transcripts.py

   1 #
   2 #  transcripts.py
   3 #  ENRAGE
   4 #
   5 #  Created by Ali Mortazavi on 1/25/08.
   6 #
   7 """ usage: python %s rpkmFile outFile [--transcriptome size] [--cells count] [--efficiency fraction]
   8            where transcriptome size is in Gbp, cell count is in arbitrary units and efficiency is a fraction
   9 """
  10
  11 import sys
  12 import optparse
  13 from commoncode import getConfigParser, getConfigFloatOption
  14
  15 def main(argv=None):
  16     if not argv:
  17         argv = sys.argv
  18
  19     print "transcripts: version 3.1"
  20     usage = "usage: python %prog rpkmFile outFile [options]"
  21
  22     parser = makeParser(usage)
  23     (options, args) = parser.parse_args(argv[1:])
  24
  25     if len(args) < 2:
  26         print usage
  27         sys.exit(1)
  28
  29     infile = args[0]
  30     outfile = args[1]
  31
  32     transcripts(infile, outfile, options.tSize, options.cellCount, options.efficiency)
  33
  34
  35 def makeParser(usage=""):
  36     parser = optparse.OptionParser(usage=usage)
  37     parser.add_option("--transcriptome", type="float", dest="tSize",
  38                       help="transcriptome size in Gbp [default 200000.0]")
  39     parser.add_option("--cells", type="float", dest="cellCount",
  40                       help="arbitrary units [default 1e6]")
  41     parser.add_option("--efficiency", type="float", dest="efficiency",
  42                       help="fraction [default 0.3]")
  43
  44     configParser = getConfigParser()
  45     section = "transcripts"
  46     tSize = getConfigFloatOption(configParser, section, "tSize", 200000.0)
  47     cellCount = getConfigFloatOption(configParser, section, "cellCount", 1e6)
  48     efficiency = getConfigFloatOption(configParser, section, "efficiency", 0.3)
  49
  50     parser.set_defaults(tSize=tSize, cellCount=cellCount, efficiency=efficiency)
  51
  52     return parser
  53
  54
  55 def transcripts(infilename, outfilename, tSize=200000, cellCount=1e6, efficiency=0.3):
  56     infile = open(infilename)
  57     outfile = open(outfilename, "w")
  58     for line in infile:
  59         fields = line.strip().split()
  60         rpkm = float(fields[-1])
  61         transcripts = rpkm * tSize
  62         transPerCell = transcripts / cellCount / efficiency
  63         outfile.write("%s\t%.1f\t%.1f\n" % (fields[0], transcripts, transPerCell))
  64
  65     infile.close()
  66     outfile.close()
  67
  68 if __name__ == "__main__":
  69     main(sys.argv)