Skip to content

ALT-aware mapping via post-extension liftover groups - #7

Draft
nh13 wants to merge 25 commits into
masterfrom
feat/alt-aware-liftgroup
Draft

ALT-aware mapping via post-extension liftover groups#7
nh13 wants to merge 25 commits into
masterfrom
feat/alt-aware-liftgroup

Conversation

@nh13

@nh13 nh13 commented Jun 25, 2026

Copy link
Copy Markdown
Owner

ALT-aware mapping via post-extension liftover groups

Draft for design feedback, opened against this fork's master. It adds optional ALT-contig
awareness: given a .alt file, reads whose confident placement is deflated by a homologous ALT
contig recover their MAPQ, while genuinely ambiguous reads (paralogs, segmental duplications) stay at
MAPQ 0. Without a .alt file the output is byte-for-byte identical to baseline minibwa — the whole
path is gated on mb_any_alt.

This is the liftover-group design — the same primitive
DRAGEN documents
("liftover groups… each scored taking properly paired alignments into account; the winning group
provides its primary representative as the primary alignment, with MAPQ from the score difference to the
second-best group") — not a port of the bwa-mem .alt / bwa-postalt.js scheme.

The main question for reviewers is at the bottom: is post-extension lifted-placement grouping the right
primitive?

Problem

GRCh38 ships ALT contigs — alternate versions of regions too variable for one linear sequence (HLA/MHC,
many segmental duplications). A read from such a region scores about equally on the main chromosome and
its ALT copy, so a normal aligner drives it to MAPQ 0 or places it on the ALT contig — lost for
primary-assembly variant calling, in exactly the regions that matter most.

Design (liftover groups, post-extension)

Two principles: (1) ALT copies are homologous, not identical — score each on its own merits and defer
the pick to the end; (2) never lift into the chainer (it breaks comput_sc's collinearity/monotonicity
gates) — all grouping is post-extension, on real DP scores.

  1. Span-lift (l2bit): load each .alt POS+CIGAR into an ALT→primary coordinate lift (holes and
    reverse strand handled; unmapped/decoy lines skipped).

  2. Per-hit lifted placement (mb_hit_place): each hit's primary footprint. Grouping key is the
    co-located lifted start (same pri_tid/strand, |Δlifted_st| ≤ tol), multi-interval so an
    SV-breakpoint-spanning ALT hit still groups via a matching sub-interval — not interval-overlap
    fraction (conflates segdup paralogs) and not a scalar k-anchor.

  3. Survival guard (mb_select_sub): diverged ALT twins survive both prune passes to reach
    reconciliation.

  4. Reconciliation (mb_reconcile_alt): group co-located hits; the non-ALT member is the
    representative; recompute group-scoped subsc/dp_max2/n_sub over representatives only, gated on
    query-span overlap with an all-cross-pairs transitive guard so genuine paralogs never merge. An
    ungrouped ALT hit whose query span overlaps a non-ALT hit inside the ALT contig's own primary region
    is folded in too (handles .alt-indel displacement and unliftable insertion holes).

  5. PE integration (pe.c): non-representative ALT twins are excluded from pair enumeration so the
    suboptimal measures the gap to the second-best group-pair; plus a post-rescue regroup and a
    cross-group demotion guard keyed on lifted placement.

    (Updated: this step used to also add --pe-pair-primary, making the SAM primary follow the
    pair-chosen endpoint. Upstream 54fb387 "r417: fixed inconsistent 0x2 flag" now does exactly that
    unconditionally — it sets h[r]->sam_pri = 1 for the pair-chosen endpoint in the paired branch — so
    the flag was redundant and has been removed. SAM output is byte-identical with and without it on
    100k HG002 pairs against hg38 with the real .alt.)

  6. ALT-seed → primary projection (seed.c): when a primary seed is dropped by max_occ subsampling
    in a segdup, inject a native primary anchor at the lifted coordinate so the primary candidate exists
    (scoped to "primary absent from the subsample"; chains normally — not a chainer-merge).

  7. Coin-flip MAPQ damp (pe.c): the residual confident-but-wrong reads are rescue-manufactured
    coin-flips — the anchor end has ≥2 co-optimal copies (a primary plus an ALT twin re-exposed as a
    distinct representative by the grouping), and the mate is placed only by rescue, fitting several of
    them equally. Mate rescue already detects this (it sets skip[r] on a rescued-pair tie) but discarded
    the signal; surface it (rescue_tie) and damp mapq_pe to ~0 for that exact triple. A no-op without
    a .alt (the second co-optimal copy only survives under lifting) — byte-identical to baseline.

CLI: --alt FILE (auto-detected as <idx>.alt), --alt-lift-tol INT (default 10), --alt-records
(emit ALT alignments).
--dbg-alt-proj / --dbg-no-alt-proj are diagnostics in the existing --dbg-* family.

Results

Ground truth is wgsim/holodeck reads with known origin; arms are bwa mem on hs38DH, minibwa --alt
on the GATK index + hs38DH.fa.alt, and minibwa no-.alt. Accuracy = primary placed within tolerance
of the true locus.

(a) Reads from the primary assembly (wgsim, ALT-homologous stratum, 66,082 reads):

Arm accuracy confident (MAPQ≥30) confident-but-wrong
bwa mem 94.0% 58,333 0
minibwa, no .alt 47.1%
minibwa --alt (pair-primary auto-on + coin-flip damp) 94.1% 58,141 8

The follow-up commits took --alt from 86.2% → 91.9% → 94.1%; the confident-and-correct gap to bwa
(reads bwa places confidently+correctly that minibwa does not) collapsed 7,162 → 731 → 210. The
coin-flip damp then halves the confident-but-wrong tail (18 → 8; ALT FDR 0.00031 → 0.00014) at a
+10-read gap (210 → 220) and no runtime cost.

(b) Reads from the ALT haplotype — the direction the feature most directly targets, with exact,
liftover-noise-free truth (holodeck from the 786 ALT/HLA contigs, 1.86M pairs; each read's origin
lifted to primary through the .alt):

Arm recovered → primary, confident stuck on ALT contig
bwa mem 82.8% 0.2%
minibwa --alt 81.6% 1.2%
minibwa, no .alt 0.0% 39.8%

Without the feature, ALT-haplotype reads are never confidently placed on primary; with it, minibwa
matches bwa. (Both arms show ~9% "wrong" here — a shared scoring artifact from single-point lift across
ALT↔primary structural indels; the minibwa−bwa delta is +0.27 pt.)

Corroboration. Against bwa mem + bwa-postalt.js on 5M PE WGS (hs38DH): 97.6% MAPQ-bucket
agreement
; the residual is bidirectional and structural, and the simulations above show minibwa is the
more conservative arm on genuine multi-mappers. Real-data DeepVariant concordance vs GIAB HG002 v4.2.1:
MHC F1 +0.39 SNP / +0.33 INDEL, segdup +0.03–0.04, control flat.

Performance & baseline. Byte-identical without .alt; no-.alt CPU equals baseline. --alt costs
≈ +5.4% CPU / +12.7% wall (1kg HG00096 5M PE WGS); the pair-primary selection and the coin-flip damp
both reuse already-computed fields — no measurable runtime change (2M-pair WGS: 283s vs 284s CPU). The
ksw2 SW kernels are untouched.

Cost and honest limitations

  • Pair-primary is now upstream's job, not this PR's. minibwa's pairing already chose the
    mate-consistent endpoint, but mb_set_sam_pri picked the primary by per-read DP/hash order — emitting
    the wrong paralog copy, which this feature's paralog guard exposed. That was fixed here behind
    --pe-pair-primary; upstream 54fb387 (r417) then fixed it unconditionally, so the flag has been
    removed and this is no longer a cost or a caveat of the feature.
  • Confident-but-wrong: 0 → 8 (FDR 0.00014). Committing the pair-chosen copy adds +18 in hard repeats
    where the truth copy scores ≤ a paralog and bwa stays at MAPQ 0; the coin-flip damp then claws that
    back to +8 by catching the rescue-manufactured ties (16 of 4M correct reads also drop below MAPQ 30, 10
    of which bwa places confidently). The base feature lowers FDR (0.00165 → 0.00027); pair-primary +
    damp lands at 0.00014 with a +10-read recall cost.
  • Residual ~220 reads are low-copy subtelomeric/segdup primary paralogs — not an ALT defect, and
    not max_occ (raising it 20× recovers ~14 for +8% runtime and worsens FDR). This is the
    segmental-duplication mapping limit DRAGEN also handles cautiously at the mapper and resolves
    downstream. A blunt n_pp-aware MAPQ cap was prototyped and refuted (it destroys 1,300–2,900 correct
    recoveries for ~6 fewer confident-wrong); the coin-flip damp succeeds where it failed because it gates
    on the rescue-tie + co-optimal-anchor conjunction, which a genuine recovery does not satisfy.
  • SV-breakpoint output (lifted, split primary-coordinate records, postalt-style) is deferred by
    design — this fixes MAPQ via grouping rather than rewriting coordinates.
  • --alt-lift-tol default 10: sweeping 10→100 changed ~50 reads out of 10M and merged no paralogs.

Validation suite

test/altlg/ — 15 end-to-end fixtures: span-lift, per-hit placement, survival guard, reconciliation
(incl. chimera/query-overlap and displaced/hole twins), PE happy/paralog/rescue, PE-pair-primary,
segdup fwd+rev, indel reject, --alt-records, --alt-lift-tol, unliftable demotion, SV-breakpoint
grouping. All green; baseline byte-identical.

For reviewers

  • Is post-extension lifted-placement grouping the primitive you'd want, or is there a cleaner design
    for ALT-awareness in minibwa?
  • Is auto-on-under-.alt the right default for pair-primary (off without .alt keeps baseline
    byte-identical), and is the rescue-tie + co-optimal-anchor conjunction the right coin-flip
    discriminator for the MAPQ damp?

Commits: (1) the feature; (2) fold displaced/unliftable ALT twins; (3) PE-pair-primary selection;
(4) chore (idiomatic --dbg-* flags, comment clarity); (5) pair-primary default-on under .alt;
(6) coin-flip MAPQ damp. Suggested reading order: l2bit.c (span-lift) → mb_hit_place /
mb_reconcile_alt in map-algo.cpe.cseed.c. The test/altlg/ fixtures double as worked
examples.

minibwa index hs38DH.fa                       # .alt auto-detected as hs38DH.fa.alt
minibwa mem hs38DH.fa r1.fq r2.fq             # pair-primary + coin-flip damp auto-on under .alt
# or explicitly against any index:
minibwa mem --alt hs38DH.fa.alt ref.fa r1.fq r2.fq

@nh13
nh13 force-pushed the feat/alt-aware-liftgroup branch from d6e59d4 to fa30992 Compare June 25, 2026 01:52
@nh13
nh13 force-pushed the feat/alt-aware-liftgroup branch from f4d71a9 to ddc1ddf Compare June 25, 2026 17:16
@nh13 nh13 changed the title Draft: ALT-aware mapping via liftover groups (DRAGEN-style) ALT-aware mapping via post-extension liftover groups Jun 25, 2026
@nh13 nh13 added the feature review vehicle for a manifest feature; never a merge vehicle label Aug 3, 2026
nh13 added 8 commits August 7, 2026 12:56
Optional ALT-contig awareness. With a .alt file (--alt FILE, or auto-detected
<idx>.alt), reads whose confident placement is deflated by a homologous ALT
contig recover their MAPQ, while genuinely ambiguous reads (paralogs, segmental
duplications) stay at MAPQ 0. Without a .alt file the aligner is byte-for-byte
identical to baseline minibwa.

Liftover-group design (not the bwa-mem .alt/postalt scheme): ALT copies are
homologous, not identical, so each is scored independently and the primary
representative is chosen post-extension.

- l2bit: load .alt POS+CIGAR into an ALT->primary span-lift (holes + reverse
  strand); unmapped/decoy lines skipped.
- mb_hit_place: per-hit lifted primary placement; grouping key = co-located
  lifted start, multi-interval (SV-breakpoint aware).
- mb_select_sub: survival guard keeps diverged ALT twins through both prune passes.
- mb_reconcile_alt: post-extension grouping by lifted placement, best-DP
  representative, group-scoped subsc/dp_max2/n_sub recompute with query-overlap
  + transitive guard (paralogs never merge).
- pe.c: PE hooks (exclude non-rep twins from pairing; rescue regroup; chimeric
  demotion guard keyed on lifted placement).
- seed.c: ALT-seed -> primary anchor projection for segdup recovery, scoped to
  primary-absent-from-subsample; per-read scratch on the kalloc arena.
- CLI: --alt, --alt-records, --alt-lift-tol.
- test/altlg: 13 end-to-end fixtures; baseline byte-identical without .alt.
An .alt CIGAR's structural indels can offset an ALT twin's lifted placement from
the read's true primary copy beyond --alt-lift-tol, or drop it into an insertion
hole (unliftable). The per-base lift then fails to co-locate the twin, so it
survives mb_reconcile_alt as a co-equal representative -- taking the SAM-primary
slot and pinning the primary's MAPQ to 0.

Fold an ungrouped ALT hit into a non-ALT hit's group when their query spans
overlap and the ALT contig's .alt-mapped primary region covers the non-ALT hit's
locus (the correspondence the .alt record establishes). Not a categorical
discard: an ALT hit with no overlapping primary stays independent, and two
primary loci never merge (paralog-safe). Within a group the non-ALT member is
the representative.

Ground truth (wgsim from GRCh38 primary, ALT-homologous reads): placement
accuracy 86.2% -> 91.9%; the confident+correct gap to bwa mem shrinks 7,162 ->
731. Byte-identical without .alt. Test: test/altlg/test-altalt.sh.
PE pairing already chooses the mate-consistent endpoint among near-equal paralog
copies, but mb_set_sam_pri marked the SAM primary by per-read DP/hash order --
which can be the other, mate-inconsistent copy. (Masked in the non-ALT path by
the chimeric demotion in pe.c; exposed where the ALT paralog-safety guard skips
that demotion for distinct lifted groups.)

--pe-pair-primary makes mb_set_sam_pri honor the pair-chosen endpoint. Off by
default (plain `minibwa mem` stays byte-identical to baseline); one O(1) index
check per pair, no measurable runtime change.

Ground truth (--alt --pe-pair-primary vs --alt): placement accuracy 91.9% ->
94.1% (bwa mem 94.0%); confident+correct gap to bwa 731 -> 210; confident-but-
wrong 11 -> 18 (FDR 0.0003) -- the cost of committing the pair-chosen copy in
hard repeats where bwa stays cautious. Test: test/altlg/test-pe-primary.sh.
Replace the MB_NO_ALT_PROJECT / MB_PROJ_TRACE environment seams in seed.c with
--dbg-no-alt-proj / --dbg-alt-proj, matching the existing --dbg-* family
(kom_dbg_flag); no getenv remains in the aligner. Rename the lettered pe.c
"hook" comments to describe what they enforce. Tests pass the flags instead of
env vars. No behaviour change; all test/altlg fixtures green, baseline
byte-identical.
PE-pair-primary recovers bwa-parity accuracy in ALT-homologous repeats but was
opt-in (off by default). Without a .alt it is a provable no-op -- the chimeric
demotion in pe.c masks it -- so defaulting it off there cost nothing AND kept
plain `minibwa mem` byte-identical to baseline. With a .alt it is the whole
accuracy win, so opt-in left it off in exactly the runs that want it.

Make it default-on when (and only when) a .alt is loaded:
- opt->pe_pair_primary is now tri-state (-1 auto / 0 off / 1 on), default auto.
- --pe-pair-primary becomes an optional-argument flag matching --long/--adap:
  bare or =yes forces on, =no forces off, absent = auto.
- Resolved to MB_F_PE_PAIR_PRI after index load: auto => on iff l2b->n_alt > 0
  (mb_idx_has_alt), so "off unless .alt" stays the resting state. l2b records the
  loaded ALT-record count; this covers both --alt FILE and an auto-detected
  sibling <idx>.alt.

Net: `minibwa mem --alt X` now delivers the recovery out of the box;
`minibwa mem` is unchanged. Ground truth (--alt, wgsim from GRCh38 primary,
ALT-homologous reads): placement accuracy 91.9% -> 94.1% (bwa mem 94.0%);
confident+correct gap to bwa 731 -> 210; confident-but-wrong 11 -> 18 (FDR
0.0003). Byte-identical without .alt. Test: test/altlg/test-pe-primary.sh
(locks =no -> copyB, auto/explicit -> copyA, chrM byte-identical).
With pair-primary on, the residual confident-but-wrong reads in ALT repeats are
"coin-flip" pairs: the anchor end has >=2 co-optimal copies (a primary plus an
ALT twin re-exposed as a distinct representative by the liftover grouping), the
mate is placed only by rescue, and rescue fits it equally beside several of the
anchor's copies. mb_matesw already detects this -- it sets skip[r] when it finds
two equal-scoring rescued pairs -- but discards the signal, so the surviving
single pair is reported confidently though it is an arbitrary pick of one copy.

Surface that signal (rescue_tie[r]) and damp mapq_pe to ~0 when all hold:
(a) one endpoint is rescue-only, (b) rescue_tie set for the anchor, and (c) the
anchor has >=2 co-optimal representatives (DP second-best within one match score
of best). (a)+(b) keep it off genuine recoveries, where the mate fits exactly
one copy. Reuses already-computed fields; no new passes.

Fires only when ALT lifting has re-exposed the second co-optimal copy as a
representative, so it is a no-op without a .alt (byte-identical to baseline;
verified 0 of 4,000,000 records change). Ground truth (--alt, on top of default
pair-primary): confident-but-wrong 18 -> 8 (ALT FDR 0.00031 -> 0.00014); 16 of
4M correct reads drop below MAPQ 30 (10 of them bwa also places confidently),
head-to-head gap to bwa 210 -> 220. No measurable runtime change. The
byte-identical-without-.alt invariant is covered by the existing test/altlg suite.
The example programs link libminibwa.a, so they need whatever extra link flags
that library was built with. $(LIBS) is empty unless something defines it, so
this is a no-op on its own and supplies the flags when a build that does set it
is in play -- for instance a parallel index build that pulls in the OpenMP
runtime, where the omission fails at link time on Linux but not on macOS.
The coin-flip MAPQ damp and the duplicate-rescued-hit guard both ran on
the baseline path. The damp carried an argument that it was a no-op
without a .alt rather than a check, and the argument was wrong: two
genuine paralogs surviving as representatives with dp_max within opt->a,
plus a rescue tie, satisfies every condition with no ALT anywhere. The
dedup argued it was 'safe' on the non-ALT path, which is true but is not
the byte-identity the manifest promises.

Measured on 100k HG002 WGS pairs against an index with no .alt: before
this commit the SAM differed from stock r421, after it the two are
byte-identical. The chrM fixture the distribution gates on has no ALT
contigs and never reaches either site, so the gate passed the whole time
the claim was false.

Both now test mb_any_alt() on the hits actually in hand, which is the
same predicate pe.c already uses for twin exclusion.
@nh13
nh13 force-pushed the feat/alt-aware-liftgroup branch from 25ff6e0 to 6e4c3d2 Compare August 7, 2026 20:03
nh13 added 16 commits August 7, 2026 15:39
…thing

Three hygiene fixes to the ALT liftover-group feature. No change to mapping
behaviour: the no-.alt path is byte-identical to stock, now demonstrably so.

Keep the lift bookkeeping off the no-.alt path. mb_select_sub allocated a
placement array and called mb_hit_place() for every kept hit whenever `l2b` was
non-null -- but l2b is the sequence index and is always present, so every read
paid for it. The survival guard it feeds can only fire on an ALT hit, and
hit->is_alt is copied from l2b->ctg[tid].is_alt, so with no .alt loaded no hit is
ever ALT and the array is written but never read. Count ALT contigs at load time
and gate on that instead; skipping the block is provably a no-op.

Give the suite a way to run. test/altlg holds 15 scripts that no target and no CI
job invoked, and several drive the api-test helpers, which `all` does not build --
so on a fresh checkout they failed for want of a binary. `make test` now builds
the helpers and runs the suite.

Replace seven vacuous byte-identity checks. Each ran the SAME command twice and
compared the output, then reported "hooks gated => inert without .alt". That only
shows the aligner is deterministic; it cannot see a hook firing. The real claim
needs an unpatched binary, so chrm_baseline() compares against $MB_STOCK and
reports itself skipped when unset rather than passing vacuously. It also indexes
inside its own tmpdir, which stops the suite leaving .mbw/.l2b beside the shared
fixture.

Unrelated pre-existing failures: with the helpers built, 11 of 15 scripts fail on
this r421 base while all 15 pass on the pre-rebase branch (ddc1ddf, r383) with the
same test files. Several failure messages anticipate exactly this ("SA ordering or
max_occ changed"), so it is not yet clear whether the rebase regressed the feature
or upstream's seeding changes invalidated the fixtures. Untouched here; the four
that pass include the new stock comparison, which confirms the no-.alt guarantee.
Eleven of the fifteen ALT scripts failed on r421 with ALT records simply absent
from the output. The cause is not the feature: `--outn` and `--alt-records` are
wired into the `map` option chain only, and the bwa-compatible `mem` chain
silently ignores unknown long options -- no diagnostic, exit status 0. Every
script asked for secondary records via `mem --outn=N`, so on r421 they asked for
nothing and then failed on the records that were never requested.

`minibwa map --outn=5` emits the twin (secondary 0x100, MAPQ 0) alongside the
MAPQ 60 primary, which is what the assertions describe. `mem` honoured `--outn`
at the branch's original base (r383) and no longer does, which is why this
surfaced only after the rebase onto r421.

Switch the 42 mapping invocations to `map`. The two subcommands differ only in
that `mem` also emits MD:Z, and nothing here asserts on MD.

With the api-test helpers built and MB_STOCK pointing at a stock r421 build, the
suite goes from 4/15 to 14/15, and the new no-.alt stock comparison passes on
every script that runs it -- the PE hooks are byte-identical to stock without a
.alt file, as the manifest claims.

test-pe-primary is the one remaining failure and is a stale fixture, not a code
defect: on r421 copyB comes out as a supplementary rather than a competing
full-length hit, and the mate-consistent copy already wins the primary slot
without the flag the test exists to exercise.
`make test` builds the api-test helper programs from tracked sources and, in two
scripts that index the shared chrM fixture in place, leaves .l2b/.mbw beside it.
Neither is source, and both showed up as untracked noise in every worktree that
ran the suite.
The flag made the SAM primary follow the PE-pair-chosen endpoint instead of the
first-by-DP/hash representative. That was a real fix when it was written: on the
pre-r417 base it moved R1's primary from the mate-inconsistent paralog copy to the
mate-consistent one, and it was verified to be a complete no-op without a .alt (0
of 400k primary records changed), because the competing copy only survives as a
representative when Hook C's r_any_alt-gated paralog-safety skip protects it.

Upstream then fixed the same thing itself. 54fb387 ("r417: fixed inconsistent 0x2
flag") sets `h[r]->sam_pri = 1` in the paired branch -- h[r] IS the pair-chosen
endpoint, &hit[r][paux.i[r]] -- and sets reset_sam_pri = 0 so the per-read pass
cannot override it. That is this flag's behaviour, unconditional and flagless.

It has also been inert since the rebase onto r421: pri_idx is assigned only inside
the paired branch, and the mb_set_sam_pri call that reads it now sits behind
upstream's `if (reset_sam_pri)`, which is false on exactly that path. Both empirical
checks agree -- --pe-pair-primary=yes vs =no differs by 0 records on its own fixture
and on 100k HG002 pairs against hg38 with the real 3,177-record .alt (ALT machinery
confirmed live there: 10,467 records differ with vs without the .alt).

So this removes redundancy, not behaviour: SAM output on that same 100k/.alt run is
byte-identical before and after. mb_set_sam_pri returns to its upstream 3-argument
signature, mb_idx_has_alt goes with it (the auto default was its only caller), and
test-pe-primary drops -- it asserted an on/off difference that upstream made
impossible. The suite is now 14/14.
minibwa.h included l2bit.h on the line directly above MB_VERSION, which
upstream bumps on every release. The two edits therefore landed in one
conflict hunk, so assembling this feature onto current master always
conflicted in minibwa.h:

  <<<<<<< master
  #define MB_VERSION "0.7-r424-dirty"
  =======
  #include "l2bit.h"

  #define MB_VERSION "0.7-r421"
  >>>>>>> feat/alt-aware-liftgroup

git rerere keys its cache on the conflict pre-image, and that pre-image
contained the version string, so a resolution recorded against one release
could never match the next one. The cache shows the result: resolved at r416
and r421, unresolved at r419 and r424, one dead entry per release. The
feature dropped out of the distribution assembly every time upstream tagged.

l2b_t is needed only by the mb_hit_place() prototype 200 lines below, so the
include moves down beside the other forward declarations. No semantic change
-- the header still compiles standalone and the aligner is untouched.
`make test` failed on any fresh checkout:

  test-output.sh: Abort trap: 6   minibwa map --outn=5 test/chrM-human.fa.gz ...

test-output.sh used the shared chrM fixture without indexing it, and the
index is a gitignored build artifact, so on a clean tree there was none.
Stock minibwa aborts in that case (`[E::main_map] failed to load the index.
ABORT!`), which is what the trap was -- not a fault in this feature.

It passed in practice only because test-projrev.sh and test-segdup.sh
indexed the shared fixture *in place* as a side effect, leaving untracked
.mbw/.l2b in test/. So the suite was order-dependent: it passed on a second
run, or whenever another script happened to go first, and failed on a clean
tree when test-output.sh led.

All three now copy chrM into TMPD and index it there, which is what
lib-baseline.sh already documents and does. The suite is hermetic and
order-independent, and test/ is left clean.
mb_idx_load() auto-detected <prefix>.alt; mb_idx_load_mmap() did not.  The same
index and the same command line therefore produced different MAPQs depending on
whether --mmap was passed -- on the ALT fixture, 60 without it and 0 with it,
because the ALT twin looked like an ordinary repeat when no lift blocks were
loaded.  Nothing warned; the only symptom was quietly worse mapping quality.

Move the resolution out of the loaders and into main_map(), via a new
mb_idx_set_alt_auto().  One place decides, so the two loaders cannot disagree,
and the precedence between --alt FILE and the adjacent file is visible where it
is chosen instead of being an override applied after the fact.  That also drops
a redundant parse: with --alt FILE, mb_idx_load() used to parse the adjacent
.alt only for mb_idx_set_alt() to reset and re-parse it immediately after.

The cost is that mb_idx_load() no longer loads .alt for library callers.  The
two api-test examples that need it already re-loaded it explicitly and now call
mb_idx_set_alt_auto(); ex-one and ex-batch inspect no ALT fields.

Covered by test/altlg/test-noalt.sh in the following commit, which asserts
--mmap parity with ALT both on and off.
ALT-aware mapping turns itself on whenever a .alt sits beside the index, which
matches bwa's convention but leaves no way to turn it back off.  That matters
more here than it does in bwa: bwa only marks ALT hits and defers the real work
to bwa-postalt.js, which you run or you don't, whereas this recomputes MAPQ
inside the aligner.  Anyone on a stock hs38DH index -- which ships a .alt -- got
different MAPQs from stock minibwa with nothing in the command line or the log
to say why, and no way to check the difference short of moving the file aside.

--no-alt loads no .alt at all.  is_alt is set only by l2b_set_alt() and is never
stored in the index, so loading nothing leaves mb_any_alt(), n_alt_ctg and
mb_hit_place() all on their stock paths; there is no new branch in the hot loop
and no ALT parse to pay for.  It beats an explicit --alt in either order, so it
is usable as a late override on a command line assembled elsewhere.

test-noalt.sh keys on the fixture's 60-vs-0 MAPQ split and checks that --no-alt
does something, that it is byte-identical to deleting the .alt outright, that it
wins over --alt in both orders, and that --mmap agrees with the normal loader
with ALT on and off -- the regression fixed in the preceding commit.
… --alt

Two gaps in the preceding two commits, both found by self-review.

Moving .alt resolution out of mb_idx_load() and into the caller updated one
caller and missed the other: main_mem() loads the index at map-main.c:709 and
never resolved .alt, so the bwa-compat `mem` subcommand silently stopped being
ALT-aware.  On the ALT fixture it went from MAPQ 60 with two supplementaries to
a single MAPQ 0 record.  The suite did not catch it because 6ca77f2 moved the
ALT tests off `mem` and onto `map` on purpose, leaving that entry point with no
coverage at all.  Of the seven mb_idx_load*() call sites, this was the only
remaining behaviour-bearing one: fastmap is SMEM-only and ex-one/ex-batch read
no ALT fields.

Separately, mb_idx_set_alt() returned void and dropped l2b_set_alt()'s -1, so a
mistyped `--alt /path/typo.alt` loaded nothing and aligned on as if ALT-awareness
were off -- no diagnostic even at -v 3, and output indistinguishable from a
correct non-ALT run.  It now returns the record count like l2b_set_alt(), and an
explicitly named file that will not load is an error.  The auto-detect branch
still ignores its status: an absent adjacent .alt is the normal case, not a
mistake, and that asymmetry is the point.

test-noalt.sh grows a `mem`-vs-`map` parity check (both with and without the
.alt, so it cannot pass for an unrelated reason) and a check that an unreadable
--alt reports the failure and emits no alignments.  Both were confirmed to fail
when their fix is reverted.  The bad-path check asserts on stderr rather than
exit status because this branch predates r422 (lh3#67), where main() still
discarded every subcommand's return value.

Also corrects l2b_set_alt()'s header doc, which said it returns the number of
ALT contigs; l2bit.c counts records there and keeps the contig count separately
in n_alt_ctg.
README told the reader to use a reference without ALT contigs and listed
"does not recognize alternate haplotypes" under Limitations, both of which this
branch makes false. On the distribution that matters more than it looks: `dist`
is the fork's default branch, so this file is the front page every clone lands
on, and update_readme() inserts the generated feature banner while deliberately
preserving upstream's prose around it. The result was a page advertising
alt-liftgroup fifty lines from a paragraph denying it.

Corrected here rather than in the renderer so the claim tracks the feature. If
alt-liftgroup ever drops from an assembly, this commit drops with it and
upstream's original wording -- which is accurate for a build without the
feature -- comes back on its own.
Three defects in loading a .alt, all silent.

The CIGAR walk treated H as consuming nothing, so a supplementary record whose
leading hard clip IS its offset into the ALT contig lifted every block to the
wrong place. H removes bases from SEQ but they still occupy the contig, so it
advances the ALT axis like S and I. Reverse records walked the reverse
complement and stored RC coordinates, which l2b_lift()'s forward binary search
then read as forward ones; they are now remapped to ALT-forward. And blocks were
appended in file order, which is non-monotonic once reverse records or
multi-record contigs exist, so the binary search was operating on unsorted data:
each contig's blocks are now sorted, with a total order on (alt_st, alt_en,
pri_tid, pri_st, rev) because qsort is unstable and ties would otherwise let
libc decide the output.

Input from a .alt is untrusted. POS < 1 would underflow to a garbage primary
coordinate; an unknown CIGAR operator or a length with no operator silently
truncated the walk; and the reverse remap subtracts from ctg->len, so a record
whose query span runs past the contig it names wrapped on uint64 to a
near-UINT64_MAX coordinate that sorted to the end of lift[] and read back as a
valid block. That last one raised MAPQ from 0 to 60 on a 300M record against a
150bp ALT contig -- silent, and in the confidence-inflating direction. All four
are now rejected per-record, rolling back that record's partial blocks.

Separately, kom_grow() uses realloc(), which does not zero new slots, so the ALT
fields of a contig were uninitialized on the plain `minibwa index` path where
l2b_set_alt() never runs. l2b_destroy() then freed a garbage ctg->lift, crashing
`minibwa index` on any multi-contig reference.

test-badalt.sh gates the four rejections. It first proves a well-formed reverse
record loads (MAPQ 60) and that a no-.alt run differs (MAPQ 0): without those
controls the rejection assertions would pass just as loudly if ALT support were
broken outright.

The suite's shared helpers (lib.sh, fixlib.sh) land here rather than with the
later dedup commit, because this is the first test that sources them and a
commit whose own test cannot run is not worth bisecting to.
mb_hit_place/mb_place_t/mb_subpl_t/MB_MAX_SUBPL lived in the public minibwa.h,
which dragged the entire private l2bit interface into the public API (via an
`#include "l2bit.h"`) and exposed a function that cannot be called through the
public surface (mb_idx_t is opaque; the api-test probe already reaches into
mbpriv.h for the argument). Move them to mbpriv.h, which already includes both
minibwa.h and l2bit.h, and drop the l2bit.h include from the public header.
MB_LIFT_TOL stays public (opt->lift_tol references it).

Also add the MB_DBG_NO_ALT_SURVIVE debug bit used by the survival-guard test
seam.
…on seam

The doc described a `pref` parameter the function does not take, and claimed the
pairing pass disambiguates the primary here. It does not: selection in
mb_set_sam_pri() is per-read, and pe.c sets sam_pri on the pair-chosen endpoint
itself (`h[r]->sam_pri = 1`). A doc describing a parameter that does not exist is
worse than none.

--dbg-no-alt-survive ablates just the mb_select_sub survival guard, so a fixture
can assert the RED state -- the co-located ALT twin dropped -- as well as the
GREEN one. A guard whose test only ever sees it working cannot tell "the guard
holds" from "the fixture never exercises it".
`mem` auto-loaded an adjacent .alt with no way to point at another one or turn
it off, so the bwa-compat entry point silently had a different ALT surface from
`map`. Both now resolve through one mb_load_alt(), and `mem` accepts --alt,
--no-alt, --alt-records and --alt-lift-tol.

The resolution is also no longer silent about a user error. An explicitly named
--alt that will not load is fatal rather than a quiet fall back to stock
behaviour; one that loads zero usable records warns; and --alt together with
--no-alt says which one won. An adjacent .alt that is simply absent stays quiet,
because that is the normal case rather than a mistake.

`mem` carries its own long-option table rather than sharing map's. Sharing it
would make `mem` PARSE every knob map has while implementing only these, so one
it has no arm for -- --meth, --mmap, --outn, --eqx, --hic -- was accepted and
dropped on the floor: `mem --meth` ran without methylation and said nothing. It
now rejects unknown options and a missing argument the way map does, and handles
--version and --help itself.

--alt-records no longer consumes the non-ALT secondary budget: it is emitted as
an extra that bypasses the out_n cap, so counting it against n_sec let it crowd
out genuine secondaries the user asked for with -h.

usage_mem now documents the ALT knobs, and annotates bwa's `-j` as the primitive
minibwa replaces rather than leaving it an unexplained ignored flag.
mb_distinct_lifted_group recomputed mb_hit_place(q) for q == h[r] (the
pair-chosen endpoint, loop-invariant) on every iteration of the chimeric-hit
loop. Add an mb_distinct_lifted_group_pb() variant that takes q's precomputed
placement and hoist it out of the loop. Only runs when the mate has ALT hits, so
the no-.alt path is unchanged.
…e guard

The suite had no fixture for a reverse .alt record or for an ALT contig carrying
more than one, which is exactly where the span-lift walked the reverse complement
and where blocks arrived out of order -- the two defects the lift fix addresses.
test-revmulti.sh covers both.

test-survive.sh asserted only that a co-located ALT twin survives sub-hit
selection. A test that never sees the guard off cannot distinguish "the guard
holds" from "the fixture never reaches it", so it now drives --dbg-no-alt-survive
and asserts the twin is dropped without the guard.

The helpers were copy-pasted across the suite, in places with incompatible
per-file variants of the same name. lib.sh and fixlib.sh hold one version of
each, with signatures that are supersets of every prior call site so the
optional trailing arguments default to "match anything".

Also drops .gitignore patterns for build artifacts the suite no longer produces.
The step-2 comment claimed the transitive close-up is O(n_hit^2), but the
all-cross-pairs merge guard makes it O(n_hit^4) in the worst case. Correct the
comment: document the real bound, why it stays near-quadratic in practice
(n_hit is small under default max_occ/out_n subsampling -- a 40-copy segdup
measured n_hit=41; the grp[j]==grp[i] short-circuit plus the fact that
co-locating hits of one read also query-overlap means they merge on first
contact rather than re-triggering the guard), and why a union-find is not a
drop-in replacement (it cannot express the all-cross-pairs guard). No behavior
change.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

feature review vehicle for a manifest feature; never a merge vehicle

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant