findMotifs.py

   1 #
   2 #  findMotifs.py
   3 #  ENRAGE
   4 #
   5 try:
   6     import psyco
   7     psyco.full()
   8 except:
   9     pass
  10
  11 import sys
  12 import os
  13 import optparse
  14 from cistematic.experiments.fasta import Fasta
  15 from cistematic.programs.meme import Meme
  16 from cistematic.programs.cisGreedy import CisGreedy
  17 from commoncode import getConfigParser, getConfigOption, getConfigIntOption, getConfigFloatOption, getConfigBoolOption
  18 #TODO: cisSampler is not supported yet!
  19 #from cistematic.programs.cisSampler import CisSampler
  20
  21 print "findMotifs: version 3.5"
  22
  23 def main(argv=None):
  24     if not argv:
  25         argv = sys.argv
  26
  27     usage = "usage: python %prog explabel regions.fsa [options]"
  28
  29     parser = getParser(usage)
  30     (options, args) = parser.parse_args(argv[1:])
  31
  32     if len(args) < 2:
  33         print usage
  34         print "\n\twhere at least one of the motif finders (meme or cisGreedy) must be specified\n"
  35         sys.exit(1)
  36
  37     expbase = args[0]
  38     fsafile = args[1]
  39
  40     doCisSampler = False
  41     if "--cisSampler" in sys.argv:
  42         print "cisSampler is not supported yet! avoid using it for now"
  43         doCisSampler = True
  44
  45     findMotifs(expbase, fsafile, options.doMeme, options.doCisGreedy, options.saveLogo,
  46                options.threshold, options.numMotifs, options.maxWidth, options.maskLower,
  47                doCisSampler)
  48
  49
  50 def getParser(usage):
  51     parser = optparse.OptionParser(usage=usage)
  52     parser.add_option("--meme", action="store_true", dest="doMeme")
  53     parser.add_option("--cisGreedy", action="store_true", dest="doCisGreedy")
  54     parser.add_option("--logo", action="store_true", dest="saveLogo")
  55     parser.add_option("--threshold", type="float", dest="threshold")
  56     parser.add_option("--prefix", dest="motifPrefix")
  57     parser.add_option("--numMotifs", type="int", dest="numMotifs")
  58     parser.add_option("--maxWidth", type="int", dest="maxWidth")
  59     parser.add_option("--maskLower", action="store_true", dest="maskLower")
  60
  61     configParser = getConfigParser()
  62     section = "findMotifs"
  63     doMeme = getConfigBoolOption(configParser, section, "doMeme", False)
  64     doCisGreedy = getConfigBoolOption(configParser, section, "doCisGreedy", False)
  65     saveLogo = getConfigBoolOption(configParser, section, "saveLogo", False)
  66     threshold = getConfigFloatOption(configParser, section, "threshold", 75.)
  67     numMotifs = getConfigOption(configParser, section, "numMotifs", 10)
  68     maxWidth = getConfigIntOption(configParser, section, "maxWidth", 28)
  69     maskLower = getConfigBoolOption(configParser, section, "maskLower", False)
  70
  71
  72     parser.set_defaults(doMeme=doMeme, doCisGreedy=doCisGreedy, saveLogo=saveLogo,
  73                         threshold=threshold, numMotifs=numMotifs, maxWidth=maxWidth, maskLower=maskLower)
  74
  75     return parser
  76
  77
  78 def findMotifs(expbase, fsafile, doMeme=False, doCisGreedy=False, saveLogo=False, threshold=75.,
  79                numMotifs=10, maxWidth=28, maskLower=False, doCisSampler=False):
  80
  81     motifPrefix = expbase
  82
  83     #TODO: cisSampler is not supported yet!
  84     #if doMeme or doCisGreedy or doCisSampler:
  85     if not (doMeme or doCisGreedy):
  86         print "error: must specify at least one motif finder - exiting"
  87         sys.exit(1)
  88
  89     exp = Fasta(expbase, "%s.db" % expbase)
  90
  91     exp.initialize()
  92     if maskLower:
  93         exp.setMaskLowerCase(True)
  94
  95     if doMeme:
  96         prog4 = Meme()
  97         prog4.setMaxWidth(maxWidth)
  98         prog4.setNumMotifs(numMotifs)
  99         prog4.setModel("zoops")
 100         exp.appendProgram(prog4)
 101
 102     if doCisGreedy:
 103         prog5 = CisGreedy()
 104         prog5.setGenExpOptions([])
 105         prog5.setMaxWidth(maxWidth)
 106         prog5.setNumMotifs(numMotifs)
 107         exp.appendProgram(prog5)
 108
 109     #TODO: cisSampler is not supported yet!
 110     #if doCisSampler:
 111     #    prog6 = CisSampler()
 112     #    prog6.setGenExpOptions([])
 113     #    prog6.setMaxWidth(maxWidth)
 114     #    prog6.setNumMotifs(numMotifs)
 115     #    exp.appendProgram(prog6)
 116
 117     exp.run(fsafile)
 118     exp.createAnalysis()
 119     exp.loadAnalysis()
 120     exp.mapMotifs(threshold, verbose=False)
 121     exp.exportMotifs(prefix = motifPrefix)
 122     if saveLogo:
 123         exp.exportLogos(prefix = motifPrefix)
 124
 125     exp.draw("%s.png" % expbase, maxOccurences=4000)
 126     print "deleting database..."
 127     del exp
 128     os.remove("%s.db" % expbase)
 129
 130
 131 if __name__ == "__main__":
 132     main(sys.argv)