integrate hyalite v0.2.0 for CellRanger4

[Psy-Fer]

--clipAdapterType CellRanger4 now matches STAR exactly, both halves of the clip.

Closes #148. Dependency discussed in #197.

What changed

5' TSO trim (new)

Previously a heuristic: clip the full 30 nt TSO if the read prefix matched it within 3 mismatches. STAR instead runs an overlap alignment of the TSO against the first 91 bases of the read (ClipCR4::opalAlign, via Opal) and applies an acceptance gate.

solo::tso_clip_len_cr4 now replicates that setup exactly — ACGTN alphabet; match +1, mismatch −2, any-vs-N −2, N-vs-N 0; gap open = gap extend = 2; overlap mode with end tracking; target N-padded to 91 (ClipCR4::opalFillOneSeq) — then STAR's gate S<20 || (S==20 && L>26) || (S==21 && L>30) (ClipMate_clipChunk.cpp:44-47) and min(clip, Lread) (ClipMate_clip.cpp:53).

3' poly-A trim (corrected)

Was "trim a trailing run of A ≥ 8". STAR's ClipCR4::polyTail3p is a scored, mismatch-tolerant walk: +1 per A, −2 per non-A, accept while score*10 >= ib*7, give up at ib-score > 27, reject unless the kept score reached 20, and never trim a read under 20 bases.

The old rule was wrong in both directions:

read	old	STAR
ends in 10 A's	trims 10	trims 0
A*15 + C + A*15	trims 15	trims 31

The Rust transcription is Benjamin Demaille's from #148 and is correct — I re-derived all four of its hand-written expectations against STAR's compiled polyTail3p, including the non-obvious 34 for ACGT*5 + A*30, and they all hold. It is credited in the doc comment.

I did not merge #148 itself, because it wires the trim into clip::clip_mate, which the solo path never calls (solo goes through params.clip5p(0) + clip_read + solo::clip_adapter_cr4) — it would have been dead code on the only path where CellRanger4 is reachable. #148 also logs a warning that the TSO trim is unimplemented, and its CHANGELOG says CellRanger4 is "rejected at parse time"; neither is true (params/mod.rs:1716 accepts it, and #148's own code applies the trim). Nothing in #148 is left unlanded.

New dependency: hyalite 0.2 (#197)

The overlap alignment comes from hyalite 0.2 — pure Rust, zero dependencies of its own, so it adds one node to the tree and no C toolchain requirement.

Reads are scanned a batch at a time through Database::scan_all, mirroring STAR, which fills a 64-read Opal "database" per chunk and aligns the adapter against all of them in one call (ClipMate::clipChunk). solo::tso_clip_lens_cr4_batch does the same over rayon chunks; align_reads_solo resolves the whole batch up front and indexes the result per read.

This is a change to a hot path, so: measured over 100k reads including the per-chunk Database build, batched scan_all is 7.5× faster on AVX2 (4.7× SSE4.1) than per-read align_pair, bit-identical. Throughput is flat from chunk 128 to 10k, so the chunk size is picked to give rayon enough work units, not to tune SIMD.

Validation

Both halves are gated against STAR's own C++, not hand-written expectations. tests/data/cr4_opal_oracle.cpp links STAR's opal.cpp and reproduces ClipCR4 + ClipMate::clipChunk verbatim; the committed cr4_opal_oracle.tsv is its output over 938 reads straddling both decision boundaries — the 5' gate (TSO with 0–14 mismatches, shifts, truncations, indels, embedded Ns, short reads) and polyTail3p's own (score-20 floor, 70% density, ib-score>27 give-up, seqLen<20, A-rich sequence upstream of the tail). The regeneration recipe is in the .cpp header; CI needs no Opal build.

cr4_tso_matches_star_opal_oracle checks all 938 against the scalar path, a forced-SIMD path, and the production batch path. It previously read RUSTAR_CR4_ORACLE and silently returned when unset, so it had been passing vacuously; it now defaults to the committed fixture, with the env var as an override.

Differential vs STARsolo — 10x mouse chr19, both tools with --clipAdapterType CellRanger4:

	before	after
CR4 clip-amount diffs	154	0
soft-clipped share	39% (STAR 41%)	34.9% (STAR 35.0%)

5' leading-clip diffs are 0. The 15 remaining 3'-tail diffs are all STAR extending a short terminal segment across a junction or indel (9S→241N9M, 12S→2I10M) — annotated-junction extension, a separate known class, not a trim disagreement. Solo Gene matrix: 383/384 shared entries, 382 with identical counts, 442 vs 441 total UMIs.

Default SE/PE unchanged — CellRanger4 is opt-in and nothing else moved:

• SE 10k yeast: 8788/8926 (98.5%)
• PE 10k yeast: 8390 both-mapped, 0 half-mapped

Both exactly the recorded baselines.

Tests

New/changed, all present in the diff and passing:

• cr4_tso_matches_star_opal_oracle — the 938-read STAR differential, both columns
• cr4_tso_scan_all_batch_matches_scalar_and_is_backend_stable — batch vs scalar, every available backend, plus the production entry point
• cr4_tso_database_is_deterministic_across_backends_and_matches_scalar
• cr4_polya_trim_follows_star_not_a_run_length_rule — the poly-A cases above, and that the scan runs on the post-5'-clip read
• cr4_tso_clip_removes_exact_tso_prefix, cr4_tso_clip_is_noop_without_adapter, cr4_tso_clip_handles_short_and_empty_reads
• test_starsolo_cr4_tso_clip_is_per_read_within_batch — end-to-end, interleaving TSO-bearing and TSO-free reads so a misindexed batch scan lands a clip on the wrong read

594 tests pass, 0 clippy warnings, cargo fmt --check clean.

Known limitation (pre-existing, not from this PR)

CellRanger4 combined with a non-zero --clip5pNbases/--clip3pNbases misplaces most reads by clip5pNbases: on the same dataset the rustar-minus-STAR POS delta is −5 on 8527 reads, +5 on 3546, and 0 on only 552. It is specific to the pairing — the fixed clips alone are fine (1 leading-clip diff, ~980 other, i.e. the usual baseline), and CellRanger4 alone is fine (0 clip diffs, above). Mapped counts stay close (14747 vs 14687), so reads are misplaced, not lost.

This predates the PR and is not made worse by it. It is filed separately, and it matters for #176, which would make CellRanger4 a default on 10x geometry — at which point users would hit this without ever opting into CR4. It should be fixed before #176 merges. Still not sure what the cause is, so will need to investigate.

One dead end worth recording: STAR does run polyTail3p after the fixed trims, so our scan window is formally wrong for that pairing. I tried narrowing it to match and agreement got worse (475 → 549 leading-clip diffs), because the position bug swamps it. Reverted; there is an NB comment in the code.

[Psy-Fer]

It would be good to get your feedback on this one @BenjaminDEMAILLE

Cheers,
James

[BenjaminDEMAILLE]

Reviewed the full diff and checked both halves against STAR's own sources (ClipCR4.cpp, ClipMate_clipChunk.cpp, ClipMate_clip.cpp, ParametersClip_initialize.cpp). This is a solid piece of work: the oracle-driven approach is the right way to land a clip rule, and catching that cr4_tso_matches_star_opal_oracle had been passing vacuously is a good save on its own.

What I verified independently

• poly_tail_3p is an exact transcription of ClipCR4.cpp:82-107, including the ib1 = seqLen-1 initialisation and the unsigned semantics of ib - score > 27. Correct.
• Scoring matrix, gapOpen = gapExt = 2, alphabetLength = 5, readLen = 91, and the N-padding all match ClipCR4() and opalFillOneSeq exactly.
• The 5' gate matches ClipMate_clipChunk.cpp:44-47, and min(clip, Lread) matches ClipMate_clip.cpp:51.
• Determinism holds: par_chunks + flat_map_iter + collect is order-preserving, the chunking is independent of --runThreadN, and there is no cross-read dependence (each target is padded to 91 in isolation), so batch composition cannot move a result.
• Database::scan_all does out.clear(); out.reserve(sequences.len()) and pushes one hit per sequence, so the chunk.iter().zip(&hits) is safe by construction rather than by luck.

Things I would like changed before merge

1. see TODO.md points at a gitignored file. src/solo/mod.rs:608 refers readers to TODO.md, which .gitignore:40 excludes. No contributor can follow that reference. Please point at the tracking issue number instead.

2. DIVERGENCE.md is not updated. CONTRIBUTING.md:54 requires an entry for any deliberate divergence, and there are two live here: the CellRanger4 + clip5pNbases position bug, and CellRanger4 being solo-only (see 3).

3. CellRanger4 is a silent no-op outside --soloType. STAR arms clip types 10/11 in ParametersClip_initialize.cpp:16-31 unconditionally and applies them in ClipMate::clip for every run, solo or not. rustar only branches on cr4_clip inside align_reads_solo (src/lib.rs:2029), so a plain alignReads --clipAdapterType CellRanger4 clips nothing and says nothing. That is pre-existing, but this PR is where it becomes a deliberate choice (#148's clip_mate wiring dropped, and #148's warning removed with it). CONTRIBUTING.md:64 is explicit that a flag which parses but does nothing must not be quietly accepted. A warning at parse time when CellRanger4 is set without --soloType would close it.

4. The NB comment about the 5' window is not right, and I think it points away from the real bug. Working from the sources: STAR does the fixed N clip first and then clippedAdN = min(clippedInfo, Lread) (ClipMate_clip.cpp:14-24 then 50-54), while rustar does the TSO clip first and the fixed clip second. Both remove N + L bases from the front, and the min only binds once the read is exhausted, in which case both paths land on an empty read. So for the 5' end the two orders are equivalent, and only the 3' window genuinely differs (STAR runs polyTail3p after clip3pNbases).

That is consistent with what you observed: narrowing the window made agreement worse because the window was never the cause. My suggestion for the follow-up investigation is the soft-clip bookkeeping instead, specifically total_clip5p = cr4_5p + clip5p where clip_read has already shifted the sequence, since a double count there would produce exactly a ±clip5pNbases POS delta while leaving mapped counts almost unchanged.

5. hyalite = "0.2" as a caret range. rustar-aligner is published to crates.io, so downstream consumers resolve their own 0.2.x without this repository's Cargo.lock. On a crate that is a week old with a single maintainer, an unvetted patch release can move alignment output. I would pin =0.2.0 (or ~0.2.0) until it has some track record. More on the dependency itself in #197.

Nits

• synthetic_cr4_reads() returns 500 reads, just under CR4_SCAN_CHUNK = 512, so the unit-level batch tests never cross a chunk boundary. The 938-read oracle does cross it, so the coverage exists, but only incidentally. Bumping the synthetic corpus above 512 would make it deliberate.
• A debug_assert_eq!(hits.len(), chunk.len()) in tso_clip_lens_cr4_batch would document the scan_all contract the zip relies on.
• The PR is currently in a conflicting state and needs a rebase (CHANGELOG).