|
33
|
1 require("minfi", quietly = TRUE)
|
|
|
2 require("ChIPseeker", quietly = TRUE)
|
|
|
3 require("ChIPpeakAnno", quietly = TRUE)
|
|
|
4 require("data.table", quietly = TRUE)
|
|
|
5
|
|
|
6
|
|
|
7 args <- commandArgs(trailingOnly = TRUE)
|
|
|
8 GSMTable = args[1]
|
|
|
9 platform = args[2]
|
|
|
10 Data_Table = args[3]
|
|
|
11 cutoff = as.numeric(args[4])
|
|
|
12 clusterSize = as.numeric(args[5])
|
|
|
13 DMR = args[6]
|
|
|
14
|
|
|
15 TAB = fread(GSMTable)
|
|
|
16
|
|
|
17 IlmnInfo = fread(platform)
|
|
|
18
|
|
|
19 gmSet = fread(Data_Table)
|
|
|
20
|
|
|
21 # bumphunter Run with processed data
|
|
|
22 designMatrix <- model.matrix( ~ TAB$Phenotype)
|
|
|
23
|
|
|
24 bumps <- bumphunter(
|
|
|
25 as.matrix(gmSet),
|
|
|
26 design = designMatrix,
|
|
|
27 pos = IlmnInfo$BP,
|
|
|
28 cutoff = cutoff,
|
|
|
29 chr = IlmnInfo$CHR
|
|
|
30 )
|
|
|
31
|
|
|
32 # choose DMR's of a certain length threshold
|
|
|
33 DMRTable <- bumps$table[which(bumps$table$L >= clusterSize), ]
|
|
|
34 DMRInfo <- data.table(DMRTable$chr, DMRTable$start, DMRTable$end)
|
|
|
35
|
|
|
36
|
|
|
37
|
|
|
38 write.table(
|
|
|
39 DMRInfo,
|
|
|
40 DMR,
|
|
|
41 quote = F,
|
|
|
42 sep = "\t",
|
|
|
43 row.names = F,
|
|
|
44 col.names = F
|
|
|
45 )
|
