Map bedMethyl intervals and strands into genomic space

[SuhasSrinivasan]

Review status: Author-reviewed and ready for ONT review.

Fixes #661.

Summary

• Always emit valid half-open single-base genomic intervals from bedmethyl map-to-genome.
• Compose transcript-relative bedMethyl strands with the GTF transcript strand.
• Add an exact coordinate/strand matrix spanning positive and negative two-exon transcripts, coordinate zero, exon boundaries, both row strands, and unstranded rows.

Severity

Severity: High — scientific correctness and output validity

Rationale: The affected command exits successfully while producing reversed intervals and incorrect genomic strand labels, which can cause downstream interval tools to reject or biologically misinterpret rows.

Root cause

The map-to-genome loop used the input bedMethyl strand to choose the coordinate endpoint, subtracting one for a negative row. Strand is an annotation and must not determine BED interval direction. The loop also changed the contig and coordinate but never composed the transcript-relative row strand with the transcript model's genomic orientation.

Implementation

• Introduce map_bedmethyl_line_to_genome, which converts the transcript start with the existing exon-aware mapper, obtains the stop with checked +1, composes strand, and updates the output contig/interval together.
• Add a crate-private TranscriptModel::strand() accessor, keeping other model details encapsulated.
• Route the CLI loop through the helper and preserve the existing per-record error/skip boundary.
• Add one table-driven unit test with nine exact cases across both transcript orientations, exon boundaries, coordinate zero, +, -, and ..

Preserved behavior

• Transcript parsing, ID/version selection, input iteration, row ordering, headers, numeric count fields, and per-record error handling are unchanged.
• Positive-transcript row strands are unchanged.
• Negative-transcript output remains in existing transcript/input order rather than being genomically resorted.

Non-goals

• No change to transcript-contig prefix matching or output sorting.
• No new validation policy for malformed/non-unit input records.
• No merge, DMR, sampling, schema, or performance change.

Behavior before and after

Case	Before	After	Expected oracle
Negative row on positive transcript at genomic start 1	[1,0), strand -	[1,2), strand -	Single-base interval; positive transcript preserves strand
Positive row at transcript position 0 on negative transcript	[11,12), strand +	[11,12), strand -	Negative transcript flips strand
Negative row at transcript position 1 on negative transcript	[10,9), strand -	[10,11), strand +	Single-base interval and flipped strand
Unstranded row on negative transcript	Valid interval, strand .	Valid interval, strand .	Unstranded remains unstranded

Testing

Test environment

• Revision tested: ab28f5f86a27259db35310308f80329af4736fb6
• Tree tested and worktree state: 788fdf37cb56da054f553393804dd53c33cd9e2b; clean
• Parent revision: 5cecc3fb3a9336068d9e3c68d5c08d678153dd2c
• Toolchain: rustc 1.90.0; cargo 1.90.0
• Platform: macOS 26.6, arm64
• Reference binary: installed modkit 0.6.4
• External tools: bgzip/tabix from HTSlib 1.23.1
• Dependency resolution identity: ignored test-only Cargo.lock SHA-256 49c08c4c51b6f4320726551146d971fa9ef2183d40c3f6c631b2005965e242c0; resolved hts-sys 2.2.0; the lockfile is not in this diff

Test layer	Exact command, fixture, or matrix	Result and evidence
Core: parent-red regression	Synthetic GTF and indexed bedMethyl commands from the linked issue on installed modkit 0.6.4	Exit 0 with [1,0) and [10,9) invalid intervals; both negative-transcript strands remain transcript-relative
Core: focused regression	cargo test --offline --locked -p mod_kit bedmethyl_util::subcommands::tests::map_to_genome_composes_coordinates_and_strands -- --exact --test-threads=1	1 passed, 0 failed; nine coordinate/strand cases asserted
Unit/library tests	cargo test --offline --locked -p mod_kit -- --test-threads=1	96 passed, 3 declared ignored, 0 failed
Core: affected CLI/integration tests	cargo test --offline --locked -p modkit --test test_bedmethyl_util -- --test-threads=1	4 passed, 0 failed
Core: applicable full workspace gate	cargo test --offline --workspace --all-targets -- --test-threads=1 on the frozen exact head	Every active workspace/all-target test passed with 0 failures in the recorded serial gate
Installed-versus-fixed comparison	Three-row public fixture from the linked issue	Positive output changed from SHA-256 faec651a64d16ac971cf19b82321c3809dfa34d28e4f7977e78f3d276916e7ed to 679f14e162ff8809d54d19b701f6d01543fa2400bcfed99316671b82d765170f; negative output changed from cc357b2588e711ecc4d84c2372e48f06abc0ce8aa8b945140b403687d8b3d7e6 to f15c5dd598081bcf68215db175cf4337a98a202db962e019acd96eb67fc74f06
Coordinate/strand matrix	Two exons per positive/negative transcript; transcript positions 0–3; row strands +, -, .	Exact contig/start/stop/strand tuples pass; every stop is start+1
Core: formatting/diff checks	Direct rustfmt --edition 2021 --check on both changed Rust files; git diff --check upstream/master...HEAD; clean-worktree check	Passed; only stable-rustfmt warnings for nightly-only repository settings were emitted; worktree clean

Tests not performed

• The 1-GB direct-RNA slice was not run; the defect has a complete exact coordinate/strand oracle and does not touch BAM processing.
• Performance, RSS, and I/O benchmarks were not run because this is a local correctness-only transformation.
• cargo clippy was not run.

Scientific validation

• Population/eligibility invariant: the same parsed input rows are processed or skipped; only genomic representation changes.
• Count/category conservation invariant: every bedMethyl count field is preserved exactly.
• Coordinate/strand/interval invariant: each mapped base is [g, g+1); genomic strand is transcript strand composed with row strand.
• Determinism invariant: the mapping helper is pure for a transcript model and row; existing row order is preserved.
• Independent oracle or specialist review: coordinates and strand algebra were independently reviewed against GTF exon orientation and exact CLI outputs; no blocker remained.

Output and compatibility

• User-visible change: affected rows now have valid genomic intervals and genomic strand labels.
• Expected output differences: only interval endpoints and negative-transcript strand labels in the affected cases.
• Byte-identical controls: count fields, modification code, coverage, score/color fields, positive-transcript strand labels, and row order remain unchanged.
• CLI/API/schema compatibility: unchanged; the new helper and accessor are crate-private.
• Partial-output or failure semantics: unchanged.

Reviewer guide

1. Review map_to_genome_composes_coordinates_and_strands for the complete nine-case oracle.
2. Review map_bedmethyl_line_to_genome for checked interval construction and three-way strand composition.
3. Review the small crate-private transcript-strand accessor.
4. Confirm the CLI loop still owns the same record parse/skip behavior.
5. Rerun the focused test and the two short public CLI commands from the issue.

Checklist

• The issue contains reproducible observed and expected behavior.
• The change is limited to the linked issue’s approved scope.
• The regression is demonstrably red on the exact parent revision.
• All tests actually performed are listed above with their results.
• Unrun or inapplicable tests are disclosed.
• Scientific counts/statistics and output compatibility are explicitly checked.
• Formatting and diff-hygiene checks pass, or unrelated findings are documented.
• No generated data, private sample identifiers, or unrelated changes are included.