feat(simulate): body error injection, template-coordinate sort fix, correctness guards (SIMU3-01/03/04/05/07)

[nh13]

Summary

W7 of the final-audit burn-down: the simulate error-model / fixture-generator cluster. simulate is fgumi-specific (no fgbio oracle), so verification is against the generator's own advertised contract — template-coordinate sorted output, two-stranded duplex families, and a truth file that agrees with fgumi correct. Stacked on #531 (nh/simulate-ss-parity), the only open PR touching mapped_reads.rs/grouped_reads.rs; retargets to main on merge.

Findings

SIMU3-01 (S2) — no body sequencing errors

The generators emitted error-free reads (only the fastq UMI prefix was perturbed), so consensus/error-correction paths never saw real discordances. Added a --error-rate flag (default 0.0) to mapped-reads, grouped-reads, and fastq-reads that injects per-base substitution errors into read bodies via a shared introduce_errors_inplace (moved to common.rs). The error_rate > 0.0 guard is important: introduce_errors_inplace draws one RNG value per base regardless of rate, so guarding keeps the default byte-identical to the error-free generator. Per the maintainer's call, this ships opt-in/off-by-default so the ~12 E2E regression baselines are unchanged; a follow-up can opt specific fixtures into a positive rate.

SIMU3-05 (S1) — sort-order break near contig ends (mapped-reads)

The template-coordinate sort-key pre-pass keyed on the pre-sampled locus, but generation re-samples a new locus when the pre-sampled one is too close to a contig end for the drawn insert size — so records were emitted out of the advertised SS:template-coordinate order. A shared effective_molecule_locus replays the molecule RNG (UMI → family size → insert size → strand coin → sequence_at/sample_sequence fallback) so the key uses the emitted coordinates. For loci that don't trigger the fallback the key is unchanged.

Scope: the audit scoped SIMU3-05 to mapped-reads, and this fixes it there. grouped-reads has a separate, pre-existing template-coordinate ordering gap (present even without any fallback), so it is deliberately left for a follow-up rather than conflated here — the grouped sort key is unchanged.

SIMU3-07 (data quality) — duplex families not two-stranded

grouped-reads --duplex used split_reads, which can assign zero reads to a strand. Switched to the existing split_reads_with_minimum(family_size, 1) (identical RNG consumption — a single sample_a_fraction draw), so a duplex family with ≥ 2 reads is genuinely two-stranded. family_size == 1 is inherently single-stranded and falls back.

SIMU3-03 (S3) — truth file unverified vs nearest-neighbor collisions

correct-reads asserted edit1/edit2 → true UMI without checking whether the observed UMI is actually nearest to a different includelist entry. expected_correction_for now computes the unique nearest includelist neighbor within the correctable radius (mirroring fgumi correct with max-mismatches 2 / min-distance 1), reporting collisions/ambiguous cases as uncorrectable so the truth agrees with correct.

SIMU3-04 (S3) — panic for short UMIs

rng.random_range(3..=umi_length.min(5)) inverts (and panics) when umi_length < 3. Guarded so the "multi" branch clamps the upper bound instead.

Verification (generator's own contract, via the built binary)

• SIMU3-05: mapped-reads output on a fallback-heavy reference (12×1100bp contigs, insert mean 400/max 800) is byte-for-byte samtools sort --template-coordinate. RED: reverting the key to the pre-sampled locus left 4238/11988 records out of order.
• SIMU3-01: --error-rate 0.05 changes 1191 sequences vs 0.0, with positions/flags/names identical (errors don't move alignments).
• SIMU3-07: grouped-reads --duplex → 290/300 families two-stranded (the 10 singletons are family_size == 1).
• SIMU3-04: correct-reads --umi-length 2 exits 0 (no panic).

Tests

• common.rs: effective_molecule_locus passthrough + contig-end fallback.
• correct_reads.rs: expected_correction_for (unique-nearest, collision, beyond-radius, exact) + a short-UMI no-panic test.
• Full suite green: cargo ci-fmt && ci-lint && ci-test (2259 default) and cargo nextest run --features simulate (2150).

Follow-ups

• grouped-reads has a separate pre-existing template-coordinate ordering gap (out of SIMU3-05's scope).
• Opting specific E2E fixtures into a non-zero --error-rate.

(This commit was pushed unsigned because 1Password signing was unavailable; it will be re-signed and force-pushed once signing is back.)

Summary by CodeRabbit

• New Features
  • Added --error-rate to fastq-reads, grouped-reads, and mapped-reads.
  • Injects per-base substitution errors into generated read bodies (beyond the UMI), while preserving UMI prefixes and read lengths.
  • Default --error-rate=0.0 keeps output byte-identical to previous versions.
• Bug Fixes
  • Improved UMI correction to use nearest-neighbor matching with a fixed radius, avoiding incorrect corrections.
  • Prevented panics for very short UMIs in multi-error generation.
• Other
  • Improved template-coordinate ordering consistency at contig boundaries and enhanced test coverage.

[coderabbitai]

Warning

Review limit reached

You’ve reached a temporary PR review limit under our Fair Usage Limits Policy.

Your recent review volume is higher than typical usage, so adaptive limits are currently applied.

Next review available in: 15 minutes

Enable usage-based reviews in Billing to review now. Otherwise, wait until the next included review is available.
You're only billed for reviews past your plan's rate limits ($0.25/file).

How can I continue?

After more reviews become available, a review can be triggered using the @coderabbitai review command as a PR comment. Alternatively, push new commits to this PR.

To avoid repeated limits, reduce automatic review volume by pausing incremental auto-reviews earlier, using label-based review opt-in, excluding WIP or generated PR titles, or requesting reviews manually when the PR is ready. If your team needs uninterrupted high-volume reviews, an organization admin can enable usage-based reviews.

How do review limits work?

CodeRabbit enforces per-developer PR review limits for each organization. Most developers receive the normal plan review availability.

For paid Pro and Pro+ PR reviews, CodeRabbit uses adaptive limits for sustained high-volume activity. When a developer's recent PR review activity reaches the 95th percentile or higher among CodeRabbit users, additional reviews become available more gradually as earlier reviews age out of the rolling window.

Please refer docs for additional details.

Review details

⚙️ Run configuration

Configuration used: Path: .coderabbit.yaml

Review profile: ASSERTIVE

Plan: Pro

Run ID: 8da63517-ab22-4a61-96ee-f0e443215d59

📥 Commits

Reviewing files that changed from the base of the PR and between d780822 and 9225518.

📒 Files selected for processing (5)

• src/lib/commands/simulate/common.rs
• src/lib/commands/simulate/correct_reads.rs
• src/lib/commands/simulate/fastq_reads.rs
• src/lib/commands/simulate/grouped_reads.rs
• src/lib/commands/simulate/mapped_reads.rs

Walkthrough

Simulation commands add configurable read-body substitution errors, preserve UMI prefixes, align mapped sort keys with effective emitted loci, constrain duplex partitioning, and improve short-UMI correction modeling.

Changes

Simulation output generation

Layer / File(s)	Summary
RNG utilities and effective loci src/lib/commands/simulate/common.rs	Adds in-place substitution and deterministic effective-locus replay, with tests for valid loci and contig-end fallback.
FASTQ error-rate propagation src/lib/commands/simulate/fastq_reads.rs	Adds --error-rate, validates and propagates it, and injects substitutions into R1/R2 bodies while preserving UMI prefixes and read lengths.
Grouped and mapped read integration src/lib/commands/simulate/grouped_reads.rs, src/lib/commands/simulate/mapped_reads.rs	Propagates body errors through grouped and mapped generation, uses effective loci for mapped sort keys, constrains duplex read splitting, and adds regression coverage.

UMI correction simulation

Layer / File(s)	Summary
Correction truth and short-UMI safety src/lib/commands/simulate/correct_reads.rs	Clamps multi-error sampling for short UMIs and computes expected correction from unique nearest includelist matches, including ambiguity and radius tests.

Estimated code review effort: 4 (Complex) | ~45 minutes

Sequence Diagram(s)

sequenceDiagram
  participant CLI
  participant GenerationParams
  participant ReadGenerator
  participant introduce_errors_inplace
  participant FASTQOutput
  CLI->>GenerationParams: validate and store error_rate
  GenerationParams->>ReadGenerator: generate molecule reads
  ReadGenerator->>introduce_errors_inplace: mutate read bodies when error_rate > 0
  introduce_errors_inplace-->>ReadGenerator: substituted sequences
  ReadGenerator-->>FASTQOutput: fixed-length R1/R2 reads

Loading

🚥 Pre-merge checks | ✅ 5

✅ Passed checks (5 passed)

Check name	Status	Explanation
Description Check	✅ Passed	Check skipped - CodeRabbit’s high-level summary is enabled.
Title check	✅ Passed	The title matches the main changes: body error injection, sorting fix, and correctness guards.
Docstring Coverage	✅ Passed	Docstring coverage is 100.00% which is sufficient. The required threshold is 80.00%.
Linked Issues check	✅ Passed	Check skipped because no linked issues were found for this pull request.
Out of Scope Changes check	✅ Passed	Check skipped because no linked issues were found for this pull request.

✨ Finishing Touches

🧪 Generate unit tests (beta)

• Create PR with unit tests
• Commit unit tests in branch nh/fix-simulate-error-model

---

_{Comment @coderabbitai help to get the list of available commands.}

[nh13]

Both reviewers run (§0 step 6).

CodeRabbit CLI (--base nh/simulate-ss-parity) — 7 findings:

• (minor ×3) validate_rate was called with "--error-rate"; the established convention passes the bare flag name ("unmapped-fraction", "conversion-rate"). Fixed all three call sites to "error-rate" so the error message is well-formed.
• (trivial) collapse the four expected_correction_for tests into a parameterized #[rstest] case table (repo convention). Done.
• (trivial) add a regression test comparing effective_molecule_locus to actual generate_molecule_reads. Done — test_effective_molecule_locus_matches_generation asserts the helper returns exactly the coordinates generation emits, across seeds and loci including a contig-end locus that triggers the fallback. This directly guards the RNG-replay correctness that the SIMU3-05 fix depends on.
• (trivial) add a fastq body-error-injection unit test. Skipped — already covered end-to-end (rate 0.05 changes bases but not positions/flags); the guard is a one-line if error_rate > 0.0.
• (major, out of this diff) grouped_reads.rs/mapped_reads.rs sort-order header uses HeaderTag::from(...) pattern matches instead of noodles' SORT_ORDER/GROUP_ORDER/SUBSORT_ORDER constants (lines 159-172). This is base fix(simulate): write SS sub-sort tag via SortOrder accessors (R2-HDR-01) #531 code, not touched by this commit → noted for fix(simulate): write SS sub-sort tag via SortOrder accessors (R2-HDR-01) #531.

Local CR-style review — no actionable findings; one nitpick (the correctable radius/min-distance in expected_correction_for are fixed at 2/1, documented and matching the generator's edit distribution).

cargo ci-fmt && ci-lint && ci-test green (2259 default); cargo nextest run --features simulate green (2151).

Commit is unsigned (1Password unavailable); will be re-signed + force-pushed once signing returns.

[nh13]

Follow-up on the grouped-reads note in the PR body: I investigated it and it's not a bug — retracting the "separate pre-existing gap" framing.

The generated grouped-reads output has 0 descending steps in the template primary key (verified both without any fallback and on a fallback-heavy reference, 12,066 records each) — i.e., it is validly template-coordinate sorted. The earlier "68/3059 records differ from samtools sort --template-coordinate" observation is benign tie-break noise among equal-key records: fgumi's own sort --order template-coordinate also differs from samtools by ~4444 for both mapped- and grouped-reads, so fgumi and samtools simply order equal keys differently — neither is wrong.

For contrast, the mapped-reads SIMU3-05 bug state (pre-sampled sort key) had 36 genuine descending steps → 0 after the fix, which is why mapped is a real S1 and grouped is not. No grouped change is needed; scoping SIMU3-05 to mapped-reads was correct, and there is no residual follow-up here.

[codecov]

Codecov Report

❌ Patch coverage is 99.67532% with 1 line in your changes missing coverage. Please review.
✅ Project coverage is 92.80%. Comparing base (e9a1ae1) to head (9225518).
⚠️ Report is 1 commits behind head on main.

Files with missing lines	Patch %	Lines
src/lib/commands/simulate/common.rs	98.18%	1 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main     #541      +/-   ##
==========================================
+ Coverage   92.73%   92.80%   +0.07%     
==========================================
  Files         166      166              
  Lines      101603   101888     +285     
==========================================
+ Hits        94218    94557     +339     
+ Misses       7385     7331      -54     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[nh13]

Overlap with #576 (simulate template-coordinate sort)

#576 takes a different, broader approach to the sort ordering: it deletes simulate's bespoke key (simulate/sort.rs, TemplateCoordKey::for_f1r2_pair) and instead emits records in molecule order to a temp BAM, then sorts into the output with the canonical fgumi-sort engine (RawExternalSorter / SortOrder::TemplateCoordinate) — the same one fgumi sort/fgumi group use.

Because it sorts the emitted records, #576:

• subsumes this PR's SIMU3-05 fix (the contig-end resampling case is handled for free — the key can't disagree with the emitted coordinates because there is no separate key), and
• also closes the grouped-reads ordering gap this PR explicitly deferred — both mapped-reads and grouped-reads are correct by construction, for every read orientation. Verified: the F2R1 (reverse-R1) case that was mis-ordered goes 44/5994 → 0 vs samtools sort --template-coordinate, and output is byte-identical to fgumi sort --order template-coordinate.

Proposal: let #576 own the simulate sort ordering, and rebase this PR to drop the SIMU3-05 sort portion — keeping its genuinely separate improvements: body error injection (SIMU3-01), two-stranded duplex families (SIMU3-07), truth-file verification (SIMU3-03), and the short-UMI panic guard (SIMU3-04). They don't overlap #576 and are worth keeping.

(These edit the same mapped_reads.rs/grouped_reads.rs lines, so whichever merges first, the other will need a rebase.)

[nh13]

@coderabbitai review

[coderabbitai]

✅ Action performed

Review finished.

Note: CodeRabbit is an incremental review system and does not re-review already reviewed commits. This command is applicable only when automatic reviews are paused.

[coderabbitai]

Actionable comments posted: 2

🤖 Prompt for all review comments with AI agents

Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

Inline comments:
In `@src/lib/commands/simulate/fastq_reads.rs`:
- Around line 539-555: Isolate body-error randomness from molecule-generation
randomness by using deterministic, independent per-read/per-mate RNG streams for
introduce_errors_inplace in src/lib/commands/simulate/fastq_reads.rs:539-555,
and apply the same change to grouped_reads.rs:559-581 and
mapped_reads.rs:473-497 so qualities and later family generation remain
unchanged. Extend the assertions in fastq_reads.rs:1074-1147 to preserve UMI
prefixes, strand/truth fields, and qualities; in grouped_reads.rs:1020-1065
preserve non-sequence fields and qualities; and in mapped_reads.rs:1089-1145
preserve positions, flags, and qualities.

In `@src/lib/commands/simulate/grouped_reads.rs`:
- Around line 244-249: Update grouped-reads ordering to account for the fallback
locus generated later in the molecule flow: re-key fallback molecules using
effective_molecule_locus before sorting, matching the mapped-reads behavior, or
pass emitted records through the canonical sorter. Remove or revise the nearby
template-coordinate ordering note so it no longer advertises behavior the
implementation does not guarantee.

🪄 Autofix (Beta)

Fix all unresolved CodeRabbit comments on this PR:

• Push a commit to this branch (recommended)
• Create a new PR with the fixes

---

ℹ️ Review info

⚙️ Run configuration

Configuration used: Path: .coderabbit.yaml

Review profile: ASSERTIVE

Plan: Pro

Run ID: 8102d586-abcc-495b-9cd2-c6d06313e732

📥 Commits

Reviewing files that changed from the base of the PR and between 9676dc1 and 5784271.

📒 Files selected for processing (5)

• src/lib/commands/simulate/common.rs
• src/lib/commands/simulate/correct_reads.rs
• src/lib/commands/simulate/fastq_reads.rs
• src/lib/commands/simulate/grouped_reads.rs
• src/lib/commands/simulate/mapped_reads.rs

[coderabbitai]

Actionable comments posted: 1

🤖 Prompt for all review comments with AI agents

Verify each finding against current code. Fix only still-valid issues, skip the
rest with a brief reason, keep changes minimal, and validate.

Inline comments:
In `@src/lib/commands/simulate/grouped_reads.rs`:
- Around line 415-421: Add a focused test for the duplex grouping path using
family_size >= 2, preferably exactly 2, and iterate across multiple RNG seeds;
assert the resulting split from split_reads_with_minimum (or the emitted
MI-tagged pairs) always contains at least one A read and one B read. Keep the
existing orientation-flip test unchanged and target the SIMU3-07
minimum-per-strand guarantee directly.

🪄 Autofix (Beta)

Fix all unresolved CodeRabbit comments on this PR:

• Push a commit to this branch (recommended)
• Create a new PR with the fixes

---

ℹ️ Review info

⚙️ Run configuration

Configuration used: Path: .coderabbit.yaml

Review profile: ASSERTIVE

Plan: Pro

Run ID: 9b16cd0e-a0d6-4131-8c60-e2fe510d5585

📥 Commits

Reviewing files that changed from the base of the PR and between 5784271 and d780822.

📒 Files selected for processing (5)

• src/lib/commands/simulate/common.rs
• src/lib/commands/simulate/correct_reads.rs
• src/lib/commands/simulate/fastq_reads.rs
• src/lib/commands/simulate/grouped_reads.rs
• src/lib/commands/simulate/mapped_reads.rs