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