Anomaly and change-point detection in time series

Before this lesson: decomposition and the remainder (lesson 01), forecasts and residuals (lesson 03), the Kalman filter and its innovations on the Nile data (lesson 07), and prediction intervals (lesson 10). From the probability-and-statistics track: the normal distribution, quantiles and the idea of a hypothesis test; from the ml-fundamentals track: precision and recall. Related: drift of feature distributions in mlops lesson 07. Time: about 90 minutes, plus the notebook.

---

Why this matters #

• Operations. Server latency, error rates, payments and sensor readings are monitored by detectors that run every minute. Twitter published a detector based on seasonal decomposition and a robust test (S-H-ESD; Hochenbaum, Vallis and Kejariwal, 2017). Every false alarm costs an engineer's time; every missed event costs money.
• ML systems. A model in production sees inputs whose distribution changes. Monitoring its inputs, outputs and errors over time is anomaly and change detection (mlops lesson 07, mlops lesson 08).
• Industry and science. Sensor faults and machine failures (predictive maintenance), unusual heartbeats in ECG recordings, breaks in climate station records after an instrument change, and changes of regime in river flows are all anomalies or change points. The worked example finds the 1899 drop of the Nile flow.
• Forecasting. Anomalies distort model fitting and evaluation. Detecting them and treating them as missing values (the Kalman filter skips the update) makes forecasts more robust.

1. Predict first #

2. The picture #

Four kinds of anomalies #

Kind	What it is	Example
point anomaly	one value far from all normal values	a spike of latency to 10 s
contextual anomaly	a value that is normal in general, but not at this time	high heating demand in July
collective anomaly	a group of values that is unusual together, though each value is normal	a missing heartbeat in an ECG; a day without the daily cycle
change point	a lasting change of level, variance, trend or season	a new dam lowers the river flow from then on

The first three are short events. A change point is different: after it, the new behaviour is normal. A detector that treats a change point as a series of outliers raises alarms for ever (or until the model adapts). The classical survey is Chandola, Banerjee and Kumar (2009); for change points see Truong, Oudre and Vayatis (2020).

Expected value, residual, threshold #

The figure shows a series with a daily cycle (24 hours) and three injected anomalies: a point anomaly at hour 30, a contextual anomaly at hours 137–139 (a peak-level value at night) and a collective anomaly at hours 190–201 (half a day without the cycle).

Interactive figure: loading…

This answers question 1 of section 1: a global quantile rule misses the July value, because 300 MW is a normal value in winter. The expected value must know the season.

A robust scale #

To turn residuals into scores, we divide by a scale. The mean and the standard deviation are a bad choice when anomalies are present: each outlier pulls the mean towards itself and makes the standard deviation larger, so the band "mean ± 3 sd" grows and the outliers hide inside it. This is masking. The median and the median absolute deviation (MAD) hardly move.

Interactive figure: loading…

3. Definitions and main results #

3.1 Detectors, scores and thresholds #

A detector computes an anomaly score $s_t$ for each time $t$ (or each window) from the data up to $t$ (online) or from the whole series (offline), and raises an alarm when $s_t > \lambda$ for a threshold $\lambda$. Online detection may only use the past: a two-sided smoother or a centred moving average (lesson 01) uses future values and cannot run in real time.

A change point of a series is a time $\tau$ at which its distribution changes, for example $y_t \sim \mathcal N(\mu_0, \sigma^2)$ for $t < \tau$ and $\mathcal N(\mu_1, \sigma^2)$ for $t \ge \tau$. Online change detection wants an alarm soon after $\tau$ and few alarms before it.

3.2 The residual recipe and robust scores #

1. Expected value $\hat y_t$. Options: the trend plus season of a decomposition (STL, lesson 01); a one-step forecast $\hat y_{t\mid t-1}$ from exponential smoothing, ARIMA or an ML model (lessons 04–09); the predicted observation of a Kalman filter (lesson 07).
2. Residual $r_t = y_t - \hat y_t$.
3. Robust score $z_t = (r_t - \operatorname{med}(r)) / (1.4826\cdot\operatorname{MAD}(r))$, where $\operatorname{MAD}(r) = \operatorname{med}_i |r_i - \operatorname{med}(r)|$. Use $|z_t|$ for anomalies in both directions.
4. Alarm when $|z_t| > \lambda$ (section 3.3 chooses $\lambda$).

With a forecaster, the score of a Kalman filter has a direct form: the standardised innovation $\nu_t / \sqrt{S_t}$ of lesson 07, which is $\mathcal N(0, 1)$ when the model is right. A prediction interval (lesson 10) is the same idea: alarm when $y_t$ falls outside the 99 % interval.

Why the factor 1.4826. For $X \sim \mathcal N(\mu, \sigma^2)$, the median of $|X - \mu|$ is the number $c$ with $P(|X - \mu| \le c) = 0.5$. Since $P(|X - \mu| \le c) = 2\Phi(c/\sigma) - 1$, we need $\Phi(c/\sigma) = 0.75$, so $c = \sigma\,\Phi^{-1}(0.75) = 0.6745\,\sigma$. Hence $\operatorname{MAD}/0.6745 = 1.4826\cdot\operatorname{MAD}$ estimates $\sigma$ for Gaussian data, and the robust score is on the same scale as an ordinary z-score.

Why robust. The breakdown point of an estimator is the share of data that can be moved to infinity before the estimate also goes to infinity. It is 0 for the mean and the standard deviation (one bad value is enough) and 50 % for the median and the MAD. Two failures follow from a non-robust scale: masking (outliers inflate the standard deviation and are no longer flagged) and swamping (outliers shift the mean, so normal values on the other side get flagged).

Practical points.

• Use the past only. For online detection, the expected value at $t$ must come from data before $t$: a one-step forecast or the filter, never the smoother of lesson 07 (it has seen $y_t$ and later values, and it partly follows the anomaly).
• Echoes. A forecaster that uses the anomalous value later (the seasonal naive forecast, an AR model) makes new residuals when it does. Replace a detected anomaly by its expected value (or treat it as missing in a Kalman filter) before you forecast on.
• Local scale. If the noise level changes with the season or the level, compute the median and MAD in a rolling window, or use the standardised innovations of a model whose variance changes.

3.3 Thresholds and the false-alarm budget #

Each value that is checked is a small hypothesis test of "this value is normal". With many tests, false alarms add up. If a detector checks $N$ normal values with a false-alarm probability $p$ each, the expected number of false alarms is $Np$. With one value per minute, $N = 1440$ per day, and $p = 0.01$ gives $1440 \cdot 0.01 = 14.4$ false alarms per day. This answers question 2 of section 1. The chance of at least one false alarm in a day is $1 - 0.99^{1440} \approx 1$.

Choose $\lambda$ from a budget, not from habit. Say you accept one false alarm per week from a detector that checks one value per minute: $N = 7 \cdot 1440 = 10\,080$ values, so $p = 1/10\,080 \approx 10^{-4}$. For a Gaussian score and a two-sided rule, $\lambda = \Phi^{-1}(1 - p/2) \approx 3.9$. Real residuals have heavier tails than the Gaussian, so it is safer to set $\lambda$ as an upper quantile of scores from a clean past period. With $n$ past scores, the $\lceil (n + 1)(1 - p)\rceil$-th smallest one gives a false-alarm rate of at most $p$ when the scores are exchangeable: this is split conformal prediction from lesson 10.

The trade-off. A higher threshold gives fewer false alarms but later or missed detections. Report both: false alarms per day (or the average time between false alarms) and the detection delay or recall.

3.4 Change points and the CUSUM #

A threshold on single values is weak for small lasting shifts. If the mean moves up by one standard deviation, a single value exceeds $3\sigma$ only with probability $P(Z > 2) = 0.023$, so the one-sided rule "$z_t > 3$" needs about $1/0.023 \approx 44$ steps on average to notice. The CUSUM (cumulative sum; Page, 1954) adds up the small evidence of many values.

Definition (upper CUSUM for standardised data). Let $z_t = (y_t - \mu_0)/\sigma$ with the in-control mean $\mu_0$ and standard deviation $\sigma$. With a reference value $k \ge 0$ and a decision threshold $h_c > 0$:

For a downward shift, use $z_t \to -z_t$: $C_t = \max(0, C_{t-1} - z_t - k)$. After an alarm, restart at $C_t = 0$. We write $h_c$ so that it is not confused with the forecast horizon $h$.

How it works. Before the change, $\mathbb E[z_t - k] = -k < 0$: the sum drifts down and the $\max(0, \cdot)$ keeps it near 0. After a shift of $\delta$ standard deviations, $\mathbb E[z_t - k] = \delta - k > 0$: the sum grows like a ramp with slope $\delta - k$ and crosses $h_c$ after about $h_c / (\delta - k)$ steps. The usual choice is $k = \delta/2$, half the shift that you want to detect.

Why $k = \delta/2$: the CUSUM is a likelihood-ratio test.

1. Compare "no change", $z_t \sim \mathcal N(0, 1)$, with "shift $\delta$", $z_t \sim \mathcal N(\delta, 1)$. The log likelihood ratio of one value is $\log\frac{\varphi(z_t - \delta)}{\varphi(z_t)} = \frac{-(z_t - \delta)^2 + z_t^2}{2} = \delta z_t - \frac{\delta^2}{2} = \delta\,(z_t - \delta/2)$.
2. For a change at time $j$, the log likelihood ratio of the values $z_j, \dots, z_t$ is the sum $\delta\sum_{i=j}^{t}(z_i - \delta/2)$.
3. We do not know $j$, so we take the largest of these sums over $j \le t$ (or 0 for "no change yet"). This maximum satisfies the recursion $M_t = \max(0, M_{t-1} + \delta(z_t - \delta/2))$: either the best start is before $t$ (then we add the new term to $M_{t-1}$) or there is no evidence (0).
4. Dividing by $\delta > 0$ gives the CUSUM with $k = \delta/2$. $\blacksquare$

The CUSUM is also optimal in a minimax sense: among all detectors with the same false-alarm rate, it has the smallest worst-case detection delay for a shift of exactly $\delta$ (Lorden, 1971; Moustakides, 1986).

Average run length. Detectors of changes are compared by the average run length (ARL): the mean number of steps until an alarm. $\mathrm{ARL}_0$ is the ARL with no change (the mean time between false alarms; large is good). $\mathrm{ARL}_1$ is the ARL after a shift of a given size (the mean detection delay; small is good). For a shift of one standard deviation:

Detector	$\mathrm{ARL}_0$ (no change)	$\mathrm{ARL}_1$ ($\delta = 1$)
single value $z_t > 3$ (one-sided)	$1/P(Z > 3) \approx 741$	$1/P(Z > 2) \approx 44$
CUSUM, $k = 0.5$, $h_c = 4$	$\approx 336$	$\approx 8.4$
CUSUM, $k = 0.5$, $h_c = 5$	$\approx 930$	$\approx 10.4$

With $h_c = 5$, the CUSUM has fewer false alarms than the 3-sigma rule and finds the shift four times faster. (The CUSUM values come from a simulation of the one-sided CUSUM. Montgomery's Introduction to Statistical Quality Control gives tables for the two-sided version, which watches both directions and therefore has about half the $\mathrm{ARL}_0$.)

Interactive figure: loading…

When did the change start? At an alarm at time $t_a$, the last time $\hat\tau - 1$ before $t_a$ at which $C_t = 0$ is a natural estimate of the change time: the evidence has been growing since $\hat\tau$.

Offline segmentation. When the whole series is available, change points can be placed by minimising a cost: the sum over segments of a within-segment cost (for example the squared deviations from the segment mean) plus a penalty $\beta$ per change point. Binary segmentation finds the best single split, then splits the parts again; it is fast but approximate. PELT (Killick, Fearnhead and Eckley, 2012) finds the exact optimum by dynamic programming and prunes candidates, with a linear expected cost under mild conditions. The Python library ruptures implements both (Truong, Oudre and Vayatis, 2020).

Bayesian online change-point detection

Adams and MacKay (2007) track the run length $r_t$, the time since the last change point. A hazard function $H$ gives the prior probability of a change at each step. The posterior $p(r_t \mid y_{1:t})$ is updated by a recursion: either the run grows ($r_t = r_{t-1} + 1$, weight $1 - H$ times the predictive density of $y_t$ given the current run) or a new run starts ($r_t = 0$, weight $H$). With conjugate models (for example a normal model with unknown mean), the predictive densities have closed forms, and the method gives the probability of a recent change at every step. It is popular for streams with several changes of mean or variance.

3.5 Windows, shapes and machine learning #

Collective anomalies are unusual shapes, not unusual values. Three families of methods look at windows $w_i = (y_i, \dots, y_{i+m-1})$ of length $m$:

• Features plus an outlier detector. Describe each window by features (mean, standard deviation, slope, seasonal strength, autocorrelations) and score the windows with an isolation forest (Liu, Ting and Zhou, 2008). It grows random trees that split on random features at random values; anomalies are isolated after few splits, so a short average path length means a high anomaly score.
• Distance to the nearest neighbour. The matrix profile (Yeh et al., 2016) stores, for every window, the z-normalised Euclidean distance to its most similar window elsewhere in the series (not overlapping it). A window whose nearest neighbour is far away has no similar shape anywhere: it is a discord, the classic collective anomaly. The z-normalisation (subtract the window mean, divide by its standard deviation) compares shapes, not levels. The distance of two z-normalised windows is $\sqrt{2m(1 - \rho)}$, with $\rho$ their correlation, so it lies between 0 and $2\sqrt m$.
• Neural models. Train a forecaster or an autoencoder on normal data and use its forecast or reconstruction error as the score (lesson 09). This is the residual recipe with a learned expected value.

Interactive figure: loading…

3.6 Evaluating detectors honestly #

Labels in time series are usually events: segments $[a, b]$ of anomalous time steps. Useful measures:

• Event recall: the share of true events with at least one alarm inside (or within a tolerance after the start).
• Precision of alarms: the share of alarms that fall inside a true event (or near one).
• Detection delay: the time from the start of an event to its first alarm.
• False alarms per day (or per 1000 points) on normal data.

The point-adjust trap. Many papers compute a point-wise F1 score after point adjustment: if at least one point of a true segment is flagged, all points of that segment count as detected. Kim et al. (2022) showed that this makes a detector with random scores reach a high F1. With one segment of 200 points out of 1000 and a random 5 % alarm rate, almost surely at least one alarm falls inside the segment, so all 200 points count as true positives while only about 40 false alarms remain: F1 $\approx 0.9$. This answers question 3 of section 1.

Interactive figure: loading…

Flawed benchmarks. Wu and Keogh (2021) showed that several popular benchmark data sets have problems: many anomalies can be found with one line of code (they are trivial), some labels are wrong, anomalies are very dense, and many anomalies sit near the end of the series. Results on such data say little.

A checklist.

1. Report event-level recall, precision of alarms, delay and false alarms per day; never only point-adjusted F1.
2. Always include a simple baseline: a robust z-score on the residuals of a seasonal decomposition or a seasonal naive forecast, and a random detector with the same alarm rate.
3. Choose thresholds on training or calibration data, not on the test labels.
4. Use benchmarks with few, clearly labelled anomalies per series, such as the UCR anomaly archive (Wu and Keogh, 2021).

4. Worked example #

The annual flow of the Nile at Aswan, 1871–1970 (100 values, the series of lesson 07; data: Durbin and Koopman, 2012, via statsmodels). Two questions: was there a lasting change, and which single years are anomalous?

(a) The level shift, found by CUSUM. We look for a drop in level.

1. Reference period. Take the first 20 years, 1871–1890, as "normal": mean $\mu_0 = 1070.85$, standard deviation $\sigma = 143.86$. Standardise all years: $z_t = (y_t - 1070.85)/143.86$.
2. Lower CUSUM with $k = 0.5$ and $h_c = 4$: $C_t = \max(0, C_{t-1} - z_t - 0.5)$. Up to 1898, $C_t$ goes back to 0 again and again (it is 0 in each year from 1895 to 1898).
3. 1899: $y = 774$, $z = (774 - 1070.85)/143.86 = -2.063$, so $C = \max(0, 0 + 2.063 - 0.5) = 1.563$.
4. 1900: $y = 840$, $z = -1.604$, so $C = 1.563 + 1.604 - 0.5 = 2.667$.
5. 1901: $y = 874$, $z = -1.368$, so $C = 2.667 + 1.368 - 0.5 = 3.535$.
6. 1902: $y = 694$, $z = -2.621$, so $C = 3.535 + 2.621 - 0.5 = 5.656 > 4$: alarm in 1902.
7. Change time. The last year with $C = 0$ before the alarm is 1898, so the estimated start of the change is 1899, with a delay of 3 years. The mean flow is 1097.75 for 1871–1898 and 849.97 for 1899–1970.

This agrees with the history: the first Aswan dam was built from 1898 to 1902, and lesson 07 saw the filtered level drop at the same time. Note that no single year here is extreme: the standardised innovation of the local level filter in 1899 is $-2.50$, below a typical threshold of 3. The CUSUM finds the change because it adds up four moderately low years.

(b) The 1913 low flow, a point anomaly. The smallest value of the series is 456, in 1913.

• Against the regime after the change (1899–1970): median $842.5$, $1.4826\cdot\mathrm{MAD} = 114.9$, so the robust score is $(456 - 842.5)/114.9 = -3.36$.
• Against the whole series (median 893.5, robust scale 179.4) it is only $-2.44$: the level shift inflates the scale. Compare regimes, not the whole mixed series.
• The local level filter of lesson 07 (with its maximum-likelihood variances $Q = 1469.1$, $R = 15099$) gives a standardised innovation of $-2.79$ in 1913, the largest in absolute value of all years (next: $+2.57$ in 1916 and $-2.50$ in 1899).

So 1913 is the clearest point anomaly, and 1899 starts a change point. A detector with $|z| > 3$ on the filter innovations would flag neither: on this short and noisy series, the change needs CUSUM, and the point anomaly needs the right regime.

5. In code #

NumPy checks the masking example and the worked example; the notebook adds the seasonal detectors, the average run length and the point-adjust experiment.

import numpy as np

NILE = np.array([1120, 1160, 963, 1210, 1160, 1160, 813, 1230, 1370, 1140, 995, 935, 1110, 994, 1020, 960, 1180, 799, 958, 1140,
                 1100, 1210, 1150, 1250, 1260, 1220, 1030, 1100, 774, 840, 874, 694, 940, 833, 701, 916, 692, 1020, 1050, 969,
                 831, 726, 456, 824, 702, 1120, 1100, 832, 764, 821, 768, 845, 864, 862, 698, 845, 744, 796, 1040, 759,
                 781, 865, 845, 944, 984, 897, 822, 1010, 771, 676, 649, 846, 812, 742, 801, 1040, 860, 874, 848, 890,
                 744, 749, 838, 1050, 918, 986, 797, 923, 975, 815, 1020, 906, 901, 1170, 912, 746, 919, 718, 714, 740], float)
years = np.arange(1871, 1971)

def robust_z(x):
    med = np.median(x)
    return (x - med) / (1.4826 * np.median(np.abs(x - med)))

# (a) masking: the classical z-score of the outlier is small, the robust one is large
x = np.array([1., 2., 2., 3., 100.])
print(np.median(x), np.median(np.abs(x - np.median(x))))   # 2.0 1.0
print(np.round((x - x.mean()) / x.std(ddof=1), 2))         # [-0.47 -0.45 -0.45 -0.42  1.79]
print(np.round(robust_z(x), 2))                            # [-0.67  0.    0.    0.67 66.1 ]

# (b) the 1913 low flow, against the regime after 1899
after = NILE[years >= 1899]
z1913 = (456 - np.median(after)) / (1.4826 * np.median(np.abs(after - np.median(after))))
print(round(z1913, 2))                                     # -3.36

# (c) lower CUSUM for a drop in level, reference period 1871-1890
ref = NILE[:20]
z = (NILE - ref.mean()) / ref.std(ddof=1)
k, h_c, C = 0.5, 4.0, 0.0
for t in range(len(z)):
    C = max(0.0, C - z[t] - k)        # adds up the evidence for a drop
    if C > h_c:
        print("alarm in", years[t], "C =", round(C, 2))    # alarm in 1902 C = 5.66
        break

6. Where it shows up #

• Monitoring of services. Alert rules on latency, error rates and traffic use seasonal baselines and robust thresholds; S-H-ESD (Hochenbaum et al., 2017) combines an STL-like decomposition with the median, the MAD and a generalised ESD test. Forecast intervals of time-series models (lesson 10) are used the same way.
• ML in production. Drift of input features, predictions and errors is a change-point problem; CUSUM-type rules and two-sample tests over windows are standard (mlops lesson 07).
• Manufacturing and quality control. CUSUM and Shewhart charts (statistical process control) watch the mean and the variance of a production process; the ARL tables come from this field.
• Health and science. The matrix profile finds unusual heartbeats in ECGs and unusual cycles in sensor data; change-point methods find breaks in climate station records (homogenisation) and regime shifts in hydrology, as in the Nile example.
• Finance and payments. Fraud detection combines per-customer expected behaviour (a contextual baseline) with robust scores and models over windows of transactions.

7. Common mistakes #

• A global threshold on seasonal data. It misses contextual anomalies and flags normal peaks. Remove the season first (decomposition, seasonal forecast, filter).
• Mean and standard deviation on contaminated data. They are pulled by the anomalies (masking, swamping). Use the median and $1.4826\cdot\mathrm{MAD}$, or estimate the scale on a clean period.
• Using future data in an online detector. Centred moving averages and smoothers look ahead. Use one-step forecasts or the filter.
• A threshold without a false-alarm budget. $|z| > 3$ checked every minute gives about 4 false alarms a day even for perfect Gaussian data ($1440 \cdot 0.0027 \approx 3.9$).
• Treating a change point as many outliers. After a lasting shift every value looks anomalous until the model adapts. Detect the change, then refit or reset the model.
• Forgetting echoes. A seasonal naive or AR forecast uses the anomaly later and creates false alarms one season later. Replace detected anomalies before forecasting on.
• Reporting point-adjusted F1. It rewards random alarms. Report event recall, alarm precision, delay and false alarms, with a random and a simple baseline.

8. Check yourself #

1. Give one example each of a point, a contextual and a collective anomaly in hourly electricity demand.
   Answer

   Point: a meter error gives 5000 MW for one hour. Contextual: a demand level that is normal at 18:00 appears at 03:00. Collective: on one day the evening peak is missing although each hourly value is in the normal range.
2. Compute the robust z-score of 13 in the data $\{10, 11, 12, 9, 10, 13\}$.
   Answer

   Sorted: $9, 10, 10, 11, 12, 13$; median $(10 + 11)/2 = 10.5$. Absolute deviations: $1.5, 0.5, 0.5, 0.5, 1.5, 2.5$; sorted $0.5, 0.5, 0.5, 1.5, 1.5, 2.5$, MAD $= (0.5 + 1.5)/2 = 1$. Robust z $= (13 - 10.5)/(1.4826 \cdot 1) = 1.69$.
3. You want at most one false alarm per day from a two-sided Gaussian rule that checks one value every 10 minutes. Which threshold $\lambda$ do you need?
   Answer

   $N = 144$ checks per day, so $p = 1/144 \approx 0.0069$. Two-sided: $\lambda = \Phi^{-1}(1 - 0.0069/2) = \Phi^{-1}(0.99653) \approx 2.70$.
4. Run the upper CUSUM with $k = 0.5$, $h_c = 2$ on $z = (0.2, 1.5, 1.2, 0.9, 1.4)$. When is the first alarm?
   Answer

   $C_1 = \max(0, 0.2 - 0.5) = 0$; $C_2 = 1.0$; $C_3 = 1.0 + 0.7 = 1.7$; $C_4 = 1.7 + 0.4 = 2.1 > 2$: alarm at $t = 4$.
5. Why does a z-normalised distance find a collective anomaly that a threshold on the values misses?
   Answer

   Each value of the abnormal window is inside the normal range, so no value-based threshold fires. The z-normalised distance compares the shape of the whole window with all other windows; an abnormal shape has no close neighbour, so its matrix profile value is the largest.
6. A random detector flags 5 % of the points. The data have one anomalous segment of 200 points among 1000. Estimate its point-wise F1 without adjustment.
   Answer

   About $0.05 \cdot 200 = 10$ true positives, $0.05 \cdot 800 = 40$ false positives and $190$ false negatives. Precision $10/50 = 0.2$, recall $10/200 = 0.05$, F1 $= 2 \cdot 0.2 \cdot 0.05/(0.25) = 0.08$.

Recall cardsAnswer in your head first, then reveal. Add the cards to your review deck to see them again just before you would forget.

1. Name the four kinds of anomalies in time series.

   Show answer

2. What are the four steps of a residual-based detector?

   Show answer

3. Why is the MAD multiplied by 1.4826?

   Show answer

4. What is masking?

   Show answer

5. A rule checks one value per minute with a false-alarm rate of 1 % per value. How many false alarms per day?

   Show answer

6. Write the upper CUSUM recursion and its alarm rule.

   Show answer

7. How should you choose the reference value $k$ of a CUSUM, and why?

   Show answer

8. What are $\mathrm{ARL}_0$ and $\mathrm{ARL}_1$?

   Show answer

9. What is a discord in the matrix profile?

   Show answer

10. Why is point-adjusted F1 misleading?

    Show answer

Add these cards to my reviewOpen the review →

9. Interview angle #

• "How would you detect anomalies in a metric with daily and weekly seasonality?" — Model the expected value with the season (STL or a seasonal forecaster, one step ahead), compute robust z-scores of the residuals (median, MAD, maybe in a rolling window), set the threshold from a false-alarm budget on a clean period, and replace flagged values before forecasting on. Add a CUSUM on the residuals for slow shifts.
• "Mean and standard deviation or median and MAD?" — Median and MAD: breakdown point 50 % against 0 %, no masking; multiply the MAD by 1.4826 to get the Gaussian scale.
• "Explain CUSUM." — Running sum of standardised deviations minus a reference $k$, floored at 0; alarm above $h_c$. It is the maximised log-likelihood ratio for a shift $2k$; tune $h_c$ for the wanted $\mathrm{ARL}_0$.
• "A vendor claims F1 = 0.95 on a public benchmark. What do you ask?" — Which protocol (point adjustment?), event-level metrics and delay, false alarms per day, the score of a random detector and of a robust residual baseline on the same data, and whether the threshold was chosen on the test labels.

10. Practice #

The notebook for this lesson runs in Colab or in the browser. It has five exercises, each with a check cell and a folded solution:

1. Write robust_zscore(r) with the median and the MAD, and test it on data with outliers.
2. Build a detector from the seasonal naive forecast and robust z-scores on a series with injected anomalies; measure event recall and precision, and find the echoes.
3. Write cusum(z, k, h) and estimate $\mathrm{ARL}_0$ and the delay for a shift of 1 by simulation, for two values of $k$.
4. Find the Nile level shift with the lower CUSUM, and the 1913 anomaly with the standardised innovations of the local level filter.
5. Write point_adjust_f1 and compare it with the point-wise F1 for a random detector.

11. Watch and read #

Anomaly Detection : Time Series Talkritvikmath · 10 min

Detecting anomalies in a time series and adjusting the series for them, with code.

Robust Anomaly Detection + Seasonal-Trend Decomposition : Time Series Talkritvikmath · 6 min

STL decomposition used for robust anomaly detection, with code.

23a CUSUMMatthias Schonlau · 5 min

A short lecture on the CUSUM chart.

Why Most Time Series Anomaly Detection Results are MeaninglessEamonn Keogh · 2 min

A two-minute summary of the problems with popular anomaly benchmarks.

Sean Law - STUMPY: Modern Time Series Analysis with Matrix Profiles | SciPy 2024SciPy · 30 min

The matrix profile with the STUMPY library: finding patterns and anomalies, with a live demo.

• Chandola, Banerjee and Kumar (2009), Anomaly detection: A survey, ACM Computing Surveys 41(3): the kinds of anomalies and the main families of methods.
• Page (1954), Continuous inspection schemes, Biometrika 41: the CUSUM. Montgomery, Introduction to Statistical Quality Control: CUSUM and Shewhart charts with ARL tables.
• Hochenbaum, Vallis and Kejariwal (2017), Automatic anomaly detection in the cloud via statistical learning (S-H-ESD).
• Truong, Oudre and Vayatis (2020), Selective review of offline change point detection methods (binary segmentation, PELT, the ruptures library).
• Adams and MacKay (2007), Bayesian online changepoint detection.
• The UCR matrix profile page: papers, code and tutorials (Yeh et al., 2016).
• Wu and Keogh (2021), Current time series anomaly detection benchmarks are flawed and are creating the illusion of progress.
• Kim et al. (2022), Towards a rigorous evaluation of time-series anomaly detection (AAAI): the point-adjust problem.
• Liu, Ting and Zhou (2008), Isolation forest, ICDM: the isolation forest.