ChIPpeakAnno "examplePattern.fa" file
1
0
Entering edit mode
Paul Shannon ▴ 750
@paul-shannon-5161
Last seen 9.6 years ago
I offer below a short example of using rGADEM with MotifDb which you may find useful. The example motivates and then demonstrates "seeded" search, in which a candidate TF is provided to rGADEM. (The rGADEM vignette illustrates the simpler task of de novo motif search. - Paul On Feb 27, 2013, at 9:03 AM, Zhu, Lihua (Julie) wrote: > Jose, > > Thanks for your kind word! > > You might want to check out rGADEM package for de nova motif discovery. > > Best regards, > > Julie > library(rGADEM) library(MotifDb) library(BSgenome.Hsapiens.UCSC.hg19) path <- system.file("extdata/Test_100.fasta", package="rGADEM") sequences <- readDNAStringSet(path, "fasta") length(sequences) # 49 # you can use the UCSC genome browser to BLAT the first 4 sequences, one at a time, # against human hg19, and see what TFs bind there # http://genome.ucsc.edu/cgi-bin/hgBlat?command=start as.character(substr(sequences[1:4], 1, 80)) # Once you blat, click through to view the blat match in the genome browser, # and enable the "Transcription Factor ChIP-seq from ENCODE" track. # The FOXA1 Tf has a high score in this region of chr1 # for most of these sequences # obtain the reported motif for this TF from the Bioc MotifDb package query(MotifDb, 'foxa1') # MotifDb object of length 1 # | Created from downloaded public sources: 2012-Nov-01 # | 1 position frequency matrices from 1 source: # | JASPAR_CORE: 1 # | 1 organism/s # | Hsapiens: 1 # Hsapiens-JASPAR_CORE-FOXA1-MA0148.1 # now use rGADEM to do a seeded search: pwm.foxa1 <- as.list(MotifDb["Hsapiens-JASPAR_CORE-FOXA1-MA0148.1"]) gadem.foxa1 <-GADEM(sequences,verbose=1,genome=Hsapiens,Spwm=pwm.foxa1, seed=TRUE) consensus(gadem.foxa1) # [1] "nTGTTTACwyw" "GmwGrrrsswGGvAGn" # how do the 49 sequences match to these two motifs? mean(sapply(gadem.foxa1 at motifList[[1]]@alignList, function(x) x at pval)) mean(sapply(gadem.foxa1 at motifList[[2]]@alignList, function(x) x at pval)) On Feb 27, 2013, at 9:03 AM, Zhu, Lihua (Julie) wrote: > Jose, > > Thanks for your kind word! > > You might want to check out rGADEM package for de nova motif discovery. > > Best regards, > > Julie >
rGADEM MotifDb rGADEM MotifDb • 1.3k views
ADD COMMENT
0
Entering edit mode
@jose-luis-lavin-5797
Last seen 9.6 years ago
Dear Julie & Paul, Please excuse my delay answering you back. Thank you very much for your help with this issue, I'll read the rGADEM vignette to try it in my analysis. Best wishes JL 2013/3/1 Paul Shannon <pshannon@fhcrc.org> > I offer below a short example of using rGADEM with MotifDb which you may > find useful. > > The example motivates and then demonstrates "seeded" search, in which a > candidate TF is provided to rGADEM. (The rGADEM vignette illustrates the > simpler task of de novo motif search. > > - Paul > > > On Feb 27, 2013, at 9:03 AM, Zhu, Lihua (Julie) wrote: > > > Jose, > > > > Thanks for your kind word! > > > > You might want to check out rGADEM package for de nova motif discovery. > > > > Best regards, > > > > Julie > > > > library(rGADEM) > library(MotifDb) > library(BSgenome.Hsapiens.UCSC.hg19) > path <- system.file("extdata/Test_100.fasta", package="rGADEM") > sequences <- readDNAStringSet(path, "fasta") > length(sequences) # 49 > > # you can use the UCSC genome browser to BLAT the first 4 sequences, > one at a time, > # against human hg19, and see what TFs bind there > # http://genome.ucsc.edu/cgi-bin/hgBlat?command=start > > as.character(substr(sequences[1:4], 1, 80)) > > # Once you blat, click through to view the blat match in the genome > browser, > # and enable the "Transcription Factor ChIP-seq from ENCODE" track. > # The FOXA1 Tf has a high score in this region of chr1 > # for most of these sequences > > # obtain the reported motif for this TF from the Bioc MotifDb package > > query(MotifDb, 'foxa1') > > # MotifDb object of length 1 > # | Created from downloaded public sources: 2012-Nov-01 > # | 1 position frequency matrices from 1 source: > # | JASPAR_CORE: 1 > # | 1 organism/s > # | Hsapiens: 1 > # Hsapiens-JASPAR_CORE-FOXA1-MA0148.1 > > # now use rGADEM to do a seeded search: > > pwm.foxa1 <- as.list(MotifDb["Hsapiens-JASPAR_CORE-FOXA1-MA0148.1"]) > gadem.foxa1 <-GADEM(sequences,verbose=1,genome=Hsapiens,Spwm=pwm.foxa1, > seed=TRUE) > consensus(gadem.foxa1) > > # [1] "nTGTTTACwyw" "GmwGrrrsswGGvAGn" > > # how do the 49 sequences match to these two motifs? > mean(sapply(gadem.foxa1@motifList[[1]]@alignList, function(x) x@pval)) > mean(sapply(gadem.foxa1@motifList[[2]]@alignList, function(x) x@pval)) > > > On Feb 27, 2013, at 9:03 AM, Zhu, Lihua (Julie) wrote: > > > Jose, > > > > Thanks for your kind word! > > > > You might want to check out rGADEM package for de nova motif discovery. > > > > Best regards, > > > > Julie > > > > -- -- José Luis Lavín Trueba, Ph.D Genome Analysis Platform CIC bioGUNE Parque Tecnológico de Bizkaia Building 502, Floor 0 48160 Derio-Spain Tel.: +34 946 572 524 Fax: +34 946 568 732 [[alternative HTML version deleted]]
ADD COMMENT

Login before adding your answer.

Traffic: 780 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6