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 >

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]]