solo: cbMinP posterior threshold, oneExact guard, adapter-anchored geometry

[BenjaminDEMAILLE]

STARsolo barcode handling: cbMinP posterior threshold, the oneExact guard, adapter-anchored geometry, featuresGeneField3.

What changed

cbMinP and oneExact. resolve_multi_cb computed the count+quality posterior for a 1MM_multi barcode, computed the normalising total, then ignored it: any candidate with positive weight won by argmax. STAR requires two more conditions, both now enforced:

• the winner must hold at least 97.5% of the total posterior mass (cbMinP), so a read whose barcode cannot be singled out is discarded rather than guessed;
• unless pseudocounts are in use, the winner must have been observed as an exact match at least once (oneExact). A candidate supported only by a pseudocount is eligible, not evidence.

Adapter-anchored geometry for CB_UMI_Complex: --soloAdapterSequence and --soloAdapterMismatchesNmax, with the four STAR anchor codes (read start, read end, adapter start, adapter end).

The adapter search follows STAR's localAlignHammingDist (SequenceFuns.cpp:341): every offset is scored and the one with the fewest mismatches wins, ties go to the leftmost of the best, and an N in the adapter is a wildcard that costs nothing. Since the barcode positions are measured from the adapter, the chosen offset decides what the read is counted as.

--soloOutFormatFeaturesGeneField3 controls the third features.tsv column; - suppresses it.

Why

Without cbMinP a barcode is assigned whenever any candidate has positive weight, including when two candidates are nearly tied — exactly the case the threshold exists to reject.

This changes default counting behaviour

An existing unit test asserted that priors of 10 against 3 resolve to the higher one: a 77% posterior share, which now correctly resolves to nothing. The test encoded the missing threshold rather than checking it. It now covers a decisive 1000-against-3 prior (accepted), a tie under pseudocounts (rejected), and a winner never seen exactly with pseudocounts off (rejected).

Fixed after review

The adapter search took the leftmost acceptable offset, not the best-scoring one — the opposite of STAR's rule, so an early approximate match beat a later exact one and shifted the whole barcode. Reported by @Psy-Fer with the upstream source. The original tests used single-candidate or exact fixtures, where the two rules agree; there is now a test where they differ, one for the tie-break, and one for N in the adapter.

--soloChemistry is removed. It is not a STAR parameter — it came from STAR-rs — and widening the CLI beyond STAR's is a positioning decision for this project, not something to arrive inside a PR about the barcode posterior. The geometry it set remains reachable through --soloCBstart, --soloCBlen, --soloUMIstart and --soloUMIlen. --soloCBwhitelistDir was described in an earlier version of this text but was never in the diff; that claim is gone too.

--soloCBtype String stays refused: it is a real STAR parameter, and refusing an unimplemented value is the honest answer.

Verification

Gate green: 568 lib + 26 integration, clippy -D warnings, cargo fmt --check, MSRV 1.89. Tests cover the anchored geometry, the three adapter-search rules, the field-3 suppression, and umi_valid_rejects_every_homopolymer (STAR's umiL = 0 quirk lets poly-C/G/T through; documented in docs-old/dev/divergences.md).

The solo differential against real STARsolo passes (test/solo_diff_docker.sh, one run, Linux container with STAR 2.7.11b):

--- 3' SE (CellRanger-style) ---
PASS: rustar-aligner matrix matches real STARsolo exactly.
--- 5' PE (--soloBarcodeMate 1) ---
PASS: Gene raw matrix byte-identical to STARsolo (3 entries)
PASS: GeneFull raw matrix byte-identical to STARsolo (3 entries)

Both runs use --soloCBmatchWLtype 1MM_multi_Nbase_pseudocounts, which is the path cbMinP and oneExact gate, so the thresholds are exercised rather than bypassed.

What this does not establish: the harness fixtures are small and hand-computed, so they show the changed path agrees with STAR on cases where the answer is known, not that counts are unchanged on a real dataset. Reads whose barcode cannot be singled out now land nowhere, so raw matrix counts will move on real data — that is the intended effect of the threshold.

[Psy-Fer]

The --soloChemistry arg doesn't exist in STAR, only in STAR-rs. This isn't made clear at all in the PR

[Psy-Fer]

--soloCBwhitelistDir also doesn't exist in the PR???

[Psy-Fer]

Also the adapter comment claims STAR takes leftmost when STAR actually takes best-scoring

STAR (SequenceFuns.cpp:341, localAlignHammingDist):

uint32 distBest = query.size();
for (ii = 0 .. text.size()-query.size()) {
    dist1 = <mismatches at offset ii>;
    if (dist1 < distBest) { distBest = dist1; pos = ii; }
}
return distBest; 

This PR (find_adapter)

for start in 0..=(seq.len() - adapter.len()) {
    ... count mm ...
    if mm <= max_mm { return Some(...); }
}

These are opposites. The tests written use exact (0-mismatch) or single-match fixtures, so leftmost and best-scoring coincide, which is why it was missed.

[BenjaminDEMAILLE]

All three fixed in 20494db and c8a5d7d.

Adapter search. You are right, and it was a real bug, not just a wrong comment: find_adapter now scores every offset and keeps the fewest mismatches, ties to the leftmost, with N in the adapter as a wildcard that costs nothing — localAlignHammingDist lines 347-358. New test find_adapter_takes_the_best_scoring_match_not_the_leftmost uses a fixture where the rules disagree (approximate match at 0, exact at 5), plus one for the tie-break and one for N.

--soloChemistry. Removed. Not a STAR parameter, so it does not belong in a PR about the barcode posterior. The geometry is still reachable through --soloCBstart, --soloCBlen, --soloUMIstart, --soloUMIlen.

--soloCBwhitelistDir. It was never in the diff — the PR text described it and the code did not have it. The text is corrected.

PR body rewritten. 568 lib + 26 integration green, clippy -D warnings, fmt, MSRV 1.89.

[Psy-Fer]

Please see #154 and the updated CONTRIBUTING.md

[BenjaminDEMAILLE]

All three findings were right, and the adapter one was a genuine correctness bug rather than a documentation problem. Rather than force-push over the review history here, the corrected work is #165. What changed:

Adapter search. Your reading of localAlignHammingDist is what the code now does: score every offset, keep the fewest-mismatch one. The old version returned the first offset inside the mismatch budget, the opposite rule. It mattered more than "where a match is reported", because CB_UMI_Complex measures barcode positions from the adapter, so the wrong offset shifts the whole barcode and changes what the read is counted as. Your point about the fixtures was also correct: single-candidate or exact-match fixtures cannot distinguish the two rules. find_adapter_takes_the_best_scoring_match_not_the_leftmost is the discriminating case, one mismatch at offset 0 against an exact match at offset 5.

Going back to that function properly turned up two more things the old code had wrong: ties go to the leftmost of the best, since STAR replaces its incumbent only on a strictly smaller distance, and an N in the adapter is a wildcard rather than a mismatch, so charging it spends budget STAR does not spend.

--soloChemistry and --soloCBwhitelistDir. Both gone, and they were one mistake rather than two. Neither exists in STAR; they came from the STAR-rs donor and were presented here as parity work without saying so. CONTRIBUTING.md is clear that a non-STAR flag needs sign-off and must not be presented as faithfulness, so reverting was the right call rather than relabelling. If chemistry presets are worth having, they should be proposed as an addition on their own merits.

#154 / CONTRIBUTING.md. Read. Two consequences beyond the above: #165 is rebased onto current main, and it adds the DIVERGENCE.md entry that was missing here. The divergence is that rustar-aligner rejects all four homopolymer UMIs while STAR under CB_UMI_Complex rejects only poly-A, and I checked the source rather than asserting it: the four packed constants are built in the SoloReadBarcode constructor from pSolo.umiL (SoloReadBarcode.cpp:16-21), but that path only learns its UMI length after extracting one (SoloReadBarcode_getCBandUMI.cpp:353-354), so all four stay 0. Poly-A packs to 0 and is caught; the other three do not and slip through.

#165 supersedes this one. Happy to close this as soon as you confirm that is how you want it tracked.