Question: How to resize a GenomicRange objects centered on a DNA motif
0
gravatar for Ming Tang
4.9 years ago by
Ming Tang2.5k
Houston/MD Anderson Cancer Center
Ming Tang2.5k wrote:

Hi there,

I have a bed file from MACS, I did the motif-enrichment by MEME-ChIP and found the primary motif. Now, I want to find the secondary motif, but I want to resize the bed file to 500bp  centered on the primary motif. Is there an easy way to do this? I know one can use resize function for GenomicRanges object, but it can only center on the middle of the ranges.

Thank you,

Ming

chip-seq R • 2.1k views
ADD COMMENTlink modified 4.9 years ago by micahgearhart40 • written 4.9 years ago by Ming Tang2.5k
4
gravatar for seidel
4.9 years ago by
seidel6.8k
United States
seidel6.8k wrote:

If you know the location of the primary motif for the features in your bed file, simply use the starts and ends of the motifs to create a new GRanges object, and resize that with a fixed center.

ADD COMMENTlink written 4.9 years ago by seidel6.8k
1
gravatar for micahgearhart
4.9 years ago by
United States
micahgearhart40 wrote:

This function returns an aligned GRanges object based on matchPattern or will also work with matchPWM.  Its a little slow because it has to retrieve the DNA sequence for each range but I wanted to make sure there was only one site under each peak.

returnAlignedGR <- function(gr,pattern,width) {
#  pattern<-"WWWGWWGAANTTCWWCWWW"
  c <-as.matrix("chr0",nrow=length(gr))
  s <- as.matrix(rep(0,width(gr[1])))
  for (i in 1:length(gr)) {
  seq<-as.character(seqnames(gr[i])@values)
  c[i]<-seq
  start<-ranges(gr[i])@start
  width2<-ranges(gr[i])@width
  peak <- DNAString(Mmusculus[[seq]],start=start,nchar=width2)
  site <-  matchPattern(pattern,peak,fixed=F)
  if(length(site) == 1) {
                          s[i]<-start+start(site)+8 #assuming motif is 15mer, 8 is center
                          } else {
                            s[i]= 0
                          }
  }
  gr_outRD<-RangedData(space=c, IRanges(start=s-width, end=s+width ), strand="*")
  return(as(gr_outRD,"GRanges"))
  }

 

ADD COMMENTlink written 4.9 years ago by micahgearhart40
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1926 users visited in the last hour