MLOps Model Monitoring

Real-World Scenario

A bank’s fraud model worked well in January. By August its precision had fallen sharply, and the team found out from customer complaints rather than from a dashboard. Two things had changed. The customers were different, which is data drift, and fraud itself had changed, which is concept drift. The first can be caught the week it starts by comparing inputs with the training data. The second shows up only when labelled outcomes arrive, and those arrive late.

Every example runs locally with NumPy, SciPy and scikit-learn, on seeded synthetic data, so the numbers are reproducible.

The Data: Training Time vs This Week

Two samples of three features: what the model was trained on, and what production sent this week. Save this as data.py.

import numpy as np

rng = np.random.default_rng(42)
REFERENCE = {                                   # what the model was trained on
    "age":     rng.normal(35, 10, 10_000),
    "income":  rng.lognormal(10.8, 0.5, 10_000),
    "balance": rng.exponential(5_000, 10_000),
}
PRODUCTION = {                                  # this week's traffic
    "age":     rng.normal(42, 12, 1_000),       # older customers
    "income":  rng.lognormal(10.5, 0.6, 1_000), # lower incomes
    "balance": rng.exponential(5_000, 1_000),   # unchanged
}

Age and income have shifted; balance hasn’t. The checks below should agree with that and nothing else.

Detecting Drift With the KS Test

The two-sample Kolmogorov-Smirnov test compares two distributions without assuming either is normal. Its statistic is the largest gap between their cumulative distributions.

from scipy import stats
from data import REFERENCE, PRODUCTION

for feature in REFERENCE:
    result = stats.ks_2samp(REFERENCE[feature], PRODUCTION[feature])
    drift = "DRIFT" if result.pvalue < 0.05 else "ok"
    print(f"{feature:8} KS statistic {result.statistic:.3f}  p-value {result.pvalue:.2e}  {drift}")

$ python ks.py
age      KS statistic 0.295  p-value 1.61e-70  DRIFT
income   KS statistic 0.254  p-value 6.66e-52  DRIFT
balance  KS statistic 0.026  p-value 5.58e-01  ok

Age and income are flagged, with KS statistics of 0.295 and 0.254 and p-values so small they are effectively zero. Balance isn’t: its statistic of 0.026 is within what two samples of one distribution produce by chance. So far the test does exactly what you’d want.

The Failure: Significant Isn’t Important

Now shift ages by 0.3 years, a change nobody would retrain a model for, and grow both samples, as they grow with real traffic.

import numpy as np
from scipy import stats

rng = np.random.default_rng(1)
for n in (1_000, 100_000, 1_000_000):
    reference = rng.normal(35, 10, n)           # both samples grow, as with real traffic
    production = rng.normal(35.3, 10, n)        # ages shifted by 0.3 years: nobody cares
    result = stats.ks_2samp(reference, production)
    flag = "DRIFT" if result.pvalue < 0.05 else "ok"
    print(f"n={n:>9,}  KS statistic {result.statistic:.3f}  p-value {result.pvalue:.1e}  {flag}")

$ python ks_trap.py
n=    1,000  KS statistic 0.056  p-value 8.7e-02  ok
n=  100,000  KS statistic 0.015  p-value 3.2e-10  DRIFT
n=1,000,000  KS statistic 0.014  p-value 9.8e-80  DRIFT

At 1,000 rows the tiny shift passes, as it should. At 100,000 and 1,000,000 rows the KS statistic settles near 0.014–0.015, the true size of this small shift, but the p-value collapses to 3.2e-10 and then 9.8e-80, and the check reports drift. A p-value measures how surprising the gap would be if nothing had changed, not how big the gap is, and with enough data even a trivial gap becomes very surprising. Production traffic is large, so a p-value-only alert fires constantly. Alert on the size of the shift instead, such as the KS statistic or PSI, and use the p-value at most as a secondary filter.

Measuring the Size of a Shift With PSI

The Population Stability Index, long used in credit scoring, measures the size of a shift directly. Save this as psi.py.

import numpy as np
from data import REFERENCE, PRODUCTION

def psi(reference, production, bins=10):
    """Population Stability Index over bins cut at the reference quantiles."""
    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)
    prod = np.histogram(production, edges)[0] / len(production)
    ref, prod = np.clip(ref, 1e-4, None), np.clip(prod, 1e-4, None)   # avoid log(0)
    return float(np.sum((prod - ref) * np.log(prod / ref)))

def band(value):
    return "stable" if value < 0.1 else "moderate shift" if value < 0.25 else "major shift"

if __name__ == "__main__":
    for feature in REFERENCE:
        value = psi(REFERENCE[feature], PRODUCTION[feature])
        print(f"{feature:8} PSI {value:.3f}  {band(value)}")

$ python psi.py
age      PSI 0.465  major shift
income   PSI 0.345  major shift
balance  PSI 0.023  stable

PSI bins the training values into deciles and compares how production fills those bins. Age (0.465) and income (0.345) are well past the 0.25 “major shift” line, and balance (0.023) is stable. Unlike the p-value, PSI doesn’t grow just because you have more rows, and the conventional bands give teams a shared vocabulary for when to act. The np.clip line matters: an empty production bin would otherwise make log(0) and an infinite PSI.

Concept Drift: When the Inputs Look Fine

Drift checks compare inputs. Here the inputs never change, but the meaning of a positive case does, month by month, and labels arrive two months after each prediction.

import numpy as np
from sklearn.linear_model import LogisticRegression

rng = np.random.default_rng(7)
def month_of_data(month, n=2_000):
    X = rng.normal(size=(n, 2))                      # inputs never change: no data drift
    weight = 1.0 - 0.35 * month                      # but what makes a case "fraud" does
    y = (X[:, 0] + weight * X[:, 1] > 0).astype(int)
    return X, y

model = LogisticRegression().fit(*month_of_data(0))
TODAY, LABEL_DELAY = 6, 2                            # labels arrive two months late
print("month  flagged  accuracy")
for month in range(1, TODAY + 1):
    X, y = month_of_data(month)
    flagged = model.predict(X).mean()                # visible immediately
    known = month <= TODAY - LABEL_DELAY
    accuracy = f"{(model.predict(X) == y).mean():.1%}" if known else "labels not in yet"
    print(f"{month:5}  {flagged:7.1%}  {accuracy}")

$ python label_delay.py
month  flagged  accuracy
    1    48.4%  92.8%
    2    50.4%  85.0%
    3    50.6%  73.8%
    4    49.8%  64.8%
    5    49.7%  labels not in yet
    6    50.4%  labels not in yet

Accuracy falls from 92.8% to 64.8% over four months, and nothing visible without labels moves. The inputs are drawn from the same distribution every month, so KS and PSI would report no drift, and the share of cases flagged stays at about 50%. The model is wrong about which half. Months 5 and 6 are probably worse still, but their labels haven’t arrived, so their accuracy is unknown.

This is why input monitoring is not enough on its own. Track accuracy for every period as its labels land, and shorten the delay where you can: for example, have analysts label a small random sample of this week’s cases, so concept drift shows up in days rather than months.

Turning Checks Into Alerts

Put the checks together, with thresholds and severities, into something a scheduler can run and a pager can read.

import json
from data import REFERENCE, PRODUCTION
from psi import psi

def run_checks(production, accuracy=None, psi_limit=0.25, accuracy_floor=0.85):
    alerts = []
    for feature, values in production.items():
        value = psi(REFERENCE[feature], values)
        if value >= psi_limit:
            alerts.append({"type": "data_drift", "feature": feature, "psi": round(value, 3),
                           "severity": "high" if value >= 0.4 else "medium"})
    if accuracy is not None and accuracy < accuracy_floor:
        alerts.append({"type": "performance", "accuracy": accuracy,
                       "severity": "critical" if accuracy < accuracy_floor - 0.1 else "high"})
    return alerts

print(json.dumps(run_checks(PRODUCTION, accuracy=0.738), indent=1))
print("baseline week:", run_checks({f: v[:1_000] for f, v in REFERENCE.items()}, accuracy=0.93))

$ python monitor.py
[
 {
  "type": "data_drift",
  "feature": "age",
  "psi": 0.465,
  "severity": "high"
 },
 {
  "type": "data_drift",
  "feature": "income",
  "psi": 0.345,
  "severity": "medium"
 },
 {
  "type": "performance",
  "accuracy": 0.738,
  "severity": "critical"
 }
]
baseline week: []

This week’s traffic raises three alerts: high-severity drift on age (PSI 0.465 is past the 0.4 severity line) and medium on income, and a critical accuracy alert for the 73.8% measured in month 3 above. The second line is a sanity check worth keeping in any monitor: a sample drawn from the training data raises no alerts. A monitor that fires on its own baseline trains people to ignore it. In production this function runs on a schedule, the alerts go to a queue or pager, and a high-severity alert opens a retraining ticket instead of retraining automatically. Retraining on drifted data without a person looking can bake in a data bug.

Practice

1. PSI depends on how many bins you use. Compute it for age and balance with 5, 10, 20 and 50 bins, and decide which setting you’d standardize on for 1,000-row weekly samples.

   Solution

   from data import REFERENCE, PRODUCTION
   from psi import psi
   
   for feature in ("age", "balance"):
       row = [f"{bins} bins {psi(REFERENCE[feature], PRODUCTION[feature], bins):.3f}"
              for bins in (5, 10, 20, 50)]
       print(f"{feature:8}", "  ".join(row))

   $ python sol_bins.py
   age      5 bins 0.405  10 bins 0.465  20 bins 0.515  50 bins 0.598
   balance  5 bins 0.010  10 bins 0.023  20 bins 0.031  50 bins 0.067

   More bins mean higher PSI for both features, including balance, which hasn’t drifted at all: at 50 bins, each bin holds about 20 production rows, and sampling noise alone pushes balance to 0.067. Ten bins is the usual choice for samples this size. Whatever you pick, keep it fixed, because PSI values computed with different bin counts aren’t comparable.

2. One odd week (a marketing campaign, a data glitch) shouldn’t page anyone. Compute weekly PSI for age and alert only when it has been over 0.25 for three weeks running.

   Solution

   import numpy as np
   from data import REFERENCE
   from psi import psi
   
   rng = np.random.default_rng(3)
   shift = {1: 0, 2: 0, 3: 9, 4: 0, 5: 6, 6: 7, 7: 7}         # week -> mean age shift in years
   streak = 0
   for week, years in shift.items():
       value = psi(REFERENCE["age"], rng.normal(35 + years, 10, 1_000))
       streak = streak + 1 if value >= 0.25 else 0
       line = f"week {week}  PSI {value:.3f}  over the limit {streak} week(s) running"
       print(line + ("  -> ALERT" if streak == 3 else ""))

   $ python sol_debounce.py
   week 1  PSI 0.018  over the limit 0 week(s) running
   week 2  PSI 0.009  over the limit 0 week(s) running
   week 3  PSI 0.761  over the limit 1 week(s) running
   week 4  PSI 0.008  over the limit 0 week(s) running
   week 5  PSI 0.324  over the limit 1 week(s) running
   week 6  PSI 0.487  over the limit 2 week(s) running
   week 7  PSI 0.488  over the limit 3 week(s) running  -> ALERT

   Week 3’s spike (PSI 0.761) resets the next week and never alerts. The sustained shift from week 5 alerts in week 7, on its third consecutive week. The cost of debouncing is two weeks of delay on real drift, so use it for medium-severity drift and keep single-week alerts for severe or safety-critical features.

3. KS and PSI as written are for numeric features. Check a categorical one, payment method, with a chi-square test, and show how the category shares moved.

   Solution

   from scipy.stats import chi2_contingency
   
   methods = ["card", "bank transfer", "wallet"]
   training = [7_000, 2_500, 500]               # payment method counts at training time
   this_week = [560, 190, 250]                  # wallets jumped from 5% to 25%
   
   result = chi2_contingency([training, this_week])
   print(f"chi-square {result.statistic:.1f}  p-value {result.pvalue:.1e}")
   for name, a, b in zip(methods, training, this_week):
       print(f"{name:14} {a / sum(training):5.1%} -> {b / sum(this_week):5.1%}")

   $ python sol_categorical.py
   chi-square 572.6  p-value 4.5e-125
   card           70.0% -> 56.0%
   bank transfer  25.0% -> 19.0%
   wallet          5.0% -> 25.0%

   The chi-square test flags the change decisively, and the share table shows what it is: wallets went from 5% to 25% of payments, mostly at the expense of cards. The same large-sample caution applies as with KS, so alert on how much the shares moved. The share table is also what the person handling the alert needs to see.