Usage¶

Online detection¶

The online detector processes the series one point at a time and keeps the posterior over the run length, the number of observations since the last change. Two helpers turn that posterior into changepoints.

from functools import partial

import torch

from bayesian_changepoint_detection import (
    StudentT,
    changepoint_probabilities,
    constant_hazard,
    get_map_changepoints,
    online_changepoint_detection,
)

torch.manual_seed(42)
data = torch.cat([
    torch.randn(50) + 0,  # first segment: mean 0
    torch.randn(50) + 3,  # second segment: mean 3
    torch.randn(50) + 0,  # third segment: mean 0
])

hazard = partial(constant_hazard, 250)  # prior: one change every ~250 points
likelihood = StudentT(alpha=0.1, beta=0.01, kappa=1, mu=0)  # unknown mean and variance

R, map_run_lengths = online_changepoint_detection(data, hazard, likelihood)

# Index of the first point of each new segment on the MAP run-length path.
# min_separation merges starts closer than that many points when the
# posterior hesitates between neighbors.
print(get_map_changepoints(R, min_separation=10))  # tensor([ 50, 100])

# Or a probability per position, judged `lag` observations later.
probs = changepoint_probabilities(R, lag=10)  # probs[t] refers to data index t
print(torch.where(probs[1:] > 0.5)[0] + 1)  # tensor([ 50, 100])

R[r, t] is P(run length = r | first t observations). Why not simply threshold R[0, :]? Under a constant hazard the posterior probability of run length 0 is the hazard rate at every step, whatever the data say; the evidence for a change at t shows up in the following columns, as mass at run length k in column t + k. changepoint_probabilities reads exactly that. viterbi_changepoints(data, hazard, likelihood) returns the single most probable run-length path instead, i.e. the MAP segmentation.

Streaming¶

online_changepoint_detection needs the whole series and returns the (T+1)² matrix R. For a stream of unknown length, feed observations one at a time to OnlineChangepointDetector; it keeps only the current run-length posterior. With max_run_length the memory and the time per observation stay bounded: run lengths above the bound are dropped and the posterior renormalized, which is exact until the bound is reached and a close approximation afterwards when segments are shorter than the bound.

from bayesian_changepoint_detection import OnlineChangepointDetector

detector = OnlineChangepointDetector(
    hazard, StudentT(alpha=0.1, beta=0.01, kappa=1, mu=0), max_run_length=500
)
starts = []
for x in data:  # any iterable: a socket, a file, a generator
    detector.update(x)
    # P(a segment started 10 observations ago), as changepoint_probabilities
    if detector.t > 10 and detector.changepoint_probability(lag=10) > 0.5:
        starts.append(detector.t - 10)
print(starts)  # [50, 100]

Without max_run_length each posterior equals the corresponding column of R from online_changepoint_detection (up to float32 rounding).

Direction and size of each change¶

segment_statistics summarizes the segments between changepoints: mean, standard deviation, and for every changepoint the change in mean and a Welch z-score (issue #42).

from bayesian_changepoint_detection import segment_statistics

stats = segment_statistics(data, get_map_changepoints(R, min_separation=10))
print(stats.means.tolist())         # about [0.096, 3.173, -0.165]
print(stats.mean_changes.tolist())  # about [3.078, -3.339]: up, then down

The z-score ignores that the changepoints were found in the same data, so it overstates significance; treat it as a rough guide. Offline positions are the last index of the old segment: pass positions + 1.

Offline detection¶

The offline detector sees the whole series and returns, for every position, the posterior probability that a segment ends there. It is usually sharper than the online detector; use it for retrospective analysis.

from bayesian_changepoint_detection import const_prior, offline_changepoint_detection
from bayesian_changepoint_detection.offline_likelihoods import StudentT as OfflineStudentT

prior = partial(const_prior, p=1 / (len(data) + 1))  # flat prior on segment length
Q, P, changepoint_log_probs = offline_changepoint_detection(data, prior, OfflineStudentT())

changepoint_probs = torch.exp(changepoint_log_probs).sum(0)  # P(a segment ends at t)
print(torch.where(changepoint_probs > 0.5)[0])  # tensor([49, 99])

The two detectors use different index conventions: online reports the first point of the new segment (50), offline the last point of the old one (49). See the FAQ.

Multivariate data¶

Pass a [T, d] tensor, one row per observation, and a multivariate likelihood; everything else is the same. All three detectors check their input first: data must be [T] or [T, d], non-empty, real and finite, and match the likelihood's dims. A transposed [d, T] tensor is rejected with a hint rather than read as d observations.

from bayesian_changepoint_detection import MultivariateT

dims = 3
mv_data = torch.cat([
    torch.randn(50, dims) + torch.tensor([0.0, 0.0, 0.0]),
    torch.randn(50, dims) + torch.tensor([2.0, -1.0, 1.0]),
    torch.randn(50, dims) + torch.tensor([0.0, 0.0, 0.0]),
])

R, _ = online_changepoint_detection(mv_data, hazard, MultivariateT(dims=dims))
print(get_map_changepoints(R, min_separation=10))  # tensor([ 48, 100])

The first start lands two points early on this draw: the lag-10 posterior puts 0.53 on 48, 0.13 on 49 and 0.26 on 50, and the MAP path takes the mode. Read changepoint_probabilities when the exact position matters.

Devices¶

Every likelihood and both detectors take a device argument. The default is the CPU; name a device on the likelihood ("cuda", "mps", or "auto" for the first available) to opt into an accelerator, and the detectors follow it. On a laptop the CPU is the faster choice for the online detector (measured: 6–30x faster than MPS), and the offline detector always runs on the CPU under MPS because it needs float64. How the argument is resolved, what has been measured, how to time your own workload and how much memory the tables need: docs/devices.md.

API at a glance¶

Function	Returns
`online_changepoint_detection(data, hazard, likelihood)`	`R` (run-length posterior, `[T+1, T+1]`) and the MAP run length after each point
`changepoint_probabilities(R, lag)`	`P(a new segment started at t)`, judged `lag` observations later
`get_map_changepoints(R, min_separation=1)`	indices where the MAP run-length path starts a new segment
`viterbi_changepoints(data, hazard, likelihood)`	the single most probable run-length path and its segment starts
`compute_run_length_posterior(data, hazard, likelihood)`	just `R`, for code that only wants the posterior
`segment_statistics(data, starts)`	per-segment mean and std, change in mean and z-score at each changepoint
`OnlineChangepointDetector(hazard, likelihood, max_run_length=None)`	streaming detector: `update(x)`, `run_length_posterior`, `map_run_length`, `changepoint_probability(lag)`
`offline_changepoint_detection(data, prior, likelihood)`	`Q` (log evidence), `P` (segment log likelihoods), `Pcp` (log probability of the j-th changepoint at t)
`constant_hazard(lam, r)`	hazard `1 / lam` for every run length
`negative_binomial_hazard(k, p, r)`	hazard of negative binomial segment lengths (mean `k / p`); the online counterpart of `negative_binomial_prior`
`const_prior`, `geometric_prior`, `negative_binomial_prior`	log prior on segment length for the offline detector
`online_likelihoods.StudentT`, `online_likelihoods.MultivariateT`	online conjugate models (Normal-Gamma, Normal-Wishart)
`offline_likelihoods.StudentT`, `MultivariateT`, `IndependentFeaturesLikelihood`, `FullCovarianceLikelihood`	offline segment marginal likelihoods
`online_likelihoods.Poisson`, `offline_likelihoods.Poisson`	count data: Gamma-Poisson, negative-binomial predictive
`online_likelihoods.NormalKnownVariance`, `offline_likelihoods.NormalKnownVariance`	mean changes with a known noise variance: Normal-Normal, Normal predictive
`get_device`, `get_device_info`, `to_tensor`	device helpers

All public functions have NumPy-style docstrings with the formulas and the paper they come from; the API reference on the documentation site is generated from them.