Filter out overlapping variants

In the initial implementation, variants overlapping a variant of interest were included in the output.
This update fixes this bug, using post-processing to filter tabix output by rsID.

rsidx/search.py

@@ -20,6 +20,16 @@ def trim_rsid(rsid):
     return rsidstr
 
 
+def filter_by_rsid(instream, rsidlist, header=False):
+    for line in instream:
+        if line.startswith('#'):
+            yield line
+            continue
+        chrom, coord, rsid, *values = line.split('\t')
+        if rsid[2:] in rsidlist:
+            yield line
+
+
 def search(rsidlist, dbconn, vcffile, header=False):
     c = dbconn.cursor()
     rsids = ', '.join(map(trim_rsid, rsidlist))
@@ -39,7 +49,7 @@ def fmt(row):
         tabixcmd.append('-h')
     tabixcmd.extend(coords)
     proc = Popen(tabixcmd, stdout=PIPE, universal_newlines=True)
-    for line in proc.stdout:
+    for line in filter_by_rsid(proc.stdout, rsids, header=header):
         yield line
 
 

rsidx/tests/test_search.py

@@ -99,3 +99,18 @@ def test_search_stdout(capsys):
     terminal = capsys.readouterr()
     outlines = terminal.out.strip().split('\n')
     assert len(outlines) == 5
+
+
+@pytest.mark.parametrize('doheader,numlines', [
+    (False, 1),
+    (True, 57 + 1),
+])
+def test_search_overlapping_variants(doheader, numlines):
+    rsidlist = ['rs8051733']
+    vcffile = data_file('overlap.vcf.gz')
+    idxfile = data_file('overlap.sqlite3')
+    conn = sqlite3.connect(idxfile)
+    outlines = list(rsidx.search.search(rsidlist, conn, vcffile, doheader))
+    assert len(outlines) == numlines
+    assert '\trs8051733\t' in outlines[-1]
+    assert '\rs967556605\t' not in outlines[-1]