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