solo: --soloUMIfiltering MultiGeneUMI_All (closes #144) - #152
Open
BenjaminDEMAILLE wants to merge 3 commits into
Open
solo: --soloUMIfiltering MultiGeneUMI_All (closes #144)#152BenjaminDEMAILLE wants to merge 3 commits into
BenjaminDEMAILLE wants to merge 3 commits into
Conversation
BenjaminDEMAILLE
changed the title
BenjaminDEMAILLE
changed the title
BenjaminDEMAILLE
changed the title
BenjaminDEMAILLE
changed the title
BenjaminDEMAILLE
added a commit
to BenjaminDEMAILLE/rustar-aligner
that referenced
this pull request
Jul 29, 2026
CONTRIBUTING.md requires the description to match the code; the entries for the other themes split out of scverse#152 belong to their own PRs.
BenjaminDEMAILLE
added a commit
to BenjaminDEMAILLE/rustar-aligner
that referenced
this pull request
Jul 29, 2026
CONTRIBUTING.md requires the description to match the code; the entries for the other themes split out of scverse#152 belong to their own PRs.
…GeneUMI `--soloUMIfiltering MultiGeneUMI_All` resolved to the same variant as `MultiGeneUMI`, which is neither what STAR does nor what the option is documented to do. Of the three available behaviours it was the only one nobody had asked for. In STAR the option is a no-op: it is parsed and stored, but its consumption site tests only the `MultiGeneUMI` flag, so selecting it leaves the filter entirely off. Documented, it removes a UMI seen in more than one gene from *all* of them, rather than from the losers only. `UmiFiltering::MultiGeneUmiAll` now exists and does the documented thing: a UMI appearing in several genes is evidence of a collision or of chimeric amplification, so it is discarded outright rather than attributed to whichever gene happened to read deepest. Single-gene UMIs are untouched, which the test checks across every mode. Raised upstream as scverse#144 before changing it, since "be faithful to STAR" and "do what the flag says" genuinely point in opposite directions here. Also adds `docs-old/dev/divergences.md`, recording this and the homopolymer-UMI rule, so deliberate differences are written down rather than rediscovered as surprises in a differential run. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
BenjaminDEMAILLE
force-pushed
the
bd/solo-umi-filter
branch
from
July 29, 2026 09:33
5b4456e to
1383b87
Compare
BenjaminDEMAILLE
changed the title
The homopolymer-UMI rule describes behaviour already on main, not something this diff changes, and its locking test does not exist yet. It belongs in its own change rather than riding along here.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
This was referenced Jul 30, 2026
Psy-Fer
added a commit
that referenced
this pull request
Aug 5, 2026
* fix: --runThreadN 1 ran on every core, not on one The rayon global pool was configured only when `--runThreadN` was greater than 1. Skipping the build at 1 does not give one thread: it leaves rayon's default, which is one worker per logical core. So `--runThreadN 1` ran the whole machine. Measured on 200k reads, before: 2.44 s wall at **1300% CPU**. After: 26.23 s wall at 100% CPU. The old figure was not a fast single-threaded run, it was a sixteen-way run wearing the wrong flag. This matters beyond the flag reading falsely. A scheduler or a container given one CPU gets sixteen worker threads; on a shared machine the run oversubscribes every other job; and with a thread-caching allocator each of those threads keeps its own heap, which is the very cost the comment above this code says the pool sizing exists to avoid. It also means the project's thread-invariance checks were weaker than they read: the `--runThreadN 1` leg was not a one-thread leg. Verified now that it is one: records are byte-identical between 1 and 8 threads on 200k real reads, and byte-identical to the previous binary's output at `--runThreadN 1`. Only the `@PG` `CL:` line differs between thread counts, because it records the command line. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> * docs(changelog): record the --runThreadN 1 fix Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> * refactor: drop Transcript::read_seq, which nothing ever reads The field was filled with a full copy of the read at every finalised transcript and assigned again in four more places, and no code in the crate reads it. Not "reads it rarely": the compiler was asked, and after deleting the field every one of the 97 resulting errors is a struct literal or an assignment. There is no read site, in `src` or in tests. Measured on 200k real reads at 8 threads: about 50 transcripts are finalised per read, so removing it takes **15.1 million allocations off a 200 million total** (7.6%) and 3.87 GB of copying. Wall clock does not move, and that is worth recording rather than hiding: six interleaved rounds at 87-92% CPU idle give medians 20.90 s against 20.65 s with the direction mixed, inside the run-to-run spread. mimalloc is fast enough that seventy-five small allocations per read do not surface. The reason to remove it is that it is dead weight, not that it is slow. That number also calibrates #168 downward: if removing 7.6% of the allocations changes nothing measurable, the rest of the allocation programme is unlikely to be worth a new dependency. `Transcript::read_seq` is `pub`, so this is an API removal and needs sign-off. Nothing outside the crate can be relying on its contents being meaningful, though, since it is only ever written. Output-neutral: SAM byte-identical on 200k real reads. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> * docs(changelog): record the Transcript::read_seq removal Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> * test: an A/B harness that refuses to average over a busy machine Every perf measurement on this project needs the same three guards, and I got each of them wrong at least once in a single session: An earlier harness sampled CPU idle once, before the first run. Three separate measurements then drifted below the threshold mid-series and it kept going, leaving me to spot the contaminated rounds by eye in the output. This samples idle before and after every run, drops a round if any of the four samples falls short, and prints how many it dropped, so a median over four surviving rounds cannot be mistaken for a median over six. The check is on CPU idle rather than load average, because load average is an exponential average over minutes: it refused to measure at 2.24 on a machine whose cores were all free. It reports the spread within each side next to the difference between the medians, and says so in as many words when the difference is smaller. Two changes I measured looked like wins on medians alone and were inside the spread. Both sides run as ./rustar-aligner with --outFileNamePrefix ./ from inside their own directory, because the @pg CL: line records argv verbatim: running ./old against ./new is enough to make the output differ, which cost me two false "output is not neutral" alarms. The header documents the fourth trap, which no script can enforce: timing a total hides the part that changed. BAM writing is 1-4% of a yeast run, so a total dominated by alignment cannot resolve a change to the writer. Run the None configuration alongside and read the difference. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> * feat(solo): --soloCellReadStats CB writes CellReads.stats One row per cell barcode with the fifteen counters STARsolo reports: how the barcode matched, whether the read mapped to one locus or several, whether it landed on a feature, where in the gene and on which strand, whether it was mitochondrial, and whether it reached the matrix. The per-cell UMI and gene totals come from the raw matrix rather than from the read counters, so they agree with what the matrix says by construction. Reads whose barcode never resolved are summed into a single `CBnotInPasslist` row instead of being dropped. That row is the reason the file is useful: it is the difference between "these cells look thin" and "most of the input never reached a cell at all". The region columns split by strand — an antisense read counts under `exonicAS` or `intronicAS`, never under `exonic` or `intronic`. `--genomeChrSetMitochondrial` names the chromosomes behind the `mito` column. Without it the column is zero throughout, which is honest: no chromosome was declared mitochondrial. D24 comes with it. STAR emits these rows by walking a libc++ `unordered_map`, which at these sizes is the reverse of each barcode's first appearance. That is reproduced, including across threads: the per-read accumulator merges in read order, so a threaded run writes the same file as a serial one. It stops being reproducible past the point where libc++ rehashes, since the order then depends on the bucket count. The values never differ, only which line they sit on. Recorded in docs-old/dev/divergences.md. * docs: record the CellReads.stats row order in DIVERGENCE.md Section 3.2, in the format CONTRIBUTING.md asks for. * docs(changelog): keep only this PR's entry CONTRIBUTING.md requires the description to match the code; the entries for the other themes split out of #152 belong to their own PRs. * refactor(solo): drop CellReadStats::merge, which nothing calls Reads are folded in under a mutex, so there are no per-thread partials to merge; the function was reachable only from its own test. CONTRIBUTING.md rules out shipping a function no production path reaches, and the PR description claimed its test as evidence of thread-safety that the mutex actually provides. * feat(solo): --runMode soloCellFiltering Cell-calls an existing raw count matrix without aligning anything, taking the raw directory and an output prefix as STAR does: `--runMode soloCellFiltering /path/to/raw/ /path/to/out/prefix`. Cell calling is a decision about a matrix, not about reads. Re-calling with different `--soloCellFilter` parameters should not mean re-aligning 400 million reads, and a matrix produced by another tool should be callable too. The matrix is streamed into the same temp-body form the align path builds, so `called_cells` and `emptydrops_called` are the identical code here and there rather than a second implementation free to drift. Counts are rounded on the way in: a multimapper matrix carries real values, and the filters work on UMI totals. `--runMode` becomes a token list, because that is what STAR's is: the mode followed by its arguments. The mode itself is now validated rather than falling back to `alignReads`, so a typo is refused instead of quietly running something else. The standalone `emptydrops` binary still exists and still carries its own copy of the algorithm, which no longer matches this one. Removing it means moving `test/solo_genefull_compare.py` and `test/solo_genefull_h5_compare.py` to the new mode first, so it is left alone here rather than broken. --------- Co-authored-by: Benjamin Demaille <benjamin.demaille@icloud.com> 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
2 participants
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.
--soloUMIfiltering MultiGeneUMI_Alldoes what it is documented to do. Closes #144.What changed
UmiFiltering::MultiGeneUmiAllbecomes its own variant instead of an alias, and removes a UMI seen in more than one gene from all of them.MultiGeneUMIflag, so selecting it leaves the filter offMultiGeneUMI: keep the deepest-read geneWhy
Reproducing the no-op ships a flag that silently does nothing to anyone who read STAR's documentation. #144 records the choice, and it was raised there before any code changed, because "be faithful to STAR" and "do what the flag says" genuinely point in opposite directions here.
Divergence
DIVERGENCE.md§1.2, with the STAR source checked (SoloFeature_collapseUMIall.cpp,ParametersSolo.cpp). Needs maintainer sign-off per CONTRIBUTING. Inverting it is a one-line change, since the test asserts the behaviour rather than the direction.Verification
multigene_umi_all_drops_the_umi_from_every_gene— a UMI seen in two genes disappears from both, whereMultiGeneUMIkeeps the deeper onemultigene_umi_all_parses_to_its_own_variant— it is no longer an aliasGate: 562 lib + 26 integration tests,
cargo clippy --all-targets -- -D warnings,cargo fmt --check, MSRV 1.89 — all green.This is a counting change on the
MultiGeneUMI_Allpath only. The default (-) and the other modes produce identical output, so no existing run changes.Not run:
test/solo_diff_docker.sh. Since the default is unaffected, the harness would compare an unchanged path; the changed path has no STAR behaviour to compare against, because in STAR it does nothing. Sign-off on the divergence is what this needs rather than an oracle run.This PR was previously bundled with four other solo themes; those are now #156 (EmptyDrops_CR statistics), #157 (CellReads.stats), #158 (Transcript3p) and #159 (soloCellFiltering), per the one-theme rule in CONTRIBUTING.md.