Question

DMRcate: Cpgids in DMR

1

Entering edit mode

yoursbassanio ▴ 10

@yoursbassanio-12717

Last seen 3.6 years ago

The below is my Dmrcate run and I do get the DMRS. but how to get the CPGids that are part of these DMRS

corfit <- duplicateCorrelation(myMs, design=designMatrix, block = biolrep)
dmrCate.<-cpg.annotate("array", myMs, what="M", arraytype = "EPIC", analysis.type="differential", design=designMatrix, coef=2, block=biolrep, correlation=corfit$consensus.correlation,fdr=0.01)
dmrcoutput_all <- dmrcate(dmrCate,lambda=500, C=5, pcutoff="fdr",min.cpgs=1,mc.cores=1)
wgbs.ranges_all.snp7 <- extractRanges(dmrcoutput_all, genome = "hg19")
write.csv (wgbs.ranges_all.snp7,file="DMRCate_DifferentailRange.txt")
write.csv(dmrcoutput_all.snp7$results,file="DMRCate_ALL.txt")

how to get the 17 Cpgids from the result shown below

"","seqnames","start","end","width","strand","no.cpgs","minfdr","Stouffer","maxbetafc","meanbetafc","overlapping.promoters"
"1","chr6",31539539,31540750,1212,"*",17,2.06654851019045e-120,2.16758469969414e-149,-0.162022961391297,-0.106211449278892,"LTA-001, LTA-002, LTA-004, LTA-003"

I tried with gRanges but unable to understand how to perform it

methylation Epic annotation DMRcate • 1.4k views

ADD COMMENT • link updated 4.5 years ago by James W. MacDonald 65k • written 4.6 years ago by yoursbassanio ▴ 10

1

Entering edit mode

 gr <- granges(OBJECT)
locs.ranges<-names(gr[seqnames(gr) == "chr6" & start(gr)>=31539539 & end(gr) <=31540750 ])

The above two lines I was able to subset the ids belonging to 1 location of interest.

But now How can I loop this for multiple locations present in a data.frame of dmrcoutput_all?

ADD REPLY • link 4.6 years ago yoursbassanio ▴ 10

1

Entering edit mode

There would be multiple ways(more efficient way ) to what do than the below script which I have used solve my own problem.

But It would be great to know If there is any efficient way to do the same.

OBJECT <- is the last object from which beta/mvalue for DMRCate is give as input

    gr <- granges(OBJECT)

 #read DMRcate outputstringsAsFactors = FALSE ;else there issue by creating levels
 mydata<-read.csv(file="DMRCate_Output.csv",header=TRUE,stringsAsFactors = FALSE)  

    for(n in 1:nrow(mydata)) 
    {

      # chr is the second column
      chr<-mydata2[n,2]   

      # DMR start is in 3rd column
       start=mydata2[n,3]   

       #DMR end is in 4th column
       end=mydata2[n,4]   

     #strand where all "*" so not included while retrieving
      #Start and end pos included 
      infos<-names(gr[seqnames(gr) == chr & start(gr)>=start & end(gr) <=end ])  

      allcpg<-toString(infos)
      mydata$cpgids<-allcpg
    }

    write.table(mydata,file="DMRCate_Model3_DifferentailRange_V3V4_withCpgids",sep="\t");

P.S: Here the output from DMRcate is written to a file

ADD REPLY • link 4.6 years ago yoursbassanio ▴ 10

0

Entering edit mode

 gr <- granges(OBJECT)
locs.ranges<-names(gr[seqnames(gr) == "chr6" & start(gr)>=31539539 & end(gr) <=31540750 ])

The above two lines I was able to subset the ids belonging to 1 location of interest.

How can I loop this for multiple locations present in a data.frame of dmrcoutput_all?

ADD REPLY • link 4.6 years ago yoursbassanio ▴ 10

score 0 · Answer 1 · 2019-10-08

0

Entering edit mode

James W. MacDonald 65k

@james-w-macdonald-5106

Last seen 1 minute ago

United States

You don't say what your starting object is, so you are asking us to read your mind. If I were to do so, I would divine that you used minfi to import and process your data, and you have a GenomicRatioSet that you got your M-values from. And since that's obviously what you did, you can just do

firstrow <- subsetByOverlaps(OBJECT, wgbs.ranges_all.snp7[1,])

Which will give you a GenomicRatioSet that only contains the data that overlap the first DMR. Which you can do things with. Also, read the vignettes for GenomicRanges.

ADD COMMENT • link 4.5 years ago James W. MacDonald 65k

0

Entering edit mode

Also, if you aren't answering your own question (which it appears you weren't?), don't use the Add your answer box, but instead ADD REPLY or ADD COMMENT.

ADD REPLY • link 4.5 years ago James W. MacDonald 65k