Variant filter

CHANGELOG.md

@@ -20,6 +20,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ### `Changed`
 
+- [#200](https://github.com/IntGenomicsLab/lrsomatic/pull/200) - In tumour-only mode with both `deepvariant` and `deepsomatic` selected, DeepVariant's germline calls are now kept only where DeepSomatic's verdict is `GERMLINE` or `PON` (`INFO/DS_VERDICT`); with `deepvariant` alone they are unfiltered and may include somatic variants (@robert-a-forsyth).
+- [#200](https://github.com/IntGenomicsLab/lrsomatic/pull/200) - The phased somatic VCF is now split from the merged germline+somatic VCF by an `INFO/SOMATIC` provenance flag (new `VCFTAG` module) instead of by position, so germline records at somatic coordinates no longer leak into it; ClairS-TO germline records keep their original `FILTER` in `INFO/ORIG_FILTER` (@robert-a-forsyth).
+- [#200](https://github.com/IntGenomicsLab/lrsomatic/pull/200) - Added `--smallvar_filter_pass` (default `true`), which restricts each small variant caller's output to `PASS` records before the caller consensus, phasing, VEP and the report; the published per-caller VCFs are unchanged (@robert-a-forsyth).
 - [#201](https://github.com/IntGenomicsLab/lrsomatic/pull/201) - `CLAIRSTO` and `CLAIRSTO_VERDICT_TAG` now pull the fork image from Docker Hub: `oras://docker.io/ljwharbers/clairs-to-sif:0.5.1-verdict-chm13-c0687e8-flat` under Singularity/Apptainer and `docker.io/ljwharbers/clairs-to:0.5.1-verdict-chm13-c0687e8-flat` otherwise, instead of `ghcr.io/ljwharbers/clairs-to`. The `-cpu` SIF on ghcr failed with `PROTOCOL_ERROR` on slow links: ghcr redirects every blob download to an Azure URL that expires at the next 5-minute mark and resets a stream still open then, and Apptainer resumes neither an `oras://` nor a `docker://` download. Docker Hub's download URLs are valid for 50 minutes and only checked when the request starts. `-flat` is the same software copied into an empty image in a few layers (3.3 GB instead of 7 GB); the software and its outputs are unchanged. `docs/usage.md` describes `pullTimeout`, Docker Hub's anonymous pull limit, pre-pulling, and how to recover the remaining `oras://ghcr.io` SIFs resumably (@ljwharbers).
 - [#199](https://github.com/IntGenomicsLab/lrsomatic/pull/199) - `CLAIRSTO` and `CLAIRSTO_VERDICT_TAG` now run the `-cpu` rebuild of the fork image (`0.5.1-verdict-chm13-c0687e8-cpu`), which swaps PyTorch's CUDA build for the CPU build of the same version. The software is otherwise unchanged, but the Apptainer SIF drops from 6.53 GB to 3.46 GB. The old image could not be pulled on a normal VSC link: Apptainer fetches an `oras://` SIF as a single unresumable stream, and the signed blob URL ghcr redirects to expires on a 15-minute wall-clock boundary, so 6.53 GB needed 7.3 MB/s sustained and was otherwise cut mid-transfer with `PROTOCOL_ERROR` (@ljwharbers).
 - [#197](https://github.com/IntGenomicsLab/lrsomatic/pull/197) - `CLAIRSTO` now runs `ghcr.io/ljwharbers/clairs-to:0.5.1-verdict-chm13-c0687e8` (ClairS-TO 0.5.1) instead of `docker.io/hkubal/clairs-to:v0.4.2`: a fork that lets Verdict read its CNA resources from `--cna_resource_dir`, fixes four places where Verdict's Python port of ASCAT departed from R, and disables Verdict with a warning when its resources cannot be read. **GRCh38 results move as well as CHM13 ones.** Revert to the upstream image once HKU-BAL/ClairS-TO carries these changes. The module also selects the SIF under `-profile apptainer` and sets explicit output prefixes (@ljwharbers).
@@ -41,6 +44,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ### `Fixed`
 
+- [#200](https://github.com/IntGenomicsLab/lrsomatic/pull/200) - `*_var_combine = 'all'` now keeps both callers' private calls; the prioritized caller's own private calls were previously dropped. Invalid `combine_method`/`prioritize_caller` values now raise an error (@robert-a-forsyth).
+- [#200](https://github.com/IntGenomicsLab/lrsomatic/pull/200) - The caller consensus now splits multi-allelic records before intersecting, so they match across callers, and rejoins them before phasing (@robert-a-forsyth).
+- [#200](https://github.com/IntGenomicsLab/lrsomatic/pull/200) - `LRSOMATICREPORT` now renders each sample as soon as its own inputs are ready instead of waiting for the whole batch, including samples without ASCAT's optional raw segments output (@robert-a-forsyth).
+- [#200](https://github.com/IntGenomicsLab/lrsomatic/pull/200) - Bumped `WAKHAN` from 0.4.3 to 0.4.4, fixing a deterministic `StatisticsError` crash in BAF binning that aborted the run (@robert-a-forsyth).
 - [#196](https://github.com/IntGenomicsLab/lrsomatic/pull/196) - `LRSOMATICREPORT` now points `XDG_CACHE_HOME` at the task directory alongside `HOME` and `TMPDIR`. Singularity/Apptainer inherit the host environment, so on sites that set it outside the bind-mounted work tree the render died with `Read-only file system (os error 30): mkdir '<...>/.cache/quarto'` (@AmberVerhasselt, @ljwharbers).
 - [#193](https://github.com/IntGenomicsLab/lrsomatic/pull/193) - `--vep_eve https://evemodel.org/api/proteins/bulk/download/` was rejected at launch because the "needs preparing" check keyed on a `.zip` suffix; it now checks whether the value is already a prepared bgzipped file (@AmberVerhasselt).
 - [#188](https://github.com/IntGenomicsLab/lrsomatic/pull/188) - `MODKIT_PILEUP` now runs a patched modkit 0.6.4 ([ljwharbers/modkit@pacbio-conflict-fix](https://github.com/ljwharbers/modkit/tree/pacbio-conflict-fix)): `ghcr.io/ljwharbers/modkit:0.6.4-pacbiofix-6e0afa2` under Docker and `oras://ghcr.io/ljwharbers/modkit-sif:0.6.4-pacbiofix-6e0afa2` under Singularity/Apptainer. Stock modkit 0.4.3-0.6.4 dropped 32-65 % of reads from recent PacBio HiFi BAMs and returned empty `--cpg` pileups ([nanoporetech/modkit#612](https://github.com/nanoporetech/modkit/issues/612), fix proposed in [nanoporetech/modkit#720](https://github.com/nanoporetech/modkit/pull/720)), and ignored `--phased`/`--modified-bases` for PacBio BAMs with 6mA calls. The image is `linux/amd64` only and Conda is not supported (use `--skip_modkit` there); return to the biocontainer once a release includes the fix (@ljwharbers).

conf/modules.config

@@ -121,13 +121,26 @@ process {
 
     withName: '.*:BCFTOOLS_NORM' {
         ext.prefix = { "${meta.id}.${meta.caller}_norm" }
+        // Split multi-allelics so BCFTOOLS_ISEC, which matches exact REF/ALT, can intersect them per ALT.
+        // BCFTOOLS_NORM_REJOIN rejoins them before phasing.
         ext.args = {
-            "-Oz"
+            "-m -any -Oz"
         }
         publishDir = [
             enabled: false
         ]
     }
+    // Rejoin split multi-allelics; LongPhase keys variants by position and cannot take two records at one POS.
+    withName: '.*:GERMLINE_CONSENSUS:BCFTOOLS_NORM_REJOIN' {
+        ext.prefix = { "${meta.id}_germline_sorted" }
+        ext.args   = { '-m +any --output-type z --write-index=tbi' }
+        publishDir = [ enabled: false ]
+    }
+    withName: '.*:SOMATIC_CONSENSUS:BCFTOOLS_NORM_REJOIN' {
+        ext.prefix = { "${meta.id}_somatic_sorted" }
+        ext.args   = { '-m +any --output-type z --write-index=tbi' }
+        publishDir = [ enabled: false ]
+    }
 
     withName: '.*:BCFTOOLS_ISEC' {
         ext.prefix = { "${meta.id}_isec" }
@@ -138,29 +151,60 @@ process {
             enabled: false
         ]
     }
-    withName: '.*STANDARDIZE_AF' {
-        ext.prefix = { "${meta.id}.${meta.caller}_standardized" }
-        ext.args   = {
-            meta.rename_to == 'VAF'
-                ? "--rename-annots <(printf 'FORMAT/AF\\tFORMAT/VAF\\n') -Oz -W=tbi"
-                : "--rename-annots <(printf 'FORMAT/VAF\\tFORMAT/AF\\n') -Oz -W=tbi"
+    withName: '.*:BCFTOOLS_ANNOTATE' {
+        ext.prefix = { "${meta.id}.${meta.caller}" }
+        // Stamp INFO/CALLER; in 'all' mode (meta.rename_to set) also unify the FORMAT allele frequency key.
+        ext.args = {
+            def rename = meta.rename_to == 'VAF'
+                ? "--rename-annots <(printf 'FORMAT/AF\\tFORMAT/VAF\\n') "
+                : meta.rename_to == 'AF'
+                    ? "--rename-annots <(printf 'FORMAT/VAF\\tFORMAT/AF\\n') "
+                    : ""
+            rename + '''-h <(echo '##INFO=<ID=CALLER,Number=1,Type=String,Description="Variant caller">') \
+                -c CHROM,POS,REF,ALT,INFO/CALLER \
+                -Oz \
+                -W=tbi'''
         }
         publishDir = [
             enabled: false
         ]
     }
-    withName: '.*:BCFTOOLS_ANNOTATE' {
-        ext.prefix = { "${meta.id}.${meta.caller}" }
+    // GERMLINE VERDICT TRANSFER (tumor-only): DeepSomatic's verdict filters DeepVariant's germline calls.
+    withName: '.*:DS_VERDICT_QUERY' {
+        ext.prefix = { "${meta.id}.ds_verdict" }
+        // Exclude PASS/RefCall rather than include GERMLINE/PON: bcftools fails on an undeclared FILTER.
         ext.args = {
-            '''-h <(echo '##INFO=<ID=CALLER,Number=1,Type=String,Description="Variant caller">') \
-                -c CHROM,POS,REF,ALT,INFO/CALLER \
+            "-e 'FILTER=\"PASS\" || FILTER=\"RefCall\"' -f '%CHROM\t%POS\t%REF\t%ALT\t%FILTER\n'"
+        }
+        publishDir = [
+            enabled: false
+        ]
+    }
+
+    withName: '.*:DS_VERDICT_ANNOTATE' {
+        ext.prefix = { "${meta.id}.deepvariant_verdict" }
+        ext.args = {
+            '''-h <(echo '##INFO=<ID=DS_VERDICT,Number=1,Type=String,Description="DeepSomatic FILTER at this site, set only where it is not PASS or RefCall (GERMLINE or PON in practice).">') \
+                -c CHROM,POS,REF,ALT,INFO/DS_VERDICT \
                 -Oz \
                 -W=tbi'''
         }
         publishDir = [
             enabled: false
         ]
     }
+
+    withName: '.*:DS_GERMLINE_SELECT' {
+        ext.prefix = { "${meta.id}.deepvariant_germline" }
+        // Keep only DeepSomatic-adjudicated germline sites; all others lack DS_VERDICT or fail this test.
+        ext.args = {
+            "-i 'INFO/DS_VERDICT=\"GERMLINE\" || INFO/DS_VERDICT=\"PON\"' --output-type z --write-index=tbi"
+        }
+        publishDir = [
+            enabled: false
+        ]
+    }
+
     withName: '.*:BCFTOOLS_QUERY' {
         ext.args = {
             "-f '%CHROM\t%POS\t%REF\t%ALT\t${meta.caller}\n'"
@@ -401,8 +445,22 @@ process {
             enabled: false
         ]
     }
+    // Provenance flags stamped before the phasing merge, so the somatic arm is selected by origin.
+    withName: '.*:TAG_SOMATIC' {
+        ext.prefix = { "${meta.id}_somatic_tagged" }
+        publishDir = [
+            enabled: false
+        ]
+    }
+    withName: '.*:TAG_GERMLINE' {
+        ext.prefix = { "${meta.id}_germline_tagged" }
+        publishDir = [
+            enabled: false
+        ]
+    }
     withName: '.*:PHASING_HAPLOTYPING:BCFTOOLS_VIEW' {
         ext.prefix = { "somatic_smallvariants" }
+        ext.args   = { "-i 'INFO/SOMATIC=1'" }
         publishDir = [
             path: { "${params.outdir}/${meta.id}/variants/phased" },
             mode: params.publish_dir_mode,
@@ -496,26 +554,26 @@ process {
         ]
     }
     withName: '.*:GERMLINE_CONSENSUS:BCFTOOLS_SORT' {
-        ext.prefix = { "${meta.id}_germline_sorted" }
+        ext.prefix = { "${meta.id}_germline_split_sorted" }
         ext.args = {'-Oz -W=tbi'}
         publishDir = [
             enabled: false
         ]
     }
     withName: '.*:SOMATIC_CONSENSUS:BCFTOOLS_SORT' {
-        ext.prefix = { "${meta.id}_somatic_sorted" }
+        ext.prefix = { "${meta.id}_somatic_split_sorted" }
         ext.args = {'-Oz -W=tbi'}
         publishDir = [
             enabled: false
         ]
     }
     withName: '.*:GERMLINE_CONSENSUS:BCFTOOLS_SORT_CONSENSUS' {
-        ext.prefix = { "${meta.id}_germline_sorted" }
+        ext.prefix = { "${meta.id}_germline_split_sorted" }
         ext.args   = { '-Oz -W=tbi' }
         publishDir = [ enabled: false ]
     }
     withName: '.*:SOMATIC_CONSENSUS:BCFTOOLS_SORT_CONSENSUS' {
-        ext.prefix = { "${meta.id}_somatic_sorted" }
+        ext.prefix = { "${meta.id}_somatic_split_sorted" }
         ext.args   = { '-Oz -W=tbi' }
         publishDir = [ enabled: false ]
     }
@@ -575,7 +633,7 @@ process {
     }
 
     withName: '.*:CLAIR3' {
-        ext.args = { "--sample_name=${meta.id}" }
+        // --sample_name is passed by the module itself, from ext.prefix (which defaults to meta.id).
         publishDir = [
             path: { "${params.outdir}/${meta.id}/variants/clair3" },
             mode: params.publish_dir_mode,
@@ -703,6 +761,20 @@ process {
         ]
     }
 
+    // PASS-only per-caller copies for downstream steps; ClairS-TO is covered by VCFSPLIT.
+    // --write-index=tbi: the module's index output is optional and the join needs it.
+    withName: '.*:(CLAIR3|CLAIRS|DEEPVARIANT|DEEPSOMATIC)_PASS_FILTER' {
+        ext.args   = '--apply-filters PASS --output-type z --write-index=tbi'
+        publishDir = [
+            enabled: false
+        ]
+    }
+
+    withName: '.*:CLAIR3_PASS_FILTER'      { ext.prefix = { "${meta.id}_clair3_pass" } }
+    withName: '.*:CLAIRS_PASS_FILTER'      { ext.prefix = { "${meta.id}_clairs_pass" } }
+    withName: '.*:DEEPVARIANT_PASS_FILTER' { ext.prefix = { "${meta.id}_deepvariant_pass" } }
+    withName: '.*:DEEPSOMATIC_PASS_FILTER' { ext.prefix = { "${meta.id}_deepsomatic_pass" } }
+
     withName : '.*:SIGPROFILER_MATRIXGENERATOR' {
         ext.args   = { params.sigprofiler_matrix_args ?: '' }
         publishDir = [