perf(align): memoize the junction-position scan inside stitchWindowAligns - #256
Open
BenjaminDEMAILLE wants to merge 1 commit into
Open
perf(align): memoize the junction-position scan inside stitchWindowAligns#256BenjaminDEMAILLE wants to merge 1 commit into
BenjaminDEMAILLE wants to merge 1 commit into
Conversation
…igns
`find_best_junction_position` is the single hottest function in the
aligner: on a 50k-pair yeast PE run it accounts for ~60% of alignment
time. The scan itself is not wasteful, but it is repeated.
`stitchWindowAligns`' include/exclude recursion reaches the same
(exon A end, seed B) pair through many different branch paths, and each
path re-runs the identical scan. The scan is a pure function of its
arguments, and within one window `read_seq`, the genome, `is_reverse`
and `n_genome` are all fixed, so six coordinates identify a scan
completely: the exon A read end and genome end, the read and genome
gaps, the previous exon length and the next seed length.
Add `JunctionScanCache`, a per-window `FxHashMap` on that key, and a
`find_best_junction_position_cached` wrapper that consults it. The
uncached function is untouched, so a hit returns exactly what a fresh
scan would have. The cache is created per window in `stitch_seeds_core`
and threaded down the recursion, which is what keeps the key complete:
it never outlives the read, genome and strand it was filled for. An
empty `HashMap` does not allocate, so windows that stitch nothing pay
nothing.
Measured on 50k yeast read pairs (Apple M4 Max, quiet machine,
best-of-6 interleaved rounds, `--outSAMtype None`):
threads before after change
1 19.51s 16.76s -14.1%
8 2.55s 2.18s -14.5%
Output is byte-identical: `Aligned.out.sam` (records and header alike,
modulo the `@PG` CL line naming the binary) and `SJ.out.tab` compare
equal against the pre-change binary on the same input. 592 tests pass,
0 clippy warnings.
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
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.
What
find_best_junction_positionis the single hottest function in the aligner. On a 50k-pair yeast PE run (sample-based profile, single thread) it is ~60% of alignment time, insidestitch_seeds_coreat ~76%.The scan itself is not wasteful, but it is repeated.
stitchWindowAligns' include/exclude recursion reaches the same (exon A end, seed B) pair through many different branch paths, and every path re-runs the identical scan.The scan is a pure function of its arguments. Within one window
read_seq, the genome,is_reverseandn_genomeare fixed, so six coordinates identify a scan completely: exon A read end, exon A genome end, read gap, genome gap, previous exon length, next seed length.This PR adds
JunctionScanCache(a per-windowFxHashMapon that key) and afind_best_junction_position_cachedwrapper. The uncached function is untouched, so a cache hit returns exactly what a fresh scan would have returned. The cache is created per window institch_seeds_coreand threaded down the recursion, which is what keeps the key complete: it never outlives the read, genome and strand it was filled for. An emptyHashMapdoes not allocate, so windows that stitch nothing pay nothing.Effect on the profile
find_best_junction_positioninclusive samples drop from 6,855 to 2,941, a 57% cut in the function's total cost.Benchmark
50k yeast read pairs, Apple M4 Max, quiet machine, best-of-6 interleaved rounds,
--outSAMtype None:Correctness
Output is byte-identical.
Aligned.out.sam(records and header, modulo the@PGCL line naming the binary) andSJ.out.tabcompare equal against the pre-change binary on the same input.--all-targets --release)cargo fmt --checkclean🤖 Generated with Claude Code