### R code from vignette source 'MotifDb.Rnw' ################################################### ### code chunk number 1: libraries ################################################### library (MotifDb) library (MotIV) library (seqLogo) ################################################### ### code chunk number 2: MotifDb.Rnw:71-90 ################################################### MotIV.toTable = function (match) { if (length (match@bestMatch) == 0) return (NA) alignments = match@bestMatch[[1]]@aligns df = data.frame (stringsAsFactors=FALSE) for (alignment in alignments) { x = alignment name = x@TF@name eVal = x@evalue sequence = x@sequence match = x@match strand = x@strand df = rbind (df, data.frame (name=name, eVal=eVal, sequence=sequence, match=match, strand=strand, stringsAsFactors=FALSE)) } # for alignment return (df) } # MotIV.toTable ################################################### ### code chunk number 3: sources ################################################### length (MotifDb) sort (table (values (MotifDb)$dataSource), decreasing=TRUE) ################################################### ### code chunk number 4: organisms ################################################### sort (table (values (MotifDb)$organism), decreasing=TRUE) ################################################### ### code chunk number 5: metadata ################################################### colnames (values (MotifDb)) ################################################### ### code chunk number 6: queryHuman ################################################### query (MotifDb, 'hsapiens') ################################################### ### code chunk number 7: querySox ################################################### query (MotifDb, 'sox') ################################################### ### code chunk number 8: queryYeastHomeo ################################################### query (query (MotifDb, 'cerevisiae'), 'homeo') ################################################### ### code chunk number 9: homeoVariety ################################################### unique (grep ('homeo', values(MotifDb)$bindingDomain, ignore.case=T, v=T)) unique (grep ('homeo', values(MotifDb)$tfFamily, ignore.case=T, v=T)) ################################################### ### code chunk number 10: grepHuman ################################################### mdb.human <- MotifDb [grep ('Hsapiens', values (MotifDb)$organism)] mdb.sox <- MotifDb [grep ('sox', values (MotifDb)$geneSymbol, ignore.case=TRUE)] yeast.indices = grepl ('scere', values (MotifDb)$organism, ignore.case=TRUE) homeo.indices.domain = grepl ('homeo', values (MotifDb)$bindingDomain, ignore.case=TRUE) homeo.indices.family = grepl ('homeo', values (MotifDb)$tfFamily, ignore.case=TRUE) yeast.homeo.indices = yeast.indices & (homeo.indices.domain | homeo.indices.family) yeast.homeoDb = MotifDb [yeast.homeo.indices] ################################################### ### code chunk number 11: withHomeo ################################################### yeast.homeo.indices <- with(values(MotifDb), grepl('scere', organism, ignore.case=TRUE) & (grepl('homeo', bindingDomain, ignore.case=TRUE) | grepl('homeo', tfFamily, ignore.case=TRUE))) ################################################### ### code chunk number 12: subsetHuman ################################################### if (interactive ()) subset (MotifDb, organism=='Hsapiens') ################################################### ### code chunk number 13: subsetSox ################################################### if (interactive ()) subset (MotifDb, tolower (geneSymbol) == 'sox4') ################################################### ### code chunk number 14: subsetYeastHomeo ################################################### if (interactive ()) subset (MotifDb, organism=='Scerevisiae' & bindingDomain=='Homeo') ################################################### ### code chunk number 15: findEgr1 ################################################### # subset is convenient: if (interactive ()) as.list (subset (MotifDb, tolower (geneSymbol) == 'egr1')) # grep returns indices which allow for more flexibility indices = grep ('egr1', values (MotifDb)$geneSymbol, ignore.case=TRUE) length (indices) ################################################### ### code chunk number 16: MotifDbViews ################################################### MotifDb [indices] ################################################### ### code chunk number 17: as.list ################################################### as.list (MotifDb [indices]) ################################################### ### code chunk number 18: as.metadata ################################################### noquote (t (as.data.frame (values (MotifDb [indices])))) ################################################### ### code chunk number 19: egr1-multi-grep ################################################### geneSymbol.rows = grep ('Egr1', values (MotifDb)$geneSymbol, ignore.case=TRUE) organism.rows = grep ('Mmusculus', values (MotifDb)$organism, ignore.case=TRUE) source.rows = grep ('JASPAR', values (MotifDb)$dataSource, ignore.case=TRUE) egr1.mouse.jaspar.rows = intersect (geneSymbol.rows, intersect (organism.rows, source.rows)) print (egr1.mouse.jaspar.rows) egr1.motif <- MotifDb [egr1.mouse.jaspar.rows] ################################################### ### code chunk number 20: MotifDbViews ################################################### ################################################### ### code chunk number 21: subsetSearchForEgr1 ################################################### if (interactive ()) { egr1.motif <- subset (MotifDb, organism=='Mmusculus' & dataSource=='JASPAR_CORE' & geneSymbol=='Egr1') } ################################################### ### code chunk number 22: examine-egr1 ################################################### egr1.motif as.list (egr1.motif) noquote (t (as.data.frame (values (egr1.motif)))) ################################################### ### code chunk number 23: egr1 ################################################### seqLogo (as.list (egr1.motif)[[1]]) ################################################### ### code chunk number 24: motifmatch ################################################### egr1.hits <- motifMatch (as.list (egr1.motif) [1], as.list (MotifDb), top=11) # 'MotIV.toTable' -- defined above (and hidden) -- will become part of MotIV in the upcoming release tbl.hits <- MotIV.toTable (egr1.hits) print (tbl.hits) ################################################### ### code chunk number 25: three.mice.metadata ################################################### if (interactive ()) noquote (t (as.data.frame (subset (values (MotifDb), geneId=='13653')))) ################################################### ### code chunk number 26: fly.Sr.metadata ################################################### noquote (t (as.data.frame (values (MotifDb)[grep ('FBgn0003499', values (MotifDb)$geneId),]))) ################################################### ### code chunk number 27: logo1 ################################################### seqLogo (MotifDb [[tbl.hits$name[1]]]) ################################################### ### code chunk number 28: logo2 ################################################### seqLogo (MotifDb [[tbl.hits$name[2]]]) ################################################### ### code chunk number 29: logo3 ################################################### seqLogo (MotifDb [[tbl.hits$name[3]]]) ################################################### ### code chunk number 30: logo4 ################################################### seqLogo (MotifDb [[tbl.hits$name[4]]]) ################################################### ### code chunk number 31: logo5 ################################################### seqLogo (MotifDb [[tbl.hits$name[5]]]) ################################################### ### code chunk number 32: logo6 ################################################### seqLogo (MotifDb [[tbl.hits$name[6]]]) ################################################### ### code chunk number 33: export ################################################### matrix.output.file = tempfile () # substitute your preferred filename here meme.text = export (MotifDb, matrix.output.file, 'meme') metadata.output.file = tempfile () # substitute your preferred filename here write.table (as.data.frame (values (MotifDb)), file=metadata.output.file, sep='\t', row.names=TRUE, col.names=TRUE, quote=FALSE)