quantwave_core/regimes/
pelt.rs

1//! Changepoint Detection (Killick et al. 2012)
2//!
3//! Source: Killick, R., Fearnhead, P., & Eckley, I. A. (2012). 
4//! "Optimal Detection of Changepoints with a Linear Computational Cost." 
5//! Journal of the American Statistical Association, 107(500), 1590-1598.
6//!
7//! Implementation of the Pruned Exact Linear Time (PELT) algorithm for exact segmentation.
8//! PELT identifies change points by minimizing a cost function over all possible partitions.
9
10use serde::{Deserialize, Serialize};
11
12/// A PELT changepoint detector.
13#[derive(Debug, Clone, Serialize, Deserialize)]
14pub struct PELT {
15    penalty: f64,
16    min_dist: usize,
17}
18
19impl PELT {
20    /// Creates a new PELT detector.
21    ///
22    /// # Arguments
23    /// * `penalty` - The penalty (beta) for adding a new changepoint (e.g., ln(n)).
24    /// * `min_dist` - Minimum distance between changepoints.
25    pub fn new(penalty: f64, min_dist: usize) -> Self {
26        Self { penalty, min_dist }
27    }
28
29    /// Normal log-likelihood cost function for change in mean.
30    /// C(y_s:t) = (t-s) * var(y_s:t)
31    fn cost(&self, data: &[f64], start: usize, end: usize) -> f64 {
32        if end <= start { return 0.0; }
33        let n = (end - start) as f64;
34        let mut sum = 0.0;
35        let mut sum_sq = 0.0;
36        for i in start..end {
37            sum += data[i];
38            sum_sq += data[i] * data[i];
39        }
40        let mean = sum / n;
41        let var = (sum_sq / n) - (mean * mean);
42        n * var.max(0.0)
43    }
44
45    /// Detect changepoints in a batch of data.
46    /// Returns indices of changepoints.
47    pub fn detect(&self, data: &[f64]) -> Vec<usize> {
48        let n = data.len();
49        if n < self.min_dist * 2 { return vec![]; }
50
51        let mut f = vec![0.0; n + 1];
52        let mut cp = vec![0; n + 1];
53        let mut r = vec![0]; // Potential last changepoints
54
55        f[0] = -self.penalty;
56
57        for t in 1..=n {
58            let mut min_val = f64::MAX;
59            let mut best_tau = 0;
60
61            for &tau in &r {
62                if t - tau < self.min_dist { continue; }
63                let val = f[tau] + self.cost(data, tau, t) + self.penalty;
64                if val < min_val {
65                    min_val = val;
66                    best_tau = tau;
67                }
68            }
69
70            f[t] = min_val;
71            cp[t] = best_tau;
72
73            // Pruning step
74            let mut next_r = vec![0];
75            for &tau in &r {
76                if f[tau] + self.cost(data, tau, t) <= f[t] {
77                    next_r.push(tau);
78                }
79            }
80            next_r.push(t);
81            r = next_r;
82        }
83
84        // Backtrack to find changepoints
85        let mut changepoints = Vec::new();
86        let mut curr = cp[n];
87        while curr > 0 {
88            changepoints.push(curr);
89            curr = cp[curr];
90        }
91        changepoints.sort();
92        changepoints
93    }
94}
quantwave_core/regimes/pelt.rs

quantwave_core/regimes/
pelt.rs
