Skip to main content

chronix_analytics/anomaly/
modified_zscore.rs

1//! Modified Z-Score detector using Median Absolute Deviation (MAD).
2
3use crate::anomaly::error::AnomalyError;
4use crate::anomaly::traits::{
5    AnomalyDetector, AnomalyScore, DetectorType, scale_floor, validate_lengths,
6};
7
8/// MAD-based detector robust to outliers in training data.
9///
10/// Modified Z-Score = 0.6745 × (value − median) / MAD
11#[derive(serde::Serialize, serde::Deserialize)]
12pub struct ModifiedZScoreDetector {
13    threshold: f64,
14    median: f64,
15    mad: f64,
16    /// Mean absolute deviation about the median, the fallback scale when the
17    /// MAD is exactly zero.
18    mean_ad: f64,
19    fitted: bool,
20}
21
22impl ModifiedZScoreDetector {
23    /// Creates a new detector. `threshold` defaults to 3.5.
24    pub fn new(threshold: Option<f64>) -> Self {
25        Self {
26            threshold: threshold.unwrap_or(3.5),
27            median: 0.0,
28            mad: 0.0,
29            mean_ad: 0.0,
30            fitted: false,
31        }
32    }
33
34    fn compute_median(data: &mut [f64]) -> f64 {
35        data.sort_unstable_by(f64::total_cmp);
36        let n = data.len();
37        if n.is_multiple_of(2) {
38            (data[n / 2 - 1] + data[n / 2]) / 2.0
39        } else {
40            data[n / 2]
41        }
42    }
43
44    /// Iglewicz & Hoaglin (1993), *How to Detect and Handle Outliers*, §3.
45    ///
46    /// The primary statistic is `M_i = 0.6745·(x − x̃)/MAD`. On **quantised**
47    /// data — a 0.5 °C thermostat, a percentage rounded to whole numbers, an
48    /// integer counter — more than half the sample can equal the median, and
49    /// then the MAD is *exactly* zero: a real and common case, not a
50    /// degenerate one. Dividing by a `1e-12` floor there gave the value one
51    /// quantisation step away a modified z-score of ~1e7, so the detector
52    /// flagged the ordinary neighbouring bucket as an extreme outlier.
53    ///
54    /// Iglewicz & Hoaglin give the alternative for exactly this case:
55    /// `M_i = (x − x̃)/(1.253314·MeanAD)`, where `MeanAD` is the *mean*
56    /// absolute deviation about the median, whose expectation is `σ/1.253314`
57    /// for normal data — hence the constant, and hence no `0.6745` factor
58    /// here. Only when both dispersion estimates are zero (a genuinely
59    /// constant sample) does the relative floor apply.
60    #[inline]
61    fn score_value(&self, value: f64) -> f64 {
62        let dev = value - self.median;
63        if self.mad > 0.0 {
64            (0.6745 * dev / self.mad).abs()
65        } else if self.mean_ad > 0.0 {
66            (dev / (1.253_314 * self.mean_ad)).abs()
67        } else {
68            (0.6745 * dev / scale_floor(0.0, self.median)).abs()
69        }
70    }
71
72    #[inline]
73    fn normalize(raw: f64, threshold: f64) -> f64 {
74        let x = (raw - threshold) * 2.0;
75        1.0 / (1.0 + (-x).exp())
76    }
77}
78
79impl AnomalyDetector for ModifiedZScoreDetector {
80    fn fit(&mut self, _timestamps: &[i64], values: &[f64]) -> Result<(), AnomalyError> {
81        let _start = std::time::Instant::now();
82        if values.len() < 3 {
83            return Err(AnomalyError::InsufficientData {
84                min: 3,
85                got: values.len(),
86            });
87        }
88        let mut sorted = values.to_vec();
89        self.median = Self::compute_median(&mut sorted);
90        let mut deviations: Vec<f64> = values.iter().map(|&v| (v - self.median).abs()).collect();
91        self.mean_ad = deviations.iter().sum::<f64>() / deviations.len() as f64;
92        self.mad = Self::compute_median(&mut deviations);
93        self.fitted = true;
94        metrics::histogram!("chronix_anomaly_fit_duration_seconds", "method" => "modified_zscore")
95            .record(_start.elapsed().as_secs_f64());
96        Ok(())
97    }
98
99    fn detect(
100        &mut self,
101        timestamps: &[i64],
102        values: &[f64],
103    ) -> Result<Vec<AnomalyScore>, AnomalyError> {
104        let _start = std::time::Instant::now();
105        validate_lengths(timestamps, values)?;
106        if !self.fitted {
107            return Err(AnomalyError::NotFitted);
108        }
109        let mut scores = Vec::with_capacity(values.len());
110        for (i, &v) in values.iter().enumerate() {
111            let raw = self.score_value(v);
112            scores.push(AnomalyScore {
113                timestamp: timestamps[i],
114                value: v,
115                score: Self::normalize(raw, self.threshold),
116                is_anomaly: raw > self.threshold,
117                method: DetectorType::ModifiedZScore,
118                threshold: self.threshold,
119                details: format!("modified_z={raw:.4}"),
120            });
121        }
122        let anomaly_count = scores.iter().filter(|s| s.is_anomaly).count();
123        metrics::counter!("chronix_anomaly_detected_total", "method" => "modified_zscore")
124            .increment(anomaly_count as u64);
125        metrics::histogram!("chronix_anomaly_detect_duration_seconds", "method" => "modified_zscore").record(_start.elapsed().as_secs_f64());
126        Ok(scores)
127    }
128
129    fn detect_point(&mut self, timestamp: i64, value: f64) -> Result<AnomalyScore, AnomalyError> {
130        if !self.fitted {
131            return Err(AnomalyError::NotFitted);
132        }
133        let raw = self.score_value(value);
134        Ok(AnomalyScore {
135            timestamp,
136            value,
137            score: Self::normalize(raw, self.threshold),
138            is_anomaly: raw > self.threshold,
139            method: DetectorType::ModifiedZScore,
140            threshold: self.threshold,
141            details: format!("modified_z={raw:.4}"),
142        })
143    }
144
145    fn detector_type(&self) -> DetectorType {
146        DetectorType::ModifiedZScore
147    }
148}
149
150#[cfg(test)]
151mod tests {
152    use super::*;
153
154    #[test]
155    fn robust_to_outliers_in_training() {
156        // Training data with outliers — MAD should be more robust than std dev
157        let ts: Vec<i64> = (0..100).map(|i| i * 1_000_000_000).collect();
158        let mut vals: Vec<f64> = vec![50.0; 100];
159        // Inject a few heavy outliers in training data
160        vals[0] = 500.0;
161        vals[1] = -400.0;
162        vals[99] = 600.0;
163
164        let mut det = ModifiedZScoreDetector::new(Some(3.5));
165        det.fit(&ts, &vals).unwrap();
166
167        // Median should be near 50.0
168        assert!((det.median - 50.0).abs() < 1.0);
169
170        // A normal value should not be flagged
171        let s = det.detect_point(999, 50.0).unwrap();
172        assert!(!s.is_anomaly);
173
174        // A genuine outlier should be flagged
175        let s = det.detect_point(999, 500.0).unwrap();
176        assert!(s.is_anomaly);
177    }
178
179    #[test]
180    fn detect_outliers() {
181        let ts: Vec<i64> = (0..200).map(|i| i * 1_000_000_000).collect();
182        let mut vals: Vec<f64> = (0..200)
183            .map(|i| (i as f64 * 0.1).sin() * 5.0 + 50.0)
184            .collect();
185        vals[75] = 200.0;
186        vals[175] = -100.0;
187
188        let mut det = ModifiedZScoreDetector::new(None);
189        det.fit(&ts, &vals).unwrap();
190        let scores = det.detect(&ts, &vals).unwrap();
191        assert!(scores[75].is_anomaly);
192        assert!(scores[175].is_anomaly);
193    }
194
195    #[test]
196    fn insufficient_data() {
197        let mut det = ModifiedZScoreDetector::new(None);
198        assert!(det.fit(&[0, 1], &[1.0, 2.0]).is_err());
199    }
200
201    #[test]
202    fn constant_data_does_not_flag_same_values() {
203        let ts: Vec<i64> = (0..50).map(|i| i * 1_000_000_000).collect();
204        let vals = vec![42.0; 50];
205
206        let mut det = ModifiedZScoreDetector::new(Some(3.5));
207        det.fit(&ts, &vals).unwrap();
208
209        // Same value should NOT be flagged
210        let scores = det.detect(&ts, &vals).unwrap();
211        let anomalies: Vec<_> = scores.iter().filter(|s| s.is_anomaly).collect();
212        assert!(
213            anomalies.is_empty(),
214            "Constant data should not flag identical values as anomalies, got {} anomalies",
215            anomalies.len(),
216        );
217    }
218
219    #[test]
220    fn constant_data_flags_genuine_deviation() {
221        let ts: Vec<i64> = (0..50).map(|i| i * 1_000_000_000).collect();
222        let vals = vec![42.0; 50];
223
224        let mut det = ModifiedZScoreDetector::new(Some(3.5));
225        det.fit(&ts, &vals).unwrap();
226
227        // A genuinely different value should be flagged
228        let s = det.detect_point(999, 100.0).unwrap();
229        assert!(
230            s.is_anomaly,
231            "Large deviation from constant median should be anomaly"
232        );
233        assert!(s.score.is_finite(), "Score should be finite, not INFINITY");
234    }
235}