qualia_core_db/solvers/statistics/
timeseries.rs

1//! Time-series kernels — autocorrelation, moving average, exponential smoothing.
2//!
3//! The canonical home for the elementary time-series transforms used by the
4//! domain libraries (`specialized_libs::statistical_computing`). Like the rest
5//! of `solvers::statistics`, these operate over caller-owned slices; the
6//! series-producing transforms write into a caller-provided `out` slice and
7//! return the number of elements written (mirroring `histogram_into`), so no
8//! allocation is imposed by this layer.
9
10use super::descriptive::mean;
11
12/// Sample autocorrelation at `lag`, using the standard biased estimator
13/// (normalised by the total sum of squares, mean-centred):
14///
15/// ```text
16/// r_k = Σ_{t=k}^{n-1} (x_t − x̄)(x_{t−k} − x̄)  /  Σ_{t=0}^{n-1} (x_t − x̄)²
17/// ```
18///
19/// `r_0` is always `1.0` for non-constant data. Returns `None` if the slice is
20/// empty, `lag >= n`, or the series is constant (zero variance → undefined).
21pub fn autocorrelation(values: &[f64], lag: usize) -> Option<f64> {
22    let n = values.len();
23    if n == 0 || lag >= n {
24        return None;
25    }
26    let m = mean(values)?;
27    let mut denom = 0.0;
28    let mut i = 0;
29    while i < n {
30        let d = values[i] - m;
31        denom += d * d;
32        i += 1;
33    }
34    if denom == 0.0 {
35        return None; // constant series — autocorrelation undefined
36    }
37    let mut num = 0.0;
38    let mut t = lag;
39    while t < n {
40        num += (values[t] - m) * (values[t - lag] - m);
41        t += 1;
42    }
43    Some(num / denom)
44}
45
46/// Simple moving average with the given `window`, written into `out`.
47///
48/// Produces `n − window + 1` values, where `out[i]` is the mean of
49/// `values[i..i+window]`. Returns the number of values written, or `None` if
50/// `window == 0`, `window > n`, or `out` is too small to hold the result.
51/// Uses a running-sum sweep (O(n), not O(n·window)).
52pub fn moving_average_into(values: &[f64], window: usize, out: &mut [f64]) -> Option<usize> {
53    let n = values.len();
54    if window == 0 || window > n {
55        return None;
56    }
57    let count = n - window + 1;
58    if out.len() < count {
59        return None;
60    }
61    let mut acc = 0.0;
62    let mut i = 0;
63    while i < window {
64        acc += values[i];
65        i += 1;
66    }
67    let w = window as f64;
68    out[0] = acc / w;
69    let mut j = 1;
70    while j < count {
71        acc += values[j + window - 1] - values[j - 1];
72        out[j] = acc / w;
73        j += 1;
74    }
75    Some(count)
76}
77
78/// Single (Brown's) exponential smoothing with factor `alpha ∈ (0, 1]`,
79/// written into `out` (same length as `values`).
80///
81/// `s_0 = x_0`; `s_t = alpha·x_t + (1 − alpha)·s_{t−1}`. Returns the number of
82/// values written, or `None` if the series is empty, `out` is too small, or
83/// `alpha` is not in `(0, 1]`.
84pub fn exponential_smoothing_into(values: &[f64], alpha: f64, out: &mut [f64]) -> Option<usize> {
85    let n = values.len();
86    if n == 0 || out.len() < n || !(alpha > 0.0 && alpha <= 1.0) {
87        return None;
88    }
89    out[0] = values[0];
90    let mut t = 1;
91    while t < n {
92        out[t] = alpha * values[t] + (1.0 - alpha) * out[t - 1];
93        t += 1;
94    }
95    Some(n)
96}
97
98#[cfg(test)]
99mod tests {
100    use super::*;
101
102    const EPS: f64 = 1e-12;
103
104    #[test]
105    fn autocorr_lag0_is_one() {
106        let v = [1.0, 2.0, 3.0, 4.0, 5.0];
107        assert!((autocorrelation(&v, 0).unwrap() - 1.0).abs() < EPS);
108    }
109
110    #[test]
111    fn autocorr_lag1_known_value() {
112        // mean=3; deviations [-2,-1,0,1,2]; num = 2+0+0+2 = 4; denom = 10 → 0.4
113        let v = [1.0, 2.0, 3.0, 4.0, 5.0];
114        assert!((autocorrelation(&v, 1).unwrap() - 0.4).abs() < EPS);
115    }
116
117    #[test]
118    fn autocorr_constant_is_none() {
119        let v = [7.0, 7.0, 7.0];
120        assert_eq!(autocorrelation(&v, 1), None);
121    }
122
123    #[test]
124    fn moving_average_window2() {
125        let v = [1.0, 2.0, 3.0, 4.0];
126        let mut out = [0.0; 3];
127        assert_eq!(moving_average_into(&v, 2, &mut out), Some(3));
128        assert!((out[0] - 1.5).abs() < EPS);
129        assert!((out[1] - 2.5).abs() < EPS);
130        assert!((out[2] - 3.5).abs() < EPS);
131    }
132
133    #[test]
134    fn moving_average_rejects_bad_window() {
135        let v = [1.0, 2.0];
136        let mut out = [0.0; 2];
137        assert_eq!(moving_average_into(&v, 0, &mut out), None);
138        assert_eq!(moving_average_into(&v, 3, &mut out), None);
139    }
140
141    #[test]
142    fn exponential_smoothing_known_value() {
143        // alpha=0.5: s0=1, s1=0.5*2+0.5*1=1.5, s2=0.5*3+0.5*1.5=2.25
144        let v = [1.0, 2.0, 3.0];
145        let mut out = [0.0; 3];
146        assert_eq!(exponential_smoothing_into(&v, 0.5, &mut out), Some(3));
147        assert!((out[0] - 1.0).abs() < EPS);
148        assert!((out[1] - 1.5).abs() < EPS);
149        assert!((out[2] - 2.25).abs() < EPS);
150    }
151
152    #[test]
153    fn exponential_smoothing_rejects_bad_alpha() {
154        let v = [1.0, 2.0];
155        let mut out = [0.0; 2];
156        assert_eq!(exponential_smoothing_into(&v, 0.0, &mut out), None);
157        assert_eq!(exponential_smoothing_into(&v, 1.5, &mut out), None);
158    }
159}
qualia_core_db/solvers/statistics/timeseries.rs

qualia_core_db/solvers/statistics/
timeseries.rs
