Monitoring drift in practice

A model is trained on one slice of the world and then used on whatever arrives next. Drift is the gap between the two, and this section is how to measure it: keep a reference sample, compare each day's inputs and predictions against it with a statistic such as PSI, and watch the outputs when the true labels will not arrive for weeks. Most of the work is deciding what deserves an alert, because a monitor that fires on every feature is a monitor nobody reads, including on the day it is right.

Compare against a reference#

Drift is always relative, so the first decision is what to compare against. Keep a fixed reference sample — usually the validation data the model was signed off on, or a recent week everyone agrees was healthy — and compare each new day's inputs and predictions with it. The tempting alternative, comparing today against yesterday, sounds more responsive and is worse: a distribution that moves one per cent a day never trips a day-on-day comparison and has moved thirty per cent by the end of the month. Compare against a fixed point, and move the fixed point deliberately, when you retrain.

Population stability index#

PSI bins the reference into equal-sized buckets and measures how much the share in each bucket moved. It is a symmetrised relative entropy, which matters only in that it is an effect size: it answers how much the distribution moved, and it does not care how many rows you have. It works for any numeric feature and for model scores, and its thresholds are conventions to tune, not laws.

import numpy as np

def psi(reference, current, bins=10):
    edges = np.quantile(reference, np.linspace(0, 1, bins + 1))
    edges[0], edges[-1] = -np.inf, np.inf
    ref = np.histogram(reference, edges)[0] / len(reference)
    cur = np.histogram(current, edges)[0] / len(current)
    ref, cur = np.clip(ref, 1e-6, None), np.clip(cur, 1e-6, None)
    return float(np.sum((cur - ref) * np.log(cur / ref)))

# common rule of thumb: below 0.1 stable, 0.1 to 0.25 worth a look, above 0.25 investigate

A detector that fires#

Reading the definition is not the same as watching it work. Below is a five-feature stream where only two features drift, gradually, and three are stationary. The reference is twenty thousand rows from launch; each day is five thousand rows of what arrived.

rng = np.random.default_rng(0)
FEATURES = ["amount", "age_days", "n_sessions", "latency_ms", "basket_size"]

def sample(n, shift=0.0, rng=rng):
    """Five features. 'shift' moves only amount and n_sessions."""
    return np.column_stack([
        rng.lognormal(3.0 + shift, 0.8, n),          # amount: mean drifts up
        rng.normal(400, 120, n),                     # age_days: stable
        rng.poisson(4 + 6 * shift, n),               # n_sessions: drifts up
        rng.normal(180, 40, n),                      # latency_ms: stable
        rng.normal(3.2, 1.1, n),                     # basket_size: stable
    ])

reference = sample(20000)

print(f"{'day':>4} {'shift':>6} " + " ".join(f"{f:>12}" for f in FEATURES))
for day, shift in [(1, 0.0), (10, 0.05), (20, 0.10), (30, 0.20), (40, 0.35)]:
    cur = sample(5000, shift)
    row = [psi(reference[:, j], cur[:, j]) for j in range(5)]
    print(f"{day:>4} {shift:>6.2f} " + " ".join(f"{v:>12.3f}" for v in row))

 day  shift       amount     age_days   n_sessions   latency_ms  basket_size
   1   0.00        0.004        0.001        0.002        0.001        0.002
  10   0.05        0.010        0.003        0.024        0.001        0.002
  20   0.10        0.022        0.003        0.084        0.001        0.003
  30   0.20        0.064        0.001        0.329        0.002        0.002
  40   0.35        0.169        0.003        0.823        0.003        0.001

That table is what a working detector looks like, and the useful part is the three columns that do nothing. They establish that 0.004 is what stationary looks like at five thousand rows a day, which is the only way to know that n_sessions at 0.084 on day twenty is real rather than sampling. Without stationary features in the same table you have a number with no scale, and the 0.1 and 0.25 rules of thumb become the only thing you have, on a problem where they may be far too loose or far too tight.

Why a p-value is the wrong tool here#

The obvious alternative is a two-sample test: the Kolmogorov-Smirnov test for numeric features, chi-squared for categorical ones. Both work, both are one line, and both are a trap at production sample sizes. Here is the same day-twenty comparison scored both ways.

from scipy.stats import ks_2samp

cur = sample(5000, 0.10)
for j, f in enumerate(FEATURES):
    stat, p = ks_2samp(reference[:, j], cur[:, j])
    print(f"{f:>12}  KS {stat:.3f}  p {p:.2e}  PSI {psi(reference[:, j], cur[:, j]):.3f}")

      amount  KS 0.062  p 7.11e-14  PSI 0.020
    age_days  KS 0.015  p 3.46e-01  PSI 0.002
  n_sessions  KS 0.122  p 4.36e-52  PSI 0.085
  latency_ms  KS 0.010  p 8.51e-01  PSI 0.002
 basket_size  KS 0.018  p 1.52e-01  PSI 0.003

Read the amount row. Its PSI is 0.020, comfortably inside anybody's stable band, and its KS p-value is seven parts in a hundred trillion. The test is not wrong: with twenty thousand reference rows and five thousand current rows, a shift that small genuinely is not chance. That is the problem. A significance test answers whether the difference is real, and at production volumes every difference is real. The question you actually have is whether the difference is large enough to matter, and that is what an effect size such as PSI, or the KS statistic itself rather than its p-value, is for.

So use the statistic and ignore the p-value. The KS statistic column works as well as PSI here and is one import: 0.010 to 0.018 for the stationary features, 0.122 for the one that moved. Where PSI earns its keep is that it is defined on categorical features too, and that the banding convention is widely enough known that a risk team will recognise the number.

When labels arrive late#

For fraud or churn, the true outcome may take weeks. Until then, watch what you can: input drift, the distribution of predicted scores, the share of decisions above the threshold, and proxies such as complaints or manual review overturns. When labels do arrive, compute performance by the week the prediction was made, not the week the label landed, or a single bad week smears across the three months of cohorts that were still open.

The failure input drift cannot see#

Everything so far watches the inputs. There is a common failure that leaves the inputs untouched and the ranking metrics untouched, and quietly destroys the probabilities: the base rate moves. A model trained when ten per cent of cases were positive is serving a population where forty per cent are. Its features still separate the classes exactly as well as they did. Its probabilities are now wrong by a quarter, and every downstream expected-value decision is wrong with them.

from sklearn.linear_model import LogisticRegression
from sklearn.metrics import roc_auc_score, brier_score_loss

def make(n, base_rate, rng):
    """Three features; only the share of positives changes between calls."""
    y = (rng.random(n) < base_rate).astype(int)
    x = rng.normal(0, 1.0, size=(n, 3)) + y[:, None] * 1.6
    return x, y

def ece(y, p, bins=10):
    """Expected calibration error: the mean gap between predicted and actual, by bin."""
    idx = np.clip(np.digitize(p, np.linspace(0, 1, bins + 1)[1:-1]), 0, bins - 1)
    return float(sum(
        (idx == b).mean() * abs(y[idx == b].mean() - p[idx == b].mean())
        for b in range(bins) if (idx == b).any()))

rng2 = np.random.default_rng(1)
Xtr, ytr = make(40000, 0.10, rng2)
clf = LogisticRegression(max_iter=1000).fit(Xtr, ytr)

print("month  true rate     AUC   Brier     ECE   mean p   flagged   precision")
for month, rate in [(0, 0.10), (2, 0.15), (4, 0.20), (6, 0.30), (8, 0.40)]:
    X, y = make(40000, rate, rng2)
    p = clf.predict_proba(X)[:, 1]
    hit = p > 0.5
    print(f"{month:>5} {rate:>10.2f} {roc_auc_score(y, p):>7.4f} "
          f"{brier_score_loss(y, p):>7.4f} {ece(y, p):>7.4f} {p.mean():>8.4f} "
          f"{hit.mean():>9.3f} {y[hit].mean():>11.3f}")

month  true rate     AUC   Brier     ECE   mean p   flagged   precision
    0       0.10  0.9741  0.0309  0.0027   0.1006     0.087       0.840
    2       0.15  0.9740  0.0410  0.0165   0.1348     0.123       0.892
    4       0.20  0.9760  0.0477  0.0310   0.1662     0.157       0.925
    6       0.30  0.9744  0.0689  0.0677   0.2327     0.229       0.953
    8       0.40  0.9745  0.0867  0.1000   0.2999     0.301       0.968

Every number a normal dashboard carries is telling you the model is fine or improving. AUC does not move at all, because it depends only on the ranking, and the ranking is untouched by a change in how many positives there are; any monotone rescaling of the scores leaves it identical, which is exactly the class of damage happening here. Precision at the 0.5 threshold rises from 0.840 to 0.968, which will be reported as a win. Meanwhile the average predicted probability is 0.2999 against a true rate of 0.3999, ECE has gone from 0.0027 to 0.1000, and Brier score has nearly tripled. Anything sized from the probability — expected fraud loss, how many reviewers to roster, which customers clear a cost threshold — is a quarter low and nobody is looking at the metric that says so.

Two consequences for what you monitor. Put a calibration metric on the dashboard, not only a ranking one: mean predicted probability against realised rate is the cheap version and needs no binning. And notice that the prediction-distribution monitor would have caught this on its own — the share flagged went from 8.7 per cent to 30.1 per cent — which is the argument for watching outputs even when you cannot watch outcomes.

The repair is cheaper than retraining and usually better. Fit one logistic regression from the model's own log-odds to the recent labels, which is Platt scaling and adjusts the intercept without touching what the model learned. On eight thousand labelled rows from the drifted month it moves Brier from 0.0858 to 0.0586 and ECE from 0.1005 to 0.0033, leaving AUC at 0.9755 exactly where it was. Recalibrate first; retrain when the ranking itself has decayed.

Alerts people do not ignore#

The statistics are the easy half. The hard half is that the detector above produces five numbers a day per model, a real system has eighty features, and a dashboard with eighty statistics on it is not monitoring — it is a place where drift goes to be unnoticed. With eighty features and a threshold tuned so a stationary feature trips one day in a hundred, roughly one feature alerts every single day, forever, and the team learns within a fortnight that these alerts mean nothing.

• Alert on a handful of signals, not on every feature: the prediction distribution, calibration against realised outcomes, and the five to ten features the model actually leans on, taken from feature importance.
• Compute the rest and store them without alerting. They are what you read during an investigation, and they are worthless as pages.
• Require persistence. Two consecutive days, or three of five, removes nearly all of the single-day noise and delays a genuine signal by one day.
• Alert on the input-data checks first. Most drift alerts that fire turn out to be an upstream unit change or a broken join, and the data-quality suite catches those with a clearer message.
• Every alert links to a dashboard and a runbook naming what to check first, in order. An alert nobody knows how to action is an alert that gets muted.

What to skip#

• Skip the p-values. At production volume they are all significant; use the effect size.
• Skip drift detection on features the model barely uses. It is real, and it does not matter, and it costs you the credibility of the alerts that do.
• Skip the elaborate multivariate detectors until the per-feature ones are in place and trusted. Most production drift is one feature or one upstream release.
• Worth memorising: the PSI bands, that AUC is blind to calibration drift, and that comparing today against yesterday cannot see a slow trend.