From 57c61c561a44addafe04f9b44dfd6de01cff5f95 Mon Sep 17 00:00:00 2001 From: SuhasSrinivasan <32346517+SuhasSrinivasan@users.noreply.github.com> Date: Tue, 4 Aug 2026 23:08:12 -0700 Subject: [PATCH] Apply localize minimum coverage before aggregation --- modkit-core/src/localise/subcommand.rs | 1 + modkit-core/src/localise/util.rs | 2 ++ modkit/tests/test_localize.rs | 48 ++++++++++++++++++++++++++ 3 files changed, 51 insertions(+) diff --git a/modkit-core/src/localise/subcommand.rs b/modkit-core/src/localise/subcommand.rs index 0a36c14..551320e 100644 --- a/modkit-core/src/localise/subcommand.rs +++ b/modkit-core/src/localise/subcommand.rs @@ -271,6 +271,7 @@ impl EntryLocalize { .map(|gr| { gr.into_localized_mod_counts( &tabix_index, + min_cov, self.stranded_features, stranded_features, self.io_threads, diff --git a/modkit-core/src/localise/util.rs b/modkit-core/src/localise/util.rs index 54fcffb..79cc8a0 100644 --- a/modkit-core/src/localise/util.rs +++ b/modkit-core/src/localise/util.rs @@ -190,6 +190,7 @@ impl GenomeRegion { pub(super) fn into_localized_mod_counts( self, index: &HtsTabixHandler, + min_coverage: u64, strand_rule: Option, stranded_features: Option, io_threads: usize, @@ -203,6 +204,7 @@ impl GenomeRegion { let anchor_point = self.midpoint(); let loc_counts = bedmethyl_records .into_par_iter() + .filter(|bm| bm.valid_coverage >= min_coverage) .filter(|bm| { stranded_features .map(|f| { diff --git a/modkit/tests/test_localize.rs b/modkit/tests/test_localize.rs index 2225773..a62dd90 100644 --- a/modkit/tests/test_localize.rs +++ b/modkit/tests/test_localize.rs @@ -1,4 +1,6 @@ use crate::common::run_modkit; +use std::fs; +use tempfile::tempdir; mod common; @@ -9,3 +11,49 @@ fn test_localise_helps() { let _ = run_modkit(&["localise", "--help"]) .expect("failed to run modkit localise help"); } + +#[test] +fn test_localize_min_coverage_filters_before_aggregation() { + let temp_dir = tempdir().unwrap(); + let regions = temp_dir.path().join("regions.bed"); + let genome_sizes = temp_dir.path().join("genome-sizes.tsv"); + fs::write(®ions, "chr20\t9681998\t9681999\nchr20\t9838537\t9838538\n") + .unwrap(); + fs::write(&genome_sizes, "chr20\t64444167\n").unwrap(); + + let run = |min_coverage: u64| { + let output = temp_dir.path().join(format!("min-{min_coverage}.tsv")); + let min_coverage = min_coverage.to_string(); + run_modkit(&[ + "localize", + "../tests/resources/lung_00733-m_adjacent-normal_5mc-5hmc_chr20_cpg_pileup.bed.gz", + "--regions", + regions.to_str().unwrap(), + "--genome-sizes", + genome_sizes.to_str().unwrap(), + "--window", + "1", + "--min-coverage", + &min_coverage, + "--threads", + "1", + "--io-threads", + "1", + "--out-file", + output.to_str().unwrap(), + ]) + .unwrap(); + fs::read_to_string(output).unwrap().replace("\r\n", "\n") + }; + + assert_eq!( + run(1), + "mod_code\toffset\tn_valid\tn_mod\tpercent_modified\n\ + C\t-1\t24\t3\t12.5\n" + ); + assert_eq!( + run(3), + "mod_code\toffset\tn_valid\tn_mod\tpercent_modified\n\ + C\t-1\t23\t2\t8.695652\n" + ); +}