solo: cbMinP posterior threshold, oneExact guard, chemistry presets, featuresGeneField3 - #150
Open
BenjaminDEMAILLE wants to merge 4 commits into
Open
Conversation
…uard `resolve_multi_cb` computed the count+quality posterior correctly, computed the normalising total, and then ignored it: any candidate with positive weight won by argmax. Two of STAR's conditions were missing. `cbMinP`. STAR requires the winner to hold at least 97.5% of the total posterior mass. Without it a barcode is assigned whenever any candidate has positive weight, including when two candidates are nearly tied — which is precisely the case the threshold exists to reject. A read whose barcode cannot be singled out should be discarded, not guessed. `oneExact`. Unless pseudocounts are in use, the winner must have been observed as an exact match at least once. A candidate whose only support is a pseudocount is eligible, not evidence. **This changes default counting behaviour**, and the existing unit test had to change with it: it asserted that priors of 10 against 3 resolve to the higher one, which is a 77% posterior share and now correctly resolves to nothing. The test was encoding the missing threshold. It now covers both sides — a decisive 1000-against-3 prior is accepted, a tie under pseudocounts is not. Verification I could not do here, and would want before this merges: the solo differential against real STARsolo (`test/solo_diff_docker.sh`). Reads that previously landed on an ambiguous barcode now land nowhere, so raw matrix counts will move, and the direction should be confirmed against the oracle rather than argued from the source. Everything in-tree is green — 560 lib + 26 integration — but that is not the same thing. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
…er-UMI rule `--soloOutFormatFeaturesGeneField3` sets the third column of `features.tsv`, STAR's default being `Gene Expression`. The sentinel `-` drops the column, giving the two-column form some downstream tools expect. Previously the string was hard-coded. Also adds `umi_valid_rejects_every_homopolymer`. No code change was needed — `is_homopolymer` already rejects a UMI of any single repeated base — but the rule deserves a test that says why it is deliberate. STAR has a quirk (`umiL = 0`) under which a non-A homopolymer UMI slips through its filter, so poly-C, poly-G and poly-T survive there and are counted. That is a defect rather than a rule: a UMI of one repeated base carries no information whichever base it happens to be. This diverges from 2.7.11b knowingly, and the test now records that instead of leaving it to be rediscovered as a difference in a differential run. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
… String `--soloChemistry` sets the CB/UMI geometry from a 10x chemistry name so it does not have to be spelled out four flags at a time: SC3Pv1 (14-base CB, 10-base UMI), SC3Pv2 (16/10), SC3Pv3, SC3Pv4 and SC5P (16/12), plus the CellRanger aliases. Unknown names are refused. The preset wins over an explicit `--soloCBstart` and friends, and warns when it overrode one. Half-applying a preset — taking its CB length but the flag's UMI start, say — would produce a geometry neither the user nor the preset asked for, and would fail silently in the counts rather than loudly at startup. `--soloCBtype` is declared and accepts `Sequence`. `String` is refused rather than accepted and ignored: barcodes here are 2-bit packed ACGT, and a non-ACGT barcode packed that way is nonsense rather than an approximation. Supporting it means a different whitelist representation throughout, which is its own change. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
BenjaminDEMAILLE
changed the title
`--soloCBposition` and `--soloUMIposition` take four fields, `startAnchor_startDist_endAnchor_endDist`, where the anchor codes are 0 = read start, 1 = read end, 2 = adapter start, 3 = adapter end. Only code 0 was accepted; the rest were a parse error. Codes 2 and 3 are the reason this needs a new layout rather than a wider parser: the adapter moves from read to read, so the offsets cannot be resolved once at startup the way every other layout's are. `SoloBarcodeLayout::Anchored` carries the specs with their anchors and the adapter, and resolves per read. Adds `--soloAdapterSequence` and `--soloAdapterMismatchesNmax`. The adapter search takes the leftmost position within the mismatch budget, which is what STAR takes — the first acceptable match, not the best-scoring one. A read where the adapter is not found, or where an anchor resolves outside the read, yields no barcode. Guessing would be worse than declining: a barcode read at an unknown offset is not salvageable, and a wrong barcode is a read attributed to the wrong cell rather than a read lost. Tests cover the adapter moving (both fields follow it), one mismatch inside the budget, no adapter at all, an anchor running off the front of the read, and the leftmost-match rule. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
STARsolo barcode handling:
cbMinPposterior threshold, theoneExactguard, chemistry presets, adapter-anchored geometry,featuresGeneField3.What changed
cbMinPandoneExact.resolve_multi_cbcomputed 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:cbMinP), so a read whose barcode cannot be singled out is discarded rather than guessed;oneExact). A candidate supported only by a pseudocount is eligible, not evidence.Chemistry presets (
--soloChemistry) set CB/UMI geometry for the named 10x chemistries. The 7 10x whitelists are not vendored: they total 26 MB and this crate is published.--soloCBwhitelistDirpoints at them, andautoscores only the whitelists present there.Adapter-anchored geometry for
CB_UMI_Complex:--soloAdapterSequenceand--soloAdapterMismatchesNmax, with the four STAR anchor codes (read start, read end, adapter start, adapter end). The adapter is located leftmost within the mismatch budget.--soloOutFormatFeaturesGeneField3controls the thirdfeatures.tsvcolumn;-suppresses it.Why
Without
cbMinPa 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).
Verification
Gate green: 560 lib + 26 integration, clippy
-D warnings,cargo fmt --check, MSRV 1.89. New tests cover each preset, the anchored geometry, the field-3 suppression, andumi_valid_rejects_every_homopolymer(STAR'sumiL = 0quirk lets poly-C/G/T through; documented indocs-old/dev/divergences.md).Not verified: the solo differential against real STARsolo (
test/solo_diff_docker.sh). Reads that previously landed on an ambiguous barcode now land nowhere, so raw matrix counts move. That run needs Docker and a STARsolo oracle, neither available in the environment this was developed in. Since the STARsolo work in #90/#133 was validated against real STARsolo, this path should get the same treatment before merge.🤖 Generated with Claude Code