comparison filter_gd_snp.py @ 22:95a05c1ef5d5

update to devshed revision aaece207bd01
author Richard Burhans <burhans@bx.psu.edu>
date Mon, 11 Mar 2013 11:28:06 -0400
parents
children 8997f2ca8c7a
comparison
equal deleted inserted replaced
21:d6b961721037 22:95a05c1ef5d5
1 #!/usr/bin/env python
2
3 import sys
4 import subprocess
5 from Population import Population
6
7 ################################################################################
8
9 def convert_non_negative_int(string_value):
10 try:
11 val = int(string_value)
12 except:
13 print >> sys.stderr, '"%s" is not an integer' % string_value
14 sys.exit(1)
15
16 if val < 0:
17 print >> sys.stderr, '"%d" is negative' % val
18 sys.exit(1)
19
20 return val
21
22
23 def convert_percent(string_value):
24 if string_value.endswith('%'):
25 val = convert_non_negative_int(string_value[:-1])
26 if val > 100:
27 print >> sys.stderr, 'percentage: "%d" > 100' % val
28 sys.exit(1)
29 val = val * -1
30 else:
31 val = convert_non_negative_int(string_value)
32
33 return str(val)
34
35 ################################################################################
36
37 if len(sys.argv) < 9:
38 print >> sys.stderr, "Usage"
39 sys.exit(1)
40
41 input, p1_input, output, lo, hi, lo_ind, lo_ind_qual = sys.argv[1:8]
42 individual_metadata = sys.argv[8:]
43
44 p_total = Population()
45 p_total.from_tag_list(individual_metadata)
46
47 p1 = Population()
48 p1.from_population_file(p1_input)
49
50 if not p_total.is_superset(p1):
51 print >> sys.stderr, 'There is an individual in the population that is not in the SNP table'
52 sys.exit(1)
53
54 lo = convert_percent(lo)
55 hi = convert_percent(hi)
56
57 ################################################################################
58
59 prog = 'filter_snps'
60
61 args = []
62 args.append(prog)
63 args.append(input)
64 args.append(lo)
65 args.append(hi)
66 args.append(lo_ind)
67 args.append(lo_ind_qual)
68
69 columns = p1.column_list()
70
71 for column in sorted(columns):
72 args.append(column)
73
74 fh = open(output, 'w')
75
76 #print "args:", ' '.join(args)
77 p = subprocess.Popen(args, bufsize=-1, stdin=None, stdout=fh, stderr=sys.stderr)
78 rc = p.wait()
79 fh.close()
80
81 sys.exit(0)
82