Mercurial > repos > iuc > medaka_consensus

annotate convert_VCF_info_fields.py @ 7:fbf007713188 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
author iuc
date Wed, 10 Feb 2021 08:28:25 +0000
parents
children 6ca28e6144fa
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
7
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
1 #!/usr/bin/env python3
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
2
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
3 # Takes in VCF file annotated with medaka tools annotate and converts
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
4 #
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
5 # Usage statement:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
6 # python convert_VCF_info_fields.py in_vcf.vcf out_vcf.vcf
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
7
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
8 # 10/21/2020 - Nathan P. Roach, natproach@gmail.com
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
9
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
10 import sys
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
11 from collections import OrderedDict
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
12 from math import log10
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
13
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
14 from scipy.stats import fisher_exact
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
15
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
16
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
17 def pval_to_phredqual(pval):
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
18 return round(-10 * log10(pval))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
19
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
20
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
21 def parseInfoField(info):
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
22 info_fields = info.split(';')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
23 info_dict = OrderedDict()
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
24 for info_field in info_fields:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
25 code, val = info_field.split('=')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
26 info_dict[code] = val
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
27 return info_dict
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
28
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
29
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
30 def annotateVCF(in_vcf_filepath, out_vcf_filepath):
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
31 in_vcf = open(in_vcf_filepath, 'r')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
32 out_vcf = open(out_vcf_filepath, 'w')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
33 to_skip = set(['SC', 'SR'])
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
34 for i, line in enumerate(in_vcf):
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
35 if i == 1:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
36 out_vcf.write("##convert_VCF_info_fields=0.1\n")
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
37 if line[0:2] == "##":
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
38 if line[0:11] == "##INFO=<ID=":
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
39 id_ = line[11:].split(',')[0]
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
40 if id_ in to_skip:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
41 continue
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
42 out_vcf.write(line)
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
43 elif line[0] == "#":
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
44 out_vcf.write('##INFO=<ID=DPSPS,Number=2,Type=Integer,Description="Spanning Reads Allele Frequency By Strand">\n')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
45 out_vcf.write('##INFO=<ID=AF,Number=1,Type=Float,Description="Spanning Reads Allele Frequency">\n')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
46 out_vcf.write('##INFO=<ID=FAF,Number=1,Type=Float,Description="Forward Spanning Reads Allele Frequency">\n')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
47 out_vcf.write('##INFO=<ID=RAF,Number=1,Type=Float,Description="Reverse Spanning Reads Allele Frequency">\n')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
48 out_vcf.write('##INFO=<ID=SB,Number=1,Type=Integer,Description="Phred-scaled strand bias of spanning reads at this position">\n')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
49 out_vcf.write('##INFO=<ID=DP4,Number=4,Type=Integer,Description="Counts for ref-forward bases, ref-reverse, alt-forward and alt-reverse bases in spanning reads">\n')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
50 out_vcf.write('##INFO=<ID=AS,Number=4,Type=Integer,Description="Total alignment score to ref and alt allele of spanning reads by strand (ref fwd, ref rev, alt fwd, alt rev) aligned with parasail match 5, mismatch -4, open 5, extend 3">\n')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
51 out_vcf.write(line)
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
52 else:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
53 fields = line.split('\t')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
54 info_dict = parseInfoField(fields[7])
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
55 sr_list = [int(x) for x in info_dict["SR"].split(',')]
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
56 sc_list = [int(x) for x in info_dict["SC"].split(',')]
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
57 if len(sr_list) == len(sc_list):
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
58 variant_list = fields[4].split(',')
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
59 dpsp = int(info_dict["DPSP"])
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
60 ref_fwd, ref_rev = 0, 1
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
61 dpspf, dpspr = (int(x) for x in info_dict["AR"].split(','))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
62 for i in range(0, len(sr_list), 2):
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
63 dpspf += sr_list[i]
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
64 dpspr += sr_list[i + 1]
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
65 for j, i in enumerate(range(2, len(sr_list), 2)):
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
66 dp4 = (sr_list[ref_fwd], sr_list[ref_rev], sr_list[i], sr_list[i + 1])
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
67 dp2x2 = [[dp4[0], dp4[1]], [dp4[2], dp4[3]]]
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
68 _, p_val = fisher_exact(dp2x2)
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
69 sb = pval_to_phredqual(p_val)
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
70
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
71 as_ = (sc_list[ref_fwd], sc_list[ref_rev], sc_list[i], sc_list[i + 1])
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
72
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
73 info = []
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
74 for code in info_dict:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
75 if code in to_skip:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
76 continue
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
77 val = info_dict[code]
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
78 info.append("%s=%s" % (code, val))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
79
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
80 info.append("DPSPS=%d,%d" % (dpspf, dpspr))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
81
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
82 if dpsp == 0:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
83 info.append("AF=NaN")
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
84 else:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
85 af = dp4[2] + dp4[3] / dpsp
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
86 info.append("AF=%.6f" % (af))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
87 if dpspf == 0:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
88 info.append("FAF=NaN")
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
89 else:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
90 faf = dp4[2] / dpspf
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
91 info.append("FAF=%.6f" % (faf))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
92 if dpspr == 0:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
93 info.append("RAF=NaN")
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
94 else:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
95 raf = dp4[3] / dpspr
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
96 info.append("RAF=%.6f" % (raf))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
97 info.append("SB=%d" % (sb))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
98 info.append("DP4=%d,%d,%d,%d" % (dp4))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
99 info.append("AS=%d,%d,%d,%d" % (as_))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
100 new_info = ';'.join(info)
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
101 fields[4] = variant_list[j]
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
102 fields[7] = new_info
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
103 out_vcf.write("%s" % ("\t".join(fields)))
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
104 else:
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
105 print("WARNING - SR and SC are different lengths, skipping variant")
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
106 print(line.strip()) # Print the line for debugging purposes
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
107 in_vcf.close()
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
108 out_vcf.close()
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
109
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
110
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
111 if __name__ == "__main__":
fbf007713188 "planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/medaka commit 9b7d28ac59ad082874670ee989836631ba8d7fb4"
iuc
parents:
diff changeset
112 annotateVCF(sys.argv[1], sys.argv[2])