feat: add --only-passing option to statSTR (#229)

Co-authored-by: Arya Massarat <[email protected]>

trtools/prancSTR/README.rst

@@ -45,12 +45,12 @@ Other general parameters:
 * :code:`--region <string>`: Restrict to the region chr:start-end. VCF file must be bgzipped and indexed to use this option.
 * :code:`--samples <string>`: Restrict to the given list of samples. Samples are comma separated.
 * :code:`--vcftype <string>`: Specify the tool which generated the vcf call file for STRs. Currently this will fail if using anything other than :code:`hipstr` VCFs.
-* :code:`--only-passing`: Filters out the VCF records with non-passing FILTER column
-* :code:`--output-all`: Force tool to output results for all loci. Overrides :code:``--only-passing``.
+* :code:`--only-passing`: Filters out the VCF records with non-passing FILTER column.
+* :code:`--output-all`: Force tool to output results for all loci. Overrides :code:`--only-passing`.
 * :code:`--readfield <string>`: Specify which VCF format field output by HipSTR to utilize for extracting read information. We recommend setting this to "MALLREADS". "ALLREADS" is also accepted but we have found that it produces unreliable results.
 * :code:`--debug`: Print helpful debug messages.
 * :code:`--quiet`: Restrict printing of any messages.
-* :code:`--version`: Print the version of the tool
+* :code:`--version`: Print the version of the tool.
 
 Notes:
 

trtools/statSTR/README.rst

@@ -30,6 +30,7 @@ Optional general parameters:
 * :code:`--sample-prefixes <string>`: The prefixes to name output for each samples group. By default uses 1, 2, 3 etc. Must be sample length as :code:`--samples`.
 * :code:`--region <string>`: Restrict to specific regions (chrom:start-end). Requires the input VCF to be bgzipped and tabix indexed.
 * :code:`--precision <int>`: How much precision to use when writing stats (default = 3)
+* :code:`--only-passing`: Filters out the VCF records with non-passing FILTER column.
 
 For specific statistics available, see below.
 

trtools/statSTR/statSTR.py

@@ -453,6 +453,8 @@ def getargs(): # pragma: no cover
     filter_group.add_argument("--region", help="Restrict to the region "
                               "chrom:start-end. Requires file to bgzipped and"
                               " tabix indexed.", type=str)
+    filter_group.add_argument("--only-passing", help="Only process records "
+                              " where FILTER==PASS", action="store_true")
     stat_group_name = "Stats group"
     stat_group = parser.add_argument_group(stat_group_name)
     stat_group.add_argument("--thresh", help="Output threshold field (max allele size, used for GangSTR strinfo).", action="store_true")
@@ -574,6 +576,10 @@ def main(args):
             nrecords += 1
 
             trrecord = trh.HarmonizeRecord(vcftype, record)
+
+            if args.only_passing and record.FILTER is not None:
+                continue
+
             if args.plot_afreq and num_plotted <= MAXPLOTS:
                 PlotAlleleFreqs(trrecord, args.out, sample_indexes=sample_indexes, sampleprefixes=sample_prefixes)
                 num_plotted += 1

trtools/statSTR/tests/test_statSTR.py

@@ -17,6 +17,7 @@ def args(tmpdir):
     args.sample_prefixes = None
     args.plot_afreq = False
     args.region = None
+    args.only_passing = False
     args.thresh = False
     args.afreq = False
     args.acount = False
@@ -49,6 +50,17 @@ def test_RightFile(args, vcfdir):
     retcode = main(args)
     assert retcode==0
 
+# Test rhe only passing option
+def test_OnlyPassing(args, vcfdir):
+    fname = os.path.join(vcfdir, "CEU_test.vcf.gz")
+    args.vcf = fname
+
+    # With only passing
+    args.only_passing = True
+    args.region = None
+    retcode = main(args)
+    assert retcode==0
+
 # Test all the statSTR options
 def test_Stats(args, vcfdir, capsys):
     fname = os.path.join(vcfdir, "few_samples_few_loci.vcf.gz")