Skip to content

Pool the empirical error table across reference individuals (SRP434573), then examine and update the paper #49

Description

@davmlaw

Summary

The SRP434573 real-data validation currently builds one empirical error table per mixture pair, each from only that pair's two reference individuals. Because a site only supplies the alt->ref direction where a reference individual is hom-alt (the minority case), most sites carry just one substitution direction, and the low-rate tail is estimated from few reads. The seven individuals (F1-F3, M1-M4) share the panel and sequencing chemistry, so the substitution background is shared and can be pooled across them.

Pooling should give near-complete both-direction per-site coverage (some individual is hom-ref, some hom-alt, at nearly every polymorphic site) and a tighter low-rate tail. Both hit exactly the low-dilution regime (0.5%, 1%) where the presence test currently falls back to the flat --error-rate default at the many hom-alt-donor markers.

Follow-up to #41 (panel-of-normals error table) and #23 (per-patient error tables).

What is already wired up (pending TAU-side regen)

Code and config are in place; the table itself needs regenerating on /tau (BAMs are not accessible in dev).

  • pipeline/Snakefile: new opt-in build_pooled_error_table (default off), enumerates unique reference individuals (deduplicated, so an individual shared by several mixtures is counted once), and four rules: pooled_targets -> pileup_refs_panel_pooled + pileup_refs_bg_pooled -> error_table_pooled producing pooled.error_table.tsv. Per-patient tables still build as fallback.
  • paper/public_data/SRP434573/config.yaml: build_pooled_error_table: true.
  • paper/scripts/run_srp434573_allomix.py: resolve_error_table() prefers the pooled table for every mixture, else the mixture's own, else the flat default. Wired into the real, semi-synthetic two-person, and semi-synthetic three-person runs.
  • paper/public_data/SRP434573/README.md: "Error tables" section documents the pooled build and copy step.

Verified in dev: Snakefile parses, all four rules register, the pooled.error_table.tsv target resolves its full rule chain (stops only at the /tau FASTQ wall), run script compiles and passes ruff. Fallback is safe: with no pooled table committed, the build uses the existing per-patient tables.

Tasks

  • Regenerate TAU-side and commit the pooled table:
    snakemake -s pipeline/Snakefile --configfile paper/public_data/SRP434573/config.yaml --cores 16
    cp output/genotypes/SRP434573/pooled.error_table.tsv paper/public_data/SRP434573/genotypes/
  • Examine the result before trusting it. Check per-direction site coverage: e_altref currently covers ~300-550 sites per per-patient table; pooled should approach the full ~1300 for both directions. If it does not jump, the union-targets or sample-naming assumption needs a second look.
  • Rebuild the affected real-data rules and compare low-dilution presence numbers (0.5%, 1% rows) before vs after pooling. This is where the change should move.
  • Sanity-check that pooling did not distort any per-site rate (spot-check a few high-rate sites against their per-patient values; the rate should be a read-weighted average, not shifted by one individual).
  • Alter the paper if the results differ significantly. If the pooled table materially changes the real-data presence LoD or the 0.5%/1% detection numbers, update the affected figures/facts and the Methods/Results text. If the change is small, note that pooling was used but the numbers are unchanged. Either way the Methods "training cohort" wording already covers pooling.

Notes

  • Off by default in the general pipeline on purpose: pooling is only sound when the patients share a run / flowcell / panel / chemistry. A pipeline run over patients from unrelated runs keeps the per-patient tables.
  • Pooling is for the stable substitution background only. Run-specific index hopping and cross-sample contamination stay in the separate contamination machinery (contamination table + per-sample floor), not the error table.
  • For SRP434573 we cannot confirm a single flowcell (SRA stripped the @RG PU tags), so "same flowcell" is framed as a shared-panel/shared-chemistry assumption, not an asserted fact.

Metadata

Metadata

Assignees

No one assigned

    Labels

    enhancementNew feature or request

    Type

    No type

    Fields

    No fields configured for issues without a type.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions