annotate tablemerger.py @ 6:e711f63e99ff draft

Uploaded
author rreumerman
date Fri, 05 Apr 2013 05:06:42 -0400
parents bd5692103d5b
children 8de0ffc2166f
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
4
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
1 '''Takes tab-delimited SNP tables from user input and merges them into one.'''
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
2
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
3 import sys
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
4 files = []
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
5 filenames = []
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
6
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
7 try:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
8 output = open(sys.argv[1], "w")
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
9 output.write('Position\tReference')
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
10 except:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
11 exit("No output file given or unable to open output file.")
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
12 for name in sys.argv[2:]:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
13 try:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
14 files.append(open(name, "rU"))
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
15 except:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
16 continue
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
17
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
18 # Fetch headers and print them to output file;
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
19 headers = [header.readline()[:-1].split('\t')[2:] for header in files]
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
20 columns = [len(strains) for strains in headers]
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
21 for strain in [a for b in headers for a in b]:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
22 output.write('\t'+strain)
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
23 output.flush()
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
24
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
25 file_active = [True]*len(files)
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
26 snps = [row.readline()[:-1].split('\t') for row in files]
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
27
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
28 while True in file_active:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
29 for h in range(0,len(snps)):
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
30 if file_active[h]:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
31 cur_pos = [h]
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
32 lowest = int(snps[h][0])
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
33 break
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
34 i = 1
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
35
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
36 # Determine lowest position
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
37 while i < len(snps):
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
38 if int(snps[i][0]) < lowest and file_active[i]:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
39 lowest = int(snps[i][0])
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
40 cur_pos = [i]
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
41 elif int(snps[i][0]) == lowest:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
42 cur_pos.append(i)
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
43 i+=1
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
44
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
45 # Check if all SNPs sharing a position have the same reference base, exit with message otherwise;
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
46 if len(cur_pos) > 1:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
47 ref_base = snps[cur_pos[0]][1].lower()
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
48 for j in cur_pos[1:]:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
49 if snps[j][1].lower() != ref_base:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
50 error = '\nError: Reference bases not the same for position %s, present in following files:' % lowest
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
51 for k in cur_pos:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
52 error += ' '+filenames[k]
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
53 exit(error+'.')
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
54
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
55 # Write line to output file containing position, ref base and snps/empty cells;
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
56 output.write('\n'+snps[cur_pos[0]][0]+'\t'+snps[cur_pos[0]][1].lower())
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
57 for l,row in enumerate(snps):
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
58 if l in cur_pos:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
59 for snp in row[2:]:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
60 output.write('\t'+snp)
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
61 else:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
62 output.write('\t'*columns[l])
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
63
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
64 # Read new line in files that had snp at current position;
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
65 for m in cur_pos:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
66 line = files[m].readline()
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
67 if line == '': file_active[m] = False
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
68 else:
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
69 snps[m] = line.split('\t')
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
70 snps[m][-1] = snps[m][-1].rstrip()# Remove newline character at end of line;
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
71
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
72 for it in files: it.close()
bd5692103d5b Uploaded
rreumerman
parents:
diff changeset
73 output.close()