Change-Point Detection (offline and online)

"For three weeks the dashboard was green and I was proud. Then someone pulled last quarter's data and discovered the mean had moved on a Tuesday in March, quietly, while I was busy congratulating myself for catching spikes. The spikes were never the problem. The new normal was."

A Regime That Changed When the Dashboard Was Not Looking

In Section 8.2 we built statistical control charts (Shewhart, EWMA, CUSUM) that monitor a stream and raise an alarm when a process drifts off target. Those charts are tuned to detect a shift quickly, but they answer a yes-or-no monitoring question and do not, by themselves, partition a recorded history into its constituent regimes. This section completes the picture. We treat the change point as a first-class object to be located, not merely flagged, and we develop the two regimes of the problem that recur throughout temporal AI: the retrospective analysis of a fixed dataset and the sequential analysis of a live stream. The persistent-shift idea connects backward to the non-stationarity of Section 3.1, where we first defined stationarity and warned that real series violate it, and forward to concept drift, the supervised-learning cousin of change-point detection that Chapter 20 treats in full. Throughout we use the notation of the unified table in Appendix A: $x_{1:n} = (x_1, \dots, x_n)$ for the observed series, $\tau$ for a change-point location, and $\mathcal{C}(\cdot)$ for a segment cost.

1. Anomalies Versus Change Points Beginner

The single most useful idea in this section is also the simplest: an anomaly is transient and a change point is permanent. An anomaly is a brief deviation, one or a few samples that depart from the local pattern and then the series returns to its previous behavior, the lone temperature spike, the one fraudulent transaction, the sensor glitch. A change point is a persistent shift in the generating distribution: after the change point the series keeps obeying a different law (a new mean, a new variance, a new correlation structure, a new trend) and never reverts on its own. The mental test is the recovery test. Mask a candidate region and ask what the series does afterward. If it snaps back to the old regime, the region was an anomaly. If it settles into a new regime, you have found a change point.

This distinction matters because it changes what you should do. An anomaly is usually something to filter, flag, or repair; a change point is something to adapt to. A forecaster that treats a regime shift as a string of anomalies will keep "correcting" data that is in fact the new truth, and a forecaster that treats an anomaly as a regime shift will throw away a perfectly good model over a single glitch. The two failures are mirror images, and both are expensive: the first floods an operations team with alerts about a new normal that is not going away (the practical example at the end of this section is exactly this story), and the second resets a working model every time a sensor hiccups. Figure 8.3.1 contrasts the two on a shared time axis.

Formally, model the series as independent draws whose distribution can change at a finite set of unknown locations $0 = \tau_0 < \tau_1 < \cdots < \tau_m < \tau_{m+1} = n$. On segment $j$, spanning samples $\tau_{j}+1$ through $\tau_{j+1}$, the data follow a single distribution $F_j$, and the change points are exactly the boundaries where $F_j \ne F_{j+1}$. A pure mean-shift model writes $x_t = \mu_j + \varepsilon_t$ with $\mu_j$ constant within a segment; a variance-change model holds the mean fixed and lets $\sigma_j^2$ jump; a correlation-change model lets the dependence structure shift while marginals stay put. The detector's job is to estimate the number $m$ and the locations $\{\tau_j\}$ from the data alone. This is harder than anomaly detection precisely because $m$ is unknown: every additional change point can only reduce the within-segment cost, so something must stop the model from declaring a change at every sample. That something is the penalty, and it is the subject of the next subsection.

2. Offline Detection: Cost Plus Penalty, PELT, and Binary Segmentation Intermediate

Offline, or retrospective, detection has the whole series in hand and asks for the best segmentation of all of it at once. The dominant framework is cost-plus-penalty. Assign each candidate segment a cost $\mathcal{C}(x_{a+1:b})$ that measures how badly a single regime fits the samples from $a+1$ to $b$ (low cost means the segment is homogeneous), then choose the set of change points $\tau_1 < \cdots < \tau_m$ that minimizes the total segment cost plus a penalty $\beta$ paid once per change point:

The two terms pull in opposite directions, and that tension is the whole story. The cost term always falls as you add change points (more segments fit the data more tightly, and with one segment per sample the cost can be driven to zero), so without the penalty the optimum is the absurd segmentation that calls every sample its own regime. The penalty $\beta m$ charges a fixed price for each split, so a change point survives only if it lowers the cost by more than $\beta$. The penalty is therefore the knob that sets the number of change points: small $\beta$ yields many, large $\beta$ yields few. The most common penalties are information criteria carried over from Section 5.7: the Bayesian Information Criterion sets $\beta = p \ln n$ where $p$ is the number of parameters per segment, penalizing splits by the log of the sample size, while the Akaike criterion uses the lighter $\beta = 2p$.

The natural cost for a mean-and-variance Gaussian segment is the negative log-likelihood of the maximum-likelihood fit on that segment. For a segment of length $L = b - a$ with sample variance $\hat\sigma^2$, the Gaussian negative log-likelihood reduces (up to constants) to

so a segment whose samples scatter tightly around their own mean is cheap and one that must straddle a hidden shift is expensive. The choice of cost function is where you tell the detector what kind of change you care about, and it is the single most consequential modelling decision in offline detection. An $\ell_2$ cost $\mathcal{C} = \sum_t (x_t - \bar x)^2$ targets mean shifts only and is blind to a regime that changes variance while keeping its mean: feed it a segment whose mean is constant but whose noise doubles and it sees nothing to split. A Gaussian mean-and-variance cost (the log-variance form above) catches both. A kernel cost, which compares segments through a feature map such as the radial-basis kernel, targets arbitrary distributional change including shifts in higher moments or in the shape of the distribution. A rank-based or nonparametric cost trades a little power for robustness to outliers, catching shifts a Gaussian model would smear. The cost-plus-penalty machinery is identical across all of them; you swap the per-segment cost to swap the kind of change you detect, which is exactly why the ruptures library exposes the cost as a one-word argument in the worked example below.

How do we actually minimize the objective? Naively there are $2^{n-1}$ possible segmentations, far too many to enumerate. Two classes of algorithm make it tractable, and they trade exactness against speed.

Exact dynamic programming and PELT. Let $F(b)$ be the optimal total cost of segmenting the prefix $x_{1:b}$. The last segment ends at $b$ and begins at some earlier change point $t$, which gives the Bellman recursion

solved left to right; the change points are recovered by tracing back the minimizing $t$ at each step. Plain dynamic programming costs $O(n^2)$ because each $b$ scans all earlier $t$. The Pruned Exact Linear Time algorithm (PELT) of Killick, Fearnhead, and Eckley (2012) keeps the exactness but prunes the candidate set: a candidate change point $t$ can be discarded forever once $F(t) + \mathcal{C}(x_{t+1:b}) + K \ge F(b)$ for a constant $K$, because such a $t$ can never be the optimal last change point for any later endpoint. With this pruning, and under mild conditions on the change-point spacing, PELT runs in expected $O(n)$ time while still returning the globally optimal segmentation, which is why it is the default exact method in practice.

Binary segmentation. The classic greedy alternative is binary segmentation. Scan the whole series for the single best split (the location that most reduces the total cost), accept it if the reduction beats the penalty, then recurse independently on the left and right pieces, stopping when no further split pays for itself. Binary segmentation is fast, $O(n \log n)$, and trivial to implement, but it is approximate: because it commits to one split before looking for the next, it can misplace change points when two real changes interact, and it is biased when segments are short. Refinements (wild binary segmentation, which scans many random sub-intervals, and the SeedBS variant) repair much of this at modest extra cost. Table 8.3.1 lays the offline methods side by side.

3. The Likelihood-Ratio and CUSUM View of a Single Change Intermediate

Before tackling many change points, it pays to understand one perfectly, because the single-change case has a clean optimal statistic that every multiple-change method inherits. Suppose the series is Gaussian with known variance $\sigma^2$, and we test the hypothesis that the mean shifts exactly once. Under the null $H_0$ there is no change: all $n$ samples share one mean $\mu_0$. Under the alternative $H_1(\tau)$ the mean is $\mu_0$ up to time $\tau$ and $\mu_1 \ne \mu_0$ afterward. For a fixed candidate $\tau$, the generalized log-likelihood ratio (with the means replaced by their sample estimates) reduces to a comparison of within-segment fit, and maximizing it over all candidate locations gives the test statistic

where $\bar x_{1:\tau}$ and $\bar x_{\tau+1:n}$ are the pre- and post-candidate means. The statistic is large when some split produces two segments whose means are far apart, and the weight $\tau(n-\tau)/n$ rewards balanced splits (a change in the middle is easier to detect than one at the very edge, where one segment has almost no data). The estimated change point is the $\tau$ that achieves the maximum, and a change is declared when $\Lambda$ exceeds a threshold chosen to control the false-alarm rate. This single-change likelihood-ratio scan is exactly what binary segmentation calls at every recursion, so understanding it here explains the inner loop of subsection two.

The same idea, recast for a stream that arrives one sample at a time, is the CUSUM statistic we met as a control chart in Section 8.2. CUSUM accumulates the log-likelihood ratio of "post-change" against "pre-change" for each new sample and resets the accumulation whenever it would go negative, so it tracks the running evidence that a change has already happened. For a shift in the mean of a Gaussian, the one-sided upward CUSUM is the recursion

where $\mu_0$ is the in-control mean and $k$ is a slack (reference value) usually set to half the smallest shift worth detecting. The statistic stays near zero while the data sit at $\mu_0$, but once the mean shifts upward the increments become persistently positive and $S_t$ climbs steadily until it crosses an alarm threshold $h$. The link to this section's framework is exact: CUSUM is the sequential, recursive form of the same likelihood-ratio test that the offline scan computes in one shot over the whole series. The offline scan asks "where is the single best split?" with all data visible; CUSUM asks "has the split already happened?" with data arriving live. The connection is summarized below before we turn to the fully Bayesian streaming detector.

Code 8.3.0 makes the equivalence tangible: it runs the recursive CUSUM on a stream that shifts halfway through and reports the alarm time, so you can watch the accumulating evidence cross the threshold in real time.

import numpy as np

def cusum_detect(x, mu0, k, h):
    """One-sided upward CUSUM. Accumulates evidence that the mean rose
    above mu0; resets to 0 when evidence goes negative; alarms at h.
    Returns (alarm_time or None, the running statistic S)."""
    S = 0.0
    track = np.empty(len(x))
    alarm = None
    for t, xt in enumerate(x):
        S = max(0.0, S + (xt - mu0) - k)     # accumulate, floored at zero
        track[t] = S
        if alarm is None and S > h:           # first crossing is the detection
            alarm = t
    return alarm, track

rng = np.random.default_rng(3)
stream = np.concatenate([rng.normal(0.0, 1.0, 200),   # in-control: mean 0
                         rng.normal(1.0, 1.0, 200)])  # shift to mean 1 at t=200
alarm, _ = cusum_detect(stream, mu0=0.0, k=0.5, h=4.0)
print("true shift at: 200   CUSUM alarm at:", alarm, " delay:", alarm - 200)

true shift at: 200   CUSUM alarm at: 208  delay: 8

4. Online Detection: Bayesian Online Change-Point Detection Advanced

CUSUM gives a single threshold-crossing alarm. Bayesian Online Change-Point Detection (BOCPD), introduced by Adams and MacKay (2007), gives something richer for the streaming setting: at every time step it maintains a full posterior distribution over the run length, the number of samples since the most recent change point. Let $r_t$ be the run length at time $t$. If $r_t = 0$ a change just occurred; if $r_t = 12$ the current regime has lasted twelve samples. BOCPD recursively updates the posterior $P(r_t \mid x_{1:t})$ as each new sample arrives, and a change is signalled when probability mass collapses onto $r_t = 0$.

The recursion has two ingredients. First, a hazard function $H(r)$ giving the prior probability that a change occurs given the current run length; a constant hazard $H(r) = 1/\lambda$ encodes a geometric prior on segment lengths with expected length $\lambda$. Second, a predictive distribution $\pi^{(r)}(x_t)$ that scores the new sample under the segment of length $r$ (for a Gaussian with a conjugate prior, this is a closed-form Student-$t$). The joint distribution over run length and data then updates by a message-passing recursion that splits the mass at each step into a "growth" branch (the run continues, $r_t = r_{t-1} + 1$) and a "change" branch (the run resets, $r_t = 0$):

Normalizing the joint at each step yields the run-length posterior $P(r_t \mid x_{1:t})$. The intuition is that each existing run-length hypothesis predicts the next sample; hypotheses that predict it well grow their mass, hypotheses that predict it badly shrink, and whenever the newest sample is poorly explained by every ongoing run, mass floods into $r_t = 0$ and a change is declared. It helps to picture the posterior as a histogram over "how long has the current regime lasted." During a stable regime that histogram develops a single tall peak that marches steadily to the right (the run-length most-probable estimate ticks up by one each sample, since the regime keeps lasting). At a true change point the new sample is improbable under every ongoing run, the change branch collects almost all the mass, the peak collapses back to $r_t = 0$, and the march begins again from zero. Reading off the location of that peak over time, as the worked example does, turns the full posterior into a simple change-detection signal: a sawtooth that climbs through each regime and drops at each boundary.

The cost is one growing distribution per step (the support of $r_t$ expands by one each sample), so practical implementations cap or prune the run length to keep memory bounded, discarding run-length hypotheses whose posterior mass has fallen below a threshold. That pruning is not a detail; it is exactly the streaming constraint that the online-learning methods of Chapter 20 wrestle with at scale, where a detector cannot store all of history and must summarize the past in bounded state. BOCPD also illustrates a clean separation that recurs in sequential AI: the hazard function encodes a prior belief about how often regimes change, the predictive distribution encodes a model of what each regime looks like, and the recursion combines them. Swap in a richer predictive (a neural one, per the research frontier) and the same skeleton detects changes in series no Gaussian could model.

Every online detector lives or dies by one trade-off: detection delay against false-alarm rate. Make the detector eager (low threshold, high hazard) and it catches real changes quickly but also fires on noise; make it cautious and it stops false alarms but reacts late, sometimes far too late for the change to be actionable. There is no free lunch: for a fixed amount of data per sample, reducing the expected detection delay necessarily raises the false-alarm probability, and the only ways to improve both at once are to gather more informative data or to accept a longer averaging window. This is the streaming-detection analogue of the bias-variance trade-off, and it governs how aggressively you can tune $H$, $k$, $h$, or the BOCPD hazard.

5. Worked Example: Detecting Regimes From Scratch and With ruptures Advanced

We now run the full pipeline on a synthetic regime-switching signal, the sensor-IoT running thread of Part II: imagine a vibration sensor on a machine whose operating mode changes a few times over a shift, each mode having its own mean amplitude. Code 8.3.1 builds a piecewise-constant signal with three true change points and implements the single-change likelihood-ratio scan of subsection three, then drives it greedily (binary segmentation) to recover all the changes from scratch.

import numpy as np

rng = np.random.default_rng(8)

# Piecewise-constant signal: four regimes, three true change points.
true_cps = [200, 450, 600]               # sample indices where the mean shifts
means    = [0.0, 3.0, 1.0, 4.0]          # mean amplitude of each regime
segs = np.split(np.arange(800), true_cps)
signal = np.concatenate([m + rng.normal(0, 0.8, len(s))
                         for m, s in zip(means, segs)])

def best_single_split(x):
    """Generalized likelihood-ratio scan for ONE mean change in x.
    Returns (best_index, statistic) maximizing tau(n-tau)/n * (mean gap)^2."""
    n = len(x)
    if n < 4:
        return None, 0.0
    csum = np.concatenate([[0.0], np.cumsum(x)])     # prefix sums for O(1) segment means
    best_stat, best_tau = 0.0, None
    for tau in range(2, n - 1):                       # need >=2 points each side
        left_mean  = csum[tau] / tau
        right_mean = (csum[n] - csum[tau]) / (n - tau)
        stat = (tau * (n - tau) / n) * (left_mean - right_mean) ** 2
        if stat > best_stat:
            best_stat, best_tau = stat, tau
    return best_tau, best_stat

def binseg(x, threshold, offset=0, found=None):
    """Binary segmentation: recursively accept the best split while it
    beats the threshold, recursing on each resulting sub-segment."""
    if found is None:
        found = []
    tau, stat = best_single_split(x)
    if tau is not None and stat > threshold:
        found.append(offset + tau)                    # record absolute location
        binseg(x[:tau], threshold, offset, found)     # recurse left
        binseg(x[tau:], threshold, offset + tau, found)  # recurse right
    return sorted(found)

detected = binseg(signal, threshold=40.0)
print("true change points    :", true_cps)
print("detected change points:", detected)

true change points    : [200, 450, 600]
detected change points: [201, 449, 601]

The from-scratch code makes the mechanism transparent, but in production you would reach for the ruptures library, which packages PELT, binary segmentation, window-based, and dynamic-programming search behind a uniform interface with a menu of cost functions. Code 8.3.2 is the library shortcut: the exact same detection, both with the optimal PELT search and with the greedy Binseg, in a handful of lines.

import ruptures as rpt

# Same signal as Code 8.3.1. ruptures expects a column vector.
sig = signal.reshape(-1, 1)

# PELT: exact, penalized, unknown number of change points.
pelt = rpt.Pelt(model="l2", min_size=20).fit(sig)
cps_pelt = pelt.predict(pen=30.0)        # 'pen' is the per-change penalty beta

# Binary segmentation: greedy, ask for the known number of changes.
binseg = rpt.Binseg(model="l2", min_size=20).fit(sig)
cps_binseg = binseg.predict(n_bkps=3)    # request exactly 3 break points

print("PELT change points  :", cps_pelt)     # trailing value is n (end marker)
print("Binseg change points:", cps_binseg)

PELT change points  : [201, 449, 601, 800]
Binseg change points: [201, 449, 601, 800]

Finally, the streaming view. Code 8.3.3 sketches Bayesian Online Change-Point Detection on the same signal fed one sample at a time, tracking the run-length posterior and flagging a change whenever the most probable run length collapses toward zero. This is the from-scratch counterpart to the recursion of subsection four, kept deliberately compact with a Gaussian predictive of fixed variance.

import numpy as np
from scipy.stats import norm

def bocpd_stream(x, hazard=1/250, mu0=2.0, kappa=0.1, sig=0.8, run_cap=400):
    """Bayesian Online Change-Point Detection with a Gaussian predictive.
    Returns, for each t, the most probable run length. A sharp drop to ~0
    marks a detected change. Run length is capped to bound memory."""
    R = np.array([1.0])                    # run-length posterior, starts at r=0
    means = np.array([mu0])                # running posterior mean per run length
    counts = np.array([kappa])             # pseudo-count per run length
    map_run = []
    for xt in x:
        pred = norm.pdf(xt, loc=means, scale=sig)       # predictive prob per run length
        growth = R * pred * (1 - hazard)                 # run continues
        change = np.array([np.sum(R * pred * hazard)])   # run resets to 0
        R = np.concatenate([change, growth])
        R /= R.sum()                                     # normalize the posterior
        # Update sufficient statistics for each (now shifted) run length.
        new_mean = (means * counts + xt) / (counts + 1)
        means = np.concatenate([[mu0], new_mean])
        counts = np.concatenate([[kappa], counts + 1])
        if len(R) > run_cap:                             # prune to bound memory
            R, means, counts = R[:run_cap], means[:run_cap], counts[:run_cap]
            R /= R.sum()
        map_run.append(int(np.argmax(R)))                # most probable run length
    return np.array(map_run)

map_run = bocpd_stream(signal)
# A change is where the most probable run length drops sharply back toward 0.
drops = [t for t in range(1, len(map_run))
         if map_run[t] < 0.4 * map_run[t - 1] and map_run[t - 1] > 30]
print("BOCPD run-length resets near:", drops)

BOCPD run-length resets near: [203, 452, 604]

The three outputs together tell the section's whole story. The offline detectors of Codes 8.3.1 and 8.3.2 pinpoint the change points almost exactly because they see the entire series and optimize globally; the online detector of Code 8.3.3 finds the same changes but a few samples late, paying the detection-delay price for the privilege of flagging changes as they happen rather than in hindsight. Choosing between them is choosing your information regime: retrospective analysis of a recorded history wants PELT, while a live monitoring stream wants CUSUM or BOCPD, and the cost-plus-penalty and run-length-posterior skeletons we built make both the same idea seen from two sides.

Step back and see what the section assembled. We separated anomalies (transient) from change points (persistent) by the recovery test, cast offline detection as the cost-plus-penalty optimization solved exactly by PELT or greedily by binary segmentation, derived the single-change likelihood-ratio statistic and connected it to the CUSUM chart of Section 8.2, built the run-length posterior of Bayesian Online Change-Point Detection for the streaming regime, and confronted the detection-delay-versus-false-alarm trade-off that governs every sequential detector. The worked example tied it together, detecting the same regimes offline and online and making the delay visible. The persistent-shift theme hands directly to the rest of the book: Section 8.4 turns to forecasting-residual and reconstruction methods that detect anomalies and changes through a model's prediction errors, while the concept-drift machinery of Chapter 20 takes the change-point idea into supervised online learning, where a detected change triggers a model to adapt rather than merely a flag to fire.