twitcher 0.6.0

use clap::{Args, ValueEnum};

/// Selects how a candidate cluster is scored and validated.
///
/// The two strategies differ only in three primitives — per-event mass, cluster span, and the set
/// of validity gates — which are encapsulated by [`ClusteringSettings::event_mass`],
/// [`ClusteringSettings::span`], and [`ClusteringSettings::is_valid_cluster`]. Callers never branch
/// on the strategy themselves.
#[derive(Debug, Clone, Copy, PartialEq, Eq, ValueEnum)]
pub enum ClusterStrategy {
    /// Original behaviour: log-scaled variant complexity, reference-only span, a minimum of two
    /// records, and no edit-mass gate. This is the default.
    Legacy,
    /// Linear edit mass, symmetric span (`max(reference extent, query extent)`), a minimum of one
    /// record, and an additional `--cluster-min-mass` gate. Lets a single large indel or MNV form a
    /// cluster on its own.
    EditMass,
}

#[derive(Debug, Clone, Args)]
pub struct ClusteringSettings {
    /// Clustering strategy. `legacy` (default) uses log-scaled complexity with a reference-only
    /// span and a two-record minimum; `edit-mass` uses linear edit mass with a symmetric span, a
    /// one-record minimum, and the `--cluster-min-mass` gate.
    #[arg(
        long = "cluster-strategy",
        value_enum,
        default_value_t = ClusterStrategy::Legacy
    )]
    pub strategy: ClusterStrategy,

    /// Minimum number of mutation events required to form a cluster. When unset, defaults to 2
    /// under the `legacy` strategy and 1 under `edit-mass` (where a single sufficiently large
    /// event can form a cluster on its own).
    #[arg(long = "cluster-min-records", value_name = "NUM")]
    pub min_records: Option<usize>,

    /// Minimum edit mass required for a cluster to be reported (`edit-mass` strategy only). Edit
    /// mass is the sum of edited bases across all events (see [`edit_len`]): a SNP contributes 1,
    /// an `n`-base MNV/indel contributes `n`. Defaults to 2.0, which excludes a lone SNP or 1-base
    /// indel while admitting any larger single edit or multi-event cluster.
    #[arg(long = "cluster-min-mass", value_name = "NUM")]
    pub min_mass: Option<f64>,

    /// Maximum gap in bases between consecutive mutations that still allows them to be in the
    /// same cluster. A gap larger than this value will break the cluster.
    #[arg(long = "cluster-max-gap", value_name = "BASES", default_value_t = 20)]
    pub max_gap: usize,

    /// Minimum density required for a cluster to be reported. Density is the cluster mass divided
    /// by its span, i.e. the edited fraction of the cluster region. Higher values require more
    /// tightly packed edits.
    #[arg(
        long = "cluster-min-density",
        value_name = "NUM",
        default_value_t = 0.2
    )]
    pub min_density: f64,
}

impl Default for ClusteringSettings {
    fn default() -> Self {
        Self {
            strategy: ClusterStrategy::Legacy,
            min_records: None,
            min_mass: None,
            max_gap: 20,
            min_density: 0.2,
        }
    }
}

impl ClusteringSettings {
    /// Validates cross-field constraints that clap's derive cannot express. Call once before the
    /// settings are used to drive clustering.
    pub fn validate(&self) -> anyhow::Result<()> {
        if self.min_mass.is_some() && self.strategy != ClusterStrategy::EditMass {
            anyhow::bail!("--cluster-min-mass is only valid with --cluster-strategy edit-mass");
        }
        Ok(())
    }

    /// Effective minimum number of events, resolving the per-strategy default when unset.
    fn min_records(&self) -> usize {
        self.min_records.unwrap_or(match self.strategy {
            ClusterStrategy::Legacy => 2,
            ClusterStrategy::EditMass => 1,
        })
    }

    /// Effective minimum edit mass (`edit-mass` strategy only), resolving the default when unset.
    fn min_mass(&self) -> f64 {
        self.min_mass.unwrap_or(2.0)
    }

    /// Mass contribution of a single variant, in the units of the active strategy.
    pub fn event_mass(&self, ref_len: usize, alt_len: usize) -> f64 {
        match self.strategy {
            ClusterStrategy::Legacy => variant_complexity(ref_len, alt_len),
            ClusterStrategy::EditMass => edit_len(ref_len, alt_len),
        }
    }

    /// Cluster span derived from its reference and query extents (both in bases). The `legacy`
    /// strategy ignores the query extent; `edit-mass` takes the larger of the two.
    pub fn span(&self, ref_extent: f64, query_extent: f64) -> f64 {
        match self.strategy {
            ClusterStrategy::Legacy => ref_extent,
            ClusterStrategy::EditMass => ref_extent.max(query_extent),
        }
    }

    /// Returns true if the accumulated cluster statistics meet the configured thresholds.
    ///
    /// * `mass`  — sum of [`Self::event_mass`] over all events in the cluster
    /// * `span`  — cluster span as produced by [`Self::span`], in bases
    /// * `count` — number of distinct mutation events
    pub fn is_valid_cluster(&self, mass: f64, span: f64, count: usize) -> bool {
        if count < 1.max(self.min_records()) {
            return false;
        }
        if self.strategy == ClusterStrategy::EditMass && mass < self.min_mass() {
            return false;
        }
        mass / span.max(1.0) >= self.min_density
    }
}

/// Log-scaled complexity score for a single variant (`legacy` strategy).
///
/// Lengths use the **VCF anchor convention**: a single-base substitution is `(1, 1)`, an insertion
/// of `n` bases is `(1, n + 1)`, and a deletion of `n` bases is `(n + 1, 1)`.
///
/// Scores: SNP → `1.0`; indel of length `n` → `1.0 + log₂(n)`.
#[allow(clippy::cast_precision_loss)]
pub fn variant_complexity(ref_len: usize, alt_len: usize) -> f64 {
    match (ref_len, alt_len) {
        (1, 1) => 1.0,
        (n, m) => 1.0 + ((n.max(m) - 1) as f64).log2(),
    }
}

/// Number of edited bases for a single variant — its contribution to a cluster's edit mass
/// (`edit-mass` strategy).
///
/// Both lengths use the **VCF anchor convention**: the anchor base shared by ref and alt is
/// included, so a single-base substitution is `(1, 1)`, an insertion of `n` bases is
/// `(1, n + 1)`, and a deletion of `n` bases is `(n + 1, 1)`.
///
/// * substitution / MNV (equal lengths) → that length (SNP → `1`, `n`-base MNV → `n`)
/// * indel (unequal lengths) → `max(ref_len, alt_len) - 1` (strip the shared anchor → `n`)
#[allow(clippy::cast_precision_loss)]
pub fn edit_len(ref_len: usize, alt_len: usize) -> f64 {
    if ref_len == alt_len {
        ref_len as f64
    } else {
        (ref_len.max(alt_len) - 1) as f64
    }
}

#[cfg(test)]
#[allow(clippy::float_cmp)]
mod tests {
    use super::*;

    #[test]
    fn snp_edit_len_is_one() {
        assert_eq!(edit_len(1, 1), 1.0);
    }

    #[test]
    fn mnv_edit_len_is_block_length() {
        // 4-base MNV: ref="ACGT", alt="TGCA" → both len 4 → 4 edited bases
        assert_eq!(edit_len(4, 4), 4.0);
    }

    #[test]
    fn insertion_edit_len() {
        // 1-base ins: ref="A", alt="AT"  → ref_len=1, alt_len=2  → 1 edited base
        assert_eq!(edit_len(1, 2), 1.0);
        // 2-base ins: alt_len=3 → 2
        assert_eq!(edit_len(1, 3), 2.0);
        // 4-base ins: alt_len=5 → 4
        assert_eq!(edit_len(1, 5), 4.0);
    }

    #[test]
    fn deletion_edit_len_mirrors_insertion() {
        assert_eq!(edit_len(2, 1), edit_len(1, 2));
        assert_eq!(edit_len(3, 1), edit_len(1, 3));
    }

    #[test]
    fn snp_variant_complexity_is_one() {
        assert_eq!(variant_complexity(1, 1), 1.0);
    }

    #[test]
    fn indel_variant_complexity_is_log_scaled() {
        // 2-base indel: 1 + log2(2) = 2.0
        assert_eq!(variant_complexity(1, 3), 2.0);
        assert_eq!(variant_complexity(3, 1), 2.0);
    }

    #[test]
    fn valid_cluster_at_density_threshold() {
        let s = ClusteringSettings {
            strategy: ClusterStrategy::EditMass,
            min_records: Some(1),
            min_mass: Some(2.0),
            max_gap: 20,
            min_density: 0.2,
        };
        // mass=2.0, span=10.0 → density=0.2, exactly at threshold
        assert!(s.is_valid_cluster(2.0, 10.0, 2));
    }

    #[test]
    fn cluster_rejected_below_density_threshold() {
        let s = ClusteringSettings {
            strategy: ClusterStrategy::EditMass,
            min_records: Some(1),
            min_mass: Some(0.0),
            max_gap: 20,
            min_density: 0.2,
        };
        // mass=1.9, span=10.0 → density=0.19 < 0.2
        assert!(!s.is_valid_cluster(1.9, 10.0, 2));
    }

    #[test]
    fn cluster_rejected_below_min_mass() {
        let s = ClusteringSettings {
            strategy: ClusterStrategy::EditMass,
            min_records: Some(1),
            min_mass: Some(2.0),
            max_gap: 20,
            min_density: 0.0,
        };
        // A lone SNP: mass=1 < min_mass=2 → rejected even though density would pass.
        assert!(!s.is_valid_cluster(1.0, 1.0, 1));
        // A single 5-base deletion: mass=5 over span=6 → accepted on its own.
        assert!(s.is_valid_cluster(5.0, 6.0, 1));
    }

    #[test]
    fn legacy_has_no_min_mass_gate() {
        // Legacy ignores min_mass entirely: a lone event passes on density alone, but the
        // default two-record minimum still applies.
        let s = ClusteringSettings {
            strategy: ClusterStrategy::Legacy,
            min_records: None,
            min_mass: None,
            max_gap: 20,
            min_density: 0.0,
        };
        // Single event: below the legacy default min_records=2 → rejected.
        assert!(!s.is_valid_cluster(1.0, 1.0, 1));
        // Two events with mass 1 (would fail an edit-mass min_mass=2 gate) → accepted under legacy.
        assert!(s.is_valid_cluster(1.0, 1.0, 2));
    }

    #[test]
    fn cluster_rejected_below_min_records() {
        let s = ClusteringSettings {
            strategy: ClusterStrategy::EditMass,
            min_records: Some(3),
            min_mass: Some(0.0),
            max_gap: 20,
            min_density: 0.0,
        };
        assert!(!s.is_valid_cluster(10.0, 5.0, 2));
        assert!(s.is_valid_cluster(10.0, 5.0, 3));
    }

    #[test]
    fn validate_rejects_min_mass_without_edit_mass() {
        let s = ClusteringSettings {
            strategy: ClusterStrategy::Legacy,
            min_mass: Some(2.0),
            ..Default::default()
        };
        assert!(s.validate().is_err());

        let s = ClusteringSettings {
            strategy: ClusterStrategy::EditMass,
            min_mass: Some(2.0),
            ..Default::default()
        };
        assert!(s.validate().is_ok());
    }
}