Changeset 5753

Timestamp:

Nov 22, 2019, 1:10:32 PM (3 years ago)

Author:

Nicklas Nordborg

Message:

References #1205: Variant calling fails with "malformed BED entry"

Fixes the gc_stat.py script so that it uses the correct sequence if the variant location is within the 50 first base pairs of a chromosome.

File:

• other/pipeline/trunk/gc_stat.py (modified) (3 diffs)

other/pipeline/trunk/gc_stat.py

@@ -9 +9 @@
 
 0: Chromosome (printed to VCF#CHROM)
-1: Start (=original location - 50; not used)
+1: Start (=original location - 50; or 0 if original location < 50)
 2: Stop (=original location + 50; not used)
 3: Original location (printed to VCF#POS)
@@ -35 +35 @@
 # Column indexes from the BED file we are reading
 COL_CHR = 0
+COL_START = 1
+COL_END = 2
 COL_POS = 3
 COL_REF = 4
@@ -79 +81 @@
     cols = line.rstrip().split('\t')
     seq = cols[COL_SEQ]
-    seqLocal = cols[COL_SEQ][NORM_LEN-LOCAL_LEN:NORM_LEN+LOCAL_LEN]
+    posInSeq = int(cols[COL_POS])-int(cols[COL_START])
+    seqLocal = cols[COL_SEQ][max(posInSeq-LOCAL_LEN, 0):posInSeq+LOCAL_LEN]
     gcNorm, strengthNorm = getSequenceInfo(seq)
     gcLocal, strengthLocal = getSequenceInfo(seqLocal)