Skip to content

Harden RSS LD sketch chromosome merging and indel event identification - #1404

Merged
gaow merged 6 commits into
StatFunGen:mainfrom
jaempawi:agent/rss-ld-sketch-afreq-event-id
Aug 19, 2026
Merged

gaow merged 6 commits into
StatFunGen:mainfrom
jaempawi:agent/rss-ld-sketch-afreq-event-id

Conversation

@jaempawi

@jaempawi jaempawi commented Aug 14, 2026 •

Copy link
Copy Markdown
Contributor

Adds canonical directional event IDs for same-position REF/ALT-swap ("mirror") indel pairs in the RSS LD-sketch pipeline, so insertion vs. deletion anchoring can't flip effect-allele orientation during association harmonization.

Scope

The change spans three parts of the rss_ld_sketch module:

  1. merge_chrom (core feature) — detects indels that form an exact REF/ALT-swap pair at the same position and assigns each a canonical event ID (chr:pos:INS:x / chr:pos:DEL:x, DEL anchored at pos+1). The event ID replaces the variant ID in .pvar and .afreq; the original↔event mapping is written to a .event_id.tsv sidecar. Only mirror-pair indels are relabeled — SNPs, equal-length substitutions, and non-mirror indels keep standard IDs. If a panel has no mirror pairs, no sidecar is written and .pvar/.afreq stay fully standard. .pgen is never modified.

  2. process_block (supporting change) — VCF shard discovery now accepts .vcf.gz in addition to .bgz. Both extensions validated end-to-end.

  3. Test fixture + coverage — the rss_ld_sketch test fixture now carries one same-position REF/ALT-swap indel pair (A/AT + AT/A at chr22:16500000) injected into the existing 60-sample cohort, converted to .vcf.gz. Adds expected/event_id.tsv, regenerates expected/afreq_deterministic.tsv (673→675 variants), and adds test_rss_merge_chrom_mirror_event_ids asserting the INS:T/DEL:T relabeling and the sidecar contents. All 5 tests pass.

Validation

  • Reproduces the deployed R5 EUR chr1 panel exactly (48,014 mirror event IDs, byte-identical ID column). Same criterion verified across all four deployed R5 panels (EUR/AFR/EAS/SAS).
  • Full pipeline (generate_w → process_block → merge_chrom) runs end-to-end via SoS on the augmented fixture, producing the mirror event IDs; the no-mirror path (0 pairs → no sidecar, standard IDs) is also exercised on the real cohort.

Notes for reviewers

  • Scope grew beyond merge_chrom to include the process_block .vcf.gz discovery change and the test fixture/coverage — flagged here for transparency.
  • Known limitation (separate issue, not addressed here): merge_chrom requires ≥2 LD blocks per chromosome — a single-block chromosome fails at plink2 --pmerge-list ("requires at least two filesets").

Anak Empawi added 3 commits August 18, 2026 15:50
Narrow event-ID assignment from all indels to only indels that form an
exact REF/ALT-swap pair at the same position (the ambiguous cases where
INS vs DEL anchoring can flip effect-allele orientation). Non-mirror
indels, SNPs, and equal-length substitutions retain standard IDs.

merge_chrom now substitutes the canonical event ID into .pvar and .afreq
for the mirror set, and writes no .event_id.tsv when a panel has no
mirror pairs (leaving .pvar/.afreq fully standard). Notebook narrative
updated to match. Validated against deployed R5 EUR chr1 (reproduces
48014 mirror event IDs exactly).
@danielnachun
danielnachun force-pushed the agent/rss-ld-sketch-afreq-event-id branch from de4c783 to 2738bca Compare August 18, 2026 22:50
Anak Empawi added 3 commits August 19, 2026 09:56
Test coverage for the mirror-pair event-ID feature:
- process_block VCF discovery now accepts .vcf.gz in addition to .bgz
  (rss_ld_sketch.R); both extensions validated.
- Augment the rss_ld_sketch test fixture with one same-position REF/ALT-swap
  indel pair (A/AT + AT/A at chr22:16500000, ~40% AF, well-called) injected
  into the existing 60-sample cohort; convert fixture to .vcf.gz.
- Regenerate expected/afreq_deterministic.tsv (673->675 variants) and add
  expected/event_id.tsv (the mirror-pair mapping).
- Update test_rss_ld_sketch.py: .bgz->.vcf.gz, single-dot cohort-id, 673->675,
  and a new test_rss_merge_chrom_mirror_event_ids asserting the INS:T/DEL:T
  relabeling and the .event_id.tsv sidecar.
- Notebook narrative: name the .vcf.gz fixture, note both extensions accepted,
  show the mirror pair in the Output example, correct stale W_B50.npy -> .rds.

All 5 tests in test_rss_ld_sketch.py pass.
merge_chrom conflated two responsibilities: assembling the per-chromosome
pgen, and canonicalizing mirror-pair indel IDs. Separate them:

- Extract canonicalize_mirror_event_ids(chrom, pos, id, ref, alt) as a pure
  function that computes the mirror-pair event-ID mapping from variant records
  and returns it (0-row data.frame if no pairs). No file I/O.
- do_merge_chrom now calls it, then owns reading/writing/ID-substitution as
  before. Output is byte-identical (5 pipeline tests unchanged; full SoS
  pipeline reproduces the same 675-variant / INS:T+DEL:T result).
- Guard the CLI with sys.nframe()==0L so the worker is sourceable for unit tests.
- Add a direct unit test (test_canonicalize_mirror_event_ids_unit + helper R
  script) exercising the rule on an in-memory fixture in ~1s: mirror pair
  relabeled, non-mirror indel and SNP excluded, empty input -> empty mapping.
  Surfaced and fixed a latent empty-input edge case.
Follow-through on the canonicalization extraction: pull the remaining inline
concerns out of do_merge_chrom so each step is a named, focused function and
the merge step reads as a sequence of calls.

- reconcile_afreq(): gather per-block .afreq, validate against .pvar, reorder,
  atomically install; returns pvar.
- apply_event_id_substitution(): write the .event_id.tsv sidecar and substitute
  event IDs into .pvar/.afreq (collision- and duplicate-guarded, atomic).
- summarize_block_filters(): tally per-block .meta counts and print the summary.
- cleanup_merge_intermediates(): remove plink2 merge temps and block dirs.
- Hoist read_tab() to a file-level helper shared by the above.

do_merge_chrom drops from ~113 lines to a ~25-line orchestrator:
reconcile -> canonicalize -> (substitute | note-standard) -> summarize -> cleanup.
Behavior byte-identical: 6 tests pass and the full SoS pipeline reproduces the
same 675-variant / INS:T+DEL:T result with identical reconcile/substitution/
filter-summary messages.
@gaow
gaow merged commit 440af8d into StatFunGen:main Aug 19, 2026
3 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants