echtvar

pandora-config.example.toml

@@ -4,9 +4,6 @@
 # General filesystem layout / I/O
 #######################################
 
-# Directory where POD / run description files are located.
-pod_dir = "/data/run_data"
-
 # Root directory where all results will be written.
 result_dir = "/mnt/beegfs02/scratch/t_steimle/data/wgs"
 
@@ -25,19 +22,9 @@ threads = 5
 # Singularity bin
 singularity_bin = "module load singularity-ce && singularity"
 
-# Temporary directory used when unarchiving input data.
-unarchive_tmp_dir = "/data/unarchived"
-
-# Maximum memory available for dockerized tools, in GiB.
-docker_max_memory_go = 400
-
-# Path to the SQLite database of processed cases.
-db_cases_path = "/data/cases.sqlite"
-
 # Path to the conda activation script.
 conda_sh = "/mnt/beegfs02/software/recherche/miniconda/25.1.1/etc/profile.d/conda.sh"
 
-
 #######################################
 # Reference genome & annotations
 #######################################
@@ -62,22 +49,11 @@ refseq_gff = "/data/ref/hs1/chm13v2.0_RefSeq_Liftoff_v5.1_sorted.gff3.gz"
 # {id}         -> case identifier
 mask_bed = "{result_dir}/{id}/diag/mask.bed"
 
-# BED file with early-replicating regions.
-early_bed = "/data/ref/hs1/replication_early_25_hs1.bed"
-
-# BED file with late-replicating regions.
-late_bed = "/data/ref/hs1/replication_late_75_hs1.bed"
-
-# BED file with CpG coordinates.
-cpg_bed = "/data/ref/hs1/hs1/hs1_CpG.bed"
-
 # Panels of interest: [ [name, bed_path], ... ]
 panels = [
-  ["OncoT",         "/data/ref/hs1/V1_V2_V3_V4_V5_intersect_targets_hs1_uniq.bed"],
-  ["variable_chips","/data/ref/hs1/top_1500_sd_pos.bed"],
+  ["Repeat",         "/home/t_steimle/ref/hs1/all_repeats_chm13_final.bed"],
 ]
 
-
 #######################################
 # Sample naming / BAM handling
 #######################################
@@ -94,13 +70,12 @@ haplotagged_bam_tag_name = "HP"
 # Minimum MAPQ for reads kept during BAM filtering.
 bam_min_mapq = 40
 
-# Threads for BAM-level operations (view/sort/index…).
-bam_n_threads = 150
+# Number of threads for hts BAM reader decrompression (should be adapted to IO speed).
+bam_n_threads = 4
 
 # Number of reads sampled for BAM composition estimation.
 bam_composition_sample_size = 20000
 
-
 #######################################
 # Coverage counting / somatic-scan
 #######################################
@@ -117,7 +92,6 @@ count_n_chunks = 1000
 # Force recomputation of counting even if outputs exist.
 somatic_scan_force = false
 
-
 #######################################
 # Somatic pipeline global settings
 #######################################

src/annotation/cosmic.rs

@@ -27,6 +27,7 @@ impl FromStr for Cosmic {
     /// - The input must contain exactly three parts, separated by semicolons (`;`).
     /// - The third part must be of the form `CNT=<number>`, where `<number>` can be parsed as a `u64`.
     /// - If the first part contains the word `"MISSING"`, parsing will fail.
+    /// encode with echtvar: [{"field":"GENOME_SCREEN_SAMPLE_COUNT", "alias": "CNT"}]
     ///
     /// # Examples
     ///

src/annotation/gnomad.rs

@@ -2,6 +2,26 @@ use bitcode::{Decode, Encode};
 use serde::{Deserialize, Serialize};
 use std::str::FromStr;
 
+/// Generated with echtvar json:
+/// ```
+///[
+///  { "field": "AC",     "alias": "gnomad_ac" },
+///  { "field": "AN",     "alias": "gnomad_an" },
+///  { "field": "AF",     "alias": "gnomad_af",     "multiplier": 2000000 },
+///  { "field": "AF_oth", "alias": "gnomad_af_oth", "multiplier": 2000000 },
+///  { "field": "AF_ami", "alias": "gnomad_af_ami", "multiplier": 2000000 },
+///  { "field": "AF_sas", "alias": "gnomad_af_sas", "multiplier": 2000000 },
+///  { "field": "AF_fin", "alias": "gnomad_af_fin", "multiplier": 2000000 },
+///  { "field": "AF_eas", "alias": "gnomad_af_eas", "multiplier": 2000000 },
+///  { "field": "AF_amr", "alias": "gnomad_af_amr", "multiplier": 2000000 },
+///  { "field": "AF_afr", "alias": "gnomad_af_afr", "multiplier": 2000000 },
+///  { "field": "AF_mid", "alias": "gnomad_af_mid", "multiplier": 2000000 },
+///  { "field": "AF_asj", "alias": "gnomad_af_asj", "multiplier": 2000000 },
+///  { "field": "AF_nfe", "alias": "gnomad_af_nfe", "multiplier": 2000000 }
+///]
+///
+/// ```
+
 #[derive(Debug, PartialEq, Serialize, Deserialize, Clone, Encode, Decode)]
 pub struct GnomAD {
     pub gnomad_ac: u64,

src/annotation/mod.rs

@@ -700,7 +700,7 @@ pub trait CallerCat {
 /// # Arguments
 ///
 /// * `anns` - A slice of `Annotation` values describing variant annotations,
-///            such as gnomAD frequency, sample types, and alternate counts.
+///   such as gnomAD frequency, sample types, and alternate counts.
 ///
 /// # Returns
 ///

src/annotation/vep.rs

@@ -1,4 +1,4 @@
-use anyhow::{anyhow, Context};
+use anyhow::anyhow;
 use bitcode::{Decode, Encode};
 use hashbrown::HashMap;
 use itertools::Itertools;
@@ -17,8 +17,6 @@ use crate::{
     commands::{Command as JobCommand, LocalBatchRunner, LocalRunner, SbatchRunner},
     config::Config,
     helpers::singularity_bind_flags,
-    run,
-    runners::Run,
 };
 
 use super::ncbi::NCBIAcc;
@@ -251,17 +249,17 @@ pub enum VepConsequence {
     /// Deletion of a regulatory region
     RegulatoryRegionAblation,
     /// Amplification of a regulatory region
-    RegulatoryRegionAmplification,
-    /// Variant causing a feature to be extended
-    FeatureElongation,
-    /// Variant in a regulatory region
-    RegulatoryRegionVariant,
-    /// Variant causing a feature to be shortened
-    FeatureTruncation,
-    /// Variant in intergenic region
-    IntergenicVariant,
-    /// General sequence variant
-    SequenceVariant,
+RegulatoryRegionAmplification,
+/// Variant causing a feature to be extended
+FeatureElongation,
+/// Variant in a regulatory region
+RegulatoryRegionVariant,
+/// Variant causing a feature to be shortened
+FeatureTruncation,
+/// Variant in intergenic region
+IntergenicVariant,
+/// General sequence variant
+SequenceVariant,
 }
 
 /// Represents the severity of a variant's impact as predicted by the
@@ -270,27 +268,27 @@ pub enum VepConsequence {
 /// The impact categories are ordered from most severe (HIGH) to least severe (MODIFIER).
 #[derive(Debug, Clone, PartialEq, Eq, Serialize, Deserialize, Encode, Decode)]
 pub enum VepImpact {
-    /// High impact variants are expected to have high (disruptive) impact in the protein,
-    /// probably causing protein truncation, loss of function or triggering nonsense mediated decay.
-    HIGH,
-    /// Moderate impact variants are non-disruptive variants that might change protein effectiveness.
-    MODERATE,
-    /// Low impact variants are mostly harmless or unlikely to change protein behavior.
-    LOW,
-    /// Modifier variants are usually non-coding variants or variants affecting non-coding genes,
-    /// where predictions are difficult or there is no evidence of impact.
-    MODIFIER,
+/// High impact variants are expected to have high (disruptive) impact in the protein,
+/// probably causing protein truncation, loss of function or triggering nonsense mediated decay.
+HIGH,
+/// Moderate impact variants are non-disruptive variants that might change protein effectiveness.
+MODERATE,
+/// Low impact variants are mostly harmless or unlikely to change protein behavior.
+LOW,
+/// Modifier variants are usually non-coding variants or variants affecting non-coding genes,
+/// where predictions are difficult or there is no evidence of impact.
+MODIFIER,
 }
 
 impl Display for VepImpact {
-    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
-        write!(
-            f,
-            "{}",
-            match self {
-                VepImpact::HIGH => "HIGH",
-                VepImpact::MODERATE => "MODERATE",
-                VepImpact::LOW => "LOW",
+fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
+    write!(
+        f,
+        "{}",
+        match self {
+            VepImpact::HIGH => "HIGH",
+            VepImpact::MODERATE => "MODERATE",
+            VepImpact::LOW => "LOW",
                 VepImpact::MODIFIER => "MODIFIER",
             }
         )
@@ -645,7 +643,7 @@ impl FromStr for VEPExtra {
 }
 
 #[derive(Debug)]
-struct VepJob {
+pub struct VepJob {
     in_vcf: PathBuf,
     out_vcf: PathBuf,
     config: Config,
@@ -870,7 +868,7 @@ pub fn get_best_vep(d: &[VEP]) -> anyhow::Result<VEP> {
 #[cfg(test)]
 mod tests {
     use super::*;
-    use crate::helpers::test_init;
+    use crate::{helpers::test_init, run};
 
     #[test]
     fn vep_run() -> anyhow::Result<()> {

src/collection/prom_run.rs

@@ -1388,7 +1388,7 @@ fn sort_and_index_chunks(
     // Index every sorted BAM.
     let index_jobs: Vec<SamtoolsIndex> = original_to_sorted
         .values()
-        .map(|sorted| SamtoolsIndex::from_config(config, &sorted.to_string_lossy().to_string()))
+        .map(|sorted| SamtoolsIndex::from_config(config, &sorted.to_string_lossy().as_ref()))
         .collect();
 
     info!("Submitting {} index jobs", index_jobs.len());

src/commands/samtools.rs

@@ -631,13 +631,12 @@ impl LocalBatchRunner for SamtoolsSort {}
 #[cfg(test)]
 mod tests {
     use super::*;
-    use log::info;
 
     use crate::{
         commands::{run_many_sbatch, SlurmRunner},
         config::Config,
         helpers::test_init,
-        run, TEST_DIR,
+        TEST_DIR,
     };
 
     // #[test]

src/config.rs

@@ -13,9 +13,6 @@ const CONFIG_TEMPLATE: &str = include_str!("../pandora-config.example.toml");
 /// `{result_dir}`, `{id}`, `{time}`, `{reference_name}`, `{haplotagged_bam_tag_name}`, `{output_dir}`.
 pub struct Config {
     // === General filesystem layout / I/O ===
-    /// Directory where POD / run description files are located.
-    pub pod_dir: String,
-
     /// Root directory where all results will be written.
     pub result_dir: String,
 
@@ -31,18 +28,9 @@ pub struct Config {
     /// Software threads
     pub threads: u8,
 
-    /// Singularity bin
+    /// Singularity/Apptainer bin
     pub singularity_bin: String,
 
-    /// Temporary directory used when unarchiving input data.
-    pub unarchive_tmp_dir: String,
-
-    /// Maximum memory available for dockerized tools, in GiB.
-    pub docker_max_memory_go: u16,
-
-    /// Path to the SQLite database of processed cases.
-    pub db_cases_path: String,
-
     /// Path to the `conda.sh` activation script (used to activate envs before running tools).
     pub conda_sh: String,
 
@@ -53,7 +41,7 @@ pub struct Config {
     /// Minimum MAPQ for reads to be kept during BAM filtering.
     pub bam_min_mapq: u8,
 
-    /// Number of threads for BAM processing steps (view, sort, index…).
+    /// Number of threads for hts BAM reader
     pub bam_n_threads: u8,
 
     /// Number of reads sampled when estimating BAM composition (e.g. tumor contamination).
@@ -82,15 +70,6 @@ pub struct Config {
     /// - `{id}`: case identifier
     pub mask_bed: String,
 
-    /// BED file with CpG coordinates in the reference.
-    pub cpg_bed: String,
-
-    /// BED file with early-replicating regions (used for replication timing–based analyses).
-    pub early_bed: String,
-
-    /// BED file with late-replicating regions.
-    pub late_bed: String,
-
     /// Panels of interest (name, BED path).
     pub panels: Vec<(String, String)>,
 

src/runners.rs

@@ -3,20 +3,16 @@ use std::{
     io::{BufRead, BufReader, Write},
     path::Path,
     process::{Child, Command, Stdio},
-    sync::{
-        mpsc::{self, TryRecvError},
-        Arc, Mutex,
-    },
-    thread,
+    sync::mpsc::{self, TryRecvError},
 };
 
 use anyhow::Context;
 use chrono::{DateTime, Utc};
-use log::{info, warn};
+use log::info;
 use serde::{Deserialize, Serialize};
 use uuid::Uuid;
 
-use crate::{config::Config, io::writers::get_gz_writer, DOCKER_ID};
+use crate::io::writers::get_gz_writer;
 
 /// Trait for running a command.
 pub trait Run {
@@ -140,187 +136,6 @@ pub trait Log {
     fn log(&self) -> String;
 }
 
-/// Represents a Docker command to be run, with facilities for execution, monitoring, and logging.
-#[derive(Debug, Default)]
-pub struct DockerRun {
-    /// Arguments for the Docker command
-    pub args: Vec<String>,
-    /// Container ID of the running Docker instance
-    pub container_id: String,
-    /// Start time of the Docker command execution
-    pub start: DateTime<Utc>,
-    /// Accumulated logs of the Docker command execution
-    pub logs: Arc<Mutex<String>>,
-}
-
-impl DockerRun {
-    /// Helper function to execute a Docker command and capture its output.
-    pub fn new(args: &[&str]) -> Self {
-        DockerRun {
-            args: args.iter().map(|e| e.to_string()).collect(),
-            start: Utc::now(),
-            ..Default::default()
-        }
-    }
-}
-
-impl Run for DockerRun {
-    /// Executes the configured Docker command inside a container.
-    ///
-    /// - Sets a Ctrl-C (SIGINT) handler to stop the container if the process is interrupted.
-    /// - Adds memory limits to the Docker arguments based on configuration.
-    /// - Captures and logs the container ID.
-    /// - Spawns a thread to follow container logs in real-time.
-    ///
-    /// # Returns
-    ///
-    /// * `Ok(())` if the container starts successfully.
-    /// * An `anyhow::Error` if setup or execution fails.
-    fn run(&mut self) -> anyhow::Result<()> {
-        // Sets up a Ctrl-C handler to stop Docker containers on interrupt.
-        // ctrlc::try_set_handler(move || {
-        //     if let Ok(container_id) = DOCKER_ID.lock() {
-        //         for id in container_id.iter() {
-        //             warn!("Stopping Docker container {id}...");
-        //             let _ = Command::new("docker").args(["stop", id]).status();
-        //         }
-        //     }
-        //     std::process::exit(1);
-        // })
-        // .context("Failed to set Ctrl-C handler")?;
-
-        if let Err(e) = ctrlc::try_set_handler(move || {
-            if let Ok(container_ids) = DOCKER_ID.lock() {
-                for id in container_ids.iter() {
-                    warn!("Stopping Docker container {id}...");
-                    let _ = Command::new("docker").args(["stop", id]).status();
-                }
-            }
-            std::process::exit(1);
-        }) {
-            // Ignore if a handler was already set; propagate all other errors
-            if !matches!(e, ctrlc::Error::MultipleHandlers) {
-                return Err(e).context("Failed to set Ctrl-C handler")?;
-            }
-        }
-
-        // Configures memory limits for the Docker container.
-        let c = Config::default();
-        self.args
-            .insert(1, format!("--memory={}g", c.docker_max_memory_go));
-
-        // Spawns the Docker command and captures its output.
-        let output = Command::new("docker")
-            .args(&self.args)
-            .stdout(Stdio::piped())
-            .stderr(Stdio::inherit())
-            .output()
-            .expect("Failed to run Docker command");
-
-        let info = format!("Running with docker: {}", &self.args.join(" "));
-        info!("{info}");
-        {
-            let mut logs_lock = self.logs.lock().unwrap();
-            logs_lock.push_str(&format!("{info}\n"));
-        }
-
-        // add id to Arc
-        let id = String::from_utf8_lossy(&output.stdout).trim().to_string();
-        {
-            DOCKER_ID.lock().unwrap().push(id.clone());
-            self.container_id = id.clone();
-        }
-
-        // Spawn a thread to follow the logs
-        let log_id = id.clone();
-        let logs_clone = Arc::clone(&self.logs);
-        let _loger_thread = thread::spawn(move || {
-            let output = Command::new("docker")
-                .args(["inspect", "--format='{{.Config.Image}}'", &log_id])
-                .output()
-                .expect("Failed to execute Docker inspect");
-            let container_name = if output.status.success() {
-                String::from_utf8(output.stdout)
-                    .expect("Failed to convert output to string")
-                    .trim() // Trim to remove any trailing newlines
-                    .to_string()
-            } else {
-                "?".to_string()
-            };
-
-            #[allow(clippy::zombie_processes)] // the process is not zombie but waited by Wait trait
-            let mut child = Command::new("docker")
-                .args(["logs", "--follow", &log_id])
-                .stdout(Stdio::piped())
-                .stderr(Stdio::inherit())
-                .spawn()
-                .expect("Failed to follow Docker logs");
-
-            if let Some(stdout) = child.stdout.take() {
-                let reader = BufReader::new(stdout);
-                for line in reader.lines().map_while(Result::ok) {
-                    info!("[{container_name}] {line}");
-                    let mut logs_lock = logs_clone.lock().unwrap();
-                    logs_lock.push_str(&line);
-                    logs_lock.push('\n');
-                }
-            }
-        });
-        Ok(())
-    }
-}
-
-impl Wait for DockerRun {
-    /// Waits for the Docker container associated with this `DockerRun` instance to finish.
-    ///
-    /// This method blocks until the container exits, using `docker wait`. It captures the
-    /// exit status and logs a warning if the container did not exit successfully.
-    ///
-    /// After the container exits, its ID is removed from the global `DOCKER_ID` list
-    /// to prevent further cleanup or signal handling on an already-finished container.
-    ///
-    /// # Returns
-    ///
-    /// * `Ok(())` if the wait completes successfully (regardless of container success).
-    /// * An error if the `docker wait` command fails to execute.
-    ///
-    /// # Errors
-    ///
-    /// This function returns an `anyhow::Error` if the underlying `docker wait` command fails to launch.
-    fn wait(&mut self) -> anyhow::Result<()> {
-        let status = Command::new("docker")
-            .args(["wait", &self.container_id])
-            .status()
-            .expect("Failed to wait on Docker container");
-
-        if !status.success() {
-            let warn = format!("Docker command failed with status: {status}");
-            warn!("{warn}");
-            {
-                let mut logs_lock = self.logs.lock().unwrap();
-                logs_lock.push_str(&format!("\n{warn}\n"));
-            }
-        } else {
-            info!(
-                "Container {} exited successfully with status: {}",
-                self.container_id, status
-            );
-        }
-
-        let mut container_id_lock = DOCKER_ID.lock().unwrap();
-        container_id_lock.retain(|entry| *entry != self.container_id);
-
-        Ok(())
-    }
-}
-
-impl Log for DockerRun {
-    fn log(&self) -> String {
-        let logs_lock = self.logs.lock().unwrap();
-        logs_lock.to_string()
-    }
-}
-
 /// Represents a command to be executed, with streaming log capture,
 /// real-time monitoring, and access to its output.
 pub struct CommandRun {

src/scan/scan.rs

@@ -241,7 +241,6 @@ impl fmt::Display for BinOutlier {
 /// - `config: &Config`: A configuration object containing the following fields:
 ///   - `count_bin_size`: The size of each bin in base pairs.
 ///   - `count_n_chunks`: The number of bins per chunk for parallel processing.
-///   - `dict_file`: Path to the dictionary file containing contig names and lengths.
 ///
 /// # Returns
 /// - `anyhow::Result<()>`: Returns `Ok(())` if successful, or an error if any operation fails.
@@ -283,7 +282,6 @@ impl fmt::Display for BinOutlier {
 /// # Errors
 /// This function will return an error if:
 /// - The output directory cannot be created.
-/// - The dictionary file cannot be read.
 /// - A `Bin` object cannot be created for a specific region.
 /// - Any I/O operation (e.g., writing results) fails.
 pub fn par_whole_scan(id: &str, time_point: &str, config: &Config) -> anyhow::Result<()> {

src/variant/variants_stats.rs

@@ -166,7 +166,7 @@ use rayon::prelude::*;
 use serde::{Deserialize, Serialize, Serializer};
 
 use crate::{
-    annotation::{vep::VepImpact, Annotation, ReplicationClass},
+    annotation::{vep::VepImpact, Annotation},
     config::Config,
     helpers::bin_data,
     io::{
@@ -235,7 +235,12 @@ where
 }
 
 impl VariantsStats {
-    pub fn new(variants: &mut Variants, _id: &str, config: &Config, high_depth_ranges: &[GenomeRange]) -> anyhow::Result<Self> {
+    pub fn new(
+        variants: &mut Variants,
+        _id: &str,
+        config: &Config,
+        high_depth_ranges: &[GenomeRange],
+    ) -> anyhow::Result<Self> {
         let n = variants.data.len() as u32;
         let alteration_categories: DashMap<String, u32> = DashMap::new();
         let vep_impact: DashMap<String, u32> = DashMap::new();
@@ -343,8 +348,6 @@ impl VariantsStats {
 
         let all_somatic_rates = somatic_rates(&variants.data, &exon_ranges, config)?;
 
-        
-
         let exon_ranges_ref: Vec<&GenomeRange> = exon_ranges.iter().collect();
         let exons_high_depth = range_intersection_par(
             &high_depth_ranges.iter().collect::<Vec<&GenomeRange>>(),
@@ -378,87 +381,6 @@ impl VariantsStats {
         );
         mutation_rates.push(("Exons HighDepths".to_string(), res));
 
-        // CpG
-        let cpg_ranges: Vec<GenomeRange> = read_bed(&config.cpg_bed)?
-            .into_iter()
-            .map(|e| e.range)
-            .collect();
-        let ann = Annotation::CpG;
-        let res = variants.annotate_with_ranges(
-            &cpg_ranges,
-            Some(ann.clone()),
-            config.min_n_callers,
-            Vec::new(),
-        );
-        mutation_rates.push((ann.to_string(), res));
-
-        // CpG HighDepths
-        let cpg_high_depth = range_intersection_par(
-            &high_depth_ranges.iter().collect::<Vec<&GenomeRange>>(),
-            &cpg_ranges.iter().collect::<Vec<&GenomeRange>>(),
-        );
-        let res = variants.annotate_with_ranges(
-            &cpg_high_depth,
-            Some(ann.clone()),
-            config.min_n_callers,
-            Vec::new(),
-        );
-        mutation_rates.push(("CpG HighDepths".to_string(), res));
-
-        // Early replication
-        let early_ranges: Vec<GenomeRange> = read_bed(&config.early_bed)?
-            .into_iter()
-            .map(|e| e.range)
-            .collect();
-        let ann = Annotation::ReplicationTiming(ReplicationClass::Early);
-        let res = variants.annotate_with_ranges(
-            &early_ranges,
-            Some(ann.clone()),
-            config.min_n_callers,
-            Vec::new(),
-        );
-        mutation_rates.push((ann.to_string(), res));
-
-        // Early replication HighDepths
-        let early_ranges_high_depth = range_intersection_par(
-            &high_depth_ranges.iter().collect::<Vec<&GenomeRange>>(),
-            &early_ranges.iter().collect::<Vec<&GenomeRange>>(),
-        );
-        let res = variants.annotate_with_ranges(
-            &early_ranges_high_depth,
-            Some(ann.clone()),
-            config.min_n_callers,
-            Vec::new(),
-        );
-        mutation_rates.push(("Early replication HighDepths".to_string(), res));
-
-        // Late replication
-        let late_ranges: Vec<GenomeRange> = read_bed(&config.late_bed)?
-            .into_iter()
-            .map(|e| e.range)
-            .collect();
-        let ann = Annotation::ReplicationTiming(ReplicationClass::Late);
-        let res = variants.annotate_with_ranges(
-            &late_ranges,
-            Some(ann.clone()),
-            config.min_n_callers,
-            Vec::new(),
-        );
-        mutation_rates.push((ann.to_string(), res));
-
-        // Late replication HighDepths
-        let late_ranges_high_depth = range_intersection_par(
-            &high_depth_ranges.iter().collect::<Vec<&GenomeRange>>(),
-            &late_ranges.iter().collect::<Vec<&GenomeRange>>(),
-        );
-        let res = variants.annotate_with_ranges(
-            &late_ranges_high_depth,
-            Some(ann.clone()),
-            config.min_n_callers,
-            Vec::new(),
-        );
-        mutation_rates.push(("Late replication HighDepths".to_string(), res));
-
         for (name, path) in config.panels.iter() {
             let panel_ranges: Vec<GenomeRange> =
                 read_bed(path)?.into_iter().map(|e| e.range).collect();
@@ -644,7 +566,7 @@ pub fn somatic_depth_quality_ranges(
         .into_par_iter()
         .map(|contig| {
             let normal_path = format!("{}/{}_count.tsv.gz", cfg.normal_dir_count(id), contig);
-            let tumor_path  = format!("{}/{}_count.tsv.gz", cfg.tumoral_dir_count(id), contig);
+            let tumor_path = format!("{}/{}_count.tsv.gz", cfg.tumoral_dir_count(id), contig);
 
             let normal_rdr = get_gz_reader(&normal_path)
                 .with_context(|| format!("Failed to open normal file: {}", normal_path))?;
@@ -663,33 +585,63 @@ pub fn somatic_depth_quality_ranges(
                 let t_next = tl.next();
                 match (n_next, t_next) {
                     (None, None) => break,
-                    (Some(Err(e)), _) => return Err(anyhow::anyhow!("{} line {}: {}", normal_path, line_no + 1, e)),
-                    (_, Some(Err(e))) => return Err(anyhow::anyhow!("{} line {}: {}", tumor_path,  line_no + 1, e)),
+                    (Some(Err(e)), _) => {
+                        return Err(anyhow::anyhow!(
+                            "{} line {}: {}",
+                            normal_path,
+                            line_no + 1,
+                            e
+                        ))
+                    }
+                    (_, Some(Err(e))) => {
+                        return Err(anyhow::anyhow!(
+                            "{} line {}: {}",
+                            tumor_path,
+                            line_no + 1,
+                            e
+                        ))
+                    }
                     (Some(Ok(n_line)), Some(Ok(t_line))) => {
                         line_no += 1;
 
-                        let n = BinCount::from_tsv_row(&n_line)
-                            .with_context(|| format!("Parse error at {}: {}", normal_path, line_no))?;
-                        let t = BinCount::from_tsv_row(&t_line)
-                            .with_context(|| format!("Parse error at {}: {}", tumor_path, line_no))?;
+                        let n = BinCount::from_tsv_row(&n_line).with_context(|| {
+                            format!("Parse error at {}: {}", normal_path, line_no)
+                        })?;
+                        let t = BinCount::from_tsv_row(&t_line).with_context(|| {
+                            format!("Parse error at {}: {}", tumor_path, line_no)
+                        })?;
 
                         if n.contig != t.contig {
-                            anyhow::bail!("Contig mismatch at line {}: {} vs {}", line_no, n.contig, t.contig);
+                            anyhow::bail!(
+                                "Contig mismatch at line {}: {} vs {}",
+                                line_no,
+                                n.contig,
+                                t.contig
+                            );
                         }
                         if n.start != t.start {
-                            anyhow::bail!("Position mismatch at line {}: {} vs {}", line_no, n.start, t.start);
+                            anyhow::bail!(
+                                "Position mismatch at line {}: {} vs {}",
+                                line_no,
+                                n.start,
+                                t.start
+                            );
                         }
                         // Ensure equal bin widths
                         if n.depths.len() != t.depths.len() {
                             anyhow::bail!(
                                 "Depth vector length mismatch at line {}: {} vs {}",
-                                line_no, n.depths.len(), t.depths.len()
+                                line_no,
+                                n.depths.len(),
+                                t.depths.len()
                             );
                         }
                         if n.low_qualities.len() != t.low_qualities.len() {
                             anyhow::bail!(
                                 "LowQ vector length mismatch at line {}: {} vs {}",
-                                line_no, n.low_qualities.len(), t.low_qualities.len()
+                                line_no,
+                                n.low_qualities.len(),
+                                t.low_qualities.len()
                             );
                         }
 
@@ -699,18 +651,38 @@ pub fn somatic_depth_quality_ranges(
                         });
 
                         // NOTE: if you intended "low-quality regions" (bad), invert predicate.
-                        let lowq_mask_iter = n.low_qualities.iter().zip(&t.low_qualities).map(|(&nq, &tq)| {
-                            nq > cfg.max_depth_low_quality && tq > cfg.max_depth_low_quality
-                        });
-
-                        high_runs.extend(ranges_from_consecutive_true_iter(high_mask_iter, n.start, &n.contig));
-                        lowq_runs.extend(ranges_from_consecutive_true_iter(lowq_mask_iter, n.start, &n.contig));
+                        let lowq_mask_iter =
+                            n.low_qualities
+                                .iter()
+                                .zip(&t.low_qualities)
+                                .map(|(&nq, &tq)| {
+                                    nq > cfg.max_depth_low_quality && tq > cfg.max_depth_low_quality
+                                });
+
+                        high_runs.extend(ranges_from_consecutive_true_iter(
+                            high_mask_iter,
+                            n.start,
+                            &n.contig,
+                        ));
+                        lowq_runs.extend(ranges_from_consecutive_true_iter(
+                            lowq_mask_iter,
+                            n.start,
+                            &n.contig,
+                        ));
                     }
                     (Some(_), None) => {
-                        anyhow::bail!("Line count mismatch: {} has extra lines after {}", normal_path, line_no);
+                        anyhow::bail!(
+                            "Line count mismatch: {} has extra lines after {}",
+                            normal_path,
+                            line_no
+                        );
                     }
                     (None, Some(_)) => {
-                        anyhow::bail!("Line count mismatch: {} has extra lines after {}", tumor_path, line_no);
+                        anyhow::bail!(
+                            "Line count mismatch: {} has extra lines after {}",
+                            tumor_path,
+                            line_no
+                        );
                     }
                 }
             }
@@ -731,14 +703,9 @@ pub fn somatic_depth_quality_ranges(
     ))
 }
 
-
 /// Iterator-based version (no temporary Vec<bool>).
 /// Produces end-exclusive ranges in the same shape as your old function.
-pub fn ranges_from_consecutive_true_iter<I>(
-    mask: I,
-    start0: u32,
-    contig: &str,
-) -> Vec<GenomeRange>
+pub fn ranges_from_consecutive_true_iter<I>(mask: I, start0: u32, contig: &str) -> Vec<GenomeRange>
 where
     I: IntoIterator<Item = bool>,
 {
@@ -753,13 +720,19 @@ where
                 current_start = Some(start0 + i);
             }
         } else if let Some(s) = current_start.take() {
-            ranges.push(GenomeRange { contig, range: s..(start0 + i) });
+            ranges.push(GenomeRange {
+                contig,
+                range: s..(start0 + i),
+            });
         }
         i += 1;
     }
 
     if let Some(s) = current_start {
-        ranges.push(GenomeRange { contig, range: s..(start0 + i) });
+        ranges.push(GenomeRange {
+            contig,
+            range: s..(start0 + i),
+        });
     }
 
     ranges
@@ -771,9 +744,7 @@ pub fn merge_adjacent_ranges(mut ranges: Vec<GenomeRange>) -> Vec<GenomeRange> {
         return ranges;
     }
 
-    ranges.sort_by(|a, b| {
-        (a.contig, a.range.start).cmp(&(b.contig, b.range.start))
-    });
+    ranges.sort_by(|a, b| (a.contig, a.range.start).cmp(&(b.contig, b.range.start)));
 
     let mut merged = Vec::with_capacity(ranges.len());
     let mut cur = ranges[0].clone();