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 {
33            return 0.0;
34        }
35        let n = (end - start) as f64;
36        let mut sum = 0.0;
37        let mut sum_sq = 0.0;
38        for i in start..end {
39            sum += data[i];
40            sum_sq += data[i] * data[i];
41        }
42        let mean = sum / n;
43        let var = (sum_sq / n) - (mean * mean);
44        n * var.max(0.0)
45    }
46
47    /// Detect changepoints in a batch of data.
48    /// Returns indices of changepoints.
49    pub fn detect(&self, data: &[f64]) -> Vec<usize> {
50        let n = data.len();
51        if n < self.min_dist * 2 {
52            return vec![];
53        }
54
55        let mut f = vec![0.0; n + 1];
56        let mut cp = vec![0; n + 1];
57        let mut r = vec![0]; // Potential last changepoints
58
59        f[0] = -self.penalty;
60
61        for t in 1..=n {
62            let mut min_val = f64::MAX;
63            let mut best_tau = 0;
64
65            for &tau in &r {
66                if t - tau < self.min_dist {
67                    continue;
68                }
69                let val = f[tau] + self.cost(data, tau, t) + self.penalty;
70                if val < min_val {
71                    min_val = val;
72                    best_tau = tau;
73                }
74            }
75
76            f[t] = min_val;
77            cp[t] = best_tau;
78
79            // Pruning step
80            let mut next_r = vec![0];
81            for &tau in &r {
82                if f[tau] + self.cost(data, tau, t) <= f[t] {
83                    next_r.push(tau);
84                }
85            }
86            next_r.push(t);
87            r = next_r;
88        }
89
90        // Backtrack to find changepoints
91        let mut changepoints = Vec::new();
92        let mut curr = cp[n];
93        while curr > 0 {
94            changepoints.push(curr);
95            curr = cp[curr];
96        }
97        changepoints.sort();
98        changepoints
99    }
100}
quantwave_core/regimes/pelt.rs

quantwave_core/regimes/
pelt.rs
