Commits

dawe  committed 953ceea

added scale option

  • Participants
  • Parent commits 97ef02f

Comments (0)

Files changed (4)

File bin/hc_bam.py

   option_parser.add_argument("-l","--level", help="Map level (<10)", type=int, default=8)
   option_parser.add_argument("-x","--exclude", help="Exclude this chromosome from the analysis", action='append')
   option_parser.add_argument("-i","--include", help="Limit the analysis to this chromosome", action='append')
-  option_parser.add_argument("-L","--linearscale", help="Plot map in linear scale", action="store_true", default=False)
+  option_parser.add_argument("-L","--logscale", help="Plot map in log scale", action="store_true", default=False)
   option_parser.add_argument("-r","--resolution", help="Resolution of resulting map", type=int, default=600)
   option_parser.add_argument("-s","--summarize", help="Function to summarize data", choices=['sum', 'mean', 'max'], default='sum')
-  option_parser.add_argument("-g","--genome", help="Plot a map for the whole genome",action="store_true", default=False)
+#  option_parser.add_argument("-g","--genome", help="Plot a map for the whole genome",action="store_true", default=False)
+  option_parser.add_argument("-S","--scale", help="Scale data by this factor",type=float, default=1.0)
   options = option_parser.parse_args()
   
   
     sys.stderr.write("Processing %s [%d bp] " % (chromosome, chromosome_size))
     
     sys.stderr.write("Creating map ")
-    hc = gilbert.create_map_from_bam(fh, chromosome, chromosome_size, options.level, options.summarize)
+    hc = gilbert.create_map_from_bam(fh, chromosome, chromosome_size, options.level, options.summarize, options.scale)
 
-    if not options.linearscale:
+    if options.logscale:
       hc = np.log1p(hc)
     sys.stderr.write("\n")
     sys.stderr.write("Dumping Map into numpy array %s.npy \n" % file_prefix)

File bin/hc_bed.py

   option_parser.add_argument("-l","--level", help="Map level (<10)", type=int, default=8)
   option_parser.add_argument("-x","--exclude", help="Exclude this chromosome from the analysis", action='append')
   option_parser.add_argument("-i","--include", help="Limit the analysis to this chromosome", action='append')
-  option_parser.add_argument("-L","--linearscale", help="Plot map in linear scale", action="store_true", default=False)
+  option_parser.add_argument("-L","--logscale", help="Plot map in log scale", action="store_true", default=False)
   option_parser.add_argument("-r","--resolution", help="Resolution of resulting map", type=int, default=600)
   option_parser.add_argument("-s","--summarize", help="Function to summarize data", choices=['sum', 'mean', 'max'], default='sum')
   option_parser.add_argument("-g","--genome", help="Plot a map for the whole genome",action="store_true", default=False)
     
     hc = gilbert.create_map_from_interval_tree(tree_chrom[chromosome], chromosome_size, options.level, options.summarize)
 
-    if not options.linearscale:
+    if options.logscale:
       hc = np.log1p(hc)
     sys.stderr.write("\n")
     sys.stderr.write("Dumping Map into numpy array %s.npy \n" % file_prefix)

File bin/hc_bigwig.py

   option_parser.add_argument("-l","--level", help="Map level (<10)", type=int, default=8)
   option_parser.add_argument("-x","--exclude", help="Exclude this chromosome from the analysis", action='append')
   option_parser.add_argument("-i","--include", help="Limit the analysis to this chromosome", action='append')
-  option_parser.add_argument("-L","--linearscale", help="Plot map in linear scale", action="store_true", default=False)
+  option_parser.add_argument("-L","--logscale", help="Plot map in log scale", action="store_true", default=False)
   option_parser.add_argument("-r","--resolution", help="Resolution of resulting map", type=int, default=600)
   option_parser.add_argument("-s","--summarize", help="Function to summarize data", choices=['sum', 'mean', 'max'], default='sum')
   option_parser.add_argument("-g","--genome", help="Plot a map for the whole genome",action="store_true", default=False)
     sys.stderr.write("Creating map ")
     hc = gilbert.create_map_from_bigwig(bw_handler, chromosome, chromosome_size, options.level, options.summarize)
 
-    if not options.linearscale:
+    if options.logscale:
       hc = np.log1p(hc)
     sys.stderr.write("\n")
     sys.stderr.write("Dumping Map into numpy array %s.npy \n" % file_prefix)

File gilbert/gilbert.pyx

   plt.savefig("%s_%s.png" % (name, chromosome), dpi=resolution)
   plt.close()
   
-def create_map_from_array(data_array, int chromosome_size, int level, summarize='sum'):
+def create_map_from_array(data_array, int chromosome_size, int level, summarize='sum', double scale):
   cdef int n = 1 << level
   hc_map = np.zeros((n, n))
   counts = np.zeros((n, n))
-  cdef int position, g_value
+  cdef int position
+  cder double g_value
   
   for position in range(chromosome_size):
-    g_value = data_array[position]
+    g_value = data_array[position] * scale
     if position % 1000000 == 0:
       sys.stderr.write('.')
     (y, x) = d2xy(n, (n**2 - 1) * position / chromosome_size)  
 
   return hc_map  
 
-def create_map_from_bigwig(bw_handler, chromosome, int chromosome_size, int level, summarize = 'sum'):
+def create_map_from_bigwig(bw_handler, chromosome, int chromosome_size, int level, summarize = 'sum', float scale):
   data = bw_handler.get_as_array(chromosome, 0, chromosome_size)
   data[np.isnan(data)] = 0
-  return create_map_from_array(data, chromosome_size, level, summarize)
+  return create_map_from_array(data, chromosome_size, level, summarize, scale)
 
 
 
 
   return hc_map  
 
-def create_map_from_bam(fh, chromosome, int chromosome_size, int level, summarize='sum'):    
+def create_map_from_bam(fh, chromosome, int chromosome_size, int level, summarize='sum', double scale):    
   cdef int n = 1 << level
   hc_map = np.zeros((n, n))
   counts = np.zeros((n, n))
-  cdef int position, g_value
+  cdef int position
+  cdef double g_value
   for i in fh.pileup(chromosome, 0, chromosome_size):
     position = i.pos
-    g_value = i.n
+    g_value = i.n * scale
     if position % 1000000 == 0:
       sys.stderr.write('.')
     (y, x) = d2xy(n, (n**2 - 1) * position / chromosome_size)