nodedb_query/ts_functions/
stddev.rs

1// SPDX-License-Identifier: Apache-2.0
2
3//! Online standard deviation (Welford's algorithm).
4
5/// Online standard deviation accumulator.
6///
7/// Single-pass, numerically stable Welford algorithm.
8#[derive(Debug)]
9pub struct TsStddevAccum {
10    n: u64,
11    mean: f64,
12    m2: f64,
13}
14
15impl TsStddevAccum {
16    pub fn new() -> Self {
17        Self {
18            n: 0,
19            mean: 0.0,
20            m2: 0.0,
21        }
22    }
23
24    pub fn update(&mut self, value: f64) {
25        self.n += 1;
26        let delta = value - self.mean;
27        self.mean += delta / self.n as f64;
28        let delta2 = value - self.mean;
29        self.m2 += delta * delta2;
30    }
31
32    pub fn update_batch(&mut self, values: &[f64]) {
33        for &v in values {
34            if !v.is_nan() {
35                self.update(v);
36            }
37        }
38    }
39
40    /// Population standard deviation (`σ`). Returns `None` if n < 2.
41    pub fn evaluate_population(&self) -> Option<f64> {
42        if self.n < 2 {
43            return None;
44        }
45        Some((self.m2 / self.n as f64).sqrt())
46    }
47
48    /// Sample standard deviation (`s`). Returns `None` if n < 2.
49    pub fn evaluate_sample(&self) -> Option<f64> {
50        if self.n < 2 {
51            return None;
52        }
53        Some((self.m2 / (self.n - 1) as f64).sqrt())
54    }
55
56    /// Serialise state for partial-aggregate merging.
57    pub fn state(&self) -> [f64; 3] {
58        [self.n as f64, self.mean, self.m2]
59    }
60
61    /// Merge another accumulator using Chan's parallel algorithm.
62    pub fn merge(&mut self, other: &Self) {
63        if other.n == 0 {
64            return;
65        }
66        if self.n == 0 {
67            self.n = other.n;
68            self.mean = other.mean;
69            self.m2 = other.m2;
70            return;
71        }
72        let n_a = self.n as f64;
73        let n_b = other.n as f64;
74        let n_ab = n_a + n_b;
75        let delta = other.mean - self.mean;
76
77        self.m2 += other.m2 + delta * delta * n_a * n_b / n_ab;
78        self.mean = (n_a * self.mean + n_b * other.mean) / n_ab;
79        self.n += other.n;
80    }
81
82    /// Merge from serialised state.
83    pub fn merge_state(&mut self, state: &[f64; 3]) {
84        let other = Self {
85            n: state[0] as u64,
86            mean: state[1],
87            m2: state[2],
88        };
89        self.merge(&other);
90    }
91
92    pub fn size(&self) -> usize {
93        std::mem::size_of::<Self>()
94    }
95}
96
97impl Default for TsStddevAccum {
98    fn default() -> Self {
99        Self::new()
100    }
101}
102
103/// Batch population standard deviation. Skips NaN values.
104pub fn ts_stddev(values: &[f64]) -> Option<f64> {
105    let mut accum = TsStddevAccum::new();
106    accum.update_batch(values);
107    accum.evaluate_population()
108}
109
110#[cfg(test)]
111mod tests {
112    use super::*;
113
114    #[test]
115    fn known_values() {
116        // population stddev of [2, 4, 4, 4, 5, 5, 7, 9] = 2.0
117        let vals = [2.0, 4.0, 4.0, 4.0, 5.0, 5.0, 7.0, 9.0];
118        assert!((ts_stddev(&vals).unwrap() - 2.0).abs() < 1e-10);
119    }
120
121    #[test]
122    fn constant_values() {
123        assert!((ts_stddev(&[5.0, 5.0, 5.0, 5.0]).unwrap()).abs() < 1e-12);
124    }
125
126    #[test]
127    fn sample_vs_population() {
128        let mut acc = TsStddevAccum::new();
129        acc.update_batch(&[2.0, 4.0]);
130        let pop = acc.evaluate_population().unwrap();
131        let samp = acc.evaluate_sample().unwrap();
132        assert!(samp > pop); // sample stddev > population stddev for small n
133    }
134
135    #[test]
136    fn too_few() {
137        assert!(ts_stddev(&[1.0]).is_none());
138        assert!(ts_stddev(&[]).is_none());
139    }
140
141    #[test]
142    fn nan_skipped() {
143        let r = ts_stddev(&[2.0, f64::NAN, 4.0, 4.0, 4.0, 5.0, 5.0, 7.0, 9.0]);
144        assert!(r.is_some());
145    }
146
147    #[test]
148    fn merge_preserves_result() {
149        let mut a = TsStddevAccum::new();
150        a.update_batch(&[2.0, 4.0, 4.0, 4.0]);
151        let mut b = TsStddevAccum::new();
152        b.update_batch(&[5.0, 5.0, 7.0, 9.0]);
153        a.merge(&b);
154        assert!((a.evaluate_population().unwrap() - 2.0).abs() < 1e-10);
155    }
156}
nodedb_query/ts_functions/stddev.rs

nodedb_query/ts_functions/
stddev.rs
