Variant filter - #200
Variant filter#200robert-a-forsyth wants to merge 28 commits into
Conversation
…ysis in Deepvariant by default (toggleable)
# Conflicts: # CHANGELOG.md # subworkflows/local/phasing_haplotyping.nf # subworkflows/local/small_variant_consensus.nf # tests/clair_only.nf.test.snap # tests/consensus.nf.test.snap # tests/default.nf.test.snap # tests/union.nf.test.snap
…l joins BCFTOOLS_NORM ran with -Oz only, so multi-allelic records were never split. BCFTOOLS_ISEC matches on exact CHROM/POS/REF/ALT, so a site one caller reports as A>G,GT and the other as A>G never intersected: in consensus mode the variant was dropped, in union mode it appeared twice. ClairS-TO emits no multi-allelic records while DeepVariant emits ~2.5% (128,115 of 5,058,527 on B1975944), so every DeepVariant multi-allelic site was systematically excluded from the consensus. Adding -m -any splits them first; verified on real DeepSomatic output (45,638 -> 45,677 records over chr1:1-5Mb, 39 split, none left multi-allelic). This also makes the AF declaration mismatch harmless: ClairS-TO declares FORMAT/AF as Number=1 while DeepVariant and DeepSomatic declare Number=A, and bcftools concat only warns before keeping the first file's definition. Once every record is biallelic the two declarations are equivalent. The comment on STANDARDIZE_AF is corrected accordingly: it claimed DeepVariant emits VAF, but all four callers emit FORMAT/AF and none emits VAF, so the rename is a no-op under the default prioritize_caller='clair' and is kept as a guarantee for WAKHAN rather than as a live conversion. Two silent-failure paths in SMALL_VARIANT_CONSENSUS are closed. The caller branch had no `other` arm, so an unrecognised meta.caller removed the sample from every downstream result with a successful exit; it now errors. The DeepVariant/Clair join was a plain join, so a sample present for one caller but not the other was dropped just as quietly; it now uses failOnMismatch and failOnDuplicate, matching PHASING_HAPLOTYPING. Record counts in the merged output change, so the snapshots need regenerating. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
With --germline_var_keep deepvariant the tumour-only germline arm runs DeepVariant
on the tumour BAM. DeepVariant is a germline caller with no somatic
discrimination -- its FILTER vocabulary is only PASS/RefCall/LowQual/NoCall -- so
those calls mix germline and clonal somatic variants and cannot be separated on
their own. Published unchanged, germline_smallvariants.vcf.gz and vep/germline/
carried most of the somatic call set.
DeepSomatic evaluates the same sites and does emit a verdict: FILTER=GERMLINE
("Non somatic variants"), PON, RefCall or PASS. That verdict is now transferred
onto the DeepVariant records as INFO/DS_VERDICT and used to select the germline
arm, mirroring how ClairS-TO adjudicates its own calls through NonSomatic and
VCFSPLIT. This is a caller's own adjudication rather than positional subtraction
of the somatic call set: on B1975944 only 1.14% (57,684) of DeepVariant's
5,058,527 PASS calls carry a positive somatic verdict, against ~87% positional
overlap with the union somatic set.
Only positively adjudicated germline sites are kept (GERMLINE or PON). RefCall and
sites DeepSomatic never evaluated are dropped rather than assumed germline, so the
arm is not overpopulated with unadjudicated calls; it keeps 83.1% of the input.
The verdict stays in INFO/DS_VERDICT so the decision is auditable in the published
VCF.
The transfer is three independent bcftools invocations, so it is three aliased
instances of the existing upstream modules (DS_VERDICT_QUERY, DS_VERDICT_ANNOTATE,
DS_GERMLINE_SELECT) rather than a new bespoke process. DeepSomatic FILTER is
single-valued in practice (RefCall/GERMLINE/PON/PASS over 13.7M records), so
transferring it as a string cannot inject the ";" that would break INFO parsing.
DEEPSOMATIC now runs before DEEPVARIANT in the subworkflow because the verdict is
built from its raw VCF, before the PASS filter discards the GERMLINE records. The
deep family must therefore be enabled as a pair: --germline_var_keep deepvariant
without deepsomatic in --somatic_var_keep is rejected at launch rather than
silently producing an unadjudicated germline arm.
Paired mode is untouched -- both Clair3 and DeepVariant already run on the normal
BAM there, so the germline arm needs no adjudication.
Also corrects the --smallvar_filter_pass docs: VCFTAG normalises FILTER to PASS,
so `false` does not restore the pre-filter behaviour as previously claimed.
Germline record counts change, so the snapshots need regenerating.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
…bites VCFTAG rewrites FILTER to PASS on both arms before phasing, which is needed so that downstream tools filtering on PASS still see the germline records ClairS-TO marks NonSomatic. But it did so destructively, and every record reaching the published phased VCFs, VEP and the report therefore read FILTER=PASS with the caller's verdict gone. That also silently disabled the signature-stage filter. SIGNATURES_BCFTOOLS_VIEW ran --apply-filters PASS on a call set whose FILTER had already been normalised to PASS, so it could never remove anything: with --smallvar_filter_pass false, RefCall/LowQual/GERMLINE/PON records reached SigProfilerMatrixGenerator that previously could not, and the matrices were computed over mostly-reference sites. The parameter was documented as restoring the previous behaviour, which it did not. VCFTAG now records the caller's FILTER in INFO/ORIG_FILTER before overwriting it, as VCFSPLIT already did for the ClairS-TO arm, and the signature filter tests that field instead. VCFSPLIT's stamp is respected rather than duplicated: both the header line and the per-record field are added only when absent, since a duplicate INFO key makes the record unparseable. Multi-valued FILTER is joined with "," because ";" separates INFO fields, and a FILTER of "." is recorded as "." rather than skipped. bcftools accepts only one of -i/-e, so the ALT="*" exclusion is folded into the same include expression. Verified against the module fixture with bcftools 1.20, including a record that already carried ORIG_FILTER and one with a multi-valued FILTER: exactly one ORIG_FILTER header line, no double stamping, and htslib parses the result. The module test is extended to cover these cases. Note that the module test cannot currently run locally: nf-test 0.9.3 generates `VCFTAG(*input)` and the available Nextflow 26.04.1 rejects the spread operator under its strict syntax, while the pipeline requires >=25.10.4 so 25.04.6 will not run either. This is pre-existing and reproduces unchanged at 5b38dcf. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The test read its output through path(...).vcf, an accessor provided by an nft-vcf plugin that nf-test.config does not load, so it failed with MissingPropertyException on every run and had never passed. It now uses the built-in linesGzip accessor instead, which needs no plugin. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
STANDARDIZE_AF and BCFTOOLS_ANNOTATE were two back-to-back bcftools annotate calls on the same file with disjoint options -- one renaming the allele frequency FORMAT key, the next stamping INFO/CALLER. --rename-annots composes with -a/-c/-h, verified against a real VCF, so a single invocation does both and the alias is removed along with its modules.config block. The rename is selected by meta.rename_to, set only when combine_method is 'all'; in 'consensus' mode every surviving record comes from one caller and needs no rename. Ordering is safe because BCFTOOLS_QUERY reads only CHROM/POS/REF/ALT and does not care whether the rename has happened. The rendered ext.args was checked for all three meta cases (VAF, AF, key absent): the escapes reach printf intact and an absent key yields no --rename-annots. The joins around the annotate step now use failOnMismatch/failOnDuplicate for the same reason as the caller join: a missing annotation table should stop the run, not quietly drop the sample. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
CLAIR3 received --sample_name twice: the module passes it from ext.prefix, which defaults to meta.id, and conf/modules.config passed the identical value again. Harmless with argparse, but it reads as if the name were configurable. Removed the config copy, matching how CLAIRSTO is already handled on dev. The comment in PHASING_HAPLOTYPING claiming the somatic VCF is passed first "(higher priority in phasing)" was false: BCFTOOLS_CONCAT sorts its input file list alphabetically, so with the _germline_tagged/_somatic_tagged prefixes the germline file is always passed first. The ordering is cosmetic under -a, which emits in coordinate order; the comment now says so rather than describing an intent the module discards. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Picks up TAG_GERMLINE and TAG_SOMATIC, and restores the clairsto and lrsomatic_report versions that the previous snapshot run had rolled back. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
# Conflicts: # tests/clair_only.nf.test.snap # tests/consensus.nf.test.snap # tests/deep_only.nf.test.snap # tests/default.nf.test.snap # tests/union.nf.test.snap
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
# Conflicts: # CHANGELOG.md
|
Nextflow's strict syntax resolves `report_empty_slot()` against declared functions, not local closure variables, so it read as undefined. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
ljwharbers
left a comment
There was a problem hiding this comment.
Thanks @robert-a-forsyth, this fixes real problems: the all union dropping the prioritised caller's private calls, the -T positional somatic recovery, and the DeepVariant/DeepSomatic record explosion. The provenance tagging is the right way to do it. A few things to address before merge, most important first.
Should fix
-
Paired-only runs are now rejected at launch.
validateDeepCallerPairingfails any run withdeepvariantin--germline_var_keepbut notdeepsomaticin--somatic_var_keep, whatever the samplesheet contains.--germline_var_keep clair,deepvariant --somatic_var_keep clairon a paired samplesheet was valid (DeepVariant runs on the normal BAM there) and now errors, with a message about tumour-only mode. The check is needed only becauseTUMORONLY_SMALLVARis always invoked and now readsDEEPSOMATIC.out.
Suggestion: inTUMORONLY_SMALLVAR, runDEEPSOMATICwhendeepsomaticis insomatic_var_keepordeepvariantis ingermline_var_keep, and add its calls to the somatic set only in the first case. Then drop the launch check. -
The verdict transfer drops multi-allelic germline sites.
DS_VERDICT_QUERY/DS_VERDICT_ANNOTATEmatch on exactCHROM,POS,REF,ALT, but neither the DeepVariant nor the DeepSomatic VCF has been split at that point. If DeepVariant reportsA>G,Tand DeepSomatic reportsA>G, the record gets noDS_VERDICTandDS_GERMLINE_SELECTremoves it. It's the same mismatch this PR fixes for isec with-m -any, and it probably accounts for part of the "4.3% unevaluated". Pleasebcftools norm -m -anyboth inputs before the transfer. -
--smallvar_filter_pass falsepublishes rejected calls asFILTER=PASS.
VCFTAGrewrites FILTER to PASS on every record, on both arms and for every caller. Withfalse, RefCall/LowQual/NonSomatic records therefore reach the phased VCFs, VEP and the report marked PASS, sofalsedoes not "restore the previous unfiltered behaviour" asusage.mdsays. Only the ClairS-TO germline arm needs the PASS rewrite, andVCFSPLITalready does it. CouldVCFTAGstamp the flag without touching FILTER? -
DS_GERMLINE_SELECTignores DeepVariant's own FILTER.
With--smallvar_filter_pass false, a DeepVariant RefCall at a site DeepSomatic callsGERMLINEenters the germline arm. AddingFILTER="PASS"to the-iexpression fixes this whatever the parameter is set to.
Should address
-
Test coverage. PR CI runs only tag
small, so neither the newVCFTAGmodule test nor theextendedunion/deep_only/consensus tests (the only ones that reach the DeepSomatic verdict path) ran here. The green checks don't cover them. Please run the extended set on the cluster and post the result. Also, commit 33bb9b2 says theVCFTAGtest covers a record that already hasORIG_FILTERand a multi-valued FILTER, buttests/fixtures/vcftag_input.vcfhas neither. Please add both. -
groupKey sizes count optional outputs.
groupKey(meta, 3)for ASCAT and Wakhan countssegments_rawandsolution_dirs, which are bothoptional: true. A sample missing one waits for the whole batch (or may be dropped), which is the clog this is meant to fix. -
Split multi-allelics downstream (please check). In consensus mode,
-m -anynow leaves two records at one position in the VCFs Longphase and Wakhan read. Can you confirm both handle that, or rejoin withnorm -m +anyafter isec? -
VCF_PASS_FILTERas a subworkflow. Why is this a subworkflow when it runs only one module? If the reason is the pass-through whensmallvar_filter_passis false, that can be done where each caller is handled: includeBCFTOOLS_VIEWunder the four aliases and choose between its output and the input with a ternary. That drops a file and a layer ofwithNamenesting (.*:DEEPVARIANT_PASS_FILTER:PASS_FILTER). -
Comment and changelog length. Many of the new inline comments are long explanations: the
VCFSPLITandVCFTAGscripts,SIGNATURES_BCFTOOLS_VIEW, the verdict transfer block intumoronly_smallvar.nf,small_variant_consensus.nf. They also repeat sample-specific numbers (B1975944, B2037620) that will go stale. Please cut them to a line or two saying what each step does and why. The detail belongs in the commit messages and the PR description. Same forCHANGELOG.md: each entry should be one or two sentences on the user-visible change, matching the other entries, not a write-up of the investigation.
Docs / nits
modules/local/bcftools/view/meta.ymlstill documents the removed-Ttargets inputs, andmodules/local/vcftaghas nometa.yml.- The
usage.mdtable still lists--germline_var_keep/--somatic_var_keepas "Default =null" (they're'clair'). - The
ORIG_FILTERrow says "FILTER value in the ClairS-TO output", butVCFTAGsets it for every caller. DS_VERDICT_QUERYwrites all ~13.7M DeepSomatic rows while onlyGERMLINE/PONare used. Filtering there would cut the I/O a lot, and the select step becomes a plain "DS_VERDICTpresent" test.- The local
BCFTOOLS_VIEWnow hard-codes-i 'INFO/SOMATIC=1'in its script. It belongs inext.args, or the nf-corebcftools/viewalready used in this PR can do it. - The
ORIG_FILTERawk block is duplicated inVCFSPLITandVCFTAG, and the two already describe the field differently. - Worth a check: the new
ALT!="*"may not exclude the same records as the old--exclude 'ALT="*"'when a record has several ALTs. Keeping it an exclude (-e 'ALT="*" || …') avoids the question.
Drop the hard-coded -i from the module script and set it in the PHASING_HAPLOTYPING:BCFTOOLS_VIEW config block. Refresh meta.yml: drop the removed targets inputs, generic description, document versions output/topic. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
VCFTAG now only stamps the provenance flag, so --smallvar_filter_pass false keeps each caller's FILTER. ORIG_FILTER is stamped only on the ClairS-TO germline split, the one arm rewritten to PASS, and SIGNATURES_BCFTOOLS_VIEW goes back to --apply-filters PASS with --exclude 'ALT="*"'. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Only GERMLINE/PON verdicts are consumed downstream, so the verdict table no longer carries the ~9M RefCall and PASS rows (13.4M -> 4.3M on B1975944). The filter excludes PASS/RefCall rather than including GERMLINE/PON, because bcftools errors on a FILTER ID the header does not declare. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
BCFTOOLS_NORM splits multi-allelics (-m -any) so isec can match per ALT,
which left two records at one POS in the VCF handed to LongPhase and
Wakhan. BCFTOOLS_NORM_REJOIN (-m +any) now rejoins them at the end of
SMALL_VARIANT_CONSENSUS on both the consensus and 'all' paths. The final
file keeps its <id>_{germline,somatic}_sorted name; the sort steps now
write *_split_sorted.
Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
groupKey(meta, 3) counted optional outputs, so a sample missing one was silently dropped at channel close. ASCAT now joins its required outputs and adds segments_raw with remainder: true. Wakhan always writes a solution_* link for every row of solutions_ranks.tsv, so solution_dirs is now required and the three outputs are joined. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Drop the one-module VCF_PASS_FILTER subworkflow; each call site now picks between the filtered copy and its input on params.smallvar_filter_pass. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Tumour-only DeepVariant germline calls are adjudicated with DeepSomatic's verdict only when deepsomatic is selected; otherwise they pass through. Drops validateDeepCallerPairing, which rejected valid paired-only runs. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
…wo lines Addresses review point 9 on #200. Comments in the consensus, verdict transfer, phasing, VCFSPLIT and modules.config blocks now say what each step does in a line or two; sample IDs and cohort counts are removed from source and docs. The #200 CHANGELOG entries are rewritten to match the length of the others, with short Fixed lines for the multi-allelic rejoin and the per-sample report. No code changes. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
PASS_FILTER is now reported per alias, BCFTOOLS_NORM_REJOIN is new, and the multi-caller modes phase from rejoined multi-allelics, so union/consensus haplotagged BAMs and WhatsHap stats change. Snapshots regenerated with Nextflow 25.10.4. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
PR checklist
nf-core pipelines lint).nextflow run . -profile test,docker --outdir <OUTDIR>).nextflow run . -profile debug,test,docker --outdir <OUTDIR>).docs/usage.mdis updated.docs/output.mdis updated.CHANGELOG.mdis updated.README.mdis updated (including new tool citations and authors/contributors).