Merge pull request #28 from MuhammedHasan/master

Bug fixes for VCFSeqExtractor: #25 , #26 and compare ref with variant file.

Author: Avsecz

kipoiseq/extractors/vcf_seq.py

@@ -30,7 +30,8 @@ def _region(self, interval):
         return '%s:%d-%d' % (interval.chrom, interval.start, interval.end)
 
     def _has_variant(self, variant, sample_id):
-        return variant.gt_types[self.sample_mapping[sample_id]] != 0
+        gt_type = variant.gt_types[self.sample_mapping[sample_id]]
+        return gt_type != 0 and gt_type != 2
 
 
 class IntervalSeqBuilder(list):
@@ -161,10 +162,10 @@ def _variant_to_sequence(self, variants):
         """
         for v in variants:
             ref = Sequence(name=v.CHROM, seq=v.REF,
-                           start=v.POS, end=v.POS + len(v.REF))
+                           start=v.start, end=v.start + len(v.REF))
             # TO DO: consider alternative alleles.
             alt = Sequence(name=v.CHROM, seq=v.ALT[0],
-                           start=v.POS, end=v.POS + len(v.ALT[0]))
+                           start=v.start, end=v.start + len(v.ALT[0]))
             yield ref, alt
 
     def _split_overlapping(self, variant_pairs, anchor):

tests/data/test.vcf.gz

@@ -18,6 +18,7 @@ def test_multi_sample_vcf_fetch_variant(multi_sample_vcf):
     interval = Interval('chr1', 3, 5)
     assert len(list(multi_sample_vcf.fetch_variants(interval))) == 2
     assert len(list(multi_sample_vcf.fetch_variants(interval, 'NA00003'))) == 1
+    assert len(list(multi_sample_vcf.fetch_variants(interval, 'NA00001'))) == 0
 
     interval = Interval('chr1', 7, 12)
     assert len(list(multi_sample_vcf.fetch_variants(interval))) == 0