6.2 Splitting the data

Every claim about how well a model predicts is only as honest as the data it was tested on. Those data must be kept apart from everything the model and its builder learned from, must resemble the data the model will meet, and must be used once.

Code: data and style used in this chapter

import matplotlib.pyplot as plt
import numpy as np
import statsmodels.api as sm
from matplotlib.axes import Axes
from sklearn import datasets
from sklearn.ensemble import GradientBoostingRegressor
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import GroupShuffleSplit, train_test_split
from sklearn.neighbors import KNeighborsClassifier, KNeighborsRegressor
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

from statly import style

style.use()
rng = np.random.default_rng(1)

diabetes = datasets.load_diabetes()  # 442 patients: 10 measurements, disease progression
sun = sm.datasets.sunspots.load_pandas().data  # sunspots per year, 1700 to 2008
PART_COLORS = {"training": style.BLUE, "validation": style.GOLD, "test": style.CORAL}


def bar_of_parts(
    ax: Axes, parts: list[tuple[str, float, str]], y: float = 0, height: float = 0.5
) -> None:
    """One horizontal bar cut into labelled, coloured parts."""
    left = 0.0
    for label, width, color in parts:
        ax.barh(y, width, left=left, color=color, alpha=0.6, height=height)
        ax.text(left + width / 2, y, label, ha="center", va="center", fontsize=8.5)
        left += width

6.2.1 Three jobs, three parts

Building a predictive model involves three different jobs, and each needs its own data:

Fitting: estimating the model’s coefficients, splits or weights. This uses the training set.
Choosing: picking the method, the predictors, the amount of smoothing or penalty, the number of trees, when to stop training. This uses the validation set, or cross-validation within the training data (Chapter 6.3).
Reporting: saying how well the final, chosen model will do on new data. This uses the test set, once, at the very end.

Code

STEPS = [
    ("1  split", style.GREY),
    ("2  explore and fit\non training", style.BLUE),
    ("3  choose on\nvalidation", style.GOLD),
    ("4  refit on\ntraining +\nvalidation", style.PURPLE),
    ("5  test once", style.CORAL),
    ("6  use", style.GREEN),
]

fig, ax = plt.subplots(figsize=(10, 3.4))
bar_of_parts(
    ax,
    [
        ("training 60%", 0.6, style.BLUE),
        ("validation 20%", 0.2, style.GOLD),
        ("test 20%", 0.2, style.CORAL),
    ],
    y=1.2,
    height=0.45,
)
for i, (label, color) in enumerate(STEPS):
    x = i / len(STEPS) + 0.005
    ax.add_patch(plt.Rectangle((x, -0.2), 1 / len(STEPS) - 0.03, 0.75, color=color, alpha=0.35))
    ax.text(x + 0.5 / len(STEPS) - 0.015, 0.17, label, ha="center", va="center", fontsize=8.5)
    if i:
        ax.annotate(
            "",
            xy=(x, 0.17),
            xytext=(x - 0.025, 0.17),
            arrowprops=dict(arrowstyle="->", color=style.INK),
        )
ax.set(xlim=(0, 1), ylim=(-0.4, 1.6))
ax.set_axis_off()
plt.show()

Figure 6.2.1: Top: a typical three-way split of the data. Bottom: the order of work. The split comes first, before any exploring; the test set is locked away until the last step. After choosing, the final model is usually refitted on training and validation data together, and only then tested.

Why not choose and report on the same data? Because choosing the best of several candidates picks, in part, the one that was luckiest on those particular cases. Its score on the data used to choose it is too good, the winner’s curse of Chapter 3.10:

Code

def winner_scores(candidates: int = 20, n: int = 60) -> tuple[float, float]:
    """Validation and test accuracy of the best of several random classifiers."""
    truth = rng.integers(0, 2, n)
    validation = rng.integers(0, 2, (candidates, n)) == truth
    best = int(np.argmax(validation.mean(axis=1)))
    test = rng.integers(0, 2, n) == rng.integers(0, 2, n)
    return validation[best].mean(), test.mean()


scores = np.array([winner_scores() for _ in range(200)])
fig, ax = plt.subplots(figsize=(8, 3.1))
bins = np.linspace(0.3, 0.9, 25)
ax.hist(
    scores[:, 0],
    bins=bins,
    color=style.CORAL,
    alpha=0.6,
    label=f"validation: {scores[:, 0].mean():.0%} on average",
)
ax.hist(
    scores[:, 1],
    bins=bins,
    color=style.BLUE,
    alpha=0.6,
    label=f"fresh test: {scores[:, 1].mean():.0%}",
)
ax.axvline(0.5, color=style.INK, ls="--", lw=1)
ax.xaxis.set_major_formatter(style.PERCENT)
ax.legend(fontsize=8)
ax.set(xlabel="accuracy of the chosen classifier (truth: 50%)", ylabel="repetitions")
plt.show()

Figure 6.2.2: Twenty useless classifiers (each guesses at random, so the true accuracy of every one is 50%) compared on a validation set of 60 cases, 200 times over. The winner’s validation accuracy (red) is far above 50%; on a fresh test set (blue) it falls back to 50%. With real models the effect is smaller, but always in the same direction.

The rule follows: every look at the test set is a choice. If the test result sends you back to change the model, the test set has become a validation set, and you need a new test set.

6.2.2 How big should each part be?

Two forces pull in opposite directions. A larger training set gives a better model; a larger test set gives a more precise measurement of it. The precision of a test accuracy p measured on n cases follows from the binomial (Chapter 2.2):

$\text{SE}(\hat p) = \sqrt{\frac{p(1-p)}{n}}$

Code

X_pool, X_hold, y_pool, y_hold = train_test_split(
    diabetes.data, diabetes.target, test_size=142, random_state=0
)
SIZES = np.array([25, 50, 100, 150, 200, 250, 300])


def learning_curve(size: int, draws: int = 20) -> tuple[float, float]:
    """Mean training and held-out error (root mean squared) for a given training size."""
    train_err, test_err = [], []
    for seed in range(draws):
        rows = np.random.default_rng(seed).choice(len(y_pool), size, replace=False)
        model = GradientBoostingRegressor(
            n_estimators=100, max_depth=2, learning_rate=0.05, random_state=0
        )
        model.fit(X_pool[rows], y_pool[rows])
        train_err.append(np.sqrt(np.mean((model.predict(X_pool[rows]) - y_pool[rows]) ** 2)))
        test_err.append(np.sqrt(np.mean((model.predict(X_hold) - y_hold) ** 2)))
    return np.mean(train_err), np.mean(test_err)


curve = np.array([learning_curve(n) for n in SIZES])
fig, axes = style.subplots(2, height=3.2, width=10)
axes[0].plot(SIZES, curve[:, 1], "o-", color=style.CORAL, lw=2, label="held-out error")
axes[0].plot(SIZES, curve[:, 0], "o-", color=style.BLUE, lw=2, label="training error")
axes[0].legend(fontsize=8)
axes[0].set(
    xlabel="training patients", ylabel="error (root mean squared)", title="More training data"
)
test_sizes = np.geomspace(25, 2000, 60)
half = 1.96 * np.sqrt(0.8 * 0.2 / test_sizes)
axes[1].fill_between(test_sizes, 0.8 - half, 0.8 + half, color=style.GOLD, alpha=0.4)
axes[1].axhline(0.8, color=style.INK, lw=1)
axes[1].set_xscale("log")
axes[1].yaxis.set_major_formatter(style.PERCENT)
axes[1].set(xlabel="test cases (log scale)", ylabel="measured accuracy", title="More test data")
fig.tight_layout()
plt.show()

Figure 6.2.3: Left: a learning curve. Error of boosted trees on 142 held-out patients, trained on growing parts of the other 300 (average of 20 random draws, with the training error for comparison). More training data still helps at 300. Right: how precisely a test set measures an accuracy of 80%: the 95% range of the measured accuracy for test sets of 25 to 2,000 cases. With 100 test cases, a measured 80% is compatible with anything from about 72% to 88%.

Common choices are 60/20/20 or 70/15/15 for training, validation and test, or 80/20 when cross-validation on the training part does the choosing. With a few hundred cases, a test set is too small to say much, and nested cross-validation (Chapter 6.3) is the better tool. With millions of cases, a test set of 1% is plenty.

6.2.3 Keep the class shares: stratified splits

When one class is rare, a random split can leave the test set with very few of them, or very many, by chance. A stratified split cuts each class separately in the same proportions, so every part has the same share of each class.

Code

labels = np.r_[np.ones(10, int), np.zeros(190, int)]


def positives_in_test(stratify: bool, seed: int) -> int:
    """Number of positive cases in a 25% test set."""
    _, test = train_test_split(
        labels, test_size=0.25, random_state=seed, stratify=labels if stratify else None
    )
    return int(test.sum())


counts = {
    name: np.array([positives_in_test(name == "stratified", s) for s in range(2000)])
    for name in ("random", "stratified")
}
fig, axes = style.subplots(2, height=3, width=10)
for name, color, offset in (("random", style.GREY, -0.2), ("stratified", style.BLUE, 0.2)):
    values, freq = np.unique(counts[name], return_counts=True)
    axes[0].bar(values + offset, freq / 2000, width=0.4, color=color, alpha=0.7, label=name)
axes[0].yaxis.set_major_formatter(style.PERCENT)
axes[0].legend(fontsize=8)
axes[0].set(
    xlabel="positives in the test set", ylabel="share of splits", title="How many positives?"
)
found = [rng.binomial(n, 0.7) / n for n in counts["random"] if n > 0]
axes[1].hist(found, bins=np.linspace(0, 1, 21), color=style.GREY, alpha=0.7)
axes[1].axvline(0.7, color=style.INK, ls="--", lw=1)
axes[1].xaxis.set_major_formatter(style.PERCENT)
axes[1].set(xlabel="measured sensitivity (truth 70%)", ylabel="splits", title="Random split")
fig.tight_layout()
plt.show()

Figure 6.2.4: 200 cases of which 10 are positive (5%), split 75/25 at random 2,000 times. Left: the number of positives that land in the 50 test cases. A random split gives anywhere from none to six or more; a stratified split always gives 2 or 3. Right: what that does to a measured sensitivity, here for a classifier that truly finds 70% of positives; with only a couple of positives, the test set can say nothing useful about it.

6.2.4 Keep groups together

Rows are often not independent: several scans of the same patient, several visits of the same customer, several measurements from the same machine. If rows of one patient land in both training and test, the model can recognise the patient instead of learning the disease, and the test score measures recognition. The fix is a grouped split: all rows of a patient go to the same side.

Code

PATIENTS, VISITS = 300, 10
local = np.random.default_rng(0)
signature = local.normal(0, 3, (PATIENTS, 5))  # what makes each patient recognisable
risk = signature[:, 0] * 0.3 + local.normal(0, 1, PATIENTS)  # partly related to the signature
patient = np.repeat(np.arange(PATIENTS), VISITS)
features = signature[patient] + local.normal(0, 1, (PATIENTS * VISITS, 5))
outcome = (risk[patient] > 0).astype(int)  # the same outcome on every visit


def new_patients_accuracy(model_rows: np.ndarray) -> float:
    """Accuracy on 2,000 freshly simulated patients, one visit each."""
    fresh = np.random.default_rng(9)
    sig = fresh.normal(0, 3, (2000, 5))
    y_new = (sig[:, 0] * 0.3 + fresh.normal(0, 1, 2000) > 0).astype(int)
    model = KNeighborsClassifier(5).fit(features[model_rows], outcome[model_rows])
    return model.score(sig + fresh.normal(0, 1, (2000, 5)), y_new)


def split_accuracy(grouped: bool) -> tuple[float, np.ndarray, np.ndarray]:
    """Test accuracy of 5 nearest neighbours under a random or a grouped 75/25 split."""
    rows = np.arange(len(outcome))
    if grouped:
        train, test = next(
            GroupShuffleSplit(test_size=0.25, random_state=0).split(rows, groups=patient)
        )
    else:
        train, test = train_test_split(rows, test_size=0.25, random_state=0)
    model = KNeighborsClassifier(5).fit(features[train], outcome[train])
    return model.score(features[test], outcome[test]), train, test


random_acc, random_train, _ = split_accuracy(False)
grouped_acc, grouped_train, _ = split_accuracy(True)
fresh_acc = new_patients_accuracy(grouped_train)


def plot_assignment(ax: Axes, train: np.ndarray, y: float, title: str) -> None:
    """For ten patients, each visit as a square: filled if training, open if test."""
    in_train = np.isin(np.arange(len(outcome)), train)
    for p in range(6):
        for v in range(VISITS):
            row = p * VISITS + v
            color = plt.cm.tab10(p)
            ax.scatter(
                p * 11 + v,
                y,
                marker="s",
                s=45,
                color=color if in_train[row] else "none",
                edgecolors=color,
            )
    ax.text(-2, y + 0.35, title, fontsize=9)


fig = plt.figure(figsize=(10.5, 3.2))
grid = fig.add_gridspec(1, 2, width_ratios=[1.6, 1])
ax = fig.add_subplot(grid[0])
plot_assignment(ax, random_train, 1, "random split (filled: training, open: test)")
plot_assignment(ax, grouped_train, 0, "grouped split")
ax.set(ylim=(-0.5, 1.7), xlim=(-3, 68))
ax.set_axis_off()
ax = fig.add_subplot(grid[1])
results = [
    ("random split", random_acc),
    ("grouped split", grouped_acc),
    ("new patients", fresh_acc),
]
ax.barh(
    [r[0] for r in results],
    [r[1] for r in results],
    color=[style.GREY, style.BLUE, style.GREEN],
    alpha=0.7,
)
for i, (_, value) in enumerate(results):
    ax.text(value + 0.01, i, f"{value:.0%}", va="center", fontsize=9)
ax.invert_yaxis()
ax.xaxis.set_major_formatter(style.PERCENT)
ax.set(xlim=(0, 1.1), xlabel="accuracy")
fig.tight_layout()
plt.show()

Figure 6.2.5: Simulated data: 300 patients with 10 visits each. Every patient has a strong personal signature in the measurements, and an outcome that depends only partly on it. Left: in a random split (top) the same patient (one colour per patient, six shown) appears on both sides; in a grouped split (bottom) each patient is on one side only. Right: accuracy of a nearest-neighbour classifier. The random split reports far higher accuracy than the model achieves on new patients; the grouped split is much closer to the truth.

The question to ask is: what will be new when the model is used? If it will meet new patients, the test set must contain new patients. If it will meet new hospitals, the test set must contain hospitals it has not seen.

6.2.5 Time: train on the past, test on the future

When data come in time order, the model will be used on the future. A random split lets it peek: it trains on the days around each test day and borrows their information. A chronological split trains on everything up to a date and tests on what comes after. A gap between the two (an embargo) stops information leaking across the boundary when the predictors look back in time.

Code

LAGS = 3
series = sun.SUNACTIVITY.to_numpy()
years = sun.YEAR.to_numpy()[LAGS:]
lagged = np.column_stack([series[i : len(series) - LAGS + i] for i in range(LAGS)])
target = series[LAGS:]
cut = int(0.75 * len(target))


def holdout_error(test: np.ndarray) -> float:
    """Root mean squared error on the held-out years."""
    train = ~test
    model = KNeighborsRegressor(5).fit(lagged[train], target[train])
    return float(np.sqrt(np.mean((model.predict(lagged[test]) - target[test]) ** 2)))


random_test = np.zeros(len(target), bool)
random_test[rng.choice(len(target), len(target) - cut, replace=False)] = True
future_test = np.arange(len(target)) >= cut
errors = {
    "random years": holdout_error(random_test),
    "last 25% of years": holdout_error(future_test),
}

fig, axes = plt.subplots(2, 1, figsize=(9, 4), gridspec_kw={"height_ratios": [1.5, 1]})
for row, (name, test) in enumerate((("random", random_test), ("chronological", future_test))):
    axes[0].scatter(years[~test], np.full((~test).sum(), row), s=4, color=style.BLUE)
    axes[0].scatter(years[test], np.full(test.sum(), row), s=4, color=style.CORAL)
axes[0].set(
    yticks=[0, 1], yticklabels=["random", "chronological"], ylim=(-0.6, 1.6), xlabel="year"
)
axes[1].barh(list(errors), list(errors.values()), color=[style.GREY, style.BLUE], alpha=0.7)
for i, value in enumerate(errors.values()):
    axes[1].text(value + 0.3, i, f"{value:.1f}", va="center", fontsize=9)
axes[1].set(xlabel="error on held-out years (sunspots)")
fig.tight_layout()
plt.show()

Figure 6.2.6: Predicting next year’s sunspots from the previous three years with 5 nearest neighbours. Top: a random 25% of the years held out (red) against the last 25% held out. Bottom: the error. The random split, with test years scattered between training years, promises a smaller error than the model achieves on the years after its training data.

6.2.6 Leakage

Leakage is any way information from the test cases, or from the outcome itself, reaches the model during training. It makes test scores look better than real use will be, and it is the most common reason a model that “worked” fails in practice. The usual forms:

Preprocessing on all the data: scaling, imputing, selecting predictors or balancing classes before the split. Every such step must be learned on the training part only; a pipeline that is fitted as a whole makes this automatic. Section 6.3.7 shows how selecting predictors on all the data turns pure noise into apparent skill.
Copies across the split: duplicate rows, augmented copies of the same image, or oversampled rare cases made before splitting put near-identical cases on both sides.
Target leakage: a predictor that is only known after the outcome, or is a consequence of it, such as “was given antibiotics” when predicting infection.

Code

def make_data(n: int, seed: int) -> tuple[np.ndarray, np.ndarray]:
    """Five noisy predictors and a yes/no outcome that depends on them."""
    local = np.random.default_rng(seed)
    X = local.normal(size=(n, 5))
    y = (X @ np.array([1.0, 0.8, 0.5, 0, 0]) + local.normal(0, 1.2, n) > 0).astype(int)
    return X, y


def honest_and_leaky(kind: str) -> tuple[float, float]:
    """Accuracy on fresh data (honest) and on a leaky test split, for one kind of leakage."""
    X, y = make_data(600, 0)
    X_new, y_new = make_data(2000, 1)
    model = KNeighborsClassifier(1) if kind != "target" else LogisticRegression()
    if kind == "duplicates":
        copies = np.random.default_rng(2).choice(len(y), 240, replace=False)
        X_all, y_all = np.r_[X, X[copies]], np.r_[y, y[copies]]
    elif kind == "oversampling":
        rare = np.flatnonzero(y == 1)[:60]
        X_all, y_all = (
            np.r_[X[y == 0][:300], np.repeat(X[rare], 5, axis=0)],
            np.r_[np.zeros(300, int), np.ones(300, int)],
        )
        X_new, y_new = (
            np.r_[X_new[y_new == 0][:900], X_new[y_new == 1][:900]],
            np.r_[np.zeros(900, int), np.ones(900, int)],
        )
    else:
        X_all, y_all = X, y
    if kind == "target":
        recorded_after = y_all + np.random.default_rng(3).normal(0, 0.3, len(y_all))
        X_all = np.c_[X_all, recorded_after]
        X_new = np.c_[X_new, np.zeros(len(y_new))]  # not known when the model is used
    X_train, X_test, y_train, y_test = train_test_split(
        X_all, y_all, test_size=0.25, random_state=0
    )
    model.fit(X_train, y_train)
    return model.score(X_new, y_new), model.score(X_test, y_test)


KINDS = {
    "duplicates": "duplicates",
    "oversampling": "oversampling first",
    "target": "target leakage",
}
results = {label: honest_and_leaky(kind) for kind, label in KINDS.items()}
fig, ax = plt.subplots(figsize=(8, 2.8))
for i, (label, (honest, leaky)) in enumerate(results.items()):
    ax.barh(
        i - 0.18,
        leaky,
        height=0.36,
        color=style.GREY,
        alpha=0.7,
        label="leaky test split" if i == 0 else None,
    )
    ax.barh(
        i + 0.18,
        honest,
        height=0.36,
        color=style.BLUE,
        alpha=0.7,
        label="real use" if i == 0 else None,
    )
    ax.text(leaky + 0.01, i - 0.18, f"{leaky:.0%}", va="center", fontsize=8.5)
    ax.text(honest + 0.01, i + 0.18, f"{honest:.0%}", va="center", fontsize=8.5)
ax.set_yticks(range(len(results)), list(results))
ax.invert_yaxis()
ax.xaxis.set_major_formatter(style.PERCENT)
ax.legend(fontsize=8, loc="lower right")
ax.set(xlim=(0, 1.12), xlabel="accuracy")
plt.show()

Figure 6.2.7: Three kinds of leakage, on simulated data. Grey: the score on a test split with the leak. Blue: the same model on fresh data, as it will really be used. Duplicates: 40% of the rows appear twice, and some copies land in the test set. Oversampling first: the rare class was copied five times before splitting. Target leakage: a predictor recorded after the outcome, which is missing when the model is really used.

A practical test for leakage: if a result looks too good, find the predictors that matter most (Chapter 6.8) and ask, for each, whether it would really be known at the moment of prediction.

6.2.7 When new data differ: external validation

A test set drawn from the same pool as the training data answers: how well does the model do on more cases like these? Often the real question is different: other hospitals, other years, other countries, other machines. When the new data differ in their predictors (covariate shift) or in how predictors relate to the outcome (concept shift), the internal test score is too optimistic. External validation, testing on data from a different source, is the strongest evidence a model can have.

Code

def true_risk(age: np.ndarray) -> np.ndarray:
    """A risk that rises with age up to about 50, then falls."""
    return 1 / (1 + np.exp((age - 50) ** 2 / 100 - 1.5))


def hospital(n: int, centre: float, seed: int) -> tuple[np.ndarray, np.ndarray]:
    """Ages around `centre` and outcomes drawn from the true risk."""
    local = np.random.default_rng(seed)
    age = local.normal(centre, 8, n)
    risk = true_risk(age)
    return age[:, None], (local.uniform(size=n) < risk).astype(int)


age_a, y_a = hospital(1500, 40, 0)
age_b, y_b = hospital(1500, 65, 1)
A_train, A_test, ya_train, ya_test = train_test_split(age_a, y_a, test_size=0.3, random_state=0)
model = LogisticRegression().fit(A_train, ya_train)
grid_age = np.linspace(15, 95, 200)[:, None]

fig, axes = style.subplots(2, height=3.1, width=10)
axes[0].hist(age_a[:, 0], bins=30, color=style.BLUE, alpha=0.5, label="hospital A")
axes[0].hist(age_b[:, 0], bins=30, color=style.CORAL, alpha=0.5, label="hospital B")
twin = axes[0].twinx()
twin.plot(
    grid_age[:, 0],
    model.predict_proba(grid_age)[:, 1],
    color=style.INK,
    lw=1.5,
    label="model from A",
)
truth_age = grid_age[:, 0]
twin.plot(truth_age, true_risk(truth_age), color=style.INK, ls="--", lw=1.2, label="true risk")
twin.set(ylabel="risk", ylim=(0, 1))
twin.legend(fontsize=8, loc="center left")
axes[0].legend(fontsize=8, loc="upper left")
axes[0].set(xlabel="age", ylabel="patients")
results = [
    ("A, held out", model.score(A_test, ya_test)),
    ("B, external", model.score(age_b, y_b)),
]
axes[1].barh(
    [r[0] for r in results], [r[1] for r in results], color=[style.BLUE, style.CORAL], alpha=0.7
)
for i, (_, value) in enumerate(results):
    axes[1].text(value + 0.01, i, f"{value:.0%}", va="center", fontsize=9)
axes[1].invert_yaxis()
axes[1].xaxis.set_major_formatter(style.PERCENT)
axes[1].set(xlim=(0, 1.1), xlabel="accuracy of the model from A")
fig.tight_layout()
plt.show()

Figure 6.2.8: Simulated patients from two hospitals. Hospital A’s patients (blue) are mostly young, hospital B’s (red) mostly old, and the true risk (dashed) rises with age up to about 50 and falls after. A logistic model trained on hospital A (solid) learns only the rise. It scores well on held-out patients from A, and badly on hospital B, whose patients lie where A had few cases.

The same idea in time is temporal validation: train on earlier years, test on later ones, and keep checking a model in use, because the world it was trained on drifts away.

6.2.8 Choosing on the validation set: tuning

Most methods have settings that are not learned from the data by fitting: the penalty λ, the depth of trees, the learning rate, the number of neighbours. These hyperparameters are chosen by trying values and comparing validation scores. With several of them, a full grid of values quickly gets expensive, and random search is usually better: the same number of tries covers more distinct values of each setting, and often only one or two settings really matter.

Code

def score(setting: np.ndarray) -> np.ndarray:
    """A validation score that depends on the first setting only, peaking at 0.73."""
    return np.exp(-((setting - 0.73) ** 2) / 0.01)


grid_points = np.array([[a, b] for a in (1 / 6, 0.5, 5 / 6) for b in (1 / 6, 0.5, 5 / 6)])
random_points = np.random.default_rng(4).uniform(0, 1, (9, 2))
line = np.linspace(0, 1, 200)

fig, axes = plt.subplots(
    2, 2, figsize=(8.5, 4.4), gridspec_kw={"height_ratios": [1, 3]}, sharex=True
)
for col, (points, title) in enumerate(
    ((grid_points, "grid search"), (random_points, "random search"))
):
    axes[0, col].plot(line, score(line), color=style.INK, lw=1.2)
    axes[0, col].scatter(points[:, 0], score(points[:, 0]), color=style.CORAL, zorder=3)
    axes[0, col].set(title=title, yticks=[])
    axes[1, col].scatter(points[:, 0], points[:, 1], color=style.CORAL, s=40)
    axes[1, col].axvline(0.73, color=style.GREEN, ls="--", lw=1)
    axes[1, col].set(
        xlim=(0, 1), ylim=(0, 1), xlabel="setting 1 (matters)", ylabel="setting 2 (does not)"
    )
fig.tight_layout()
plt.show()

Figure 6.2.9: Nine tries to tune two settings, when only the first one matters (the validation score, shown as the curve above each panel, depends on it alone). A 3 × 3 grid (left) tests only three distinct values of the important setting; nine random tries (right) test nine, and land closer to the best value.

Early stopping (Chapter 8.1) is tuning too: the number of training rounds is chosen on the validation set. All of it belongs to the “choosing” job, so none of it may touch the test set.

6.2.9 A checklist

Decide what “new” means in use (new cases, patients, sites, times), and split so that the test set is new in that way.
Split first, before exploring the data.
Stratify rare classes; keep groups together; respect time.
Fit every preprocessing step inside the training data, in one pipeline.
Choose on validation data or cross-validation; never on the test set.
Test once; report the test score with its uncertainty.
If possible, validate externally, and keep monitoring a model in use.

6.2.10 In one picture

Three jobs, three parts: training fits, validation chooses, test reports, once.
Choosing on a set makes its score too good: the winner’s curse.
Test precision: $\text{SE} = \sqrt{p(1-p)/n}$ ; small test sets measure little.
Stratify rare classes, group related rows, split chronologically for time.
Beware leakage: preprocessing on all data, copies across the split, predictors known only after the outcome.
Internal tests answer “more of the same”; external validation answers “somewhere else”.
Hyperparameters are chosen on validation data; random search beats grid search.

Next: cross-validation, which lets every case take a turn as unseen data.