Skip to content

--alignEndsType on pairs: extend a mate's two ends while stitching (100M, AS:182 as STAR) - #368

Open
BenjaminDEMAILLE wants to merge 2 commits into
perf/hotpathsfrom
fix/end-to-end-pairs
Open

BenjaminDEMAILLE wants to merge 2 commits into
perf/hotpathsfrom
fix/end-to-end-pairs

Conversation

@BenjaminDEMAILLE

@BenjaminDEMAILLE BenjaminDEMAILLE commented Nov 3, 2026 •

Copy link
Copy Markdown
Contributor

Part of #354 (finding 3).

Moved from BenjaminDEMAILLE#2. The base, perf/hotpaths, is a copy of #346's branch at aeba641, pushed here so this PR can target it. If #346 gets new commits, the copy needs them too. Once #346 merges, retarget this PR to main and delete perf/hotpaths.

Stacked on #346 (perf/hotpaths, tip aeba641). Only the one commit is this PR. It merges without a conflict with the PR for finding 2, so neither is stacked on the other. Pairs only reach zero differences with both, see the tables.

Problem

--alignEndsType EndToEnd left the inner ends of a pair's mates soft-clipped. With 4 mismatches at the 3' end of each mate:

CIGARs AS
STAR 100M / 100M 182
perf/hotpaths 96M4S / 4S96M 190
This branch 100M / 100M 182

The end-to-end decision was made in finalize_transcript, on one mate's slice, after the pair had been chosen and scored. Its mismatches never counted in the pair's score or limits, and a mate that could not reach its end could still be clipped into a pair. Extend5pOfRead1, Extend5pOfReads12 and Extend3pOfRead1 went through the same code.

What STAR does

Observed output of Homebrew STAR 2.7.11b; source not read.

  • Both orientations of the pair, and mismatches at the outer ends instead, give 100M / 100M, AS:182, nM:8.
  • Extend3pOfRead1 extends only read 1's 3' end: 100M and 4S96M, AS:186, nM:4.
  • A mate that runs past the chromosome start maps to nothing in a pair. With filters at 0.3, only the other mate maps (flag 153, 100M, AS:98).
  • A single-end read past the chromosome start is unmapped "other".

Change

All in src/align/stitch.rs, plus the removed argument at four calls in src/align/read_align.rs.

  • The two extensions at the mate boundary and the two at a transcript's leaf take their flag from ext[iFrag][iEnd]. ext_frag maps a seed's mate_id to that fragment; single-end seeds carry mate_id 2 and use fragment 0.
  • An extension that hits a chromosome edge rejects the branch in the stitcher: return false at the boundary, break 'leaf at the leaf.
  • finalize_transcript loses its imate parameter and no longer forces any end. Its check for that case moved with it.
  • wt_final_score is not floored at 0 any more. An end-to-end extension can take a transcript below zero, and the window's score gate has to see that.

In Local mode every extension flag stays false, as before. Removing the floor does reach Local, for a transcript whose genomic-length penalty takes it below zero; the numbers below show no change there, on these reads.

Measured against STAR

Fields 1-11 plus NH, HI, AS, nM.

Crafted reads (reads with mismatches at the ends, with both mate orders, and reads that cannot reach an end), this branch:

Run Records STAR / here SAM and Log.final.out counts
Audit fixture e1.fq/e2.fq 2 / 2 identical, same order
Pairs, EndToEnd 12 / 12 identical, same order
Pairs, EndToEnd, --outSAMunmapped Within 22 / 22 identical, same order
Pairs, Extend5pOfRead1 18 / 18 identical, same order
Pairs, Extend5pOfReads12 16 / 16 identical, same order
Pairs, Extend3pOfRead1 18 / 18 identical, same order
Pairs, Local 20 / 20 identical, same order
Single-end, EndToEnd 18 / 18 identical, same order
Single-end, Extend5pOfRead1 20 / 20 identical, same order
Single-end, Extend3pOfRead1 20 / 20 identical, same order
Single-end, EndToEnd, no mismatch limit 19 / 19 identical, same order
Pairs, EndToEnd, filters 0.3 15 / 18 counts identical; 3 extra unmapped-mate records, which finding 2 removes (15 / 15 with both)

On perf/hotpaths the audit fixture differs, pairs under EndToEnd have 16 records against 12, and Extend3pOfRead1 has 20 against 18.

Fuzz: 4000 repeat-rich pairs, and their first mates for single-end. Reads differing from STAR:

Options perf/hotpaths This branch This branch merged with finding 2
Pairs, EndToEnd 1026 651 0
Pairs, EndToEnd, filters 0.3 2353 2012 0
Pairs, Extend5pOfRead1 662 533 0
Pairs, Extend5pOfReads12 840 566 0
Pairs, Extend3pOfRead1 769 545 0
Pairs, Local 523 523 0
Single-end, EndToEnd 1 0 0
Single-end, Local 0 0 0

So this branch alone does not bring pairs to zero. The differences left are the half-mapped reads of finding 2, which is why Local, which this PR does not touch, stays at 523. The merge's Log.final.out counts are identical to STAR's for every row but single-end EndToEnd, below.

The merge also passes cargo test --release (831 passed).

Not the same as STAR

On the single-end EndToEnd fuzz, the alignments are identical but 7 of the 335 unmapped junk reads carry another label:

Unmapped reason STAR Here
too short 117 124
other 218 211

I do not know why. I tried treating a negative score as no alignment: that gave 121 / 214, and three crafted reads scoring 0, -1 and -2 said STAR calls those "too short", so the idea was wrong and I removed it. The crafted fixtures and the audit's reads do not show it. On perf/hotpaths, a read past the chromosome start was "too short" where STAR says "other"; that case is fixed and tested.

Deliberate divergences

None, so no DIVERGENCE.md entry and no sign-off needed.

Tests

tests/end_to_end_pairs.rs, all five failing on the old src/, with expected values from STAR's output on the same fixture:

  • inner_ends_are_extended_and_charged
  • outer_ends_are_charged_to_the_pair
  • a_mate_past_the_chromosome_start_cannot_be_clipped_into_a_pair
  • one_sided_modes_extend_only_their_inner_end
  • single_end_read_past_the_chromosome_start_has_no_alignment

cargo test --release (824 passed), cargo clippy --all-targets and cargo fmt --check are clean on macOS.

Not checked

  • test/parity_suite.sh was not run.
  • No real dataset: all numbers are from synthetic reads.
  • Speed was not measured. EndToEnd now extends during stitching, on every leaf, rather than once at the end. Local mode takes the same path with the flag off.
  • --peOverlapNbasesMin was not run together with --alignEndsType.
  • Chimeric detection was not run together with --alignEndsType.

Second commit, after review

An independent review fuzzed 4000 pairs against STAR 2.7.11b under 17 option sets. With this PR and #367 merged, every set had 0 differing reads except three EndToEnd sets with 1 or 3. The single read came from the order of two checks in the end-to-end extension: STAR's extendAlign.cpp:24-31 checks the chromosome edge before the mate spacer, this code the other way round, so a mate whose 3' end sits on a chromosome edge was extended into a pair. f7da8ee swaps them.

Mate 1 = last 100 bases of the chromosome STAR before f7da8ee
mate 2 its reverse complement, EndToEnd or Extend3pOfRead1 153 100M AS:98 pair, AS:198 153 100M AS:98
mate 2 random, EndToEnd unmapped 73 100M AS:98 unmapped
mate 2 random, Extend3pOfRead1 unmapped pair, 73S6M21S unmapped

Filters at 0.3 in all rows. Locked by a_mate_ending_on_the_chromosome_end_cannot_be_extended_into_a_pair, which fails without the change. Without #367 this branch also writes the unmapped mate's record (flag 101), which STAR omits; #367 handles that.

The same review timed 200k pairs at --runThreadN 1: 12.0-12.2 s of user CPU on the base and on this branch, Local and EndToEnd alike, so no measurable cost. Remaining: EndToEnd with --peOverlapNbasesMin 10 leaves 2 reads with mates at a chromosome end unmapped where STAR maps them with 2S/1S.

🤖 Generated with Claude Code

Both ends of a mate are now extended in the stitcher, with ext[iFrag][iEnd]
applied at the mate boundary and at the transcript's outer ends. The
extension's mismatches count in the pair's score and limits. finalize_transcript
worked on one mate's slice after the pair was already chosen, so inner ends
stayed clipped and a mate that cannot reach its end could be clipped into a pair.

Transcript scores are not floored at 0 before the window's score gate, so a
negative end-to-end transcript loses to an existing one.

Part of #354 (finding 3).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
STAR's extendAlign checks the chromosome boundary first and the fragment
spacer second, so a mate whose 3' end sits on a chromosome edge is killed
under EndToEnd / Extend3pOfRead1 instead of stopping at the spacer. With
mate 1 the last 100 bases of the chromosome and mate 2 its reverse
complement, STAR reports mate 2 alone (153 100M AS:98); this branch reported
a pair with AS:198. It now matches STAR, and with an unrelated mate 2 both
leave the read unmapped.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>

This branch has not been deployed

No deployments
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant