buildMatrix.py

   1 #
   2 #  buildMatrix.py
   3 #  ENRAGE
   4 #
   5 #  Created by Ali Mortazavi on 3/6/09.
   6 #
   7 import sys
   8 import string
   9 import optparse
  10 from commoncode import writeLog, getConfigParser, getConfigOption, getConfigBoolOption
  11
  12 versionString = "buildMatrix: version 1.5"
  13 print versionString
  14
  15
  16 def main(argv=None):
  17     if not argv:
  18         argv = sys.argv
  19
  20     usage = "usage: python %prog matrix.step.N-1 data.part matrix.step.N [--rescale] [--truncate maxRPKM] [--log altlogfile]"
  21
  22     parser = makeParser(usage)
  23     (options, args) = parser.parse_args(argv[1:])
  24
  25     if len(args) < 3:
  26         print usage
  27         sys.exit(0)
  28
  29     infile = args[0]
  30     colfilename = args[1]
  31     outfilename = args[2]
  32
  33     if options.maxRPKM is not None:
  34         truncateRPKM = True
  35         maxRPKM = options.maxRPKM
  36     else:
  37         truncateRPKM = False
  38         maxRPKM = 100000000
  39
  40     buildMatrix(infile, colfilename, outfilename, truncateRPKM, maxRPKM,
  41                 options.rescale, options.logfilename)
  42
  43
  44 def makeParser(usage=""):
  45     parser = optparse.OptionParser(usage=usage)
  46     parser.add_option("--rescale", action="store_true", dest="rescale")
  47     parser.add_option("--truncate", type="int", dest="maxRPKM")
  48     parser.add_option("--log", dest="logfilename")
  49
  50     configParser = getConfigParser()
  51     section = "buildMatrix"
  52     rescale = getConfigBoolOption(configParser, section, "rescale", False)
  53     maxRPKM = getConfigOption(configParser, section, "maxRPKM", None)
  54     logfilename = getConfigOption(configParser, section, "logfilename", "buildMatrix.log")
  55
  56     parser.set_defaults(rescale=rescale, maxRPKM=maxRPKM, logfilename=logfilename)
  57
  58     return parser
  59
  60
  61 def buildMatrix(inFileName, colfilename, outfilename, truncateRPKM,
  62                 maxRPKM=100000000, rescale=False, logfilename="buildMatrix.log"):
  63
  64     writeLog(logfilename, versionString, string.join(sys.argv[1:]))
  65
  66     if "/" in colfilename:
  67         colname = colfilename.split("/")[-1]
  68     else:
  69         colname = colfilename
  70
  71     fileParts = colname.split(".")
  72     colID =  fileParts[0]
  73
  74     infile = open(inFileName)
  75     colfile = open(colfilename)
  76     outfile = open(outfilename, "w")
  77     header = infile.readline()[:-1]
  78     if header.strip() == "":
  79         header = "#\t"
  80
  81     outfile.write("%s\t%s\n" % (header, colID))
  82
  83     values = []
  84     min = 20000000000.
  85     max = -1.
  86     untruncatedMax = -1.
  87     for line in colfile:
  88         if doNotProcessLine(line):
  89             continue
  90
  91         fields = line.strip().split()
  92         val = float(fields[-1])
  93         if truncateRPKM and val > maxRPKM:
  94             if val > untruncatedMax:
  95                 untruncatedMax = val
  96
  97             val = maxRPKM
  98
  99         values.append(val)
 100         if val < min:
 101             min = val
 102
 103         if val > max:
 104             max = val
 105
 106     range = max - min
 107     if rescale:
 108         finalValues = [(val - min)/range for val in values]
 109     else:
 110         finalValues = values
 111
 112     for val in finalValues:
 113         line = infile.readline().strip()
 114         line += "\t%1.3f\n" % val
 115         outfile.write(line)
 116
 117     outfile.close()
 118
 119     if untruncatedMax > 0:
 120         max = untruncatedMax
 121
 122     message = "max value in %s was %.2f" % (colname, max)
 123     if untruncatedMax > 0:
 124         message += " but was truncated to %d" % maxRPKM
 125
 126     print message
 127     writeLog(logfilename, versionString, message)
 128
 129
 130 def doNotProcessLine(line):
 131     return line[0] == "#"
 132
 133
 134 if __name__ == "__main__":
 135     main(sys.argv)