-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathfilter_cnv_sites.py
26 lines (23 loc) · 1.02 KB
/
filter_cnv_sites.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
import sys
cnv_dict={}
with open(sys.argv[2], 'r') as cnvfile:
for line in cnvfile:
if not line.startswith('chrm'):
chrm,snppos,rd = line.strip().split('\t')
cnv_dict[chrm+'_'+snppos] = rd
with open(sys.argv[1], 'w') as rm_hetSNV_f:
#rm_hetSNV_f = open(sys.argv[1], 'w')
rm_hetSNV_f.write('#chr\tref_coord\tcA_cC_cG_cT__rd\n')
for line in sys.stdin:
if not line.startswith('#'):
(chrm,ref_coord,hap1_coord,hap2_coord,ref_allele,hap1_allele,hap2_allele,
cA,cC,cG,cT,cN,ref_allele_ratio,sum_ref_n_alt_cnts,p_binom,
warning_hap1,warning_hap2) = line.split()
rd = float(cnv_dict[chrm+"_"+ref_coord])
if 0.5 <= rd <= 1.5:
sys.stdout.write(line.strip()+'\t'+str(rd)+'\n')
else:
rm_hetSNV_f.write('\t'.join([chrm, ref_coord,
'_'.join([cA,cC,cG,cT,'',str(rd)])])+'\n')
else: sys.stdout.write(line.strip()+'\tcnv\n')
#rm_hetSNV_f.close()