From 1a54768c6e9388aae228a223772141a0d75d9e85 Mon Sep 17 00:00:00 2001 From: SuhasSrinivasan <32346517+SuhasSrinivasan@users.noreply.github.com> Date: Wed, 5 Aug 2026 16:15:59 -0700 Subject: [PATCH 1/4] test: expose dropped DMR segmentation sites --- modkit-core/src/dmr/single_site.rs | 72 ++++++++++++++++ modkit-core/src/hmm.rs | 131 ++++++++++++++++++++++++++++- modkit/tests/test_dmr.rs | 74 ++++++++++++++++ 3 files changed, 276 insertions(+), 1 deletion(-) diff --git a/modkit-core/src/dmr/single_site.rs b/modkit-core/src/dmr/single_site.rs index 6e0f6e62..4aa1a7c2 100644 --- a/modkit-core/src/dmr/single_site.rs +++ b/modkit-core/src/dmr/single_site.rs @@ -1345,3 +1345,75 @@ fn path_to_region_labels( agg } } + +#[cfg(test)] +mod segmentation_path_tests { + use super::{path_to_region_labels, States}; + + fn assert_every_position_covered_once( + path: &[States], + positions: &[u64], + regions: &[(u64, u64, States)], + ) { + for (&position, &state) in positions.iter().zip(path) { + let covering_regions = regions + .iter() + .filter(|(start, end, _)| *start <= position && position < *end) + .collect::>(); + assert_eq!( + covering_regions.len(), + 1, + "position {position} in state {state:?} was not covered exactly once" + ); + assert_eq!(covering_regions[0].2, state); + } + } + + fn check_case( + path: &[States], + positions: &[u64], + expected: &[(u64, u64, States)], + ) { + let regions = path_to_region_labels(path, positions); + assert_eq!(regions, expected); + assert_every_position_covered_once(path, positions, ®ions); + } + + #[test] + fn singleton_position_becomes_a_single_base_region() { + check_case(&[States::Same], &[10], &[(10, 11, States::Same)]); + } + + #[test] + fn two_positions_are_both_integrated() { + check_case( + &[States::Same, States::Same], + &[10, 20], + &[(10, 21, States::Same)], + ); + check_case( + &[States::Same, States::Different], + &[10, 20], + &[(10, 11, States::Same), (20, 21, States::Different)], + ); + } + + #[test] + fn multi_position_transitions_cover_each_position_once() { + check_case( + &[ + States::Same, + States::Same, + States::Different, + States::Different, + States::Same, + ], + &[10, 20, 30, 40, 50], + &[ + (10, 21, States::Same), + (30, 41, States::Different), + (50, 51, States::Same), + ], + ); + } +} diff --git a/modkit-core/src/hmm.rs b/modkit-core/src/hmm.rs index 5a86669d..f3f6f6fc 100644 --- a/modkit-core/src/hmm.rs +++ b/modkit-core/src/hmm.rs @@ -413,7 +413,106 @@ impl Projection { #[cfg(test)] mod hmm_tests { - use crate::hmm::HmmModel; + use crate::hmm::{HmmModel, States}; + + fn test_model() -> HmmModel { + HmmModel::new(0.1, 0.9, 0.3, -0.1, 0.01, 500, true).unwrap() + } + + fn emission_score(model: &HmmModel, p: f64, state: States) -> f64 { + let p = if p == 0f64 { 1e-5 } else { p }; + let (factor, log_probability) = match state { + States::Same => (model.same_state_factor, p.ln()), + States::Different => { + (model.diff_state_factor, (1f64 - p + 1e-5).ln()) + } + }; + factor * (log_probability - model.significance_factor) + } + + fn transition_score( + model: &HmmModel, + previous: States, + current: States, + diff_stay: f64, + ) -> f64 { + match (previous, current) { + (States::Same, States::Same) => model.same_to_same, + (States::Same, States::Different) => model.same_to_diff, + (States::Different, States::Same) => (1f64 - diff_stay).ln(), + (States::Different, States::Different) => diff_stay.ln(), + } + } + + /// Exhaustively score the hidden start state and every emitted state. + /// This deliberately does not use the dynamic-programming matrix or its + /// back-pointers, so it independently checks both state order and length. + fn brute_force_path( + model: &HmmModel, + scores: &[f64], + positions: &[u64], + ) -> Vec { + assert!(!scores.is_empty()); + assert_eq!(scores.len(), positions.len()); + + let probabilities = scores + .iter() + .map(|&score| (-score.max(0f64)).exp()) + .collect::>(); + let diff_stays = positions.windows(2).fold( + vec![model.dmr_prior], + |mut transitions, window| { + let gap = (window[1] - window[0]) as f64; + transitions.push(if model.linear_proj { + model.projection.linear_project_prob(gap) + } else { + model.projection.ln_project_prob(gap) + }); + transitions + }, + ); + + let state_count = scores.len() + 1; + let mut best: Option<(f64, Vec)> = None; + for encoded_path in 0..(1usize << state_count) { + let states = (0..state_count) + .map(|i| { + if encoded_path & (1usize << i) == 0 { + States::Same + } else { + States::Different + } + }) + .collect::>(); + let initial_score = match states[0] { + States::Same => model.same_to_same, + States::Different => model.same_to_diff, + }; + let total_score = probabilities + .iter() + .zip(diff_stays.iter()) + .enumerate() + .fold(initial_score, |total, (i, (&p, &diff_stay))| { + total + + transition_score( + model, + states[i], + states[i + 1], + diff_stay, + ) + + emission_score(model, p, states[i + 1]) + }); + + if best + .as_ref() + .map(|(best_score, _)| total_score > *best_score) + .unwrap_or(true) + { + best = Some((total_score, states[1..].to_vec())); + } + } + best.unwrap().1 + } #[test] fn test_prob_to_factor() { @@ -421,4 +520,34 @@ mod hmm_tests { let fact = HmmModel::prob_to_factor(sig_fact).unwrap(); dbg!(fact); } + + #[test] + fn viterbi_path_matches_independent_oracle_for_tiny_sequences() { + let model = test_model(); + let cases = [ + (vec![0.0], vec![10]), + (vec![12.0], vec![10]), + (vec![0.0, 12.0], vec![10, 20]), + (vec![12.0, 0.0], vec![10, 20]), + (vec![0.0, 0.0, 12.0, 12.0, 0.0], vec![10, 20, 30, 40, 50]), + ]; + + for (scores, positions) in cases { + let expected = brute_force_path(&model, &scores, &positions); + let actual = model.viterbi_path(&scores, &positions); + assert_eq!(actual.len(), scores.len()); + assert_eq!(actual, expected, "scores: {scores:?}"); + } + } + + #[test] + fn viterbi_path_preserves_state_transitions_in_order() { + let model = test_model(); + let scores = vec![0.0, 0.0, 12.0, 12.0, 0.0]; + let positions = vec![10, 20, 30, 40, 50]; + let expected = brute_force_path(&model, &scores, &positions); + assert!(expected.windows(2).any(|states| states[0] != states[1])); + + assert_eq!(model.viterbi_path(&scores, &positions), expected); + } } diff --git a/modkit/tests/test_dmr.rs b/modkit/tests/test_dmr.rs index 5a2991d3..dbfa61f6 100644 --- a/modkit/tests/test_dmr.rs +++ b/modkit/tests/test_dmr.rs @@ -1,6 +1,8 @@ use crate::common::{ check_against_expected_text_file, check_legal_csv, run_modkit, }; +use std::fs; +use std::path::Path; mod common; @@ -75,6 +77,78 @@ fn test_dmr_regression() { ); } +fn run_segmented_dmr( + output: &Path, + segments: &Path, + threads: &str, + io_threads: &str, +) { + run_modkit(&[ + "dmr", + "pair", + "-a", + "../tests/resources/\ + lung_00733-m_adjacent-normal_5mc-5hmc_chr20_cpg_pileup.bed.gz", + "-b", + "../tests/resources/\ + lung_00733-m_primary-tumour_5mc-5hmc_chr20_cpg_pileup.bed.gz", + "-o", + output.to_str().unwrap(), + "--segment", + segments.to_str().unwrap(), + "--ref", + "../tests/resources/GRCh38_chr20.fa", + "--header", + "--base", + "C", + "--max-coverages", + "100", + "100", + "--threads", + threads, + "--io-threads", + io_threads, + "--suppress-progress", + "--force", + ]) + .expect("segmented DMR run should succeed"); +} + +#[test] +fn segmentation_includes_last_site_and_is_thread_deterministic() { + let temp_dir = tempfile::tempdir().unwrap(); + let sites_one = temp_dir.path().join("sites-threads-1.bed"); + let segments_one = temp_dir.path().join("segments-threads-1.bed"); + let sites_four = temp_dir.path().join("sites-threads-4.bed"); + let segments_four = temp_dir.path().join("segments-threads-4.bed"); + + run_segmented_dmr(&sites_one, &segments_one, "1", "1"); + run_segmented_dmr(&sites_four, &segments_four, "4", "2"); + + let sites_one = fs::read_to_string(sites_one).unwrap(); + let sites_four = fs::read_to_string(sites_four).unwrap(); + let segments_one = fs::read_to_string(segments_one).unwrap(); + let segments_four = fs::read_to_string(segments_four).unwrap(); + assert_eq!(sites_one.as_bytes(), sites_four.as_bytes()); + assert_eq!(segments_one.as_bytes(), segments_four.as_bytes()); + + let site_rows = sites_one.lines().skip(1).collect::>(); + let segment_rows = segments_one.lines().skip(1).collect::>(); + assert_eq!(site_rows.len(), 17_271); + + let segment_site_total = segment_rows + .iter() + .map(|row| row.split('\t').nth(5).unwrap().parse::().unwrap()) + .sum::(); + assert_eq!(segment_site_total, site_rows.len()); + + let last_site_end = site_rows.last().unwrap().split('\t').nth(2).unwrap(); + let last_segment_end = + segment_rows.last().unwrap().split('\t').nth(2).unwrap(); + assert_eq!(last_site_end, "10804378"); + assert_eq!(last_segment_end, last_site_end); +} + // todo // test pair with explicit index // test multi From ba63c29527fde8c25e6f916ebac1fd62e9adb028 Mon Sep 17 00:00:00 2001 From: SuhasSrinivasan <32346517+SuhasSrinivasan@users.noreply.github.com> Date: Wed, 5 Aug 2026 16:20:18 -0700 Subject: [PATCH 2/4] fix: include every site in DMR segmentation --- modkit-core/src/dmr/single_site.rs | 2 +- modkit-core/src/hmm.rs | 21 ++++++++------------- 2 files changed, 9 insertions(+), 14 deletions(-) diff --git a/modkit-core/src/dmr/single_site.rs b/modkit-core/src/dmr/single_site.rs index 4aa1a7c2..fea3d05d 100644 --- a/modkit-core/src/dmr/single_site.rs +++ b/modkit-core/src/dmr/single_site.rs @@ -1319,7 +1319,7 @@ fn path_to_region_labels( path: &[States], positions: &[u64], ) -> Vec<(u64, u64, States)> { - assert_eq!(path.len(), positions.len() - 1); + assert_eq!(path.len(), positions.len()); if path.is_empty() { return Vec::new(); } else { diff --git a/modkit-core/src/hmm.rs b/modkit-core/src/hmm.rs index f3f6f6fc..9d91d789 100644 --- a/modkit-core/src/hmm.rs +++ b/modkit-core/src/hmm.rs @@ -187,7 +187,7 @@ impl HmmModel { assert_eq!(probs.len(), transitions.len()); let (dp_matrix, pointers) = self.viterbi_forward(&probs, &transitions); let path = self.viterbi_decode(&dp_matrix, &pointers); - assert_eq!(path.len(), scores.len() - 1); + assert_eq!(path.len(), scores.len()); path } @@ -197,21 +197,16 @@ impl HmmModel { pointers: &[PointerCell], ) -> Vec { let final_state = dp_matrix.last().unwrap().argmax(); - // dbg!(final_state); let mut path = vec![final_state]; - let mut curr_pointer = - pointers.last().unwrap().get_value(final_state).unwrap(); - for pointers in pointers.iter().rev().skip(1) { - let pointer = pointers.get_value(curr_pointer); - if let Some(pointer) = pointer { - path.push(pointer); - curr_pointer = pointer; - } else { - break; - } + let mut current_state = final_state; + // The first pointer cell is the empty start cell, and the second + // points back to the un-emitted start state. Decode only the emitted + // states, from the final score back through the second score. + for pointers in pointers.iter().skip(2).rev() { + current_state = pointers.get_value(current_state).unwrap(); + path.push(current_state); } - path.pop(); path.reverse(); path } From d80b35c7b59cb269f9b03519035bf0b94d7fed75 Mon Sep 17 00:00:00 2001 From: SuhasSrinivasan <32346517+SuhasSrinivasan@users.noreply.github.com> Date: Thu, 6 Aug 2026 08:18:59 -0700 Subject: [PATCH 3/4] test: expose DMR contig re-entry across batches --- modkit-core/src/dmr/single_site.rs | 28 +++- modkit/tests/test_dmr.rs | 146 ++++++++++++++++++ tests/resources/dmr_contig_order.fa | 6 + tests/resources/dmr_contig_order.fa.fai | 3 + tests/resources/dmr_contig_order_a.bed | 27 ++++ tests/resources/dmr_contig_order_a.bed.gz | Bin 0 -> 237 bytes tests/resources/dmr_contig_order_a.bed.gz.tbi | Bin 0 -> 142 bytes tests/resources/dmr_contig_order_b.bed | 27 ++++ tests/resources/dmr_contig_order_b.bed.gz | Bin 0 -> 234 bytes tests/resources/dmr_contig_order_b.bed.gz.tbi | Bin 0 -> 142 bytes 10 files changed, 236 insertions(+), 1 deletion(-) create mode 100644 tests/resources/dmr_contig_order.fa create mode 100644 tests/resources/dmr_contig_order.fa.fai create mode 100644 tests/resources/dmr_contig_order_a.bed create mode 100644 tests/resources/dmr_contig_order_a.bed.gz create mode 100644 tests/resources/dmr_contig_order_a.bed.gz.tbi create mode 100644 tests/resources/dmr_contig_order_b.bed create mode 100644 tests/resources/dmr_contig_order_b.bed.gz create mode 100644 tests/resources/dmr_contig_order_b.bed.gz.tbi diff --git a/modkit-core/src/dmr/single_site.rs b/modkit-core/src/dmr/single_site.rs index fea3d05d..ffac1d74 100644 --- a/modkit-core/src/dmr/single_site.rs +++ b/modkit-core/src/dmr/single_site.rs @@ -808,6 +808,9 @@ fn collapse_counts( } type ChromToSingleScores = (String, Vec>); + +fn sort_chrom_to_site_scores(_: &mut [ChromToSingleScores]) {} + fn process_batch_of_positions( batch: DmrBatchOfPositions, sample_index: Arc, @@ -816,7 +819,7 @@ fn process_batch_of_positions( let (a_lines, b_lines) = sample_index.read_bedmethyl_lines_organized_by_position(batch)?; - let chrom_to_site_scores = a_lines + let mut chrom_to_site_scores = a_lines .into_iter() // intersect a_lines and b_lines on contig/chrom, there should be a // filter upstream of this to make sure that this is not ever a miss @@ -852,9 +855,32 @@ fn process_batch_of_positions( }) .collect::>(); + sort_chrom_to_site_scores(&mut chrom_to_site_scores); + Ok(chrom_to_site_scores) } +#[cfg(test)] +mod chrom_score_order_tests { + use super::{sort_chrom_to_site_scores, ChromToSingleScores}; + + #[test] + fn chrom_scores_are_sorted_lexically() { + let mut chrom_to_site_scores = ["cc", "aa", "bb"] + .into_iter() + .map(|chrom| (chrom.to_string(), Vec::new())) + .collect::>(); + + sort_chrom_to_site_scores(&mut chrom_to_site_scores); + + let chroms = chrom_to_site_scores + .iter() + .map(|(chrom, _)| chrom.as_str()) + .collect::>(); + assert_eq!(chroms, vec!["aa", "bb", "cc"]); + } +} + struct Coverages { a_coverages: Vec, b_coverages: Vec, diff --git a/modkit/tests/test_dmr.rs b/modkit/tests/test_dmr.rs index dbfa61f6..8392da57 100644 --- a/modkit/tests/test_dmr.rs +++ b/modkit/tests/test_dmr.rs @@ -149,6 +149,152 @@ fn segmentation_includes_last_site_and_is_thread_deterministic() { assert_eq!(last_segment_end, last_site_end); } +struct DmrOutput { + sites: String, + segments: String, +} + +fn run_multi_contig_dmr( + output_dir: &Path, + label: &str, + interval_size: &str, + threads: &str, + io_threads: &str, +) -> DmrOutput { + let sites = output_dir.join(format!("sites-{label}.bed")); + let segments = output_dir.join(format!("segments-{label}.bed")); + + run_modkit(&[ + "dmr", + "pair", + "-a", + "../tests/resources/dmr_contig_order_a.bed.gz", + "-b", + "../tests/resources/dmr_contig_order_b.bed.gz", + "-o", + sites.to_str().unwrap(), + "--segment", + segments.to_str().unwrap(), + "--ref", + "../tests/resources/dmr_contig_order.fa", + "--header", + "--base", + "C", + "--max-coverages", + "100", + "100", + "--interval-size", + interval_size, + "--batch-size", + "1", + "--threads", + threads, + "--io-threads", + io_threads, + "--suppress-progress", + "--force", + ]) + .expect("multi-contig segmented DMR run should succeed"); + + DmrOutput { + sites: fs::read_to_string(sites).unwrap(), + segments: fs::read_to_string(segments).unwrap(), + } +} + +fn parse_site_keys(output: &str) -> Vec<(String, u64, u64)> { + output + .lines() + .skip(1) + .map(|row| { + let fields = row.split('\t').collect::>(); + ( + fields[0].to_string(), + fields[1].parse::().unwrap(), + fields[2].parse::().unwrap(), + ) + }) + .collect() +} + +fn parse_segments(output: &str) -> Vec<(String, u64, u64, String, usize)> { + output + .lines() + .skip(1) + .map(|row| { + let fields = row.split('\t').collect::>(); + ( + fields[0].to_string(), + fields[1].parse::().unwrap(), + fields[2].parse::().unwrap(), + fields[3].to_string(), + fields[5].parse::().unwrap(), + ) + }) + .collect() +} + +#[test] +fn multi_contig_segmentation_is_invariant_to_batch_geometry_and_threads() { + let temp_dir = tempfile::tempdir().unwrap(); + let interval_ten_one = + run_multi_contig_dmr(temp_dir.path(), "i10-t1", "10", "1", "1"); + let interval_ten_four = + run_multi_contig_dmr(temp_dir.path(), "i10-t4", "10", "4", "2"); + let interval_three_one = + run_multi_contig_dmr(temp_dir.path(), "i3-t1", "3", "1", "1"); + let interval_three_four = + run_multi_contig_dmr(temp_dir.path(), "i3-t4", "3", "4", "2"); + + for (label, other) in [ + ("interval 10, four threads", &interval_ten_four), + ("interval 3, one thread", &interval_three_one), + ("interval 3, four threads", &interval_three_four), + ] { + assert!( + interval_ten_one.sites.as_bytes() == other.sites.as_bytes(), + "site output changed for {label}" + ); + assert!( + interval_ten_one.segments.as_bytes() == other.segments.as_bytes(), + "segment output changed for {label}" + ); + } + + let site_keys = parse_site_keys(&interval_ten_one.sites); + let expected_site_keys = [("aa", 6u64), ("bb", 15u64), ("cc", 6u64)] + .into_iter() + .flat_map(|(chrom, size)| { + (0..size).map(move |position| { + (chrom.to_string(), position, position + 1) + }) + }) + .collect::>(); + assert_eq!(site_keys, expected_site_keys); + + let segments = parse_segments(&interval_ten_one.segments); + assert_eq!( + segments, + vec![ + ("aa".to_string(), 0, 6, "different".to_string(), 6), + ("bb".to_string(), 0, 15, "different".to_string(), 15), + ("cc".to_string(), 0, 6, "different".to_string(), 6), + ] + ); + + for (chrom, start, end) in site_keys { + let covering_segments = segments + .iter() + .filter(|(segment_chrom, segment_start, segment_end, _, _)| { + segment_chrom == &chrom + && *segment_start <= start + && end <= *segment_end + }) + .count(); + assert_eq!(covering_segments, 1, "site {chrom}:{start}-{end}"); + } +} + // todo // test pair with explicit index // test multi diff --git a/tests/resources/dmr_contig_order.fa b/tests/resources/dmr_contig_order.fa new file mode 100644 index 00000000..ac83f386 --- /dev/null +++ b/tests/resources/dmr_contig_order.fa @@ -0,0 +1,6 @@ +>aa +CCCCCC +>bb +CCCCCCCCCCCCCCC +>cc +CCCCCC diff --git a/tests/resources/dmr_contig_order.fa.fai b/tests/resources/dmr_contig_order.fa.fai new file mode 100644 index 00000000..84bdbae4 --- /dev/null +++ b/tests/resources/dmr_contig_order.fa.fai @@ -0,0 +1,3 @@ +aa 6 4 6 7 +bb 15 15 15 16 +cc 6 35 6 7 diff --git a/tests/resources/dmr_contig_order_a.bed b/tests/resources/dmr_contig_order_a.bed new file mode 100644 index 00000000..653af2d1 --- /dev/null +++ b/tests/resources/dmr_contig_order_a.bed @@ -0,0 +1,27 @@ +aa 0 1 C 10 + 0 1 255,0,0 10 100.00 10 0 0 0 0 0 0 +aa 1 2 C 10 + 1 2 255,0,0 10 100.00 10 0 0 0 0 0 0 +aa 2 3 C 10 + 2 3 255,0,0 10 100.00 10 0 0 0 0 0 0 +aa 3 4 C 10 + 3 4 255,0,0 10 100.00 10 0 0 0 0 0 0 +aa 4 5 C 10 + 4 5 255,0,0 10 100.00 10 0 0 0 0 0 0 +aa 5 6 C 10 + 5 6 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 0 1 C 10 + 0 1 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 1 2 C 10 + 1 2 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 2 3 C 10 + 2 3 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 3 4 C 10 + 3 4 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 4 5 C 10 + 4 5 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 5 6 C 10 + 5 6 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 6 7 C 10 + 6 7 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 7 8 C 10 + 7 8 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 8 9 C 10 + 8 9 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 9 10 C 10 + 9 10 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 10 11 C 10 + 10 11 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 11 12 C 10 + 11 12 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 12 13 C 10 + 12 13 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 13 14 C 10 + 13 14 255,0,0 10 100.00 10 0 0 0 0 0 0 +bb 14 15 C 10 + 14 15 255,0,0 10 100.00 10 0 0 0 0 0 0 +cc 0 1 C 10 + 0 1 255,0,0 10 100.00 10 0 0 0 0 0 0 +cc 1 2 C 10 + 1 2 255,0,0 10 100.00 10 0 0 0 0 0 0 +cc 2 3 C 10 + 2 3 255,0,0 10 100.00 10 0 0 0 0 0 0 +cc 3 4 C 10 + 3 4 255,0,0 10 100.00 10 0 0 0 0 0 0 +cc 4 5 C 10 + 4 5 255,0,0 10 100.00 10 0 0 0 0 0 0 +cc 5 6 C 10 + 5 6 255,0,0 10 100.00 10 0 0 0 0 0 0 diff --git a/tests/resources/dmr_contig_order_a.bed.gz b/tests/resources/dmr_contig_order_a.bed.gz new file mode 100644 index 0000000000000000000000000000000000000000..4e4c92aaee65e457c0e29e5ff2c48b5df369f829 GIT binary patch literal 237 zcmb2|=3rp}f&Xj_PR>jW7Z~PV-pza1fQRM6uT5S(4Oi}$stNA=&&Mbd;F*4@;IMXF z@wKy@OpkfhY_8ADvVEao(AXKbW!?VwhUyY0Ih1#LIGIhF5Ua9S;-rRNP(ZIpi1^A3 zH?3Qyf}G}RLU}$acV>9(Fq~MBtW++Q{lZOZS9b}Ivg~13Ha}(N<=t}9G)D&>+fW0!+G_C3WW2B`QVfl U`^s4v8067BB+bAK_8*7<024b^%m4rY literal 0 HcmV?d00001 diff --git a/tests/resources/dmr_contig_order_a.bed.gz.tbi b/tests/resources/dmr_contig_order_a.bed.gz.tbi new file mode 100644 index 0000000000000000000000000000000000000000..241c105f60749e0db51bfcfb6b41c932c19b89b7 GIT binary patch literal 142 zcmb2|=3rp}f&Xj_PR>jWg$&$jWXBnnm+Re*s$isYLR-|gX!>9YDY8$uw=VPpy5Im2?Jnf@Y zS9arpA5EO^|N3r?D-oRGAg^@CzIy%E10HM=CmDq{HMYEAkUZ+(R@l&#wP1mjgR61C z#d!@L8&oGIcq<*8=F+6w$r3Hpx&6uoUXv2`7ZPk%ZOlAA%*>Y=)n@p-+5fnN#qRpP zW@f`5^&e$rG#WoUT3ZCLO>ZtMSirQ7{q2W_kmkkVynS(|7FAatya)wz{;@B(_veu# RD+7Z(ns1~Tn8BU{5dfCkR$c%A literal 0 HcmV?d00001 diff --git a/tests/resources/dmr_contig_order_b.bed.gz.tbi b/tests/resources/dmr_contig_order_b.bed.gz.tbi new file mode 100644 index 0000000000000000000000000000000000000000..64ed8bc0414e8fe543be0dc850882ad129c919f5 GIT binary patch literal 142 zcmb2|=3rp}f&Xj_PR>jWg$&$ Date: Thu, 6 Aug 2026 08:19:55 -0700 Subject: [PATCH 4/4] fix: preserve DMR contig order across batches --- modkit-core/src/dmr/single_site.rs | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/modkit-core/src/dmr/single_site.rs b/modkit-core/src/dmr/single_site.rs index ffac1d74..b6da118b 100644 --- a/modkit-core/src/dmr/single_site.rs +++ b/modkit-core/src/dmr/single_site.rs @@ -809,7 +809,12 @@ fn collapse_counts( type ChromToSingleScores = (String, Vec>); -fn sort_chrom_to_site_scores(_: &mut [ChromToSingleScores]) {} +fn sort_chrom_to_site_scores(chrom_to_site_scores: &mut [ChromToSingleScores]) { + // The batch maps use hash iteration order, while the writer and stateful + // segmenter require the lexical contig order used by SingleSiteBatches. + chrom_to_site_scores + .sort_unstable_by(|(a_chrom, _), (b_chrom, _)| a_chrom.cmp(b_chrom)); +} fn process_batch_of_positions( batch: DmrBatchOfPositions,