Mercurial > repos > testtool > find_dmr
comparison findDMR/findDMR.R @ 12:083895dbc289 draft
Uploaded
author | testtool |
---|---|
date | Mon, 12 Jun 2017 11:05:04 -0400 |
parents | 24ac6f93cc3e |
children |
comparison
equal
deleted
inserted
replaced
11:24ac6f93cc3e | 12:083895dbc289 |
---|---|
1 require("BiocGenerics", quietly = TRUE) | 1 require("minfi", quietly = TRUE) |
2 require("data.table", quietly = TRUE) | 2 |
3 require("bumphunter", quietly = TRUE) | 3 options(warn = -1) |
4 options("download.file.method"="wget") | |
4 | 5 |
5 args <- commandArgs(trailingOnly = TRUE) | 6 args <- commandArgs(trailingOnly = TRUE) |
6 GSMTable = args[1] | 7 input1 = args[1] |
7 platform = args[2] | 8 input2 = args[2] |
8 Data_Table = args[3] | 9 cutoff = as.numeric(args[3]) |
9 cutoff = as.numeric(args[4]) | 10 B = as.numeric(args[4]) |
10 clusterSize = as.numeric(args[5]) | 11 pickCutoffQ = as.numeric(args[5]) |
11 DMR = args[6] | 12 output1 = args[6] |
13 output2 = args[7] | |
12 | 14 |
13 TAB = fread(GSMTable) | 15 GRset <- get(load(input1)) |
14 | 16 |
15 IlmnInfo = fread(platform) | 17 pheno <- fread(input2) |
16 | 18 |
17 gmSet = fread(Data_Table) | 19 designMatrix <- model.matrix(~ pheno$Phenotype) |
18 | 20 |
19 # bumphunter Run with processed data | 21 dmrs <- bumphunter(GRset, design = designMatrix, |
20 designMatrix <- model.matrix( ~ TAB$Phenotype) | 22 cutoff =cutoff, B=B, type="Beta", pickCutoff=TRUE, |
23 pickCutoffQ=pickCutoffQ) | |
21 | 24 |
22 bumps <- bumphunter( | 25 DMRTable <- dmrs$table |
23 as.matrix(gmSet), | |
24 design = designMatrix, | |
25 pos = IlmnInfo$BP, | |
26 cutoff = cutoff, | |
27 chr = IlmnInfo$CHR | |
28 ) | |
29 | 26 |
30 # choose DMR's of a certain length threshold | 27 write.table(DMRTable, output1) |
31 DMRTable <- bumps$table[which(bumps$table$L >= clusterSize), ] | |
32 DMRInfo <- data.table(DMRTable$chr, DMRTable$start, DMRTable$end) | |
33 | 28 |
29 sign=sign(DMRTable$value) | |
34 | 30 |
35 write.table( | 31 sign[sign==-1]="-" |
36 DMRInfo, | 32 sign[sign==1]="+" |
37 DMR, | 33 |
38 quote = F, | 34 dmr_track<-cbind(as.character(DMRTable$chr),DMRTable$start,DMRTable$end,sign) |
39 sep = "\t", | 35 |
40 row.names = F, | 36 write.table(dmr_track,output2,quote = FALSE, sep = "\t\t\t\t",row.names = FALSE, col.names = FALSE,append=TRUE) |
41 col.names = F | |
42 ) |