Mercurial > repos > devteam > dna_filtering
comparison histogram.py @ 0:9c0d844f4e48 draft
Imported from capsule None
author | devteam |
---|---|
date | Mon, 28 Jul 2014 11:30:12 -0400 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:9c0d844f4e48 |
---|---|
1 #!/usr/bin/env python | |
2 #Greg Von Kuster | |
3 | |
4 import sys | |
5 from rpy import * | |
6 | |
7 assert sys.version_info[:2] >= ( 2, 4 ) | |
8 | |
9 def stop_err(msg): | |
10 sys.stderr.write(msg) | |
11 sys.exit() | |
12 | |
13 def main(): | |
14 | |
15 # Handle input params | |
16 in_fname = sys.argv[1] | |
17 out_fname = sys.argv[2] | |
18 try: | |
19 column = int( sys.argv[3] ) - 1 | |
20 except: | |
21 stop_err( "Column not specified, your query does not contain a column of numerical data." ) | |
22 title = sys.argv[4] | |
23 xlab = sys.argv[5] | |
24 breaks = int( sys.argv[6] ) | |
25 if breaks == 0: | |
26 breaks = "Sturges" | |
27 if sys.argv[7] == "true": | |
28 density = True | |
29 else: density = False | |
30 if len( sys.argv ) >= 9 and sys.argv[8] == "true": | |
31 frequency = True | |
32 else: frequency = False | |
33 | |
34 matrix = [] | |
35 skipped_lines = 0 | |
36 first_invalid_line = 0 | |
37 invalid_value = '' | |
38 i = 0 | |
39 for i, line in enumerate( file( in_fname ) ): | |
40 valid = True | |
41 line = line.rstrip('\r\n') | |
42 # Skip comments | |
43 if line and not line.startswith( '#' ): | |
44 # Extract values and convert to floats | |
45 row = [] | |
46 try: | |
47 fields = line.split( "\t" ) | |
48 val = fields[column] | |
49 if val.lower() == "na": | |
50 row.append( float( "nan" ) ) | |
51 except: | |
52 valid = False | |
53 skipped_lines += 1 | |
54 if not first_invalid_line: | |
55 first_invalid_line = i+1 | |
56 else: | |
57 try: | |
58 row.append( float( val ) ) | |
59 except ValueError: | |
60 valid = False | |
61 skipped_lines += 1 | |
62 if not first_invalid_line: | |
63 first_invalid_line = i+1 | |
64 invalid_value = fields[column] | |
65 else: | |
66 valid = False | |
67 skipped_lines += 1 | |
68 if not first_invalid_line: | |
69 first_invalid_line = i+1 | |
70 | |
71 if valid: | |
72 matrix += row | |
73 | |
74 if skipped_lines < i: | |
75 try: | |
76 a = r.array( matrix ) | |
77 r.pdf( out_fname, 8, 8 ) | |
78 histogram = r.hist( a, probability=not frequency, main=title, xlab=xlab, breaks=breaks ) | |
79 if density: | |
80 density = r.density( a ) | |
81 if frequency: | |
82 scale_factor = len( matrix ) * ( histogram['mids'][1] - histogram['mids'][0] ) #uniform bandwidth taken from first 2 midpoints | |
83 density[ 'y' ] = map( lambda x: x * scale_factor, density[ 'y' ] ) | |
84 r.lines( density ) | |
85 r.dev_off() | |
86 except Exception, exc: | |
87 stop_err( "%s" %str( exc ) ) | |
88 else: | |
89 if i == 0: | |
90 stop_err("Input dataset is empty.") | |
91 else: | |
92 stop_err( "All values in column %s are non-numeric." %sys.argv[3] ) | |
93 | |
94 print "Histogram of column %s. " %sys.argv[3] | |
95 if skipped_lines > 0: | |
96 print "Skipped %d invalid lines starting with line #%d, '%s'." % ( skipped_lines, first_invalid_line, invalid_value ) | |
97 | |
98 r.quit( save="no" ) | |
99 | |
100 if __name__ == "__main__": | |
101 main() |