12 - ML Deployment in Production: Serving and Monitoring

A model is finished when someone else can call it, and it keeps working while the world changes. The notebook that gave you 97% accuracy answers neither question: it does not say how the model reaches whoever needs it, nor how you notice that it has stopped working. Deployment is the work that takes the model out of the notebook; serving is the way it answers requests (online, one at a time, or in blocks); monitoring is the set of measurements that tell you, after release, whether it still answers well.

This article closes the ML Fundamentals series and uses a single small model to carry out each of these steps: saving an artifact, exposing it with a service that validates the inputs, measuring latency, producing a data-preparation error on purpose, observing two kinds of change in the data in production, and computing how many examples are needed to trust a new model.

Everything runs in a single process, with no listening server and no open ports: the service is tested with FastAPI's test client, which calls the application directly. There is no cloud and no cluster, and you do not need them to understand where a model breaks in production.

The model is the logistic regression of the third article on handwritten digits, and the metrics are those of the previous article. The infrastructure part (containers, model registry, release pipelines) has a separate place in the Machine Learning Engineer series: Dockerfile for an ML service, model registry and where the model runs.

From the notebook to the artifact

What you take to production is a file, the artifact, plus the information to use it without surprises. For a scikit-learn model the most common format is joblib, a format based on pickle. The scikit-learn documentation lists the risks: the file requires the same environment as the training one (a model saved with one version of scikit-learn is not guaranteed to load with another) and loading a pickle can execute arbitrary code.

As alternatives the guide cites ONNX, which allows serving without Python and which it defines as the safest option, but does not cover all models, and skops, which allows partly validating the content before loading it.

The listing saves a pipeline (the StandardScaler and the model together, for a reason you will see shortly), writes a metadata file with the library versions and the SHA-256 fingerprint of the artifact, reloads the model after checking the fingerprint and checks that the predictions coincide. The last part shows what «loading is executing» means: a harmless object, with a __reduce__ method that calls print, prints a line while it is being deserialized.

# A model to take to production: artifact, metadata, and what "loading" a pickle file means
import hashlib, json, tempfile
from pathlib import Path
import joblib
import numpy as np
import scipy, sklearn
from sklearn.datasets import load_digits
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

X, y = load_digits(return_X_y=True)
Xtr, Xte, ytr, yte = train_test_split(X / 16.0, y, test_size=0.3, random_state=0, stratify=y)
model = make_pipeline(StandardScaler(), LogisticRegression(max_iter=1000)).fit(Xtr, ytr)   # scaler AND model together

with tempfile.TemporaryDirectory() as d:
    path = Path(d) / "digits-v1.joblib"
    joblib.dump(model, path)
    meta = {"name": "digits", "version": "1", "sklearn": sklearn.__version__, "numpy": np.__version__,
            "scipy": scipy.__version__, "n_features": 64, "sha256": hashlib.sha256(path.read_bytes()).hexdigest()}
    (Path(d) / "digits-v1.json").write_text(json.dumps(meta))
    # reload as the service would, after checking that the fingerprint matches
    expected = json.loads((Path(d) / "digits-v1.json").read_text())["sha256"]
    assert hashlib.sha256(path.read_bytes()).hexdigest() == expected
    loaded = joblib.load(path)
    print("artifact size:", path.stat().st_size, "bytes")
print("saved metadata:", {k: v for k, v in meta.items() if k != "sha256"}, "+ sha256 (", len(meta["sha256"]), "characters )")
print("identical predictions before and after saving:", bool((loaded.predict(Xte) == model.predict(Xte)).all()))
print(f"test accuracy: {loaded.score(Xte, yte):.4f}")

# Loading a pickle means EXECUTING code: this harmless object prints a line while it is deserialized
import pickle
class Harmless:
    def __reduce__(self):
        return (print, ("  >>> code executed during pickle.loads",))
blob = pickle.dumps(Harmless())
print("about to call pickle.loads on a file I did not write:")
pickle.loads(blob)

# Output:
# artifact size: 8137 bytes
# saved metadata: {'name': 'digits', 'version': '1', 'sklearn': '1.9.0', 'numpy': '2.4.6', 'scipy': '1.17.1', 'n_features': 64} + sha256 ( 64 characters )
# identical predictions before and after saving: True
# test accuracy: 0.9722
# about to call pickle.loads on a file I did not write:
#   >>> code executed during pickle.loads

The artifact weighs 8,137 bytes and was reloaded with the same predictions (accuracy 0.9722 on the test set). The metadata say which versions it was created with: if in production you find others, scikit-learn warns with an InconsistentVersionWarning that you can turn into an error. The SHA-256 fingerprint proves that the file has not changed since you saved it, but not that it is safe: it only proves that it is the same.

The last printed line, «code executed during pickle.loads», is the reason why a pickle file of unknown origin is never loaded, as the Python documentation warns: the artifact your service loads must come only from a place you control.

A service that validates the input

A prediction service has three tasks: load the model only once, check that the input is admissible and answer. The listing writes it with FastAPI. The model is loaded in a lifecycle handler (lifespan), at startup, and not on every request, which would be expensive.

The input is described by a Pydantic model: a list of exactly 64 numbers, each between 0 and 1, with an explicit check on the range. FastAPI's TestClient client, used as a context manager with with, runs the lifecycle and calls the application in the process, without opening a port (the installed version warns that the use of httpx by TestClient is deprecated: it does not change the results). The JSON field names and the messages are in English in this version of the listing, so what the service answers is what you read in the output.

# A prediction service in FastAPI, tested inside the process with TestClient (no server, no port)
import tempfile
from contextlib import asynccontextmanager
from pathlib import Path
import joblib
from fastapi import FastAPI
from fastapi.testclient import TestClient
from pydantic import BaseModel, Field, field_validator
from sklearn.datasets import load_digits
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

X, y = load_digits(return_X_y=True)
Xtr, Xte, ytr, yte = train_test_split(X / 16.0, y, test_size=0.3, random_state=0, stratify=y)
tmp = tempfile.TemporaryDirectory()
ARTIFACT = Path(tmp.name) / "digits-v1.joblib"
joblib.dump(make_pipeline(StandardScaler(), LogisticRegression(max_iter=1000)).fit(Xtr, ytr), ARTIFACT)

state = {}

@asynccontextmanager
async def lifecycle(app):
    state["model"] = joblib.load(ARTIFACT)       # loaded ONCE at startup, not on every request
    yield
    state.clear()

app = FastAPI(title="digits", lifespan=lifecycle)

class Request(BaseModel):
    pixel: list[float] = Field(min_length=64, max_length=64)

    @field_validator("pixel")
    @classmethod
    def in_range(cls, v):
        if any(not 0.0 <= p <= 1.0 for p in v):
            raise ValueError("pixels must be in [0, 1]")
        return v

@app.get("/health")
def health():
    return {"status": "ok", "model": "digits-v1"}

@app.post("/predict")
def predict(r: Request):
    p = state["model"].predict_proba([r.pixel])[0]
    return {"class": int(p.argmax()), "probability": round(float(p.max()), 4), "model": "digits-v1"}

with TestClient(app) as client:                      # "with" runs the lifecycle: the model gets loaded
    print("GET /health ->", client.get("/health").status_code, client.get("/health").json())
    r = client.post("/predict", json={"pixel": Xte[0].tolist()})
    print("valid POST /predict ->", r.status_code, r.json(), " (true label:", int(yte[0]), ")")
    r = client.post("/predict", json={"pixel": [0.5] * 63})
    print("63 values instead of 64 ->", r.status_code, r.json()["detail"][0]["type"])
    r = client.post("/predict", json={"pixel": [2.0] + [0.0] * 63})
    print("one pixel outside [0, 1] ->", r.status_code, r.json()["detail"][0]["msg"])
    r = client.post("/predict", json={"pixel": ["a"] * 64})
    print("text instead of numbers ->", r.status_code, r.json()["detail"][0]["type"])
    correct = sum(client.post("/predict", json={"pixel": x.tolist()}).json()["class"] == t for x, t in zip(Xte[:100], yte[:100]))
    print("first 100 test images passed through the service: correct classes =", correct)

# Output:
# GET /health -> 200 {'status': 'ok', 'model': 'digits-v1'}
# valid POST /predict -> 200 {'class': 1, 'probability': 0.9525, 'model': 'digits-v1'}  (true label: 1 )
# 63 values instead of 64 -> 422 too_short
# one pixel outside [0, 1] -> 422 Value error, pixels must be in [0, 1]
# text instead of numbers -> 422 float_parsing
# first 100 test images passed through the service: correct classes = 98

The valid request receives 200 and class 1 with probability 0.9525, which is the true label. The three invalid requests do not reach the model: 63 values instead of 64 give 422 with the error too_short, a pixel at 2.0 gives 422 with the message of your check, and text instead of numbers gives 422 with float_parsing.

Without this layer, a list of 63 values would reach the pipeline and raise an internal exception, a 500 error that hides the real problem. Passing 100 test images through the service, 98 are classified correctly. This is a minimal service: it lacks authentication, it lacks rate limiting and it lacks log handling. To really expose it you need an ASGI server such as uvicorn, which I did not start here so as not to open ports.

Latency: one request at a time against a block

Latency is the time between the request and the response; for an online service one usually looks at the high percentile (the 99th percentile, that is the time under which 99 requests out of 100 fall) more than at the mean. In a service that receives requests in blocks (batch), what counts instead is how much work is done per second.

The listing measures the same model in three ways: 500 calls of one row, one call of 500 rows and blocks of different sizes. The times are measurements and change from one run to the next, because they depend on the machine and on the load: I measured them in five runs on an AMD Ryzen 5 5600H with 2 threads and I report the ranges. Since this listing measures times it was not rerun after translating, so the labels it prints are still in Italian: «stesse probabilita'» is «same probabilities», «chiamate singole» «single calls», «in tutto» «in total», «a richiesta» «per request», «a riga» «per row», «rapporto fra i costi per riga» «ratio between the costs per row» and «blocco da N righe» «block of N rows».

# Latency: one prediction at a time against blocks (TIMES: they vary from one run to the next)
import time, tempfile
from pathlib import Path
import joblib
import numpy as np
from sklearn.datasets import load_digits
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

X, y = load_digits(return_X_y=True)
Xtr, Xte, ytr, yte = train_test_split(X / 16.0, y, test_size=0.3, random_state=0, stratify=y)
m = make_pipeline(StandardScaler(), LogisticRegression(max_iter=1000)).fit(Xtr, ytr)
righe = Xte[:500]
m.predict_proba(righe[:5])                                           # warm-up

t0 = time.perf_counter()
singole = [m.predict_proba(r.reshape(1, -1))[0] for r in righe]      # 500 calls of one row each
t_singole = time.perf_counter() - t0
t0 = time.perf_counter()
blocco = m.predict_proba(righe)                                      # one call of 500 rows
t_blocco = time.perf_counter() - t0
print("stesse probabilita':", np.allclose(np.array(singole), blocco))
print(f"500 chiamate singole: {t_singole * 1000:7.1f} ms in tutto, {t_singole / 500 * 1000:.3f} ms a richiesta")
print(f"1 chiamata da 500:    {t_blocco * 1000:7.1f} ms in tutto, {t_blocco / 500 * 1000:.3f} ms a riga")
print(f"rapporto fra i costi per riga: {t_singole / t_blocco:.0f}x")

# Block size: time per row
for n in (1, 10, 100, 1000):
    reps = max(1, 2000 // n)
    t0 = time.perf_counter()
    for _ in range(reps):
        m.predict_proba(Xte[:n])
    print(f"blocco da {n:4d} righe: {(time.perf_counter() - t0) / reps / n * 1000:.4f} ms a riga")

# Output:
# stesse probabilita': True
# 500 chiamate singole:   127.5 ms in tutto, 0.255 ms a richiesta
# 1 chiamata da 500:        0.5 ms in tutto, 0.001 ms a riga
# rapporto fra i costi per riga: 262x
# blocco da    1 righe: 0.2351 ms a riga
# blocco da   10 righe: 0.0204 ms a riga
# blocco da  100 righe: 0.0024 ms a riga
# blocco da 1000 righe: 0.0004 ms a riga

The probabilities are the same in the two ways. Five hundred single calls took 100-135 ms in total, that is about 0.20-0.27 ms per request, while one call of 500 rows took 0.4-0.5: the cost per row is 240-300 times lower. It is not the computation that is slow, it is the fixed cost of every call (validation, array conversion, Python calls), which for a block is paid only once.

The cost per row falls from 0.20-0.22 ms with one row to 0.0004-0.0005 ms with a thousand. In an online service gathering several requests into a block (micro-batching) can greatly increase the number of requests per second, but it adds the wait needed to fill the block: it is a compromise between the latency of each request and the overall throughput. This measurement concerns a linear model on a laptop; for another model the numbers change and must be measured again.

The silent bug: training-serving skew

Training-serving skew is the difference between how the data are prepared in training and how they are in production. Google's «Rules of Machine Learning» guide defines it as a difference between performance during training and during serving, which can arise from a difference in how data are handled in the two pipelines, from a change in the data between the two moments or from a feedback loop between model and algorithm, and recommends measuring it explicitly.

A classic preparation error is to reimplement the standardization in the service: forgetting it, or computing it on the block of requests instead of with the training statistics. The listing tries it on the breast cancer dataset (30 features with very different scales), with three variants and a model that expects standardized inputs.

# Training-serving skew: the same model, three ways of preparing the inputs in production
import numpy as np
from sklearn.datasets import load_breast_cancer
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

X, y = load_breast_cancer(return_X_y=True)
Xtr, Xte, ytr, yte = train_test_split(X, y, test_size=0.3, random_state=0, stratify=y)

scaler = StandardScaler().fit(Xtr)
lr = LogisticRegression(max_iter=2000).fit(scaler.transform(Xtr), ytr)        # the model sees standardized data
pipeline = make_pipeline(scaler, lr)                                          # the scaler travels with the model
print(f"positive class in the test set: {yte.mean():.3f}   test samples: {len(yte)}")

print(f"A) whole pipeline saved:                         accuracy {pipeline.score(Xte, yte):.3f}")
print(f"B) the service forgets the scaling:              accuracy {lr.score(Xte, yte):.3f}")

def with_batch_statistics(batch):
    """Typical error: standardizing the requests with the mean and deviation of the received BATCH, not of the training set."""
    return (batch - batch.mean(axis=0)) / (batch.std(axis=0) + 1e-9)

for n in (1, 5, 50):
    correct = 0
    for i in range(0, len(Xte) - len(Xte) % n, n):
        correct += (lr.predict(with_batch_statistics(Xte[i:i + n])) == yte[i:i + n]).sum()
    tot = len(Xte) - len(Xte) % n
    print(f"C) batch statistics, requests in batches of {n:2d}: accuracy {correct / tot:.3f}")

# Output:
# positive class in the test set: 0.626   test samples: 171
# A) whole pipeline saved:                         accuracy 0.959
# B) the service forgets the scaling:              accuracy 0.374
# C) batch statistics, requests in batches of  1: accuracy 0.626
# C) batch statistics, requests in batches of  5: accuracy 0.847
# C) batch statistics, requests in batches of 50: accuracy 0.940

With the whole pipeline saved the accuracy is 0.959. If the service forgets the scaling it falls to 0.374, less than chance (62.6% of the samples are of the positive class: always predicting «positive» would give 0.626, always «negative» 0.374): the model sees numbers on the wrong scale. Error C is subtle: standardizing with the block statistics works more or less if the block is large (0.940 with 50 requests), worse with 5 (0.847) and fails with a single request, because the deviation is zero and all the values become zero, so the model always answers the same class, with accuracy 0.626.

No exception, no 500 error: the service answers well, with wrong numbers. The defense is the rule of this section: save the data preparation with the model, in the same pipeline, and verify with a test that the service gives the same predictions as the training environment on the same examples.

Monitoring: inputs, outputs and outcomes

After release the data change. Sculley and colleagues, in «Hidden Technical Debt in Machine Learning Systems», argue that tests are not enough in a changing world and that real-time monitoring of the system's behavior is needed; among the starting points they indicate the comparison between the distribution of the predicted labels and that of the observed labels (prediction bias).

In common language the changes are of two kinds: covariate shift, in which the distribution of the inputs changes, and concept drift, in which the link between inputs and result changes. To compare the inputs you use the two-sample Kolmogorov-Smirnov test (in SciPy, ks_2samp) or the PSI (Population Stability Index), $\sum_i (a_i - e_i)\ln(a_i/e_i)$ on the shares $e_i$ of the reference and $a_i$ of the current data in quantile bins: the PSI is a commonly used measure and not a standard of any regulation.

The listing trains a logistic regression on synthetic data with six features, and then «puts it into production» on three batches of 5,000 rows: one with no changes, one in which feature 0 shifts by one standard deviation (covariate shift), and one with the same inputs but with 30% of the outcomes flipped (concept drift).

For each batch it prints the KS statistic and its p-value, the PSI of feature 0, the accuracy and the shares of predicted and observed class 1. In production the accuracy requires labels, which arrive late or never: here we know them because it is a simulation.

# Monitoring: what the inputs see and what the accuracy sees (with the labels, which in production arrive late)
import numpy as np
from scipy.stats import ks_2samp
from sklearn.datasets import make_classification
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split

X, y = make_classification(n_samples=30000, n_features=6, n_informative=4, n_redundant=0, class_sep=1.2,
                           flip_y=0.02, random_state=0)
Xtr, Xrest, ytr, yrest = train_test_split(X, y, test_size=0.8, random_state=0)
model = LogisticRegression(max_iter=1000).fit(Xtr, ytr)
ref = Xtr[:, 0]                                                     # reference distribution (training), feature 0

def psi(reference, current, bins=10):
    """Population Stability Index: sum of (a - e) * ln(a / e) over quantile bins of the reference."""
    cuts = np.quantile(reference, np.linspace(0, 1, bins + 1)[1:-1])
    e = np.bincount(np.searchsorted(cuts, reference), minlength=bins) / len(reference)
    a = np.bincount(np.searchsorted(cuts, current), minlength=bins) / len(current)
    e, a = np.clip(e, 1e-6, None), np.clip(a, 1e-6, None)
    return float(((a - e) * np.log(a / e)).sum())

batches = {}
batches["no change"] = (Xrest[:5000], yrest[:5000])
Xs = Xrest[5000:10000].copy(); Xs[:, 0] += 1.0                      # feature 0 shifts by +1 standard deviation
batches["covariate shift (feature 0 + 1)"] = (Xs, yrest[5000:10000])
Xc, yc = Xrest[10000:15000], yrest[10000:15000].copy()              # same inputs, but in 30% of cases the outcome flips
rng = np.random.default_rng(0); flip = rng.random(len(yc)) < 0.30
yc[flip] = 1 - yc[flip]
batches["concept drift (30% of outcomes flipped)"] = (Xc, yc)

print(f"{'production batch':42s} {'KS feat.0':>9s} {'p-value':>9s} {'PSI feat.0':>10s} {'accuracy':>12s} {'class 1 predicted':>18s} {'class 1 observed':>19s}")
for name, (Xl, yl) in batches.items():
    ks = ks_2samp(ref, Xl[:, 0])
    print(f"{name:42s} {ks.statistic:9.3f} {ks.pvalue:9.2e} {psi(ref, Xl[:, 0]):10.3f} {model.score(Xl, yl):12.3f}"
          f" {model.predict(Xl).mean():18.3f} {yl.mean():19.3f}")
print(f"accuracy on the reference (training test set, 5000 rows): {model.score(Xrest[:5000], yrest[:5000]):.3f}")

# Output:
# production batch                           KS feat.0   p-value PSI feat.0     accuracy  class 1 predicted    class 1 observed
# no change                                      0.015  5.98e-01      0.002        0.747              0.453               0.489
# covariate shift (feature 0 + 1)                0.213 2.63e-108      0.389        0.731              0.592               0.499
# concept drift (30% of outcomes flipped)        0.014  6.12e-01      0.003        0.603              0.461               0.504
# accuracy on the reference (training test set, 5000 rows): 0.747

The batch with no changes gives KS 0.015 (p-value 0.60), PSI 0.002 and accuracy 0.747. The covariate shift is very evident in the inputs: KS 0.213 with a p-value of 1e-108 and PSI 0.389, but the accuracy only drops to 0.731, because the linear model copes well with the shift. The concept drift is invisible in the inputs (KS 0.014, p-value 0.61, PSI 0.003, identical to the stable batch) but makes the accuracy collapse to 0.603.

Two opposite lessons, from the same experiment. First: an alarm on the inputs is not an alarm on the quality of the model, and with thousands of rows the p-value becomes tiny even for irrelevant shifts, so look at the magnitude (KS, PSI) as well as the p-value. Second: no check on the inputs sees the concept drift, which is discovered only with the labels; the share of predicted class 1 against the observed one gives 0.461 against 0.504, a gap of 4.3 points, almost equal to that of the stable batch (3.6), and it is not a useful signal here. You need both levels, and a way to obtain labels, even partial ones.

Releasing a new model: how many examples are needed

A new model, which on the offline test is one percentage point better, must be compared with the one in production on real traffic: first in shadow (it receives the same requests but its answers are not used), then with a share of traffic. But one point of accuracy out of 90 is a small difference, and an experiment with too few examples does not see it.

For a comparison between two proportions, with a significance level of 5% and power of 80%, the number of requests per arm is $n = \left(z_{\alpha/2}\sqrt{2\overline{p}\overline{q}} + z_{\beta}\sqrt{p_1 q_1 + p_2 q_2}\right)^2 / (p_2 - p_1)^2$. The listing computes $n$ for 0.90 against 0.91 and checks the formula with a simulation: 2,000 experiments for each size, and the fraction in which a two-proportion z test recognizes the difference (the power). It also checks the false alarms with two identical models.

# How many examples are needed to notice that model B is better than A by one point of accuracy?
import numpy as np
from scipy.stats import norm

p1, p2, alpha, power = 0.90, 0.91, 0.05, 0.80
za, zb = norm.ppf(1 - alpha / 2), norm.ppf(power)
pm = (p1 + p2) / 2
n = (za * np.sqrt(2 * pm * (1 - pm)) + zb * np.sqrt(p1 * (1 - p1) + p2 * (1 - p2))) ** 2 / (p2 - p1) ** 2
n = int(np.ceil(n))
print(f"formula: about {n} requests PER ARM ({2 * n} in total) are needed to see {p1:.2f} -> {p2:.2f} (alpha {alpha}, power {power})")

def two_proportion_z(a, b, n):
    """a and b are the successes of the two arms (n requests each)."""
    pp = (a + b) / (2 * n)
    return (b / n - a / n) / np.sqrt(2 * pp * (1 - pp) / n)

def significant_fraction(n, p_a, p_b, trials=2000, seed=0):
    rng = np.random.default_rng(seed)
    z = two_proportion_z(rng.binomial(n, p_a, trials), rng.binomial(n, p_b, trials), n)
    return float((np.abs(z) > za).mean())

for k in (n // 4, n // 2, n, 2 * n):
    print(f"n per arm = {k:6d}: the experiment detects the difference in {significant_fraction(k, p1, p2):.1%} of cases")
print(f"if the two models were equal (0.90 against 0.90), false alarms with n = {n}: {significant_fraction(n, 0.90, 0.90, seed=1):.1%}")

# Output:
# formula: about 13495 requests PER ARM (26990 in total) are needed to see 0.90 -> 0.91 (alpha 0.05, power 0.8)
# n per arm =   3373: the experiment detects the difference in 28.3% of cases
# n per arm =   6747: the experiment detects the difference in 50.2% of cases
# n per arm =  13495: the experiment detects the difference in 80.4% of cases
# n per arm =  26990: the experiment detects the difference in 98.1% of cases
# if the two models were equal (0.90 against 0.90), false alarms with n = 13495: 5.1%

About 13,495 requests per arm are needed, 26,990 in all. The simulation confirms the formula: with this $n$ the experiment recognizes the difference in 80.4% of cases, with a quarter of the $n$ only in 28.3%, with half in 50.2%. If the two models were identical, the test gives a false alarm in 5.1% of cases, like the declared significance level. If the service receives, hypothetically, 5,000 requests a day and you send 10% to the new model, in a month you have about 15,000 requests on the new arm, the same order of magnitude as the number per arm computed; with a 10/90 split the formula must be redone, but the order of magnitude remains, and you need the labels.

This is the reason why increases of one point of accuracy are hard to demonstrate in production, and why it is worth deciding in advance which difference you care about and with what traffic you can detect it.

Worked exercise

On a single feature the reference (training) has shares 0.5, 0.3 and 0.2 in three bins; in production the shares are 0.4, 0.3 and 0.3. Compute the PSI by hand. Then say whether the value, by itself, is enough for you to decide whether to retrain the model.

# Exercise: PSI by hand on three bins
import numpy as np
expected = np.array([0.5, 0.3, 0.2])
actual = np.array([0.4, 0.3, 0.3])
terms = (actual - expected) * np.log(actual / expected)
print("terms per bin:", np.round(terms, 5).tolist())
print("PSI =", round(float(terms.sum()), 5))

# Output:
# terms per bin: [0.02231, 0.0, 0.04055]
# PSI = 0.06286

What you can do now, and what you have covered

With this article you can do five things, with a measurement run for each:

• save a model with its metadata and its fingerprint, and know why you do not load a pickle of unknown origin;
• expose it with a service that validates the input and test it without opening a port;
• measure latency and understand when to gather requests into blocks;
• recognize training-serving skew and prevent it by saving the data preparation with the model;
• monitor inputs and outcomes, distinguishing covariate shift and concept drift, and compute how many examples are needed for a comparison between models.

The five abilities have in common that none gives a result alone: the artifact without metadata breaks at the first version change, the service without validation hides errors, the monitoring without labels does not see concept drift. The distance between the notebook and production is made of these small checks, each of which you can run in a minute, as you did here.

The ML Fundamentals series ends here. The twelve stages, in the order in which they build on one another:

1. Introduction to Machine Learning: taxonomy of problems and workflow;
2. Linear regression from scratch: the first model, written in Python;
3. Logistic regression: from the line to binary classification;
4. Bias-variance tradeoff: why a more flexible model is not always better;
5. Cross-validation: estimating performance and choosing the model;
6. Decision trees and Random Forest: non-linear rule-based models;
7. Support Vector Machines: maximum margin and kernels;
8. K-Means and hierarchical clustering: grouping without labels;
9. PCA: reducing dimensions;
10. Feedforward neural networks: from the perceptron to backpropagation;
11. Metrics and evaluation: judging a model with the right question;
12. this article: taking the model to production and continuing to measure it.

The order is not random: the first five stages teach how to build and judge a simple model (regression, classification, error, validation), stages six to ten widen the family of models (trees, margins, groups, principal axes, networks), and the last two close the cycle with evaluation and production. If you want to continue outside the series, the Machine Learning Engineer series retraces many of the same topics with more detail on the engineering work, including the stages on experiment tracking and on testing models.