Skip to content

Commit 4a8ae65

Browse files
committed
Update annotate assumptions on chromosome change.
1 parent 3709df7 commit 4a8ae65

1 file changed

Lines changed: 20 additions & 16 deletions

File tree

rhocall/run_annotate.py

Lines changed: 20 additions & 16 deletions
Original file line numberDiff line numberDiff line change
@@ -53,12 +53,16 @@ def run_annotate(proband_vcf, bed, quality_threshold, flag_upd_at_fraction, outp
5353
aztype = 'ND'
5454
azlength = end - start + 1
5555

56+
logger.debug("Looking at bed window {0} {1}={2}".format(chrom,start,end))
57+
5658
# print("looking for chrom %s %d" % (chrom, pos))
5759
passed_win = False
60+
found_var = False
61+
5862
while not passed_win:
5963
# print("testing var chrom %s %d" % (var.CHROM, var.start))
6064
if var.CHROM == chrom and var.end >= start and var.end <= end:
61-
65+
# variant in window - print, and go to next var
6266
if az == 1:
6367
var.INFO['AZ'] = True
6468
else:
@@ -67,35 +71,35 @@ def run_annotate(proband_vcf, bed, quality_threshold, flag_upd_at_fraction, outp
6771
var.INFO['AZQUAL']=str(qual)
6872
var.INFO['AZLENGTH']=str(azlength)
6973
var.INFO['AZTYPE']=str(aztype)
70-
74+
7175
output.write(str(var))
7276
var = next(proband_vcf)
77+
found_var = True
7378
elif var.CHROM == chrom and var.start < start:
74-
# before next win (and not in win) - write and pull new var
79+
# before next win (and on same chr) - write and pull new var
7580
output.write(str(var))
7681
var = next(proband_vcf)
7782
elif var.CHROM == chrom and var.end > end:
7883
# var is after last window position, same chr
7984
# pull new window, but no new variant and don't write var yet
8085
passed_win = True
8186
elif var.CHROM != chrom:
82-
# we ask that chromosomes are sorted the same way in both files
83-
if var.end > start:
84-
# likely, we have pulled a new win on a new chr.
85-
# dump var without label, and pull new
86-
output.write(str(var))
87-
var = next(proband_vcf)
88-
elif var.end <= start:
89-
# if we have pulled a new var on a new chrom
90-
# say, skipping one/first chr entirely, or last win extended
91-
# to chromosome end
92-
# pull new window, but no new variant and don't write var yet
87+
# first time around, assume there are many variants, at least one
88+
# per bed interval
89+
90+
# then we either just exited win by drawing var from new chr
91+
if found_var:
92+
# pull next win, without deciding the fate of this var yet
9393
passed_win = True
94-
logger.debug("")
94+
else:
95+
# or window is on new chr, and we need to draw new vars to
96+
# get there
97+
var = next(proband_vcf)
98+
output.write(str(var))
9599
else:
96100
# not found, but passed the due position?!
97101
# var = next(proband_vcf)
98-
logger.warning("Oops? Dunno what to do with this win and var sequence!")
102+
logger.error("Oops? Unexpected window/variant set!")
99103

100104

101105

0 commit comments

Comments
 (0)