Skip to content

fix: APA re-analysis P0 fixes (Fixes #19 Fixes #20 Fixes #21 Fixes #22 Fixes #23 Fixes #24 Fixes #25 ) - #26

Merged
berntpopp merged 10 commits into
mainfrom
fix/apa-p0
Oct 15, 2026
Merged

berntpopp merged 10 commits into
mainfrom
fix/apa-p0

Conversation

@berntpopp

Copy link
Copy Markdown
Member

Implements the APA re-analysis P0 fixes (berntpopp/apa-sequencing#34).

Fixes #19 Fixes #20 Fixes #21 Fixes #22 Fixes #23 Fixes #24 Fixes #25

Commits

  • 00b2b4d feat(mutect2): support semicolon-delimited multi-region tumor BAMs per sample
  • 790fcdf feat(mutect2): build cohort Panel of Normals with hotspot masking
  • f1e3e86 fix(mutect2): split multiallelics with bcftools norm after FilterMutectCalls and add extra passthrough
  • 387ad80 fix(profiles): fix BIH SLURM account/memory conflict, preserve retry scaling, and pass HMM threads
  • 00d6008 fix(scatter): add target BED restriction, wire interval inputs into DAG, and fix 1-based BED conversion
  • cce549b fix(samples): enforce sample uniqueness, string dtypes, and analysis_type-normal_bam consistency
  • 9d25729 fix(mutect2): resolve -normal sample name from BAM @rg SM tag with collision guard

Implementer report

fig.v8+.yaml), mutect2.smk, freebayes.smk

Tests & Linters

  • pytest (/home/bernt-popp/miniforge3/bin/pytest tests/ -v): 166 passed in 0.53s (including unit, schema, DAG input resolution, and shell execution tests in tests/test_helpers.py).
  • ruff check & ruff format --check (workflow/rules/helpers.py scripts/generate_config.py): passed.
  • mypy (workflow/rules/helpers.py scripts/generate_config.py): passed (0 issues in 2 source files).
  • snakemake, snakefmt, and shellcheck: binaries are not installed in this local environment; rule syntax/DAG wiring and shell blocks were verified via Python exec() and bash -euo pipefail subprocess tests.

Risks & Notes

  • Cohort PoN normal selection (#24): When gatk_resources.panel_of_normals is "", pon_genomicsdb requires at least one PoN-eligible normal BAM from get_normal_bams() (analysis_type == "tumor_normal", filtered by in_pon when that column is present); normal-only (pon_only) libraries not paired with a tumor should be included via tumor_normal rows or a pre-built cohort PoN if no tumor_normal row references them.
  • Multi-region tumor @RG SM tags (#25): gatk Mutect2 requires distinct @RG SM tags across multi-region tumor BAMs (T1;T2) if per-region genotypes are to remain separate in the output VCF; if two tumor BAMs share the same SM, GATK pools their reads into a single tumor sample column.

🤖 Implemented by Gemini 4 Argon via agy; adversarially reviewed by Claude Opus 5.5 before merge.

berntpopp and others added 7 commits October 15, 2026 19:39
…llision guard

Extract @rg SM tags via samtools in mutect2_call instead of using the normal_bam filename stem, and fail if SM tags are missing, multi-sample, or equal between tumor and normal.

Fixes #19

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
…type-normal_bam consistency

Read samples.tsv with dtype=str and verify_integrity=True, validate sample uniqueness and analysis_type-normal_bam consistency, and guard normal helpers on analysis_type == tumor_normal.

Fixes #20

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
…AG, and fix 1-based BED conversion

Add scatter.target_bed with INTERSECTION interval rules across callers and GetPileupSummaries, wire scattered interval/BED files into mutect2_call and freebayes_call inputs, and convert interval_list to BED via gatk IntervalListToBed.

Fixes #21

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
…scaling, and pass HMM threads

Use hpc-ag-scholl account and single mem_mb resource across SLURM profiles, remove static set-resources overrides so attempt-scaling lambdas apply within the 7-day limit, and pass --native-pair-hmm-threads to Mutect2.

Fixes #22

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
…ctCalls and add extra passthrough

Add params.filter_mutect_calls.extra passthrough, pin bcftools and htslib in gatk.yaml, and normalize FilterMutectCalls output with bcftools norm -m-any without atomizing complex indels.

Fixes #23

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
Add cohort PoN construction (Mutect2 tumor-only -> GenomicsDBImport -> CreateSomaticPanelOfNormals) when panel_of_normals is empty, optional in_pon normal filtering, and bcftools hotspot masking via pon_mask_bed.

Fixes #24

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
…r sample

Split semicolon-delimited tumor_bam entries across helpers and common.smk, pass multiple -I tumor arguments to Mutect2, compute contamination/segmentation per (sample, bam_basename), pass all tables to FilterMutectCalls, and skip multi-BAM samples in PureCN.

Fixes #25

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
@berntpopp

Copy link
Copy Markdown
Member Author

Adversarial review round 1 (Claude Opus 5.5): output could not be parsed

@berntpopp

Copy link
Copy Markdown
Member Author

Adversarial review round 1 (Codex GPT-6 Astra, ponytail): request_changes

Request changes for ambiguous contamination paths, duplicate tumor SM handling, omitted standalone PoN normals, and ignored FreeBayes targets. All workflow files parsed; 134 read-only tests passed. No tracked files changed. Not checked: 32 temporary-file tests, full DAG execution, or GATK/SLURM execution under the read-only sandbox. Risks include failed analyses, merged region identities, and incomplete PoN construction.

  • major workflow/rules/mutect2.smk:393: The new {sample}_{bam_basename} pattern is ambiguous because both wildcards allow underscores. Snakemake parses APA1_APA1_T.contamination.table as sample=APA1_APA1 and bam_basename=T. The normal lookup then fails or selects another analysis. This breaks existing single-BAM analyses too.. Fix: Use separate path components, such as {sample}/{bam_basename}.contamination.table, consistently for contamination, segmentation, logs, benchmarks and helpers. Add a DAG regression test with underscores in both identifiers.
  • major workflow/rules/mutect2.smk:208: Multiple tumor BAMs may share the same SM tag; the new guard only checks tumor versus normal. T1;T2 with identical tumor SM values passes validation, collapsing regions into one sample. Their separate contamination tables subsequently collide in GATK's sample-keyed map. See GATK's contamination filter.. Fix: Reject duplicate tumor SM values before calling Mutect2 and test two tumor BAMs sharing an SM.
  • major workflow/rules/helpers.py:181: The cohort PoN cannot include the standalone eligible normals required by SPEC.md. A valid germline row with tumor_bam=N1, normal_bam='.' and in_pon=yes contributes nothing because selection requires tumor_normal. Putting N1 in normal_bam instead fails the new consistency validation. Clean normals from excluded tumors are therefore silently omitted.. Fix: Support explicitly PoN-eligible standalone normals using the existing germline rows and their primary BAM column. Preserve matched-normal selection and test a cohort containing both kinds.
  • major workflow/rules/common.smk:195: With the default chromosome scatter mode, FreeBayes ignores scatter.target_bed. The helper returns only '-r chr1', and get_freebayes_interval_input returns no target dependency. Setting capture targets therefore still calls whole chromosomes and emits off-target variants.. Fix: Apply the target BED restriction together with the chromosome restriction and declare the BED input. Verify that the resulting intervals are their intersection.
  • minor tests/test_helpers.py:567: The tests named for PoN DAG construction and multi-region calling never construct a Snakemake DAG. They execute helpers and inspect strings or mocked command arguments. No test invokes Snakemake, so the explicit dry-run acceptance criteria for fix(profiles): fix BIH SLURM account/memory conflict, preserve retry scaling, and pass HMM threads #22, feat(mutect2): build cohort Panel of Normals with hotspot masking #24 and feat(mutect2): support semicolon-delimited multi-region tumor BAMs per sample #25 are unmet, and the wildcard regression remains undetected.. Fix: Add small Snakemake dry-run fixtures covering both profiles, interval scattering, generated masked PoN, and multi-region inputs with underscores. Assert scheduled rules and rendered contamination inputs.
  • minor tests/test_helpers.py:3: Ruff reports five I001 import-order violations in the changed tests, and its format check fails. Snakefmt also requests formatting changes for the new Snakefile sample-loading expression.. Fix: Sort and format the changed imports and expressions using the repository's configured tools.

#19: ok, #20: ok, #21: incomplete, #22: incomplete, #23: ok, #24: incomplete, #25: incomplete

Separate calculate_contamination wildcards with slash, reject duplicate tumor SM tags across multi-region BAMs, include standalone germline in_pon normals in PoN/NormalDB, intersect FreeBayes chromosome scatter with target_bed, and add Snakemake dry-run and linter fixes.

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
@berntpopp

Copy link
Copy Markdown
Member Author

Adversarial review round 2 (Codex GPT-6 Astra, ponytail): request_changes

Request changes: normalization loses allele-filter information, target restriction fails on empty chromosome intersections, and PoN construction misses an interval fallback. Selected tests: 49 passed; changed-file lint, formatting and diff checks passed. Tracked files remain unchanged. Full GATK/SLURM execution and file-writing integration tests were not run in this read-only environment; silent false-positive admission remains the principal risk.

  • major workflow/rules/mutect2.smk:462: Normalization corrupts allele-specific filtering. GATK emits pipe-separated AS_FilterStatus values, while bcftools splits Number=A fields on commas. With ALT=C,G, FILTER=PASS and AS_FilterStatus=SITE|weak_evidence, bcftools 1.21 produces two PASS records and AS_FilterStatus=. for G. I reproduced this locally. Downstream PASS filtering therefore admits the rejected allele. See GATK filtering implementation.. Fix: Preserve GATK's per-allele filter statuses when splitting and apply each status to the corresponding biallelic FILTER, retaining site-wide filters. Add a real multiallelic VCF regression test; the current command-recording mock cannot detect this.
  • major workflow/rules/helpers.py:126: Target restriction still schedules chromosomes absent from the BED. For example, the default chromosome list includes chrM, but an exome BED may omit it. Mutect2 then receives an empty intersection and aborts; the corresponding FreeBayes job also fails on its empty target stream. See GATK's empty-intersection handling and FreeBayes target loading.. Fix: When target_bed is supplied, restrict chromosome scatter units to configured chromosomes represented in the BED, preserving their order. Use that same list for callers and gathering, and test a BED missing a configured chromosome.
  • major workflow/rules/mutect2.smk:65: Cohort PoN construction incorrectly depends on scatter.chromosomes outside chromosome mode. A valid configuration with mode=none, chromosomes=[], target_bed='' and panel_of_normals='' generates no -L arguments, so GenomicsDBImport fails after normal calling finishes. GenomicsDBImport requires intervals for a new workspace.. Fix: Use reference-derived intervals when neither a target BED nor applicable chromosome restriction exists. Add a PoN test for the supported configuration without chromosome scattering or targets.
  • minor tests/test_helpers.py:820: The new Snakemake integration test lacks the dryrun marker. Consequently, make test-dryrun selects no tests, while make test-unit includes this integration test whenever Snakemake is available. I confirmed pytest -m dryrun collects zero tests.. Fix: Mark TestSnakemakeDryRun with pytest.mark.dryrun so the existing test commands select it correctly.

#19: ok, #20: ok, #21: incomplete, #22: ok, #23: wrong, #24: incomplete, #25: ok

Preserve GATK per-allele AS_FilterStatus and apply to biallelic FILTER across bcftools norm in filter_mutect_calls; restrict chromosome scatter units to chromosomes present in target_bed; fall back to reference-derived scattered.interval_list in pon_genomicsdb when no target_bed or chromosome scatter applies; mark TestSnakemakeDryRun with pytest.mark.dryrun.

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
@berntpopp

Copy link
Copy Markdown
Member Author

Final adjudication of the Codex review (Gemini 4 Argon, ponytail): approve

All four review findings are already resolved in HEAD (5db8871): multiallelic normalization preserves GATK AS_FilterStatus and site-wide filters, chromosome scatter restricts to chromosomes present in target_bed, cohort PoN GenomicsDBImport falls back to reference N-scattered intervals when needed, and TestSnakemakeDryRun carries the dryrun marker.

Finding Decision Reason
workflow/rules/mutect2.smk:462 AS_FilterStatus corruption during bcftools norm -m-any valid-fixed Commit 5db8871 preserves pipe-delimited AS_FilterStatus across bcftools norm -m-any and restores allele-specific and site-level FILTER flags, verified by real bcftools regression tests.
workflow/rules/helpers.py:126 chromosome scatter units not restricted to chromosomes present in target_bed valid-fixed Commit 5db8871 filters configured chromosomes via get_bed_chromosomes(target_bed) in get_scatter_units, used consistently across callers and gather rules.
workflow/rules/mutect2.smk:65 pon_genomicsdb missing -L intervals when mode=none, chromosomes=[], target_bed='' valid-fixed Commit 5db8871 adds get_pon_interval_input/build_pon_interval_arg falling back to scatter_intervals_by_ns when neither target_bed nor chromosome units apply.
tests/test_helpers.py:820 TestSnakemakeDryRun missing pytest.mark.dryrun valid-fixed Commit 5db8871 decorates TestSnakemakeDryRun with @pytest.mark.dryrun so pytest -m dryrun and pytest -m 'not dryrun' select tests properly.

Tests: 171 passed (170 unit, 1 dryrun); ruff check and ruff format --check passed

Co-Authored-By: Gemini 4 Argon <noreply@google.com>
@berntpopp
berntpopp merged commit 2955fb6 into main Oct 15, 2026
2 checks passed
@berntpopp
berntpopp deleted the fix/apa-p0 branch October 18, 2026 23:33
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment