12
|
1 options(warn = -1)
|
|
2 args <- commandArgs(trailingOnly = TRUE)
|
|
3 DMR = args[1]
|
|
4 annoPeakTable = args[2]
|
|
5
|
13
|
6 DMR <- ("test-data/DMR.bed")
|
12
|
7 DMRInfo = read.table(
|
|
8 DMR,
|
|
9 header = FALSE,
|
|
10 sep = "\t",
|
|
11 stringsAsFactors = FALSE,
|
|
12 quote = ""
|
|
13 )
|
|
14
|
14
|
15 peaks <- GRanges(seqnames = DMRInfo[, 1],
|
12
|
16 ranges = IRanges
|
|
17 (start = DMRInfo[, 2], end = DMRInfo[, 3]))
|
14
|
18 peaks[1:2]
|
|
19 ## try http:// if https:// URLs are not supported
|
|
20 source("https://bioconductor.org/biocLite.R")
|
|
21 biocLite("EnsDb.Hsapiens.v75")
|
|
22 library(EnsDb.Hsapiens.v75)
|
|
23 annoData <- toGRanges(EnsDb.Hsapiens.v75)
|
|
24 annoData[1:2]
|
12
|
25
|
14
|
26 ## keep the seqnames in the same style
|
|
27 seqlevelsStyle(peaks) <- seqlevelsStyle(annoData)
|
|
28 ## do annotation by nearest TSS
|
|
29 anno <- annotatePeakInBatch(peaks, AnnotationData=annoData)
|
|
30 anno[1:2]
|
|
31 annoPeakTable <- ('test-data/ChIPPeak.csv')
|
|
32 write.csv(anno, annoPeakTable, row.names = FALSE)
|
12
|
33
|
14
|
34
|