Hands-on lab — Day 137: Thinking in Features

"""Shared fixtures.

pytest finds this file by itself -- nothing imports it. Every experiment
is wrapped in a session-scoped fixture so the whole suite runs each one
exactly once; several of them average over a hundred or more train/test
splits and would otherwise be repeated for every assertion.

Nothing here writes to disk, binds a port or reaches the network.
"""

from __future__ import annotations

import pytest

import experiments as E


@pytest.fixture(scope="session")
def leakage():
    return E.target_leakage()


@pytest.fixture(scope="session")
def scaler_contamination():
    return E.scaling_contamination()


@pytest.fixture(scope="session")
def imputer_contamination():
    return E.imputer_contamination()


@pytest.fixture(scope="session")
def encoding():
    return E.target_encoding()


@pytest.fixture(scope="session")
def temporal():
    return E.temporal_leakage()


@pytest.fixture(scope="session")
def cyclical():
    return E.cyclical_distances()


@pytest.fixture(scope="session")
def colours():
    return E.ordinal_versus_one_hot()


@pytest.fixture(scope="session")
def interaction():
    return E.interaction()


@pytest.fixture(scope="session")
def vocabulary():
    return E.vocabulary_contamination()


@pytest.fixture(scope="session")
def audit():
    return E.audit_result()

"""The datasets for Day 137 -- "Features That Do Not Cheat".

Every generator here is seeded and returns exactly the same rows on every
machine, so every number this lab asserts is reproducible rather than
sampled. Nothing is downloaded; nothing touches the network.

Six generators, one per experiment:

* `signups`          -- a conversion table carrying one planted target leak
* `pricing`          -- a heavy-tailed feature for the scaling experiment
* `city_signups`     -- 40 city codes for the target-encoding experiment
* `sensor_log`       -- a time-ordered table whose rule changes by regime
* `paint_orders`     -- unordered colours for the ordinal-encoding trap
* `credit_lines`     -- spend and income that only separate as a ratio
* `tickets`          -- short documents for the bag-of-words experiment

Read the docstrings before the code: each one states what the honest
signal is and, where a leak is planted, exactly where it was planted.
"""

from __future__ import annotations

import numpy as np
import pandas as pd

SEED = 137


def _sigmoid(z: np.ndarray) -> np.ndarray:
    return 1.0 / (1.0 + np.exp(-z))


def signups(n: int = 400, seed: int = SEED) -> pd.DataFrame:
    """A conversion table with three honest features and one planted leak.

    Honest: `visits`, `minutes_on_site`, `discount_pct`. The label
    `converted` is drawn from a logistic model of those three, so the
    honest signal is real but far from perfect.

    Planted leak: `days_to_first_invoice`. An invoice only exists once a
    visitor has converted, so the column is a positive number for every
    converted row and the sentinel -1 for every unconverted one. At
    prediction time -- before the visitor has decided -- this column
    cannot exist at all. It is the outcome, wearing a different name.
    """
    rng = np.random.default_rng(seed)
    visits = rng.integers(1, 30, n)
    minutes = rng.gamma(2.0, 6.0, n)
    discount = rng.choice([0, 5, 10, 20], n)
    logit = -2.4 + 0.08 * visits + 0.045 * minutes + 0.035 * discount
    converted = (rng.random(n) < _sigmoid(logit)).astype(int)
    days = np.where(converted == 1, rng.integers(1, 15, n), -1)
    return pd.DataFrame(
        {
            "visits": visits.astype(float),
            "minutes_on_site": minutes,
            "discount_pct": discount.astype(float),
            "days_to_first_invoice": days.astype(float),
            "converted": converted,
        }
    )


def pricing(n: int = 500, seed: int = SEED + 1) -> pd.DataFrame:
    """A two-feature table whose first feature is heavy-tailed.

    `order_value` is lognormal, so a handful of rows sit two orders of
    magnitude above the median. A scaler fitted on all the data sees
    those rows wherever they fall; a scaler fitted on the training half
    only sees the ones that landed in training. That difference is the
    whole point of the scaling experiment, and it is why the feature is
    heavy-tailed rather than Gaussian: with a well-behaved feature the
    two scalers agree to three decimals and there is nothing to see.
    """
    rng = np.random.default_rng(seed)
    order_value = rng.lognormal(mean=3.0, sigma=1.4, size=n)
    tenure_days = rng.gamma(3.0, 90.0, size=n)
    logit = -1.1 + 0.55 * np.log(order_value) - 0.004 * tenure_days
    renewed = (rng.random(n) < _sigmoid(logit)).astype(int)
    return pd.DataFrame(
        {
            "order_value": order_value,
            "tenure_days": tenure_days,
            "renewed": renewed,
        }
    )


def panel_readings(
    n: int = 600,
    n_panels: int = 120,
    missing_rate: float = 0.60,
    seed: int = SEED + 7,
) -> pd.DataFrame:
    """Readings with a lot of gaps, grouped into small panels.

    Each solar panel has its own characteristic output, and `reading` is
    that output plus noise. `fault` depends on the reading. Sixty per
    cent of the readings are missing, and with 120 panels over 600 rows
    a panel holds five rows on average -- so a group mean fitted on a
    small training set is often estimated from one row, or from none.

    That is what makes group-mean imputation the sharpest available
    demonstration of a statistic fitted before the split. Impute over the
    whole table and a test row's gap is filled from the observed readings
    of its own panel, including the ones sitting in the test set.
    """
    rng = np.random.default_rng(seed)
    panel = rng.integers(0, n_panels, n)
    panel_output = rng.normal(0.0, 2.0, n_panels)
    reading = panel_output[panel] + rng.normal(0.0, 0.5, n)
    fault = (rng.random(n) < _sigmoid(1.2 * reading)).astype(int)
    observed = rng.random(n) >= missing_rate
    return pd.DataFrame(
        {
            "panel": pd.Series([f"P{p:03d}" for p in panel], dtype="str"),
            "reading": np.where(observed, reading, np.nan),
            "fault": fault,
        }
    )


def city_signups(n: int = 600, n_cities: int = 40, seed: int = SEED + 2) -> pd.DataFrame:
    """High-cardinality categories for the target-encoding experiment.

    40 city codes over 600 rows is 15 rows per city on average, which is
    exactly the regime where a per-category mean of the target is mostly
    noise. The cities do carry a small real effect, so an honest encoding
    is not useless -- it is just far less impressive than the naive one
    pretends.
    """
    rng = np.random.default_rng(seed)
    city = rng.integers(0, n_cities, n)
    city_effect = rng.normal(0.0, 0.35, n_cities)
    visits = rng.integers(1, 20, n)
    logit = -0.3 + 0.05 * (visits - 10) + city_effect[city]
    converted = (rng.random(n) < _sigmoid(logit)).astype(int)
    return pd.DataFrame(
        {
            "city": pd.Series([f"C{c:02d}" for c in city], dtype="str"),
            "visits": visits.astype(float),
            "converted": converted,
        }
    )


#: Alarm rate per calibration batch, in batch order. A batch is a period
#: of time, and the rate changes when the hardware is recalibrated.
BATCH_ALARM_RATE = [0.85, 0.15, 0.88, 0.12, 0.90, 0.10]


def sensor_log(per_batch: int = 60, seed: int = SEED + 3) -> pd.DataFrame:
    """A time-ordered table whose batches are periods of time.

    Rows come back in time order. `batch` is the calibration batch the
    sensor was running under, so a batch is not scattered through the
    table -- it occupies one contiguous stretch of it. The alarm rate
    changes sharply from batch to batch (a recalibration, a firmware
    change, a new supplier), while `reading` carries a mild, stable
    effect that holds across all of them.

    That combination is what a random split hides. Split at random and
    every batch has rows in training, so the model learns each batch's
    alarm rate from rows recorded at the same time as the ones it is
    scored on -- information from after the prediction moment. Split by
    time and the last batch is one the model has never seen, which is the
    situation every deployed model is actually in.
    """
    rng = np.random.default_rng(seed)
    frames = []
    for batch, rate in enumerate(BATCH_ALARM_RATE):
        reading = rng.normal(50.0, 6.0, per_batch)
        humidity = rng.normal(40.0, 9.0, per_batch)
        logit = np.log(rate / (1 - rate)) + 0.05 * (reading - 50.0)
        alarm = (rng.random(per_batch) < _sigmoid(logit)).astype(int)
        frames.append(
            pd.DataFrame(
                {
                    "batch": pd.Series([f"B{batch}"] * per_batch, dtype="str"),
                    "reading": reading,
                    "humidity": humidity,
                    "alarm": alarm,
                }
            )
        )
    out = pd.concat(frames, ignore_index=True)
    out.insert(0, "t", np.arange(len(out), dtype=int))
    return out


PAINT_COLOURS = ["amber", "cobalt", "ivory", "olive", "rose", "slate"]

#: The true return rate per colour, in the same order as PAINT_COLOURS.
#: Deliberately NOT monotone in the list position: ivory (code 2) is the
#: worst and cobalt (code 1) the best, so any model that can only move
#: monotonically with the code is guaranteed to be wrong somewhere.
PAINT_RETURN_RATE = [0.20, 0.08, 0.72, 0.15, 0.62, 0.25]


def paint_orders(per_colour: int = 200, seed: int = SEED + 4) -> pd.DataFrame:
    """Unordered categories with a deliberately non-monotone outcome.

    Colour names have no order. Alphabetical position is not a quantity,
    which is exactly why an ordinal code invites a model to interpolate
    between categories that have no midpoint.
    """
    rng = np.random.default_rng(seed)
    colours: list[str] = []
    returned: list[int] = []
    for name, rate in zip(PAINT_COLOURS, PAINT_RETURN_RATE):
        colours.extend([name] * per_colour)
        returned.extend((rng.random(per_colour) < rate).astype(int).tolist())
    return pd.DataFrame(
        {
            "colour": pd.Series(colours, dtype="str"),
            "returned": np.array(returned, dtype=int),
        }
    )


def credit_lines(n: int = 500, seed: int = SEED + 5) -> pd.DataFrame:
    """Spend and income whose marginals overlap but whose ratio does not.

    Income is drawn from the same distribution for both classes. Spend is
    income times a ratio drawn per class, and the two ratio bands do not
    overlap. So neither column separates the classes on its own, while
    `spend / income` separates them completely. This is an interaction in
    the plainest possible form.
    """
    rng = np.random.default_rng(seed)
    label = (rng.random(n) < 0.5).astype(int)
    income = rng.uniform(30_000, 140_000, n)
    ratio = np.where(label == 1, rng.uniform(0.52, 0.68, n), rng.uniform(0.32, 0.48, n))
    spend = income * ratio
    return pd.DataFrame(
        {
            "income": income,
            "spend": spend,
            "stressed": label,
        }
    )


#: The corpus vocabulary and each word's pull towards "urgent". No word
#: is decisive: the largest effect only doubles a word's odds of turning
#: up, so a document is classified by the whole bag or not at all.
TICKET_WORDS = [
    "outage", "refund", "broken", "charged", "cancel", "waiting",
    "thanks", "question", "hello", "curious", "manual", "hours",
    "order", "account", "team", "please", "about", "with", "the", "and",
]
TICKET_EFFECT = [
    0.40, 0.34, 0.38, 0.30, 0.36, 0.26,
    -0.40, -0.32, -0.36, -0.28, -0.34, -0.24,
    0.06, -0.05, 0.04, -0.03, 0.02, -0.02, 0.01, -0.01,
]


def tickets(n: int = 300, seed: int = SEED + 6) -> pd.DataFrame:
    """Short support tickets, each a bag of lower-case words.

    Every document draws its words from one shared pool; the label only
    tilts the odds. Twelve words carry a real but modest pull and eight
    carry essentially none, so no single word settles a document and the
    classifier has to add evidence up.

    On top of that sits a long tail of reference codes, each appearing in
    only a handful of documents corpus-wide. Those are the words that
    make the choice of vocabulary consequential: a code that turns up
    three times and happens to land twice on an urgent ticket looks like
    a strong feature to anything that scores words against the label.
    """
    rng = np.random.default_rng(seed)
    effect = np.array(TICKET_EFFECT, dtype=float)
    rows: list[str] = []
    labels: list[int] = []
    for _ in range(n):
        urgent = int(rng.random() < 0.5)
        sign = 1.0 if urgent else -1.0
        weights = np.exp(sign * effect)
        weights = weights / weights.sum()
        length = int(rng.integers(8, 15))
        words = rng.choice(TICKET_WORDS, size=length, replace=True, p=weights).tolist()
        if rng.random() < 0.35:
            words.append(f"ref{int(rng.integers(0, 80)):02d}")
        rows.append(" ".join(words))
        labels.append(urgent)
    return pd.DataFrame(
        {
            "text": pd.Series(rows, dtype="str"),
            "urgent": np.array(labels, dtype=int),
        }
    )

"""The nine experiments, each one a function that returns numbers.

The tests assert on what these return, and `tests/run_tests.sh` prints
them. Keeping the experiments here rather than inside the tests means the
same code produces the numbers in the lesson, the numbers in
`expected-output/` and the numbers the assertions check, so the three
cannot drift apart.

Every function is deterministic. Where an experiment averages over many
random splits it says so in its name and its docstring, and the seeds are
generated from a fixed base seed.
"""

from __future__ import annotations

import numpy as np
import pandas as pd

import data
import features as F
import models as M

HONEST_COLUMNS = ["visits", "minutes_on_site", "discount_pct"]
LEAKY_COLUMNS = HONEST_COLUMNS + ["days_to_first_invoice"]


def _fit_score(X_train, y_train, X_test, y_test, model=None) -> float:
    model = model or M.LogisticRegression()
    model.fit(X_train, y_train)
    return model.score(X_test, y_test)


# --- 1. target leakage -----------------------------------------------------


def target_leakage(test_size: int = 100, seed: int = 137) -> dict[str, float]:
    """Score the same task twice: with the planted leak and without it.

    Both runs use the same split, the same model and the same number of
    gradient steps. The only difference is one column.
    """
    frame = data.signups()
    y = frame["converted"].to_numpy()
    train_idx, test_idx = M.random_split(len(frame), test_size, seed)
    out: dict[str, float] = {}
    for name, columns in (("with_leak", LEAKY_COLUMNS), ("without_leak", HONEST_COLUMNS)):
        X = frame[columns].to_numpy(dtype=float)
        scaler = F.Standardiser().fit(X[train_idx])
        out[name] = _fit_score(
            scaler.transform(X[train_idx]),
            y[train_idx],
            scaler.transform(X[test_idx]),
            y[test_idx],
        )
    out["gap_points"] = 100.0 * (out["with_leak"] - out["without_leak"])
    return out


# --- 2. train/test contamination through a scaler --------------------------


def scaling_contamination(
    trials: int = 200, train_size: int = 60, test_size: int = 25, seed: int = 137
) -> dict[str, float]:
    """Fit the scaler on everything, or on the training rows only.

    Averaged over `trials` random splits, because with a 25-row test set
    one row is four accuracy points and a single split says nothing. The
    model is the nearest-centroid classifier, which is distance-based and
    therefore cares about scale; a logistic regression run to convergence
    is very nearly invariant to an affine change of features and shows
    almost nothing here.

    The mechanism is worth stating plainly: the contaminated scaler is a
    BETTER estimate of the population mean and spread than a scaler
    fitted on 60 rows. It got that way by reading rows it was not allowed
    to read, and the score it produces is therefore not a score you will
    ever see in production.
    """
    frame = data.pricing()
    X_all = frame[["order_value", "tenure_days"]].to_numpy(dtype=float)
    y_all = frame["renewed"].to_numpy()
    contaminated_scores: list[float] = []
    correct_scores: list[float] = []
    for trial in range(trials):
        rng = np.random.default_rng(seed + trial)
        order = rng.permutation(len(frame))
        train_idx = order[:train_size]
        test_idx = order[train_size : train_size + test_size]

        wrong = F.Standardiser().fit(X_all)  # every row, including the test rows
        right = F.Standardiser().fit(X_all[train_idx])  # training rows only

        contaminated_scores.append(
            _fit_score(
                wrong.transform(X_all[train_idx]),
                y_all[train_idx],
                wrong.transform(X_all[test_idx]),
                y_all[test_idx],
                M.NearestCentroid(),
            )
        )
        correct_scores.append(
            _fit_score(
                right.transform(X_all[train_idx]),
                y_all[train_idx],
                right.transform(X_all[test_idx]),
                y_all[test_idx],
                M.NearestCentroid(),
            )
        )
    contaminated = float(np.mean(contaminated_scores))
    correct = float(np.mean(correct_scores))
    return {
        "contaminated": contaminated,
        "correct": correct,
        "optimism_points": 100.0 * (contaminated - correct),
        "trials": float(trials),
        "train_size": float(train_size),
        "test_size": float(test_size),
        "contaminated_mean_order_value": float(F.Standardiser().fit(X_all).mean_[0]),
    }


def imputer_contamination(
    trials: int = 150, train_size: int = 120, test_size: int = 200, seed: int = 137
) -> dict[str, float]:
    """The same question asked of an imputer instead of a scaler.

    A group-mean imputer is not an affine map. Each missing value is
    filled with a number computed from the rows that happen to be
    observed in its own group -- so if the fit sees the test rows, a
    test row's gap is filled using the readings of its own panel that
    live in the test set. That is information the deployed system will
    not have, and here it is worth several accuracy points rather than
    the fraction of a point the scaler was worth.
    """
    frame = data.panel_readings()
    panel = frame["panel"].to_numpy()
    reading = frame["reading"].to_numpy(dtype=float)
    y = frame["fault"].to_numpy()
    n = len(frame)
    missing = np.isnan(reading)

    everything = F.GroupMeanImputer().fit(panel, reading)

    def build(rows, imputer) -> np.ndarray:
        filled = imputer.transform(panel[rows], reading[rows])
        return np.column_stack([filled, missing[rows].astype(float)])

    contaminated_scores: list[float] = []
    correct_scores: list[float] = []
    for trial in range(trials):
        rng = np.random.default_rng(seed + trial)
        order = rng.permutation(n)
        train_idx = order[:train_size]
        test_idx = order[train_size : train_size + test_size]
        train_only = F.GroupMeanImputer().fit(panel[train_idx], reading[train_idx])
        for imputer, bucket in ((everything, contaminated_scores), (train_only, correct_scores)):
            X_train = build(train_idx, imputer)
            X_test = build(test_idx, imputer)
            scaler = F.Standardiser().fit(X_train)
            bucket.append(
                _fit_score(
                    scaler.transform(X_train), y[train_idx], scaler.transform(X_test), y[test_idx]
                )
            )
    contaminated = float(np.mean(contaminated_scores))
    correct = float(np.mean(correct_scores))
    return {
        "contaminated": contaminated,
        "correct": correct,
        "optimism_points": 100.0 * (contaminated - correct),
        "trials": float(trials),
        "train_size": float(train_size),
        "test_size": float(test_size),
        "missing_fraction": float(missing.mean()),
        "groups_known_to_all_data_fit": float(len(everything.means_)),
    }


# --- 3. target encoding ----------------------------------------------------


def target_encoding(
    test_size: int = 150, seed: int = 137, n_folds: int = 5, trials: int = 40
) -> dict[str, float]:
    """Three ways to encode a city, scored on the same held-out rows.

    * `naive_all_data` computes the per-city mean of the target over the
      whole table, then splits. Every test row's feature was computed
      partly from its own answer.
    * `naive_train_only` computes it on the training rows only, which
      removes the test rows' contribution but still lets each training
      row see its own target.
    * `out_of_fold` encodes each training row from the folds that do not
      contain it.
    """
    frame = data.city_signups()
    y = frame["converted"].to_numpy()
    visits = frame["visits"].to_numpy(dtype=float)
    city = frame["city"].to_numpy()

    all_scores: list[float] = []
    train_scores: list[float] = []
    oof_scores: list[float] = []
    for trial in range(trials):
        train_idx, test_idx = M.random_split(len(frame), test_size, seed + trial)

        def score(train_codes, test_codes) -> float:
            X_train = np.column_stack([visits[train_idx], train_codes])
            X_test = np.column_stack([visits[test_idx], test_codes])
            scaler = F.Standardiser().fit(X_train)
            return _fit_score(
                scaler.transform(X_train), y[train_idx], scaler.transform(X_test), y[test_idx]
            )

        all_map, all_default = F.target_encode_fit(city, y)
        all_scores.append(
            score(
                F.target_encode_transform(city[train_idx], all_map, all_default),
                F.target_encode_transform(city[test_idx], all_map, all_default),
            )
        )

        train_map, train_default = F.target_encode_fit(city[train_idx], y[train_idx])
        train_scores.append(
            score(
                F.target_encode_transform(city[train_idx], train_map, train_default),
                F.target_encode_transform(city[test_idx], train_map, train_default),
            )
        )

        oof = F.target_encode_out_of_fold(
            city[train_idx], y[train_idx], n_folds=n_folds, seed=seed + trial
        )
        oof_scores.append(
            score(oof, F.target_encode_transform(city[test_idx], train_map, train_default))
        )

    naive_all = float(np.mean(all_scores))
    naive_train = float(np.mean(train_scores))
    out_of_fold = float(np.mean(oof_scores))
    return {
        "naive_all_data": naive_all,
        "naive_train_only": naive_train,
        "out_of_fold": out_of_fold,
        "gap_all_vs_oof_points": 100.0 * (naive_all - out_of_fold),
        "gap_train_vs_oof_points": 100.0 * (naive_train - out_of_fold),
        "n_folds": float(n_folds),
        "trials": float(trials),
    }


# --- 4. temporal leakage ---------------------------------------------------


def temporal_leakage(test_size: int = 60, seed: int = 137) -> dict[str, float]:
    """The same table, split two ways: at random, and by time.

    The rows are already in time order. The random split lets rows from
    the last regime into training, so the model is scored on a rule it
    has already seen. The time-ordered split holds the last regime out
    entirely, which is the only one of the two that answers the question
    anybody actually asked.
    """
    frame = data.sensor_log()
    y = frame["alarm"].to_numpy()
    reading = frame["reading"].to_numpy(dtype=float)
    batch = frame["batch"].to_numpy()

    def score(train_idx, test_idx) -> float:
        # The batch encoding is fitted on the training rows only, which is
        # the correct thing to do and is exactly what exposes the problem:
        # a batch the training rows never contained becomes an all-zero
        # block, and the model has nothing to say about it.
        seen = sorted(set(batch[train_idx].tolist()))
        train_block, _ = F.one_hot(batch[train_idx], seen)
        test_block, _ = F.one_hot(batch[test_idx], seen)
        X_train = np.column_stack([reading[train_idx], train_block])
        X_test = np.column_stack([reading[test_idx], test_block])
        scaler = F.Standardiser().fit(X_train)
        return _fit_score(
            scaler.transform(X_train), y[train_idx], scaler.transform(X_test), y[test_idx]
        )

    random_train, random_test = M.random_split(len(frame), test_size, seed)
    ordered_train, ordered_test = M.time_ordered_split(len(frame), test_size)
    random_score = score(random_train, random_test)
    ordered_score = score(ordered_train, ordered_test)
    return {
        "random_split": random_score,
        "time_ordered_split": ordered_score,
        "gap_points": 100.0 * (random_score - ordered_score),
        "batches_seen_by_random_train": float(len(set(batch[random_train].tolist()))),
        "batches_seen_by_ordered_train": float(len(set(batch[ordered_train].tolist()))),
        "test_batches_unseen_by_ordered_train": float(
            len(set(batch[ordered_test].tolist()) - set(batch[ordered_train].tolist()))
        ),
        "majority_rate_in_ordered_test": float(max(y[ordered_test].mean(), 1 - y[ordered_test].mean())),
    }


# --- 5. cyclical encoding --------------------------------------------------


def cyclical_distances(period: int = 24) -> dict[str, float]:
    """Distances between hours, raw and on the circle.

    Raw integer hours: 23 and 0 sit 23 units apart while 3 and 4 sit 1
    apart. Sine-cosine: every adjacent pair, wrap included, sits exactly
    2*sin(pi/24) apart.
    """
    hours = np.arange(period)
    raw = hours.reshape(-1, 1).astype(float)
    circle = F.cyclical_encode(hours, period)

    def distance(matrix, a: int, b: int) -> float:
        return float(np.linalg.norm(matrix[a] - matrix[b]))

    adjacent_raw = [distance(raw, h, (h + 1) % period) for h in range(period)]
    adjacent_circle = [distance(circle, h, (h + 1) % period) for h in range(period)]
    return {
        "raw_23_to_0": distance(raw, 23, 0),
        "raw_3_to_4": distance(raw, 3, 4),
        "cyclical_23_to_0": distance(circle, 23, 0),
        "cyclical_3_to_4": distance(circle, 3, 4),
        "cyclical_0_to_12": distance(circle, 0, 12),
        "cyclical_adjacent_spread": float(max(adjacent_circle) - min(adjacent_circle)),
        "raw_adjacent_spread": float(max(adjacent_raw) - min(adjacent_raw)),
        "expected_adjacent": float(2.0 * np.sin(np.pi / period)),
    }


# --- 6. ordinal versus one-hot --------------------------------------------


def ordinal_versus_one_hot(seed: int = 137) -> dict[str, object]:
    """Fit the same model on an ordinal code and on a one-hot block.

    The colours have no order, and their return rates are deliberately
    not monotone in the alphabetical code. A model reading the code as a
    number can only produce predictions that rise or fall with it; the
    one-hot model is free to give each colour its own answer.
    """
    frame = data.paint_orders()
    y = frame["returned"].to_numpy()
    colours = frame["colour"].to_numpy()
    order = list(data.PAINT_COLOURS)

    codes = F.ordinal_encode(colours, order).reshape(-1, 1)
    ordinal_model = M.LogisticRegression(learning_rate=0.3, steps=4000).fit(codes, y)
    ordinal_rates = ordinal_model.predict_proba(
        np.arange(len(order), dtype=float).reshape(-1, 1)
    )

    block, categories = F.one_hot(colours, order)
    one_hot_model = M.LogisticRegression(learning_rate=0.3, steps=4000).fit(block, y)
    one_hot_rates = one_hot_model.predict_proba(np.eye(len(order)))

    observed = np.array(
        [float(y[colours == name].mean()) for name in order], dtype=float
    )

    def is_monotone(values: np.ndarray) -> bool:
        diffs = np.diff(values)
        return bool(np.all(diffs >= -1e-12) or np.all(diffs <= 1e-12))

    return {
        "categories": categories,
        "observed_rates": observed.tolist(),
        "ordinal_predictions": ordinal_rates.tolist(),
        "one_hot_predictions": one_hot_rates.tolist(),
        "ordinal_is_monotone": is_monotone(ordinal_rates),
        "one_hot_is_monotone": is_monotone(one_hot_rates),
        "ordinal_max_error": float(np.max(np.abs(ordinal_rates - observed))),
        "one_hot_max_error": float(np.max(np.abs(one_hot_rates - observed))),
        "ordinal_accuracy": float(ordinal_model.score(codes, y)),
        "one_hot_accuracy": float(one_hot_model.score(block, y)),
    }


# --- 7. an interaction -----------------------------------------------------


def interaction(test_size: int = 150, seed: int = 137) -> dict[str, float]:
    """Score income alone, spend alone, and the ratio of the two."""
    frame = data.credit_lines()
    y = frame["stressed"].to_numpy()
    train_idx, test_idx = M.random_split(len(frame), test_size, seed)
    income = frame["income"].to_numpy(dtype=float)
    spend = frame["spend"].to_numpy(dtype=float)
    ratio = F.ratio_feature(spend, income)

    def score(column: np.ndarray) -> float:
        X = column.reshape(-1, 1)
        scaler = F.Standardiser().fit(X[train_idx])
        return _fit_score(
            scaler.transform(X[train_idx]), y[train_idx], scaler.transform(X[test_idx]), y[test_idx]
        )

    both = np.column_stack([income, spend])
    scaler = F.Standardiser().fit(both[train_idx])
    both_score = _fit_score(
        scaler.transform(both[train_idx]),
        y[train_idx],
        scaler.transform(both[test_idx]),
        y[test_idx],
    )
    return {
        "income_only": score(income),
        "spend_only": score(spend),
        "ratio_only": score(ratio),
        "income_and_spend": both_score,
    }


# --- 8. vocabulary fitted on the wrong rows --------------------------------


def vocabulary_contamination(
    top_k: int = 30, test_size: int = 100, seed: int = 137, min_docs: int = 2, trials: int = 40
) -> dict[str, float]:
    """Choose the words on everything, or on the training documents only.

    The vocabulary is selected by association with the label, so choosing
    it on all the documents means the test labels helped decide which
    features exist. That is feature selection fitted before the split --
    the second kind of leakage, wearing a text-shaped costume.
    """
    frame = data.tickets()
    documents = frame["text"].tolist()
    y = frame["urgent"].to_numpy()
    all_data_vocab = F.Vocabulary(top_k=top_k, min_docs=min_docs).fit(documents, y)

    all_scores: list[float] = []
    train_scores: list[float] = []
    shared: list[int] = []
    for trial in range(trials):
        train_idx, test_idx = M.random_split(len(frame), test_size, seed + trial)
        train_docs = [documents[i] for i in train_idx]
        test_docs = [documents[i] for i in test_idx]
        train_only_vocab = F.Vocabulary(top_k=top_k, min_docs=min_docs).fit(
            train_docs, y[train_idx]
        )
        shared.append(len(set(all_data_vocab.words) & set(train_only_vocab.words)))
        for vocabulary, bucket in ((all_data_vocab, all_scores), (train_only_vocab, train_scores)):
            bucket.append(
                _fit_score(
                    vocabulary.transform(train_docs),
                    y[train_idx],
                    vocabulary.transform(test_docs),
                    y[test_idx],
                )
            )

    # One split, examined in detail, for the unseen-word behaviour.
    train_idx, test_idx = M.random_split(len(frame), test_size, seed)
    train_docs = [documents[i] for i in train_idx]
    test_docs = [documents[i] for i in test_idx]
    train_only_vocab = F.Vocabulary(top_k=top_k, min_docs=min_docs).fit(train_docs, y[train_idx])
    unseen = train_only_vocab.unseen_words(test_docs)
    test_matrix = train_only_vocab.transform(test_docs)

    fitted_all = float(np.mean(all_scores))
    fitted_train = float(np.mean(train_scores))
    return {
        "fitted_on_all_data": fitted_all,
        "fitted_on_train_only": fitted_train,
        "gap_points": 100.0 * (fitted_all - fitted_train),
        "unseen_test_words": float(len(unseen)),
        "test_matrix_columns": float(test_matrix.shape[1]),
        "test_matrix_rows": float(test_matrix.shape[0]),
        "shared_words": float(np.mean(shared)),
        "top_k": float(top_k),
        "trials": float(trials),
    }


# --- extra: a bin boundary is a decision (Day 130's bin width again) -------


def binning_decision(bins: int = 3) -> dict[str, object]:
    """Bin the same column two ways and read two different stories off it.

    Equal-width edges come from NumPy's own edge calculator. Equal-count
    edges come from the quantiles of the same column. Neither is wrong.
    They put wildly different numbers of rows in the top bin, and the
    renewal rate you would quote for "high value orders" depends entirely
    on which one you picked.
    """
    frame = data.pricing()
    values = frame["order_value"].to_numpy(dtype=float)
    y = frame["renewed"].to_numpy()

    width_edges = F.equal_width_bins(values, bins)
    count_edges = np.quantile(values, np.linspace(0.0, 1.0, bins + 1))
    width_index = F.bin_index(values, width_edges)
    count_index = F.bin_index(values, count_edges)

    def top_bin(index: np.ndarray) -> tuple[int, float]:
        top = index == index.max()
        return int(top.sum()), float(y[top].mean())

    width_rows, width_rate = top_bin(width_index)
    count_rows, count_rate = top_bin(count_index)
    return {
        "bins": bins,
        "equal_width_edges": [round(float(e), 2) for e in width_edges],
        "equal_count_edges": [round(float(e), 2) for e in count_edges],
        "equal_width_top_bin_rows": width_rows,
        "equal_count_top_bin_rows": count_rows,
        "equal_width_top_bin_rate": width_rate,
        "equal_count_top_bin_rate": count_rate,
        "rows": int(len(values)),
    }


# --- 9. the audit ----------------------------------------------------------


def audit_table() -> pd.DataFrame:
    """The signups table plus one categorical leak and one honest column."""
    frame = data.signups().copy()
    rng = np.random.default_rng(data.SEED + 9)
    converted = frame["converted"].to_numpy()
    frame["email_template"] = pd.Series(
        np.where(converted == 1, "welcome_pack", "abandoned_cart"), dtype="str"
    )
    frame["channel"] = pd.Series(
        rng.choice(["search", "social", "direct"], len(frame)), dtype="str"
    )
    return frame


def audit_result(corr_threshold: float = 0.90) -> dict[str, object]:
    frame = audit_table()
    flags = F.leakage_audit(frame, "converted", corr_threshold=corr_threshold)
    return {
        "flagged": [f.column for f in flags],
        "rules": {f.column: f.rule for f in flags},
        "details": {f.column: f.detail for f in flags},
        "columns_checked": [c for c in frame.columns if c != "converted"],
    }

"""Feature builders, every one of them split into `fit` and `transform`.

The split is the whole discipline of this lab. `fit` looks at data and
learns a statistic -- a mean, a standard deviation, a category's average
outcome, a vocabulary. `transform` applies a statistic that has already
been learned and looks at nothing. Once the two are separate functions
you can point `fit` at the training rows only, and the question "did the
test set influence this number?" has an answer you can read off the call
site instead of guessing at it.

Nothing here imports scikit-learn -- it is not installed in this lab.
Everything is NumPy and pandas, small enough to read.
"""

from __future__ import annotations

import re
from collections import Counter
from dataclasses import dataclass, field

import numpy as np
import pandas as pd

TOKEN = re.compile(r"[a-z0-9]+")


# ---------------------------------------------------------------------------
# Scaling
# ---------------------------------------------------------------------------


class Standardiser:
    """Subtract a mean and divide by a standard deviation (Day 107).

    `fit` records the two statistics. `transform` applies them. Fitting on
    rows you will later score is the second kind of leakage, and the only
    thing stopping you is which rows you hand to `fit`.
    """

    def __init__(self) -> None:
        self.mean_: np.ndarray | None = None
        self.scale_: np.ndarray | None = None

    def fit(self, X: np.ndarray) -> "Standardiser":
        X = np.asarray(X, dtype=float)
        self.mean_ = X.mean(axis=0)
        scale = X.std(axis=0)
        # A constant column has zero spread; dividing by it would produce
        # infinities, so it is left alone rather than exploded.
        self.scale_ = np.where(scale == 0.0, 1.0, scale)
        return self

    def transform(self, X: np.ndarray) -> np.ndarray:
        if self.mean_ is None or self.scale_ is None:
            raise RuntimeError("fit before transform")
        return (np.asarray(X, dtype=float) - self.mean_) / self.scale_

    def fit_transform(self, X: np.ndarray) -> np.ndarray:
        return self.fit(X).transform(X)


class GroupMeanImputer:
    """Fill a missing value with the mean of its group (Day 125).

    `fit` records one mean per group, plus a global fallback for groups
    it never saw. Like every other statistic in this file it reads rows,
    so which rows it reads is a decision -- and with small groups it is
    the most consequential decision in this lab.
    """

    def __init__(self) -> None:
        self.means_: dict[str, float] = {}
        self.default_: float = 0.0

    def fit(self, groups, values) -> "GroupMeanImputer":
        frame = pd.DataFrame(
            {"group": pd.Series(groups).astype("str"), "value": np.asarray(values, dtype=float)}
        ).dropna(subset=["value"])
        self.means_ = {
            str(k): float(v) for k, v in frame.groupby("group")["value"].mean().items()
        }
        self.default_ = float(frame["value"].mean()) if len(frame) else 0.0
        return self

    def transform(self, groups, values) -> np.ndarray:
        values = np.asarray(values, dtype=float)
        filled = values.copy()
        for i, (group, value) in enumerate(zip(pd.Series(groups).astype("str"), values)):
            if np.isnan(value):
                filled[i] = self.means_.get(str(group), self.default_)
        return filled


# ---------------------------------------------------------------------------
# Categorical encodings
# ---------------------------------------------------------------------------


def one_hot(values, categories: list[str] | None = None) -> tuple[np.ndarray, list[str]]:
    """One column of 0/1 per category, in a fixed, explicit category order.

    The category list is returned so the caller can reuse it at transform
    time. A test row carrying a category the training data never held
    becomes an all-zero row rather than a crash or a new column.
    """
    values = pd.Series(values).astype("str")
    if categories is None:
        categories = sorted(values.unique().tolist())
    matrix = np.zeros((len(values), len(categories)), dtype=float)
    position = {name: i for i, name in enumerate(categories)}
    for row, name in enumerate(values):
        column = position.get(name)
        if column is not None:
            matrix[row, column] = 1.0
    return matrix, list(categories)


def ordinal_encode(values, order: list[str]) -> np.ndarray:
    """Replace each category with its position in `order`.

    Honest when the order is real (small, medium, large). A trap when it
    is not: the code is a number, and a model reading it as a number can
    and will interpolate between categories that have no midpoint.
    """
    position = {name: float(i) for i, name in enumerate(order)}
    return np.array([position[str(v)] for v in values], dtype=float)


def target_encode_fit(categories, y) -> tuple[dict[str, float], float]:
    """Learn the mean outcome per category, plus the global mean.

    Returns the map and the fallback. Call this on TRAINING rows only:
    it reads the target, so fitting it on everything hands each test row
    a number computed partly from its own answer.
    """
    frame = pd.DataFrame({"category": pd.Series(categories).astype("str"), "y": np.asarray(y, dtype=float)})
    means = frame.groupby("category")["y"].mean()
    return {str(k): float(v) for k, v in means.items()}, float(frame["y"].mean())


def target_encode_transform(categories, mapping: dict[str, float], default: float) -> np.ndarray:
    """Apply a learned target encoding; unseen categories get the prior."""
    return np.array([mapping.get(str(c), default) for c in categories], dtype=float)


def target_encode_out_of_fold(categories, y, n_folds: int = 5, seed: int = 0) -> np.ndarray:
    """Encode each training row from the folds that do NOT contain it.

    This is the fix, and it is worth stating precisely what it fixes: a
    row's own target never contributes to its own feature value. The
    encoding is still built from training data only; out-of-fold does not
    excuse fitting on the test set.
    """
    categories = pd.Series(categories).astype("str").to_numpy()
    y = np.asarray(y, dtype=float)
    n = len(y)
    rng = np.random.default_rng(seed)
    fold_of = rng.permutation(n) % n_folds
    encoded = np.empty(n, dtype=float)
    for fold in range(n_folds):
        held_out = fold_of == fold
        mapping, default = target_encode_fit(categories[~held_out], y[~held_out])
        encoded[held_out] = target_encode_transform(categories[held_out], mapping, default)
    return encoded


# ---------------------------------------------------------------------------
# Datetime, cyclical, binning, interactions
# ---------------------------------------------------------------------------


def calendar_features(timestamps: pd.Series) -> pd.DataFrame:
    """The obvious datetime parts, pulled out as separate columns."""
    ts = pd.to_datetime(timestamps)
    return pd.DataFrame(
        {
            "hour": ts.dt.hour.astype(float),
            "day_of_week": ts.dt.dayofweek.astype(float),
            "day_of_month": ts.dt.day.astype(float),
            "month": ts.dt.month.astype(float),
            "is_weekend": (ts.dt.dayofweek >= 5).astype(float),
        }
    )


def cyclical_encode(values, period: float) -> np.ndarray:
    """Map a wrapping quantity onto a circle: two columns, sine and cosine.

    Hour 23 and hour 0 are one hour apart in the world and 23 apart as
    integers. On the circle they are neighbours again, and every pair of
    adjacent hours sits exactly the same distance apart.
    """
    angle = 2.0 * np.pi * np.asarray(values, dtype=float) / float(period)
    return np.column_stack([np.sin(angle), np.cos(angle)])


def equal_width_bins(values, bins: int) -> np.ndarray:
    """Bin edges, from NumPy's own edge calculator (Day 130's bin width)."""
    return np.histogram_bin_edges(np.asarray(values, dtype=float), bins=bins)


def bin_index(values, edges: np.ndarray) -> np.ndarray:
    """Which bin each value falls in, given edges learned elsewhere.

    Values below the first edge or above the last are clamped into the
    end bins rather than dropped: a bin boundary is a decision, and the
    decision has to cover values the training data never contained.
    """
    idx = np.digitize(np.asarray(values, dtype=float), edges[1:-1], right=False)
    return idx.astype(float)


def ratio_feature(numerator, denominator, epsilon: float = 1e-12) -> np.ndarray:
    """One column divided by another, with a guard against a zero divisor."""
    denominator = np.asarray(denominator, dtype=float)
    return np.asarray(numerator, dtype=float) / np.where(
        np.abs(denominator) < epsilon, epsilon, denominator
    )


# ---------------------------------------------------------------------------
# Text
# ---------------------------------------------------------------------------


def tokenize(document: str) -> list[str]:
    return TOKEN.findall(str(document).lower())


@dataclass
class Vocabulary:
    """A bag-of-words vocabulary, chosen by association with the target.

    `top_k` words are kept: the ones whose presence correlates most
    strongly with the label. That selection reads the target, so which
    rows you fit it on is exactly as consequential as it is for a target
    encoding. `min_docs` drops words too rare to estimate anything from.
    """

    top_k: int = 20
    min_docs: int = 3
    words: list[str] = field(default_factory=list)
    scores: dict[str, float] = field(default_factory=dict)

    def fit(self, documents, y) -> "Vocabulary":
        y = np.asarray(y, dtype=float)
        tokenized = [set(tokenize(d)) for d in documents]
        counts: Counter[str] = Counter()
        for bag in tokenized:
            counts.update(bag)
        candidates = [w for w, c in counts.items() if c >= self.min_docs]
        scored: list[tuple[float, str]] = []
        for word in candidates:
            present = np.array([1.0 if word in bag else 0.0 for bag in tokenized])
            if present.std() == 0.0 or y.std() == 0.0:
                continue
            correlation = float(np.corrcoef(present, y)[0, 1])
            scored.append((abs(correlation), word))
        scored.sort(key=lambda pair: (-pair[0], pair[1]))
        self.words = [word for _, word in scored[: self.top_k]]
        self.scores = {word: score for score, word in scored[: self.top_k]}
        return self

    def transform(self, documents) -> np.ndarray:
        position = {word: i for i, word in enumerate(self.words)}
        matrix = np.zeros((len(documents), len(self.words)), dtype=float)
        for row, document in enumerate(documents):
            for token in tokenize(document):
                column = position.get(token)
                if column is not None:
                    matrix[row, column] += 1.0
        return matrix

    def unseen_words(self, documents) -> set[str]:
        """Tokens in these documents that the vocabulary does not carry."""
        known = set(self.words)
        seen: set[str] = set()
        for document in documents:
            seen.update(tokenize(document))
        return seen - known


# ---------------------------------------------------------------------------
# The audit
# ---------------------------------------------------------------------------


@dataclass(frozen=True)
class Flag:
    """One suspicious feature, with the reason it was flagged."""

    column: str
    rule: str
    detail: str


def leakage_audit(
    frame: pd.DataFrame,
    target: str,
    corr_threshold: float = 0.90,
) -> list[Flag]:
    """Flag columns that look like they already know the answer.

    Three rules, and each one has a name so a flag can be argued with:

    * `correlation` -- a numeric column whose absolute Pearson
      correlation with the target is at or above `corr_threshold`.
    * `separable` -- a numeric column where a single threshold splits the
      classes perfectly, which catches leaks a linear correlation misses.
    * `pure_category` -- a non-numeric column where every category occurs
      with exactly one target value, so the column is the target wearing
      different words.

    What this cannot do is the point of the exercise. It reads one table
    at one moment, so it cannot see a scaler fitted on the wrong rows, it
    cannot see that a column will be unavailable at prediction time, and
    it cannot see a value that was backfilled from the future. It catches
    the loud leaks. The quiet ones are still your job.
    """
    y = frame[target].to_numpy()
    flags: list[Flag] = []
    for column in frame.columns:
        if column == target:
            continue
        series = frame[column]
        if pd.api.types.is_numeric_dtype(series):
            values = series.to_numpy(dtype=float)
            if values.std() == 0.0:
                continue
            correlation = float(np.corrcoef(values, y.astype(float))[0, 1])
            if abs(correlation) >= corr_threshold:
                flags.append(
                    Flag(column, "correlation", f"|r| = {abs(correlation):.3f} with {target!r}")
                )
                continue
            if _perfectly_separable(values, y):
                flags.append(
                    Flag(column, "separable", f"one threshold splits {target!r} without error")
                )
        else:
            grouped = frame.groupby(series.astype("str"))[target].nunique()
            if len(grouped) > 1 and int(grouped.max()) == 1:
                flags.append(
                    Flag(column, "pure_category", f"every category maps to one {target!r} value")
                )
    return flags


def _perfectly_separable(values: np.ndarray, y: np.ndarray) -> bool:
    """True when some threshold on `values` classifies `y` with no error."""
    classes = np.unique(y)
    if len(classes) != 2:
        return False
    order = np.argsort(values)
    sorted_y = y[order]
    positives = np.cumsum(sorted_y == classes[1])
    negatives = np.cumsum(sorted_y == classes[0])
    total_positive = positives[-1]
    total_negative = negatives[-1]
    for i in range(len(sorted_y) - 1):
        if values[order][i] == values[order][i + 1]:
            continue
        left_correct = negatives[i] + (total_positive - positives[i])
        right_correct = positives[i] + (total_negative - negatives[i])
        if left_correct == len(sorted_y) or right_correct == len(sorted_y):
            return True
    return False

"""Two classifiers, written from scratch in NumPy.

scikit-learn is not installed in this lab, and that is deliberate: this
course has not taught modelling yet, so every model here is small enough
to read in one sitting. Both are deterministic -- no shuffling, no random
restarts, no early stopping on a random subset -- so a score is a fact
about the features rather than a fact about the run.

`LogisticRegression` is Day 111's gradient descent applied to the
log-loss. `NearestCentroid` is Day 107's distance, applied to two class
means. Neither is state of the art and neither needs to be: the whole
point of the day is that the feature table decides the score long before
the model does.
"""

from __future__ import annotations

import numpy as np


def sigmoid(z: np.ndarray) -> np.ndarray:
    """A numerically stable logistic function."""
    out = np.empty_like(z, dtype=float)
    positive = z >= 0
    out[positive] = 1.0 / (1.0 + np.exp(-z[positive]))
    exp_z = np.exp(z[~positive])
    out[~positive] = exp_z / (1.0 + exp_z)
    return out


def accuracy(y_true: np.ndarray, y_pred: np.ndarray) -> float:
    """The fraction of predictions that match the label."""
    y_true = np.asarray(y_true).ravel()
    y_pred = np.asarray(y_pred).ravel()
    return float(np.mean(y_true == y_pred))


class LogisticRegression:
    """Binary logistic regression trained by full-batch gradient descent.

    The update is exactly Day 111's: subtract the learning rate times the
    gradient of the mean log-loss, `X.T @ (p - y) / n`, and repeat for a
    fixed number of steps. Weights start at zero, so two fits on the same
    data give bit-identical coefficients.
    """

    def __init__(self, learning_rate: float = 0.25, steps: int = 3000, l2: float = 0.0) -> None:
        self.learning_rate = learning_rate
        self.steps = steps
        self.l2 = l2
        self.weights: np.ndarray | None = None
        self.bias: float = 0.0

    def fit(self, X: np.ndarray, y: np.ndarray) -> "LogisticRegression":
        X = np.asarray(X, dtype=float)
        y = np.asarray(y, dtype=float).ravel()
        n, d = X.shape
        self.weights = np.zeros(d)
        self.bias = 0.0
        for _ in range(self.steps):
            p = sigmoid(X @ self.weights + self.bias)
            error = p - y
            grad_w = X.T @ error / n + self.l2 * self.weights
            grad_b = float(np.mean(error))
            self.weights -= self.learning_rate * grad_w
            self.bias -= self.learning_rate * grad_b
        return self

    def predict_proba(self, X: np.ndarray) -> np.ndarray:
        if self.weights is None:
            raise RuntimeError("fit before predict")
        return sigmoid(np.asarray(X, dtype=float) @ self.weights + self.bias)

    def predict(self, X: np.ndarray) -> np.ndarray:
        return (self.predict_proba(X) >= 0.5).astype(int)

    def score(self, X: np.ndarray, y: np.ndarray) -> float:
        return accuracy(y, self.predict(X))


class NearestCentroid:
    """Assign each row to the class whose mean it is closest to.

    Distance is the Euclidean norm of Day 107, which is why this model
    cares about scale: a feature measured in thousands dominates the sum
    of squares no matter how little it says about the label.
    """

    def __init__(self) -> None:
        self.classes_: np.ndarray | None = None
        self.centroids_: np.ndarray | None = None

    def fit(self, X: np.ndarray, y: np.ndarray) -> "NearestCentroid":
        X = np.asarray(X, dtype=float)
        y = np.asarray(y).ravel()
        self.classes_ = np.unique(y)
        self.centroids_ = np.vstack([X[y == c].mean(axis=0) for c in self.classes_])
        return self

    def predict(self, X: np.ndarray) -> np.ndarray:
        if self.centroids_ is None or self.classes_ is None:
            raise RuntimeError("fit before predict")
        X = np.asarray(X, dtype=float)
        distances = np.linalg.norm(X[:, None, :] - self.centroids_[None, :, :], axis=2)
        return self.classes_[np.argmin(distances, axis=1)]

    def score(self, X: np.ndarray, y: np.ndarray) -> float:
        return accuracy(y, self.predict(X))


def random_split(n: int, test_size: int, seed: int) -> tuple[np.ndarray, np.ndarray]:
    """Row indices for a seeded random train/test split."""
    rng = np.random.default_rng(seed)
    order = rng.permutation(n)
    return order[test_size:], order[:test_size]


def time_ordered_split(n: int, test_size: int) -> tuple[np.ndarray, np.ndarray]:
    """Row indices for a split that puts the LAST rows in the test set.

    The rows must already be in time order. This is the split that tells
    you what happens when the model meets a day it has never seen.
    """
    index = np.arange(n)
    return index[: n - test_size], index[n - test_size :]

"""The reference answers for Day 137 -- "Features That Do Not Cheat".

Nine exercises, each one a measurement rather than an opinion. Read
`starter/test_features.py` and try them yourself before reading this.

Every band in here was chosen after running the experiment, not before,
and the bands are wide enough to describe the result honestly rather than
narrow enough to look impressive. Where a number is exact arithmetic --
exercise 5 -- the assertion is exact.
"""

from __future__ import annotations

import math

import numpy as np
import pytest

import data
import experiments as E
import features as F
import models as M


# --- 1 ---------------------------------------------------------------------


def test_target_leakage_is_implausibly_good_and_removing_it_is_honest(leakage):
    """A feature derived from the outcome scores perfectly. That is the bug."""
    assert leakage["with_leak"] >= 0.99, "the planted leak should be near-perfect"
    assert leakage["with_leak"] == 1.0

    # The honest band. Three weak behavioural features on a task with a
    # 44% base rate cannot do much better than this, and a result far
    # above it would itself be a reason to go looking for another leak.
    assert 0.55 <= leakage["without_leak"] <= 0.80
    assert leakage["gap_points"] >= 25.0

    # And the reason the leak is a leak: the column does not exist until
    # after the outcome it is predicting.
    frame = data.signups()
    unconverted = frame.loc[frame["converted"] == 0, "days_to_first_invoice"]
    converted = frame.loc[frame["converted"] == 1, "days_to_first_invoice"]
    assert set(unconverted.unique()) == {-1.0}
    assert converted.min() >= 1.0


# --- 2 ---------------------------------------------------------------------


def test_fitting_before_the_split_costs_nothing_for_a_scaler_and_a_lot_for_an_imputer(
    scaler_contamination, imputer_contamination
):
    """Contamination is only worth what the contaminated statistic is worth.

    The scaler result is the surprise, and it is measured over 200
    splits rather than argued: standardisation applies ONE affine map to
    both halves, so contaminating it can only change the relative
    weighting of the features, and that is worth almost nothing here.
    The group-mean imputer is not affine -- it fills each gap from its
    own group -- and contaminating it is worth several points.
    """
    assert scaler_contamination["trials"] == 200
    assert abs(scaler_contamination["optimism_points"]) < 1.0

    assert imputer_contamination["contaminated"] > imputer_contamination["correct"]
    assert imputer_contamination["optimism_points"] >= 5.0
    assert imputer_contamination["test_size"] == 200

    # The scaler really was contaminated -- its statistics differ plainly.
    # The score did not move; the numbers did.
    frame = data.pricing()
    X = frame[["order_value", "tenure_days"]].to_numpy(dtype=float)
    train_idx, _ = M.random_split(len(frame), 25, 137)
    everything = F.Standardiser().fit(X)
    train_only = F.Standardiser().fit(X[train_idx[:60]])
    assert abs(everything.mean_[0] - train_only.mean_[0]) > 5.0


# --- 3 ---------------------------------------------------------------------


def test_target_encoding_leaks_when_it_is_fitted_before_the_split(encoding):
    """Replace a category with the mean outcome and you have used the target."""
    assert encoding["naive_all_data"] > encoding["out_of_fold"]
    assert encoding["gap_all_vs_oof_points"] >= 4.0

    # Restricting the encoding to the training rows removes the whole
    # inflation. Out-of-fold then buys a little more on top, because a
    # training row's own target no longer feeds its own feature.
    assert encoding["out_of_fold"] >= encoding["naive_train_only"]
    assert encoding["naive_all_data"] - encoding["naive_train_only"] >= 0.04

    # The direct evidence, without a model in the way: the naive encoding
    # is far more correlated with the target than the honest one, and the
    # difference is the part it copied from the answer.
    frame = data.city_signups()
    city = frame["city"].to_numpy()
    y = frame["converted"].to_numpy()
    mapping, default = F.target_encode_fit(city, y)
    naive = F.target_encode_transform(city, mapping, default)
    out_of_fold = F.target_encode_out_of_fold(city, y, n_folds=5, seed=137)
    assert float(np.corrcoef(naive, y)[0, 1]) > 0.25
    assert float(np.corrcoef(out_of_fold, y)[0, 1]) < 0.10


# --- 4 ---------------------------------------------------------------------


def test_a_random_split_hides_what_a_time_ordered_split_shows(temporal):
    """The number you can trust is the smaller one."""
    assert temporal["random_split"] >= 0.80
    assert temporal["time_ordered_split"] <= 0.20
    assert temporal["gap_points"] >= 50.0

    # The mechanism, stated as an assertion: the time-ordered training
    # rows never contain the final batch, and the random ones contain
    # every batch there is.
    assert temporal["batches_seen_by_random_train"] == 6.0
    assert temporal["batches_seen_by_ordered_train"] == 5.0
    assert temporal["test_batches_unseen_by_ordered_train"] == 1.0

    # Worse than uninformed: the model is confidently wrong, scoring far
    # below the majority-class baseline for that period.
    assert temporal["time_ordered_split"] < temporal["majority_rate_in_ordered_test"]


# --- 5 ---------------------------------------------------------------------


def test_cyclical_encoding_restores_the_adjacency_of_hour_23_and_hour_0(cyclical):
    """Exact arithmetic, so these assertions are exact."""
    assert cyclical["raw_23_to_0"] == 23.0
    assert cyclical["raw_3_to_4"] == 1.0
    assert cyclical["raw_adjacent_spread"] == 22.0

    expected = 2.0 * math.sin(math.pi / 24.0)
    assert cyclical["cyclical_23_to_0"] == pytest.approx(expected, abs=1e-12)
    assert cyclical["cyclical_3_to_4"] == pytest.approx(expected, abs=1e-12)
    assert cyclical["cyclical_adjacent_spread"] < 1e-12

    # Opposite hours stay opposite: the circle has diameter 2.
    assert cyclical["cyclical_0_to_12"] == pytest.approx(2.0, abs=1e-12)

    # Every one of the 24 adjacent pairs, wrap included, is the same
    # distance apart -- which is the property the raw integer lacks.
    circle = F.cyclical_encode(np.arange(24), 24)
    distances = [
        float(np.linalg.norm(circle[h] - circle[(h + 1) % 24])) for h in range(24)
    ]
    assert max(distances) - min(distances) < 1e-12


# --- 6 ---------------------------------------------------------------------


def test_an_ordinal_code_forces_an_order_that_one_hot_does_not(colours):
    """Six colours, no order, and a model that has to invent one."""
    assert colours["ordinal_is_monotone"] is True
    assert colours["one_hot_is_monotone"] is False

    # One-hot reproduces each colour's observed rate essentially exactly.
    assert colours["one_hot_max_error"] < 1e-5
    assert colours["one_hot_predictions"] == pytest.approx(
        colours["observed_rates"], abs=1e-5
    )

    # The ordinal model cannot: it is out by more than a third somewhere.
    assert colours["ordinal_max_error"] > 0.30
    assert colours["ordinal_accuracy"] < colours["one_hot_accuracy"]

    # And the reason: the true rates are not monotone in the code.
    observed = colours["observed_rates"]
    assert observed[2] > observed[1]
    assert observed[3] < observed[2]


# --- 7 ---------------------------------------------------------------------


def test_a_ratio_separates_what_neither_column_separates(interaction):
    """Spend and income overlap. Spend over income does not."""
    assert interaction["ratio_only"] == 1.0
    assert interaction["income_only"] < 0.60
    assert interaction["spend_only"] < 0.75
    assert interaction["ratio_only"] - max(
        interaction["income_only"], interaction["spend_only"]
    ) >= 0.25

    # Honest footnote, asserted rather than hidden: because the boundary
    # here is a straight line through the origin, a linear model given
    # both raw columns can reach the same score. The ratio is still the
    # feature that makes the rule visible, and a distance-based model
    # gets nothing from the two raw columns.
    assert interaction["income_and_spend"] >= 0.95


# --- 8 ---------------------------------------------------------------------


def test_the_vocabulary_must_be_chosen_on_training_documents_only(vocabulary):
    """Which words become features is itself a fitted statistic."""
    assert vocabulary["fitted_on_all_data"] > vocabulary["fitted_on_train_only"]
    assert vocabulary["gap_points"] >= 1.5
    assert vocabulary["trials"] == 40

    # Unseen words are handled, not crashed on: the matrix keeps exactly
    # the training vocabulary's width and the test-only tokens are simply
    # not counted.
    assert vocabulary["unseen_test_words"] > 0
    assert vocabulary["test_matrix_columns"] == vocabulary["top_k"]
    assert vocabulary["test_matrix_rows"] == 100.0

    # Transforming a document made only of words the vocabulary has never
    # seen gives a row of zeros rather than an exception.
    trained = F.Vocabulary(top_k=5, min_docs=1).fit(
        ["outage refund now", "thanks hello there"], np.array([1, 0])
    )
    matrix = trained.transform(["quetzal marzipan"])
    assert matrix.shape == (1, len(trained.words))
    assert matrix.sum() == 0.0


# --- 9 ---------------------------------------------------------------------


def test_the_audit_catches_the_planted_leaks_and_leaves_honest_columns_alone(audit):
    """A reusable check, and an honest account of what it cannot see."""
    assert audit["flagged"] == ["days_to_first_invoice", "email_template"]
    assert audit["rules"]["days_to_first_invoice"] == "separable"
    assert audit["rules"]["email_template"] == "pure_category"

    for honest in ("visits", "minutes_on_site", "discount_pct", "channel"):
        assert honest in audit["columns_checked"]
        assert honest not in audit["flagged"]

    # The correlation rule alone would have missed the numeric leak: its
    # absolute correlation with the target is 0.85, under the 0.90
    # threshold. The separability rule is what earns its place.
    frame = E.audit_table()
    y = frame["converted"].to_numpy(dtype=float)
    leak = frame["days_to_first_invoice"].to_numpy(dtype=float)
    correlation = abs(float(np.corrcoef(leak, y)[0, 1]))
    assert 0.80 < correlation < 0.90
    assert F.leakage_audit(frame, "converted", corr_threshold=0.99) != []

    # Raising the threshold above 1 disables the correlation rule
    # entirely, and the separability rule still catches it.
    strict = F.leakage_audit(frame, "converted", corr_threshold=1.01)
    assert [f.column for f in strict] == ["days_to_first_invoice", "email_template"]

lesson_id: D137
day: 137
kind: guided-build
languages:
  - python
  - bash
setup_commands:
  - cd labs/sections/math-statistics-and-data/day-137-thinking-in-features
  - python3 -m venv .venv
  - .venv/bin/pip install -r requirements/requirements.txt
  - '.venv/bin/python3 -c "import pandas, numpy; print(pandas.__version__, numpy.__version__)"'
run_commands:
  - .venv/bin/pytest examples
  - .venv/bin/pytest starter
test_commands:
  - bash tests/run_tests.sh
cleanup_commands:
  - >-
    find . -path ./.venv -prune -o -type d -name '__pycache__' -print -exec rm -rf -- {}
    +
  - rm -rf .pytest_cache
  - 'rm -rf .venv  # optional: removes the lab virtual environment'
  - 'git checkout -- starter/  # optional: reset your work'
requires_network: true
requires_api_key: false
estimated_minutes: 55
last_executed: '2026-08-20'
executed_on: >-
  macOS 26.5.2 (Apple Silicon, arm64), Python 3.14.0, pandas 3.0.5, numpy 2.5.2,
  pytest 9.1.1, bash 3.2.57 -- bash tests/run_tests.sh -> 55 checks, 0 failure(s), exit
  0. pytest examples -> 9 passed in 16.34s. pytest starter -> 9 skipped (untouched
  checkout). Everything was run through a real lab-local .venv created by the documented
  setup commands. scikit-learn is NOT installed and the harness asserts its absence, so
  every model here is written from scratch in NumPy. Section 6 of the harness solves
  every exercise in a scratch copy (9 passed), deliberately breaks exercise 1s leaking
  score assertion (with_leak == 1.0 -> == 0.5), confirms a non-zero exit with a printed
  failure, restores the file and confirms 9 passed again. Separately, exercise 1s gap
  assertion was broken in examples/ itself (>= 25.0 -> >= 99.0) and the whole harness
  was re-run: it reported 55 checks, 6 failure(s) and exited 1; the file was restored
  and the harness returned to 55 checks, 0 failure(s), exit 0. Section 5 confirms
  directly that `pytest examples starter` in one invocation aborts collection with
  `import file mismatch` (both directories define a module named test_features.py).
  MEASURED BEFORE/AFTER PAIRS. Target leakage 1.0000 with the planted column against
  0.6400 without it, a 36.0 point gap. Group-mean imputer fitted on all data 0.7484
  against 0.6662 fitted on training rows, 8.22 points, averaged over 150 splits. Target
  encoding computed over the whole table 0.6215 against 0.5535 out-of-fold, 6.80 points,
  averaged over 40 splits. Random split 0.8833 against time-ordered split 0.0667, an
  81.67 point gap, with the time-ordered score BELOW the 0.9333 majority-class baseline
  for that period. Bag-of-words vocabulary chosen on all documents 0.8137 against 0.7893
  chosen on training documents, 2.45 points, averaged over 40 splits. Cyclical encoding
  is exact arithmetic: raw hours put 23 and 0 exactly 23.0 apart, the circle puts every
  adjacent pair 2*sin(pi/24) = 0.26105238444010315 apart with a spread below 1e-12. The
  audit flags exactly days_to_first_invoice (separable) and email_template
  (pure_category) and none of the four honest columns. FOUR HONESTY CALLS. FIRST and
  most important: the brief expected a scaler fitted on all the data to score HIGHER
  than a correctly fitted one; measured over 200 random splits with a 25-row test set it
  scored 0.06 points LOWER, 0.6218 against 0.6224. Standardisation applies one affine
  map to both halves, so contaminating it can only change the relative weighting of the
  features, and a logistic regression run to convergence is nearly invariant to that.
  Equal-width binning and a rank transform were tried as alternatives over 300 splits at
  three training-set sizes; binning came out 1.9 to 2.7 points in favour of the
  correctly fitted version and the rank transform within half a point either way. The
  lab therefore asserts what is true -- a contaminated scaler is worth under one point
  in either direction -- and demonstrates real contamination with a group-mean imputer,
  which is not an affine map. SECOND: the 81.67 point temporal gap comes from a
  deliberately sharp construction where the alarm rate flips between roughly 0.88 and
  0.08 between calibration batches; real regime changes are milder, and the durable part
  of the finding is that the time-ordered score falls below the majority baseline.
  THIRD: the vocabulary gap is small and configuration-sensitive -- eight combinations
  of top_k and min_docs were measured, six favoured the contaminated vocabulary by 0.6
  to 2.5 points and two favoured the honest one by 0.27, so the assertion floor is 1.5
  points rather than the 2.45 measured. FOURTH: in exercise 7 a logistic regression
  given both raw columns also reaches 1.0000, because the boundary in this construction
  is a straight line through the origin; the lab asserts that rather than hiding it.

# Requirements

`requirements.txt` pins the exact versions this lab was written and run
against on 2026-08-20. Everything else it uses — `re`, `math`,
`collections`, `dataclasses` — is in the Python standard library.

Install into a lab-local virtual environment so the pins cannot collide
with anything else on your machine:

```bash
cd labs/sections/math-statistics-and-data/day-137-thinking-in-features
python3 -m venv .venv
.venv/bin/pip install -r requirements/requirements.txt
```

Only that install step needs the network. Everything after it runs
offline: the lab reads no URL, opens no socket, and needs no API key.

## What is deliberately absent

**scikit-learn is not installed, and the harness asserts that it is
not.** That is not an oversight. This lab sits four days before the
course reaches machine learning, so every model in it is written out in
NumPy — a logistic regression trained by gradient descent and a
nearest-centroid classifier — and both are short enough to read. The
lesson describes scikit-learn's `Pipeline` and `ColumnTransformer` from
their published documentation and reproduces no output from either.

If you install scikit-learn into this environment the version check in
section 1 of the harness will fail, and so will the check that asserts
its absence. That is the harness telling the truth about what it ran
against, which is the whole point of the pins.

## Why the pins are exact

Section 1 of `tests/run_tests.sh` compares every installed version
against this file and fails on a mismatch. Every result in this lab is
generated from a seeded `numpy.random.default_rng`, so the numbers are
reproducible to the last decimal — but only against the same generator
stream. Pinning NumPy is what makes "0.6218 against 0.6224" a fact
rather than an anecdote.

pandas is pinned at 3.0.5 because two of its 3.0 behaviours show up in
this lab's code: Copy-on-Write is always on, and the default string
dtype is `str` rather than `object`. The `pd.Series(..., dtype="str")`
calls in `data.py` are explicit for that reason.

## If a pin will not install

Any recent pandas 2.2+ or NumPy 2.x will almost certainly run the lab.
The version check in section 1 will complain; the nine exercises should
still pass, because every assertion is a band rather than an equality —
except exercise 5, which is exact arithmetic and holds anywhere.
`expected-output/FIELDS.md` records exactly which captured values are
version-sensitive and which are exact everywhere.

pandas==3.0.5
numpy==2.5.2
pytest==9.1.1

# Features That Do Not Cheat — the exercise brief

A feature encodes a hypothesis about what matters. The dangerous ones
quietly encode the answer, and they announce themselves by making your
results look excellent. **A result that looks too good is a bug report.**

Nine exercises. Each one is a measurement: you assert on numbers the code
actually produced, not on numbers you hoped for. Everything is seeded, so
the same numbers come back on every machine.

## Before you start

```bash
cd labs/sections/math-statistics-and-data/day-137-thinking-in-features
python3 -m venv .venv
.venv/bin/pip install -r requirements/requirements.txt
.venv/bin/pytest starter -v          # 9 skipped, on an untouched checkout
```

Run `pytest starter` and `pytest examples` as **two separate commands**.
Both directories hold a module called `test_features.py`, and pytest
collects by dotted module name, so a combined invocation aborts.

## What is where

| File | What it holds |
| --- | --- |
| `data.py` | Seven seeded generators, one per experiment. Each docstring says where a leak was planted. |
| `features.py` | Every encoder, split into `fit` and `transform`, plus the leakage audit |
| `models.py` | A logistic regression trained by gradient descent and a nearest-centroid classifier, both in NumPy |
| `experiments.py` | The nine measurements; each returns a dictionary of numbers |
| `conftest.py` | One session-scoped fixture per experiment, so each runs once |
| `test_features.py` | Your nine exercises |

scikit-learn is **not** installed, which is why every model here is
written out. You have not met a model API yet, and you do not need one:
the whole point of the day is that the feature table decides the score
long before the model does.

## The nine exercises

1. **Target leakage, measured.** Score the same task with and without a
   column derived from the outcome. Report both numbers. The leaking one
   is 1.00, which is your bug report.
2. **Fitting before the split.** Two statistics fitted on all the data:
   a scaler and a group-mean imputer. One of them buys almost nothing
   and one buys a lot. Measure both and work out why they differ.
3. **Target encoding.** Replace a city with the mean outcome for that
   city, three ways: over everything, over training rows, and
   out-of-fold. Measure the gap the first one buys.
4. **Temporal leakage.** Split time-ordered rows at random, then by
   time. One number is trustworthy. Report both.
5. **Cyclical encoding.** Hour 23 and hour 0 are one hour apart in the
   world and 23 apart as integers. Prove that sine and cosine fix it,
   with exact distances.
6. **An ordinal code imposes an order.** Six paint colours with no
   order and a deliberately non-monotone return rate. Show that the
   model given a code can only move monotonically with it.
7. **An interaction.** Neither spend nor income separates the classes.
   Their ratio separates them perfectly. Report all three.
8. **Vocabulary fitted on training only.** Choosing which words become
   features by their association with the label is a fitted statistic,
   so it belongs on the training side of the line.
9. **A leakage audit.** Write the check once and run it on any table.
   Then be honest in your own words about what it cannot see.

## How to work

Delete the `pytest.skip(...)` line, write the assertions the docstring
describes, and run `pytest starter -v` again. When an assertion fails,
print the dictionary the fixture handed you and look at the number before
you change the assertion — the measurement is the thing you are learning
from, and moving a band to make a test pass throws it away.

"""Shared fixtures.

pytest finds this file by itself -- nothing imports it. Every experiment
is wrapped in a session-scoped fixture so the whole suite runs each one
exactly once; several of them average over a hundred or more train/test
splits and would otherwise be repeated for every assertion.

Nothing here writes to disk, binds a port or reaches the network.
"""

from __future__ import annotations

import pytest

import experiments as E


@pytest.fixture(scope="session")
def leakage():
    return E.target_leakage()


@pytest.fixture(scope="session")
def scaler_contamination():
    return E.scaling_contamination()


@pytest.fixture(scope="session")
def imputer_contamination():
    return E.imputer_contamination()


@pytest.fixture(scope="session")
def encoding():
    return E.target_encoding()


@pytest.fixture(scope="session")
def temporal():
    return E.temporal_leakage()


@pytest.fixture(scope="session")
def cyclical():
    return E.cyclical_distances()


@pytest.fixture(scope="session")
def colours():
    return E.ordinal_versus_one_hot()


@pytest.fixture(scope="session")
def interaction():
    return E.interaction()


@pytest.fixture(scope="session")
def vocabulary():
    return E.vocabulary_contamination()


@pytest.fixture(scope="session")
def audit():
    return E.audit_result()

"""The datasets for Day 137 -- "Features That Do Not Cheat".

Every generator here is seeded and returns exactly the same rows on every
machine, so every number this lab asserts is reproducible rather than
sampled. Nothing is downloaded; nothing touches the network.

Six generators, one per experiment:

* `signups`          -- a conversion table carrying one planted target leak
* `pricing`          -- a heavy-tailed feature for the scaling experiment
* `city_signups`     -- 40 city codes for the target-encoding experiment
* `sensor_log`       -- a time-ordered table whose rule changes by regime
* `paint_orders`     -- unordered colours for the ordinal-encoding trap
* `credit_lines`     -- spend and income that only separate as a ratio
* `tickets`          -- short documents for the bag-of-words experiment

Read the docstrings before the code: each one states what the honest
signal is and, where a leak is planted, exactly where it was planted.
"""

from __future__ import annotations

import numpy as np
import pandas as pd

SEED = 137


def _sigmoid(z: np.ndarray) -> np.ndarray:
    return 1.0 / (1.0 + np.exp(-z))


def signups(n: int = 400, seed: int = SEED) -> pd.DataFrame:
    """A conversion table with three honest features and one planted leak.

    Honest: `visits`, `minutes_on_site`, `discount_pct`. The label
    `converted` is drawn from a logistic model of those three, so the
    honest signal is real but far from perfect.

    Planted leak: `days_to_first_invoice`. An invoice only exists once a
    visitor has converted, so the column is a positive number for every
    converted row and the sentinel -1 for every unconverted one. At
    prediction time -- before the visitor has decided -- this column
    cannot exist at all. It is the outcome, wearing a different name.
    """
    rng = np.random.default_rng(seed)
    visits = rng.integers(1, 30, n)
    minutes = rng.gamma(2.0, 6.0, n)
    discount = rng.choice([0, 5, 10, 20], n)
    logit = -2.4 + 0.08 * visits + 0.045 * minutes + 0.035 * discount
    converted = (rng.random(n) < _sigmoid(logit)).astype(int)
    days = np.where(converted == 1, rng.integers(1, 15, n), -1)
    return pd.DataFrame(
        {
            "visits": visits.astype(float),
            "minutes_on_site": minutes,
            "discount_pct": discount.astype(float),
            "days_to_first_invoice": days.astype(float),
            "converted": converted,
        }
    )


def pricing(n: int = 500, seed: int = SEED + 1) -> pd.DataFrame:
    """A two-feature table whose first feature is heavy-tailed.

    `order_value` is lognormal, so a handful of rows sit two orders of
    magnitude above the median. A scaler fitted on all the data sees
    those rows wherever they fall; a scaler fitted on the training half
    only sees the ones that landed in training. That difference is the
    whole point of the scaling experiment, and it is why the feature is
    heavy-tailed rather than Gaussian: with a well-behaved feature the
    two scalers agree to three decimals and there is nothing to see.
    """
    rng = np.random.default_rng(seed)
    order_value = rng.lognormal(mean=3.0, sigma=1.4, size=n)
    tenure_days = rng.gamma(3.0, 90.0, size=n)
    logit = -1.1 + 0.55 * np.log(order_value) - 0.004 * tenure_days
    renewed = (rng.random(n) < _sigmoid(logit)).astype(int)
    return pd.DataFrame(
        {
            "order_value": order_value,
            "tenure_days": tenure_days,
            "renewed": renewed,
        }
    )


def panel_readings(
    n: int = 600,
    n_panels: int = 120,
    missing_rate: float = 0.60,
    seed: int = SEED + 7,
) -> pd.DataFrame:
    """Readings with a lot of gaps, grouped into small panels.

    Each solar panel has its own characteristic output, and `reading` is
    that output plus noise. `fault` depends on the reading. Sixty per
    cent of the readings are missing, and with 120 panels over 600 rows
    a panel holds five rows on average -- so a group mean fitted on a
    small training set is often estimated from one row, or from none.

    That is what makes group-mean imputation the sharpest available
    demonstration of a statistic fitted before the split. Impute over the
    whole table and a test row's gap is filled from the observed readings
    of its own panel, including the ones sitting in the test set.
    """
    rng = np.random.default_rng(seed)
    panel = rng.integers(0, n_panels, n)
    panel_output = rng.normal(0.0, 2.0, n_panels)
    reading = panel_output[panel] + rng.normal(0.0, 0.5, n)
    fault = (rng.random(n) < _sigmoid(1.2 * reading)).astype(int)
    observed = rng.random(n) >= missing_rate
    return pd.DataFrame(
        {
            "panel": pd.Series([f"P{p:03d}" for p in panel], dtype="str"),
            "reading": np.where(observed, reading, np.nan),
            "fault": fault,
        }
    )


def city_signups(n: int = 600, n_cities: int = 40, seed: int = SEED + 2) -> pd.DataFrame:
    """High-cardinality categories for the target-encoding experiment.

    40 city codes over 600 rows is 15 rows per city on average, which is
    exactly the regime where a per-category mean of the target is mostly
    noise. The cities do carry a small real effect, so an honest encoding
    is not useless -- it is just far less impressive than the naive one
    pretends.
    """
    rng = np.random.default_rng(seed)
    city = rng.integers(0, n_cities, n)
    city_effect = rng.normal(0.0, 0.35, n_cities)
    visits = rng.integers(1, 20, n)
    logit = -0.3 + 0.05 * (visits - 10) + city_effect[city]
    converted = (rng.random(n) < _sigmoid(logit)).astype(int)
    return pd.DataFrame(
        {
            "city": pd.Series([f"C{c:02d}" for c in city], dtype="str"),
            "visits": visits.astype(float),
            "converted": converted,
        }
    )


#: Alarm rate per calibration batch, in batch order. A batch is a period
#: of time, and the rate changes when the hardware is recalibrated.
BATCH_ALARM_RATE = [0.85, 0.15, 0.88, 0.12, 0.90, 0.10]


def sensor_log(per_batch: int = 60, seed: int = SEED + 3) -> pd.DataFrame:
    """A time-ordered table whose batches are periods of time.

    Rows come back in time order. `batch` is the calibration batch the
    sensor was running under, so a batch is not scattered through the
    table -- it occupies one contiguous stretch of it. The alarm rate
    changes sharply from batch to batch (a recalibration, a firmware
    change, a new supplier), while `reading` carries a mild, stable
    effect that holds across all of them.

    That combination is what a random split hides. Split at random and
    every batch has rows in training, so the model learns each batch's
    alarm rate from rows recorded at the same time as the ones it is
    scored on -- information from after the prediction moment. Split by
    time and the last batch is one the model has never seen, which is the
    situation every deployed model is actually in.
    """
    rng = np.random.default_rng(seed)
    frames = []
    for batch, rate in enumerate(BATCH_ALARM_RATE):
        reading = rng.normal(50.0, 6.0, per_batch)
        humidity = rng.normal(40.0, 9.0, per_batch)
        logit = np.log(rate / (1 - rate)) + 0.05 * (reading - 50.0)
        alarm = (rng.random(per_batch) < _sigmoid(logit)).astype(int)
        frames.append(
            pd.DataFrame(
                {
                    "batch": pd.Series([f"B{batch}"] * per_batch, dtype="str"),
                    "reading": reading,
                    "humidity": humidity,
                    "alarm": alarm,
                }
            )
        )
    out = pd.concat(frames, ignore_index=True)
    out.insert(0, "t", np.arange(len(out), dtype=int))
    return out


PAINT_COLOURS = ["amber", "cobalt", "ivory", "olive", "rose", "slate"]

#: The true return rate per colour, in the same order as PAINT_COLOURS.
#: Deliberately NOT monotone in the list position: ivory (code 2) is the
#: worst and cobalt (code 1) the best, so any model that can only move
#: monotonically with the code is guaranteed to be wrong somewhere.
PAINT_RETURN_RATE = [0.20, 0.08, 0.72, 0.15, 0.62, 0.25]


def paint_orders(per_colour: int = 200, seed: int = SEED + 4) -> pd.DataFrame:
    """Unordered categories with a deliberately non-monotone outcome.

    Colour names have no order. Alphabetical position is not a quantity,
    which is exactly why an ordinal code invites a model to interpolate
    between categories that have no midpoint.
    """
    rng = np.random.default_rng(seed)
    colours: list[str] = []
    returned: list[int] = []
    for name, rate in zip(PAINT_COLOURS, PAINT_RETURN_RATE):
        colours.extend([name] * per_colour)
        returned.extend((rng.random(per_colour) < rate).astype(int).tolist())
    return pd.DataFrame(
        {
            "colour": pd.Series(colours, dtype="str"),
            "returned": np.array(returned, dtype=int),
        }
    )


def credit_lines(n: int = 500, seed: int = SEED + 5) -> pd.DataFrame:
    """Spend and income whose marginals overlap but whose ratio does not.

    Income is drawn from the same distribution for both classes. Spend is
    income times a ratio drawn per class, and the two ratio bands do not
    overlap. So neither column separates the classes on its own, while
    `spend / income` separates them completely. This is an interaction in
    the plainest possible form.
    """
    rng = np.random.default_rng(seed)
    label = (rng.random(n) < 0.5).astype(int)
    income = rng.uniform(30_000, 140_000, n)
    ratio = np.where(label == 1, rng.uniform(0.52, 0.68, n), rng.uniform(0.32, 0.48, n))
    spend = income * ratio
    return pd.DataFrame(
        {
            "income": income,
            "spend": spend,
            "stressed": label,
        }
    )


#: The corpus vocabulary and each word's pull towards "urgent". No word
#: is decisive: the largest effect only doubles a word's odds of turning
#: up, so a document is classified by the whole bag or not at all.
TICKET_WORDS = [
    "outage", "refund", "broken", "charged", "cancel", "waiting",
    "thanks", "question", "hello", "curious", "manual", "hours",
    "order", "account", "team", "please", "about", "with", "the", "and",
]
TICKET_EFFECT = [
    0.40, 0.34, 0.38, 0.30, 0.36, 0.26,
    -0.40, -0.32, -0.36, -0.28, -0.34, -0.24,
    0.06, -0.05, 0.04, -0.03, 0.02, -0.02, 0.01, -0.01,
]


def tickets(n: int = 300, seed: int = SEED + 6) -> pd.DataFrame:
    """Short support tickets, each a bag of lower-case words.

    Every document draws its words from one shared pool; the label only
    tilts the odds. Twelve words carry a real but modest pull and eight
    carry essentially none, so no single word settles a document and the
    classifier has to add evidence up.

    On top of that sits a long tail of reference codes, each appearing in
    only a handful of documents corpus-wide. Those are the words that
    make the choice of vocabulary consequential: a code that turns up
    three times and happens to land twice on an urgent ticket looks like
    a strong feature to anything that scores words against the label.
    """
    rng = np.random.default_rng(seed)
    effect = np.array(TICKET_EFFECT, dtype=float)
    rows: list[str] = []
    labels: list[int] = []
    for _ in range(n):
        urgent = int(rng.random() < 0.5)
        sign = 1.0 if urgent else -1.0
        weights = np.exp(sign * effect)
        weights = weights / weights.sum()
        length = int(rng.integers(8, 15))
        words = rng.choice(TICKET_WORDS, size=length, replace=True, p=weights).tolist()
        if rng.random() < 0.35:
            words.append(f"ref{int(rng.integers(0, 80)):02d}")
        rows.append(" ".join(words))
        labels.append(urgent)
    return pd.DataFrame(
        {
            "text": pd.Series(rows, dtype="str"),
            "urgent": np.array(labels, dtype=int),
        }
    )

"""The nine experiments, each one a function that returns numbers.

The tests assert on what these return, and `tests/run_tests.sh` prints
them. Keeping the experiments here rather than inside the tests means the
same code produces the numbers in the lesson, the numbers in
`expected-output/` and the numbers the assertions check, so the three
cannot drift apart.

Every function is deterministic. Where an experiment averages over many
random splits it says so in its name and its docstring, and the seeds are
generated from a fixed base seed.
"""

from __future__ import annotations

import numpy as np
import pandas as pd

import data
import features as F
import models as M

HONEST_COLUMNS = ["visits", "minutes_on_site", "discount_pct"]
LEAKY_COLUMNS = HONEST_COLUMNS + ["days_to_first_invoice"]


def _fit_score(X_train, y_train, X_test, y_test, model=None) -> float:
    model = model or M.LogisticRegression()
    model.fit(X_train, y_train)
    return model.score(X_test, y_test)


# --- 1. target leakage -----------------------------------------------------


def target_leakage(test_size: int = 100, seed: int = 137) -> dict[str, float]:
    """Score the same task twice: with the planted leak and without it.

    Both runs use the same split, the same model and the same number of
    gradient steps. The only difference is one column.
    """
    frame = data.signups()
    y = frame["converted"].to_numpy()
    train_idx, test_idx = M.random_split(len(frame), test_size, seed)
    out: dict[str, float] = {}
    for name, columns in (("with_leak", LEAKY_COLUMNS), ("without_leak", HONEST_COLUMNS)):
        X = frame[columns].to_numpy(dtype=float)
        scaler = F.Standardiser().fit(X[train_idx])
        out[name] = _fit_score(
            scaler.transform(X[train_idx]),
            y[train_idx],
            scaler.transform(X[test_idx]),
            y[test_idx],
        )
    out["gap_points"] = 100.0 * (out["with_leak"] - out["without_leak"])
    return out


# --- 2. train/test contamination through a scaler --------------------------


def scaling_contamination(
    trials: int = 200, train_size: int = 60, test_size: int = 25, seed: int = 137
) -> dict[str, float]:
    """Fit the scaler on everything, or on the training rows only.

    Averaged over `trials` random splits, because with a 25-row test set
    one row is four accuracy points and a single split says nothing. The
    model is the nearest-centroid classifier, which is distance-based and
    therefore cares about scale; a logistic regression run to convergence
    is very nearly invariant to an affine change of features and shows
    almost nothing here.

    The mechanism is worth stating plainly: the contaminated scaler is a
    BETTER estimate of the population mean and spread than a scaler
    fitted on 60 rows. It got that way by reading rows it was not allowed
    to read, and the score it produces is therefore not a score you will
    ever see in production.
    """
    frame = data.pricing()
    X_all = frame[["order_value", "tenure_days"]].to_numpy(dtype=float)
    y_all = frame["renewed"].to_numpy()
    contaminated_scores: list[float] = []
    correct_scores: list[float] = []
    for trial in range(trials):
        rng = np.random.default_rng(seed + trial)
        order = rng.permutation(len(frame))
        train_idx = order[:train_size]
        test_idx = order[train_size : train_size + test_size]

        wrong = F.Standardiser().fit(X_all)  # every row, including the test rows
        right = F.Standardiser().fit(X_all[train_idx])  # training rows only

        contaminated_scores.append(
            _fit_score(
                wrong.transform(X_all[train_idx]),
                y_all[train_idx],
                wrong.transform(X_all[test_idx]),
                y_all[test_idx],
                M.NearestCentroid(),
            )
        )
        correct_scores.append(
            _fit_score(
                right.transform(X_all[train_idx]),
                y_all[train_idx],
                right.transform(X_all[test_idx]),
                y_all[test_idx],
                M.NearestCentroid(),
            )
        )
    contaminated = float(np.mean(contaminated_scores))
    correct = float(np.mean(correct_scores))
    return {
        "contaminated": contaminated,
        "correct": correct,
        "optimism_points": 100.0 * (contaminated - correct),
        "trials": float(trials),
        "train_size": float(train_size),
        "test_size": float(test_size),
        "contaminated_mean_order_value": float(F.Standardiser().fit(X_all).mean_[0]),
    }


def imputer_contamination(
    trials: int = 150, train_size: int = 120, test_size: int = 200, seed: int = 137
) -> dict[str, float]:
    """The same question asked of an imputer instead of a scaler.

    A group-mean imputer is not an affine map. Each missing value is
    filled with a number computed from the rows that happen to be
    observed in its own group -- so if the fit sees the test rows, a
    test row's gap is filled using the readings of its own panel that
    live in the test set. That is information the deployed system will
    not have, and here it is worth several accuracy points rather than
    the fraction of a point the scaler was worth.
    """
    frame = data.panel_readings()
    panel = frame["panel"].to_numpy()
    reading = frame["reading"].to_numpy(dtype=float)
    y = frame["fault"].to_numpy()
    n = len(frame)
    missing = np.isnan(reading)

    everything = F.GroupMeanImputer().fit(panel, reading)

    def build(rows, imputer) -> np.ndarray:
        filled = imputer.transform(panel[rows], reading[rows])
        return np.column_stack([filled, missing[rows].astype(float)])

    contaminated_scores: list[float] = []
    correct_scores: list[float] = []
    for trial in range(trials):
        rng = np.random.default_rng(seed + trial)
        order = rng.permutation(n)
        train_idx = order[:train_size]
        test_idx = order[train_size : train_size + test_size]
        train_only = F.GroupMeanImputer().fit(panel[train_idx], reading[train_idx])
        for imputer, bucket in ((everything, contaminated_scores), (train_only, correct_scores)):
            X_train = build(train_idx, imputer)
            X_test = build(test_idx, imputer)
            scaler = F.Standardiser().fit(X_train)
            bucket.append(
                _fit_score(
                    scaler.transform(X_train), y[train_idx], scaler.transform(X_test), y[test_idx]
                )
            )
    contaminated = float(np.mean(contaminated_scores))
    correct = float(np.mean(correct_scores))
    return {
        "contaminated": contaminated,
        "correct": correct,
        "optimism_points": 100.0 * (contaminated - correct),
        "trials": float(trials),
        "train_size": float(train_size),
        "test_size": float(test_size),
        "missing_fraction": float(missing.mean()),
        "groups_known_to_all_data_fit": float(len(everything.means_)),
    }


# --- 3. target encoding ----------------------------------------------------


def target_encoding(
    test_size: int = 150, seed: int = 137, n_folds: int = 5, trials: int = 40
) -> dict[str, float]:
    """Three ways to encode a city, scored on the same held-out rows.

    * `naive_all_data` computes the per-city mean of the target over the
      whole table, then splits. Every test row's feature was computed
      partly from its own answer.
    * `naive_train_only` computes it on the training rows only, which
      removes the test rows' contribution but still lets each training
      row see its own target.
    * `out_of_fold` encodes each training row from the folds that do not
      contain it.
    """
    frame = data.city_signups()
    y = frame["converted"].to_numpy()
    visits = frame["visits"].to_numpy(dtype=float)
    city = frame["city"].to_numpy()

    all_scores: list[float] = []
    train_scores: list[float] = []
    oof_scores: list[float] = []
    for trial in range(trials):
        train_idx, test_idx = M.random_split(len(frame), test_size, seed + trial)

        def score(train_codes, test_codes) -> float:
            X_train = np.column_stack([visits[train_idx], train_codes])
            X_test = np.column_stack([visits[test_idx], test_codes])
            scaler = F.Standardiser().fit(X_train)
            return _fit_score(
                scaler.transform(X_train), y[train_idx], scaler.transform(X_test), y[test_idx]
            )

        all_map, all_default = F.target_encode_fit(city, y)
        all_scores.append(
            score(
                F.target_encode_transform(city[train_idx], all_map, all_default),
                F.target_encode_transform(city[test_idx], all_map, all_default),
            )
        )

        train_map, train_default = F.target_encode_fit(city[train_idx], y[train_idx])
        train_scores.append(
            score(
                F.target_encode_transform(city[train_idx], train_map, train_default),
                F.target_encode_transform(city[test_idx], train_map, train_default),
            )
        )

        oof = F.target_encode_out_of_fold(
            city[train_idx], y[train_idx], n_folds=n_folds, seed=seed + trial
        )
        oof_scores.append(
            score(oof, F.target_encode_transform(city[test_idx], train_map, train_default))
        )

    naive_all = float(np.mean(all_scores))
    naive_train = float(np.mean(train_scores))
    out_of_fold = float(np.mean(oof_scores))
    return {
        "naive_all_data": naive_all,
        "naive_train_only": naive_train,
        "out_of_fold": out_of_fold,
        "gap_all_vs_oof_points": 100.0 * (naive_all - out_of_fold),
        "gap_train_vs_oof_points": 100.0 * (naive_train - out_of_fold),
        "n_folds": float(n_folds),
        "trials": float(trials),
    }


# --- 4. temporal leakage ---------------------------------------------------


def temporal_leakage(test_size: int = 60, seed: int = 137) -> dict[str, float]:
    """The same table, split two ways: at random, and by time.

    The rows are already in time order. The random split lets rows from
    the last regime into training, so the model is scored on a rule it
    has already seen. The time-ordered split holds the last regime out
    entirely, which is the only one of the two that answers the question
    anybody actually asked.
    """
    frame = data.sensor_log()
    y = frame["alarm"].to_numpy()
    reading = frame["reading"].to_numpy(dtype=float)
    batch = frame["batch"].to_numpy()

    def score(train_idx, test_idx) -> float:
        # The batch encoding is fitted on the training rows only, which is
        # the correct thing to do and is exactly what exposes the problem:
        # a batch the training rows never contained becomes an all-zero
        # block, and the model has nothing to say about it.
        seen = sorted(set(batch[train_idx].tolist()))
        train_block, _ = F.one_hot(batch[train_idx], seen)
        test_block, _ = F.one_hot(batch[test_idx], seen)
        X_train = np.column_stack([reading[train_idx], train_block])
        X_test = np.column_stack([reading[test_idx], test_block])
        scaler = F.Standardiser().fit(X_train)
        return _fit_score(
            scaler.transform(X_train), y[train_idx], scaler.transform(X_test), y[test_idx]
        )

    random_train, random_test = M.random_split(len(frame), test_size, seed)
    ordered_train, ordered_test = M.time_ordered_split(len(frame), test_size)
    random_score = score(random_train, random_test)
    ordered_score = score(ordered_train, ordered_test)
    return {
        "random_split": random_score,
        "time_ordered_split": ordered_score,
        "gap_points": 100.0 * (random_score - ordered_score),
        "batches_seen_by_random_train": float(len(set(batch[random_train].tolist()))),
        "batches_seen_by_ordered_train": float(len(set(batch[ordered_train].tolist()))),
        "test_batches_unseen_by_ordered_train": float(
            len(set(batch[ordered_test].tolist()) - set(batch[ordered_train].tolist()))
        ),
        "majority_rate_in_ordered_test": float(max(y[ordered_test].mean(), 1 - y[ordered_test].mean())),
    }


# --- 5. cyclical encoding --------------------------------------------------


def cyclical_distances(period: int = 24) -> dict[str, float]:
    """Distances between hours, raw and on the circle.

    Raw integer hours: 23 and 0 sit 23 units apart while 3 and 4 sit 1
    apart. Sine-cosine: every adjacent pair, wrap included, sits exactly
    2*sin(pi/24) apart.
    """
    hours = np.arange(period)
    raw = hours.reshape(-1, 1).astype(float)
    circle = F.cyclical_encode(hours, period)

    def distance(matrix, a: int, b: int) -> float:
        return float(np.linalg.norm(matrix[a] - matrix[b]))

    adjacent_raw = [distance(raw, h, (h + 1) % period) for h in range(period)]
    adjacent_circle = [distance(circle, h, (h + 1) % period) for h in range(period)]
    return {
        "raw_23_to_0": distance(raw, 23, 0),
        "raw_3_to_4": distance(raw, 3, 4),
        "cyclical_23_to_0": distance(circle, 23, 0),
        "cyclical_3_to_4": distance(circle, 3, 4),
        "cyclical_0_to_12": distance(circle, 0, 12),
        "cyclical_adjacent_spread": float(max(adjacent_circle) - min(adjacent_circle)),
        "raw_adjacent_spread": float(max(adjacent_raw) - min(adjacent_raw)),
        "expected_adjacent": float(2.0 * np.sin(np.pi / period)),
    }


# --- 6. ordinal versus one-hot --------------------------------------------


def ordinal_versus_one_hot(seed: int = 137) -> dict[str, object]:
    """Fit the same model on an ordinal code and on a one-hot block.

    The colours have no order, and their return rates are deliberately
    not monotone in the alphabetical code. A model reading the code as a
    number can only produce predictions that rise or fall with it; the
    one-hot model is free to give each colour its own answer.
    """
    frame = data.paint_orders()
    y = frame["returned"].to_numpy()
    colours = frame["colour"].to_numpy()
    order = list(data.PAINT_COLOURS)

    codes = F.ordinal_encode(colours, order).reshape(-1, 1)
    ordinal_model = M.LogisticRegression(learning_rate=0.3, steps=4000).fit(codes, y)
    ordinal_rates = ordinal_model.predict_proba(
        np.arange(len(order), dtype=float).reshape(-1, 1)
    )

    block, categories = F.one_hot(colours, order)
    one_hot_model = M.LogisticRegression(learning_rate=0.3, steps=4000).fit(block, y)
    one_hot_rates = one_hot_model.predict_proba(np.eye(len(order)))

    observed = np.array(
        [float(y[colours == name].mean()) for name in order], dtype=float
    )

    def is_monotone(values: np.ndarray) -> bool:
        diffs = np.diff(values)
        return bool(np.all(diffs >= -1e-12) or np.all(diffs <= 1e-12))

    return {
        "categories": categories,
        "observed_rates": observed.tolist(),
        "ordinal_predictions": ordinal_rates.tolist(),
        "one_hot_predictions": one_hot_rates.tolist(),
        "ordinal_is_monotone": is_monotone(ordinal_rates),
        "one_hot_is_monotone": is_monotone(one_hot_rates),
        "ordinal_max_error": float(np.max(np.abs(ordinal_rates - observed))),
        "one_hot_max_error": float(np.max(np.abs(one_hot_rates - observed))),
        "ordinal_accuracy": float(ordinal_model.score(codes, y)),
        "one_hot_accuracy": float(one_hot_model.score(block, y)),
    }


# --- 7. an interaction -----------------------------------------------------


def interaction(test_size: int = 150, seed: int = 137) -> dict[str, float]:
    """Score income alone, spend alone, and the ratio of the two."""
    frame = data.credit_lines()
    y = frame["stressed"].to_numpy()
    train_idx, test_idx = M.random_split(len(frame), test_size, seed)
    income = frame["income"].to_numpy(dtype=float)
    spend = frame["spend"].to_numpy(dtype=float)
    ratio = F.ratio_feature(spend, income)

    def score(column: np.ndarray) -> float:
        X = column.reshape(-1, 1)
        scaler = F.Standardiser().fit(X[train_idx])
        return _fit_score(
            scaler.transform(X[train_idx]), y[train_idx], scaler.transform(X[test_idx]), y[test_idx]
        )

    both = np.column_stack([income, spend])
    scaler = F.Standardiser().fit(both[train_idx])
    both_score = _fit_score(
        scaler.transform(both[train_idx]),
        y[train_idx],
        scaler.transform(both[test_idx]),
        y[test_idx],
    )
    return {
        "income_only": score(income),
        "spend_only": score(spend),
        "ratio_only": score(ratio),
        "income_and_spend": both_score,
    }


# --- 8. vocabulary fitted on the wrong rows --------------------------------


def vocabulary_contamination(
    top_k: int = 30, test_size: int = 100, seed: int = 137, min_docs: int = 2, trials: int = 40
) -> dict[str, float]:
    """Choose the words on everything, or on the training documents only.

    The vocabulary is selected by association with the label, so choosing
    it on all the documents means the test labels helped decide which
    features exist. That is feature selection fitted before the split --
    the second kind of leakage, wearing a text-shaped costume.
    """
    frame = data.tickets()
    documents = frame["text"].tolist()
    y = frame["urgent"].to_numpy()
    all_data_vocab = F.Vocabulary(top_k=top_k, min_docs=min_docs).fit(documents, y)

    all_scores: list[float] = []
    train_scores: list[float] = []
    shared: list[int] = []
    for trial in range(trials):
        train_idx, test_idx = M.random_split(len(frame), test_size, seed + trial)
        train_docs = [documents[i] for i in train_idx]
        test_docs = [documents[i] for i in test_idx]
        train_only_vocab = F.Vocabulary(top_k=top_k, min_docs=min_docs).fit(
            train_docs, y[train_idx]
        )
        shared.append(len(set(all_data_vocab.words) & set(train_only_vocab.words)))
        for vocabulary, bucket in ((all_data_vocab, all_scores), (train_only_vocab, train_scores)):
            bucket.append(
                _fit_score(
                    vocabulary.transform(train_docs),
                    y[train_idx],
                    vocabulary.transform(test_docs),
                    y[test_idx],
                )
            )

    # One split, examined in detail, for the unseen-word behaviour.
    train_idx, test_idx = M.random_split(len(frame), test_size, seed)
    train_docs = [documents[i] for i in train_idx]
    test_docs = [documents[i] for i in test_idx]
    train_only_vocab = F.Vocabulary(top_k=top_k, min_docs=min_docs).fit(train_docs, y[train_idx])
    unseen = train_only_vocab.unseen_words(test_docs)
    test_matrix = train_only_vocab.transform(test_docs)

    fitted_all = float(np.mean(all_scores))
    fitted_train = float(np.mean(train_scores))
    return {
        "fitted_on_all_data": fitted_all,
        "fitted_on_train_only": fitted_train,
        "gap_points": 100.0 * (fitted_all - fitted_train),
        "unseen_test_words": float(len(unseen)),
        "test_matrix_columns": float(test_matrix.shape[1]),
        "test_matrix_rows": float(test_matrix.shape[0]),
        "shared_words": float(np.mean(shared)),
        "top_k": float(top_k),
        "trials": float(trials),
    }


# --- extra: a bin boundary is a decision (Day 130's bin width again) -------


def binning_decision(bins: int = 3) -> dict[str, object]:
    """Bin the same column two ways and read two different stories off it.

    Equal-width edges come from NumPy's own edge calculator. Equal-count
    edges come from the quantiles of the same column. Neither is wrong.
    They put wildly different numbers of rows in the top bin, and the
    renewal rate you would quote for "high value orders" depends entirely
    on which one you picked.
    """
    frame = data.pricing()
    values = frame["order_value"].to_numpy(dtype=float)
    y = frame["renewed"].to_numpy()

    width_edges = F.equal_width_bins(values, bins)
    count_edges = np.quantile(values, np.linspace(0.0, 1.0, bins + 1))
    width_index = F.bin_index(values, width_edges)
    count_index = F.bin_index(values, count_edges)

    def top_bin(index: np.ndarray) -> tuple[int, float]:
        top = index == index.max()
        return int(top.sum()), float(y[top].mean())

    width_rows, width_rate = top_bin(width_index)
    count_rows, count_rate = top_bin(count_index)
    return {
        "bins": bins,
        "equal_width_edges": [round(float(e), 2) for e in width_edges],
        "equal_count_edges": [round(float(e), 2) for e in count_edges],
        "equal_width_top_bin_rows": width_rows,
        "equal_count_top_bin_rows": count_rows,
        "equal_width_top_bin_rate": width_rate,
        "equal_count_top_bin_rate": count_rate,
        "rows": int(len(values)),
    }


# --- 9. the audit ----------------------------------------------------------


def audit_table() -> pd.DataFrame:
    """The signups table plus one categorical leak and one honest column."""
    frame = data.signups().copy()
    rng = np.random.default_rng(data.SEED + 9)
    converted = frame["converted"].to_numpy()
    frame["email_template"] = pd.Series(
        np.where(converted == 1, "welcome_pack", "abandoned_cart"), dtype="str"
    )
    frame["channel"] = pd.Series(
        rng.choice(["search", "social", "direct"], len(frame)), dtype="str"
    )
    return frame


def audit_result(corr_threshold: float = 0.90) -> dict[str, object]:
    frame = audit_table()
    flags = F.leakage_audit(frame, "converted", corr_threshold=corr_threshold)
    return {
        "flagged": [f.column for f in flags],
        "rules": {f.column: f.rule for f in flags},
        "details": {f.column: f.detail for f in flags},
        "columns_checked": [c for c in frame.columns if c != "converted"],
    }

"""Feature builders, every one of them split into `fit` and `transform`.

The split is the whole discipline of this lab. `fit` looks at data and
learns a statistic -- a mean, a standard deviation, a category's average
outcome, a vocabulary. `transform` applies a statistic that has already
been learned and looks at nothing. Once the two are separate functions
you can point `fit` at the training rows only, and the question "did the
test set influence this number?" has an answer you can read off the call
site instead of guessing at it.

Nothing here imports scikit-learn -- it is not installed in this lab.
Everything is NumPy and pandas, small enough to read.
"""

from __future__ import annotations

import re
from collections import Counter
from dataclasses import dataclass, field

import numpy as np
import pandas as pd

TOKEN = re.compile(r"[a-z0-9]+")


# ---------------------------------------------------------------------------
# Scaling
# ---------------------------------------------------------------------------


class Standardiser:
    """Subtract a mean and divide by a standard deviation (Day 107).

    `fit` records the two statistics. `transform` applies them. Fitting on
    rows you will later score is the second kind of leakage, and the only
    thing stopping you is which rows you hand to `fit`.
    """

    def __init__(self) -> None:
        self.mean_: np.ndarray | None = None
        self.scale_: np.ndarray | None = None

    def fit(self, X: np.ndarray) -> "Standardiser":
        X = np.asarray(X, dtype=float)
        self.mean_ = X.mean(axis=0)
        scale = X.std(axis=0)
        # A constant column has zero spread; dividing by it would produce
        # infinities, so it is left alone rather than exploded.
        self.scale_ = np.where(scale == 0.0, 1.0, scale)
        return self

    def transform(self, X: np.ndarray) -> np.ndarray:
        if self.mean_ is None or self.scale_ is None:
            raise RuntimeError("fit before transform")
        return (np.asarray(X, dtype=float) - self.mean_) / self.scale_

    def fit_transform(self, X: np.ndarray) -> np.ndarray:
        return self.fit(X).transform(X)


class GroupMeanImputer:
    """Fill a missing value with the mean of its group (Day 125).

    `fit` records one mean per group, plus a global fallback for groups
    it never saw. Like every other statistic in this file it reads rows,
    so which rows it reads is a decision -- and with small groups it is
    the most consequential decision in this lab.
    """

    def __init__(self) -> None:
        self.means_: dict[str, float] = {}
        self.default_: float = 0.0

    def fit(self, groups, values) -> "GroupMeanImputer":
        frame = pd.DataFrame(
            {"group": pd.Series(groups).astype("str"), "value": np.asarray(values, dtype=float)}
        ).dropna(subset=["value"])
        self.means_ = {
            str(k): float(v) for k, v in frame.groupby("group")["value"].mean().items()
        }
        self.default_ = float(frame["value"].mean()) if len(frame) else 0.0
        return self

    def transform(self, groups, values) -> np.ndarray:
        values = np.asarray(values, dtype=float)
        filled = values.copy()
        for i, (group, value) in enumerate(zip(pd.Series(groups).astype("str"), values)):
            if np.isnan(value):
                filled[i] = self.means_.get(str(group), self.default_)
        return filled


# ---------------------------------------------------------------------------
# Categorical encodings
# ---------------------------------------------------------------------------


def one_hot(values, categories: list[str] | None = None) -> tuple[np.ndarray, list[str]]:
    """One column of 0/1 per category, in a fixed, explicit category order.

    The category list is returned so the caller can reuse it at transform
    time. A test row carrying a category the training data never held
    becomes an all-zero row rather than a crash or a new column.
    """
    values = pd.Series(values).astype("str")
    if categories is None:
        categories = sorted(values.unique().tolist())
    matrix = np.zeros((len(values), len(categories)), dtype=float)
    position = {name: i for i, name in enumerate(categories)}
    for row, name in enumerate(values):
        column = position.get(name)
        if column is not None:
            matrix[row, column] = 1.0
    return matrix, list(categories)


def ordinal_encode(values, order: list[str]) -> np.ndarray:
    """Replace each category with its position in `order`.

    Honest when the order is real (small, medium, large). A trap when it
    is not: the code is a number, and a model reading it as a number can
    and will interpolate between categories that have no midpoint.
    """
    position = {name: float(i) for i, name in enumerate(order)}
    return np.array([position[str(v)] for v in values], dtype=float)


def target_encode_fit(categories, y) -> tuple[dict[str, float], float]:
    """Learn the mean outcome per category, plus the global mean.

    Returns the map and the fallback. Call this on TRAINING rows only:
    it reads the target, so fitting it on everything hands each test row
    a number computed partly from its own answer.
    """
    frame = pd.DataFrame({"category": pd.Series(categories).astype("str"), "y": np.asarray(y, dtype=float)})
    means = frame.groupby("category")["y"].mean()
    return {str(k): float(v) for k, v in means.items()}, float(frame["y"].mean())


def target_encode_transform(categories, mapping: dict[str, float], default: float) -> np.ndarray:
    """Apply a learned target encoding; unseen categories get the prior."""
    return np.array([mapping.get(str(c), default) for c in categories], dtype=float)


def target_encode_out_of_fold(categories, y, n_folds: int = 5, seed: int = 0) -> np.ndarray:
    """Encode each training row from the folds that do NOT contain it.

    This is the fix, and it is worth stating precisely what it fixes: a
    row's own target never contributes to its own feature value. The
    encoding is still built from training data only; out-of-fold does not
    excuse fitting on the test set.
    """
    categories = pd.Series(categories).astype("str").to_numpy()
    y = np.asarray(y, dtype=float)
    n = len(y)
    rng = np.random.default_rng(seed)
    fold_of = rng.permutation(n) % n_folds
    encoded = np.empty(n, dtype=float)
    for fold in range(n_folds):
        held_out = fold_of == fold
        mapping, default = target_encode_fit(categories[~held_out], y[~held_out])
        encoded[held_out] = target_encode_transform(categories[held_out], mapping, default)
    return encoded


# ---------------------------------------------------------------------------
# Datetime, cyclical, binning, interactions
# ---------------------------------------------------------------------------


def calendar_features(timestamps: pd.Series) -> pd.DataFrame:
    """The obvious datetime parts, pulled out as separate columns."""
    ts = pd.to_datetime(timestamps)
    return pd.DataFrame(
        {
            "hour": ts.dt.hour.astype(float),
            "day_of_week": ts.dt.dayofweek.astype(float),
            "day_of_month": ts.dt.day.astype(float),
            "month": ts.dt.month.astype(float),
            "is_weekend": (ts.dt.dayofweek >= 5).astype(float),
        }
    )


def cyclical_encode(values, period: float) -> np.ndarray:
    """Map a wrapping quantity onto a circle: two columns, sine and cosine.

    Hour 23 and hour 0 are one hour apart in the world and 23 apart as
    integers. On the circle they are neighbours again, and every pair of
    adjacent hours sits exactly the same distance apart.
    """
    angle = 2.0 * np.pi * np.asarray(values, dtype=float) / float(period)
    return np.column_stack([np.sin(angle), np.cos(angle)])


def equal_width_bins(values, bins: int) -> np.ndarray:
    """Bin edges, from NumPy's own edge calculator (Day 130's bin width)."""
    return np.histogram_bin_edges(np.asarray(values, dtype=float), bins=bins)


def bin_index(values, edges: np.ndarray) -> np.ndarray:
    """Which bin each value falls in, given edges learned elsewhere.

    Values below the first edge or above the last are clamped into the
    end bins rather than dropped: a bin boundary is a decision, and the
    decision has to cover values the training data never contained.
    """
    idx = np.digitize(np.asarray(values, dtype=float), edges[1:-1], right=False)
    return idx.astype(float)


def ratio_feature(numerator, denominator, epsilon: float = 1e-12) -> np.ndarray:
    """One column divided by another, with a guard against a zero divisor."""
    denominator = np.asarray(denominator, dtype=float)
    return np.asarray(numerator, dtype=float) / np.where(
        np.abs(denominator) < epsilon, epsilon, denominator
    )


# ---------------------------------------------------------------------------
# Text
# ---------------------------------------------------------------------------


def tokenize(document: str) -> list[str]:
    return TOKEN.findall(str(document).lower())


@dataclass
class Vocabulary:
    """A bag-of-words vocabulary, chosen by association with the target.

    `top_k` words are kept: the ones whose presence correlates most
    strongly with the label. That selection reads the target, so which
    rows you fit it on is exactly as consequential as it is for a target
    encoding. `min_docs` drops words too rare to estimate anything from.
    """

    top_k: int = 20
    min_docs: int = 3
    words: list[str] = field(default_factory=list)
    scores: dict[str, float] = field(default_factory=dict)

    def fit(self, documents, y) -> "Vocabulary":
        y = np.asarray(y, dtype=float)
        tokenized = [set(tokenize(d)) for d in documents]
        counts: Counter[str] = Counter()
        for bag in tokenized:
            counts.update(bag)
        candidates = [w for w, c in counts.items() if c >= self.min_docs]
        scored: list[tuple[float, str]] = []
        for word in candidates:
            present = np.array([1.0 if word in bag else 0.0 for bag in tokenized])
            if present.std() == 0.0 or y.std() == 0.0:
                continue
            correlation = float(np.corrcoef(present, y)[0, 1])
            scored.append((abs(correlation), word))
        scored.sort(key=lambda pair: (-pair[0], pair[1]))
        self.words = [word for _, word in scored[: self.top_k]]
        self.scores = {word: score for score, word in scored[: self.top_k]}
        return self

    def transform(self, documents) -> np.ndarray:
        position = {word: i for i, word in enumerate(self.words)}
        matrix = np.zeros((len(documents), len(self.words)), dtype=float)
        for row, document in enumerate(documents):
            for token in tokenize(document):
                column = position.get(token)
                if column is not None:
                    matrix[row, column] += 1.0
        return matrix

    def unseen_words(self, documents) -> set[str]:
        """Tokens in these documents that the vocabulary does not carry."""
        known = set(self.words)
        seen: set[str] = set()
        for document in documents:
            seen.update(tokenize(document))
        return seen - known


# ---------------------------------------------------------------------------
# The audit
# ---------------------------------------------------------------------------


@dataclass(frozen=True)
class Flag:
    """One suspicious feature, with the reason it was flagged."""

    column: str
    rule: str
    detail: str


def leakage_audit(
    frame: pd.DataFrame,
    target: str,
    corr_threshold: float = 0.90,
) -> list[Flag]:
    """Flag columns that look like they already know the answer.

    Three rules, and each one has a name so a flag can be argued with:

    * `correlation` -- a numeric column whose absolute Pearson
      correlation with the target is at or above `corr_threshold`.
    * `separable` -- a numeric column where a single threshold splits the
      classes perfectly, which catches leaks a linear correlation misses.
    * `pure_category` -- a non-numeric column where every category occurs
      with exactly one target value, so the column is the target wearing
      different words.

    What this cannot do is the point of the exercise. It reads one table
    at one moment, so it cannot see a scaler fitted on the wrong rows, it
    cannot see that a column will be unavailable at prediction time, and
    it cannot see a value that was backfilled from the future. It catches
    the loud leaks. The quiet ones are still your job.
    """
    y = frame[target].to_numpy()
    flags: list[Flag] = []
    for column in frame.columns:
        if column == target:
            continue
        series = frame[column]
        if pd.api.types.is_numeric_dtype(series):
            values = series.to_numpy(dtype=float)
            if values.std() == 0.0:
                continue
            correlation = float(np.corrcoef(values, y.astype(float))[0, 1])
            if abs(correlation) >= corr_threshold:
                flags.append(
                    Flag(column, "correlation", f"|r| = {abs(correlation):.3f} with {target!r}")
                )
                continue
            if _perfectly_separable(values, y):
                flags.append(
                    Flag(column, "separable", f"one threshold splits {target!r} without error")
                )
        else:
            grouped = frame.groupby(series.astype("str"))[target].nunique()
            if len(grouped) > 1 and int(grouped.max()) == 1:
                flags.append(
                    Flag(column, "pure_category", f"every category maps to one {target!r} value")
                )
    return flags


def _perfectly_separable(values: np.ndarray, y: np.ndarray) -> bool:
    """True when some threshold on `values` classifies `y` with no error."""
    classes = np.unique(y)
    if len(classes) != 2:
        return False
    order = np.argsort(values)
    sorted_y = y[order]
    positives = np.cumsum(sorted_y == classes[1])
    negatives = np.cumsum(sorted_y == classes[0])
    total_positive = positives[-1]
    total_negative = negatives[-1]
    for i in range(len(sorted_y) - 1):
        if values[order][i] == values[order][i + 1]:
            continue
        left_correct = negatives[i] + (total_positive - positives[i])
        right_correct = positives[i] + (total_negative - negatives[i])
        if left_correct == len(sorted_y) or right_correct == len(sorted_y):
            return True
    return False

"""Two classifiers, written from scratch in NumPy.

scikit-learn is not installed in this lab, and that is deliberate: this
course has not taught modelling yet, so every model here is small enough
to read in one sitting. Both are deterministic -- no shuffling, no random
restarts, no early stopping on a random subset -- so a score is a fact
about the features rather than a fact about the run.

`LogisticRegression` is Day 111's gradient descent applied to the
log-loss. `NearestCentroid` is Day 107's distance, applied to two class
means. Neither is state of the art and neither needs to be: the whole
point of the day is that the feature table decides the score long before
the model does.
"""

from __future__ import annotations

import numpy as np


def sigmoid(z: np.ndarray) -> np.ndarray:
    """A numerically stable logistic function."""
    out = np.empty_like(z, dtype=float)
    positive = z >= 0
    out[positive] = 1.0 / (1.0 + np.exp(-z[positive]))
    exp_z = np.exp(z[~positive])
    out[~positive] = exp_z / (1.0 + exp_z)
    return out


def accuracy(y_true: np.ndarray, y_pred: np.ndarray) -> float:
    """The fraction of predictions that match the label."""
    y_true = np.asarray(y_true).ravel()
    y_pred = np.asarray(y_pred).ravel()
    return float(np.mean(y_true == y_pred))


class LogisticRegression:
    """Binary logistic regression trained by full-batch gradient descent.

    The update is exactly Day 111's: subtract the learning rate times the
    gradient of the mean log-loss, `X.T @ (p - y) / n`, and repeat for a
    fixed number of steps. Weights start at zero, so two fits on the same
    data give bit-identical coefficients.
    """

    def __init__(self, learning_rate: float = 0.25, steps: int = 3000, l2: float = 0.0) -> None:
        self.learning_rate = learning_rate
        self.steps = steps
        self.l2 = l2
        self.weights: np.ndarray | None = None
        self.bias: float = 0.0

    def fit(self, X: np.ndarray, y: np.ndarray) -> "LogisticRegression":
        X = np.asarray(X, dtype=float)
        y = np.asarray(y, dtype=float).ravel()
        n, d = X.shape
        self.weights = np.zeros(d)
        self.bias = 0.0
        for _ in range(self.steps):
            p = sigmoid(X @ self.weights + self.bias)
            error = p - y
            grad_w = X.T @ error / n + self.l2 * self.weights
            grad_b = float(np.mean(error))
            self.weights -= self.learning_rate * grad_w
            self.bias -= self.learning_rate * grad_b
        return self

    def predict_proba(self, X: np.ndarray) -> np.ndarray:
        if self.weights is None:
            raise RuntimeError("fit before predict")
        return sigmoid(np.asarray(X, dtype=float) @ self.weights + self.bias)

    def predict(self, X: np.ndarray) -> np.ndarray:
        return (self.predict_proba(X) >= 0.5).astype(int)

    def score(self, X: np.ndarray, y: np.ndarray) -> float:
        return accuracy(y, self.predict(X))


class NearestCentroid:
    """Assign each row to the class whose mean it is closest to.

    Distance is the Euclidean norm of Day 107, which is why this model
    cares about scale: a feature measured in thousands dominates the sum
    of squares no matter how little it says about the label.
    """

    def __init__(self) -> None:
        self.classes_: np.ndarray | None = None
        self.centroids_: np.ndarray | None = None

    def fit(self, X: np.ndarray, y: np.ndarray) -> "NearestCentroid":
        X = np.asarray(X, dtype=float)
        y = np.asarray(y).ravel()
        self.classes_ = np.unique(y)
        self.centroids_ = np.vstack([X[y == c].mean(axis=0) for c in self.classes_])
        return self

    def predict(self, X: np.ndarray) -> np.ndarray:
        if self.centroids_ is None or self.classes_ is None:
            raise RuntimeError("fit before predict")
        X = np.asarray(X, dtype=float)
        distances = np.linalg.norm(X[:, None, :] - self.centroids_[None, :, :], axis=2)
        return self.classes_[np.argmin(distances, axis=1)]

    def score(self, X: np.ndarray, y: np.ndarray) -> float:
        return accuracy(y, self.predict(X))


def random_split(n: int, test_size: int, seed: int) -> tuple[np.ndarray, np.ndarray]:
    """Row indices for a seeded random train/test split."""
    rng = np.random.default_rng(seed)
    order = rng.permutation(n)
    return order[test_size:], order[:test_size]


def time_ordered_split(n: int, test_size: int) -> tuple[np.ndarray, np.ndarray]:
    """Row indices for a split that puts the LAST rows in the test set.

    The rows must already be in time order. This is the split that tells
    you what happens when the model meets a day it has never seen.
    """
    index = np.arange(n)
    return index[: n - test_size], index[n - test_size :]

"""Your exercises for Day 137 -- "Features That Do Not Cheat".

Nine exercises. Every test below currently calls `pytest.skip(...)` --
delete the skip line and replace it with real assertions. Read
`00_brief.md` for what each exercise is asking, `experiments.py` for the
measurement each one is about, and `features.py` for the encoders under
test.

Check yourself at any point:

    pytest starter -v

Never run `pytest starter examples` in one command -- both directories
hold a module named `test_features.py` and pytest collects by dotted
module name, so the two collide. Run them as two separate commands.

The reference answers live in `examples/test_features.py`. Read them
AFTER you have tried, never before.
"""

from __future__ import annotations

import math

import numpy as np
import pytest

import data
import experiments as E
import features as F
import models as M


# --- 1 ---------------------------------------------------------------------


def test_target_leakage_is_implausibly_good_and_removing_it_is_honest(leakage):
    """Exercise 1. `leakage` holds two scores for the same task.

    The table has one column, `days_to_first_invoice`, that cannot exist
    at prediction time: an invoice only happens after a visitor converts.

    Assert that:
      * `leakage["with_leak"]` is at least 0.99 -- implausibly good;
      * `leakage["without_leak"]` lands in an honest band, 0.55 to 0.80;
      * `leakage["gap_points"]` is at least 25.
    Then prove WHY it is a leak, straight from `data.signups()`: every
    unconverted row carries the sentinel -1 and every converted row
    carries a positive number.
    """
    pytest.skip("Exercise 1: assert the leaking score, the honest score and the gap.")


# --- 2 ---------------------------------------------------------------------


def test_fitting_before_the_split_costs_nothing_for_a_scaler_and_a_lot_for_an_imputer(
    scaler_contamination, imputer_contamination
):
    """Exercise 2. Two statistics fitted on all the data before splitting.

    Both fixtures average over many random splits, because with a small
    test set a single split says nothing.

    Assert that:
      * the scaler's `optimism_points` is smaller than 1 point in
        absolute value -- contaminating it buys essentially nothing;
      * the imputer's `contaminated` score beats its `correct` score,
        and `optimism_points` is at least 5;
      * the scaler's statistics really did differ: fit a `Standardiser`
        on all of `data.pricing()` and on 60 training rows, and assert
        the two `mean_[0]` values differ by more than 5.

    Ask yourself why the two differ so much. One of them applies a single
    affine map to both halves. The other does not.
    """
    pytest.skip("Exercise 2: measure the optimism from a contaminated scaler and imputer.")


# --- 3 ---------------------------------------------------------------------


def test_target_encoding_leaks_when_it_is_fitted_before_the_split(encoding):
    """Exercise 3. Replacing a category with the mean outcome uses the target.

    Assert that:
      * `naive_all_data` beats `out_of_fold`, by at least 4 points;
      * `out_of_fold` is at least as good as `naive_train_only`;
      * `naive_all_data` beats `naive_train_only` by at least 0.04.
    Then show the leak without a model in the way: build both encodings
    over the whole of `data.city_signups()` and assert the naive one
    correlates with the target above 0.25 while the out-of-fold one
    correlates below 0.10.
    """
    pytest.skip("Exercise 3: measure the naive and out-of-fold target encodings.")


# --- 4 ---------------------------------------------------------------------


def test_a_random_split_hides_what_a_time_ordered_split_shows(temporal):
    """Exercise 4. The same table, split at random and split by time.

    Assert that:
      * `random_split` is at least 0.80 and `time_ordered_split` is at
        most 0.20, a gap of at least 50 points;
      * the random training rows cover all 6 batches and the time-ordered
        ones cover 5, with exactly 1 test batch unseen;
      * the time-ordered score is BELOW the majority-class baseline for
        that period -- the model is not uninformed, it is confidently
        wrong.
    """
    pytest.skip("Exercise 4: compare the random and time-ordered splits.")


# --- 5 ---------------------------------------------------------------------


def test_cyclical_encoding_restores_the_adjacency_of_hour_23_and_hour_0(cyclical):
    """Exercise 5. Hour 23 and hour 0 are one hour apart. Prove it.

    This one is exact arithmetic, so assert exact values:
      * raw distance from 23 to 0 is 23.0, and from 3 to 4 is 1.0;
      * on the circle both are 2*sin(pi/24), to 1e-12;
      * the spread across all 24 adjacent pairs is under 1e-12;
      * hours 0 and 12 sit exactly 2.0 apart -- the circle's diameter.
    Then recompute the 24 adjacent distances yourself from
    `F.cyclical_encode(np.arange(24), 24)` rather than trusting the
    fixture, and assert they are all equal.
    """
    pytest.skip("Exercise 5: assert the raw and cyclical distances exactly.")


# --- 6 ---------------------------------------------------------------------


def test_an_ordinal_code_forces_an_order_that_one_hot_does_not(colours):
    """Exercise 6. Six paint colours, no order, one non-monotone outcome.

    Assert that:
      * `ordinal_is_monotone` is True and `one_hot_is_monotone` is False;
      * one-hot reproduces the observed per-colour rates to within 1e-5;
      * the ordinal model's largest error is above 0.30, and its accuracy
        is worse than one-hot's;
      * the observed rates really are not monotone in the code.
    """
    pytest.skip("Exercise 6: show that the ordinal model can only move monotonically.")


# --- 7 ---------------------------------------------------------------------


def test_a_ratio_separates_what_neither_column_separates(interaction):
    """Exercise 7. Spend and income overlap; spend over income does not.

    Assert the separation for all three: `ratio_only` is 1.0,
    `income_only` is under 0.60, `spend_only` is under 0.75, and the
    ratio beats the better component by at least 0.25.

    Then assert the honest footnote too: `income_and_spend` is at least
    0.95, because this particular boundary is a straight line through the
    origin and a linear model can find it from the raw columns.
    """
    pytest.skip("Exercise 7: report the separation for income, spend and the ratio.")


# --- 8 ---------------------------------------------------------------------


def test_the_vocabulary_must_be_chosen_on_training_documents_only(vocabulary):
    """Exercise 8. Which words become features is a fitted statistic too.

    Assert that:
      * `fitted_on_all_data` beats `fitted_on_train_only`, by at least
        1.5 points, averaged over 40 splits;
      * `unseen_test_words` is above zero, the matrix is exactly `top_k`
        columns wide and 100 rows tall -- unseen words are dropped, not
        crashed on.
    Then build a two-document vocabulary of your own, transform a
    document made entirely of words it has never seen, and assert you get
    a row of zeros rather than an exception.
    """
    pytest.skip("Exercise 8: measure the vocabulary leak and the unseen-word handling.")


# --- 9 ---------------------------------------------------------------------


def test_the_audit_catches_the_planted_leaks_and_leaves_honest_columns_alone(audit):
    """Exercise 9. Turn the day into a check you can run on any table.

    `E.audit_table()` is the signups table plus two extra columns: a
    categorical leak (`email_template`) and an honest one (`channel`).

    Assert that:
      * exactly `days_to_first_invoice` and `email_template` are flagged,
        by the `separable` and `pure_category` rules respectively;
      * none of `visits`, `minutes_on_site`, `discount_pct` or `channel`
        is flagged, and all four were actually checked;
      * the numeric leak's absolute correlation with the target is
        between 0.80 and 0.90 -- UNDER the default threshold, so the
        correlation rule alone would have missed it;
      * with `corr_threshold=1.01`, which disables the correlation rule
        completely, both leaks are still caught.
    """
    pytest.skip("Exercise 9: assert the audit catches the planted leaks and nothing else.")

#!/usr/bin/env bash
# Tests for the Day 137 lab. Run from the lab directory:
#   bash tests/run_tests.sh
#
# The harness proves the day's claims by running the real experiments and
# reading the numbers they produce -- never by reading source:
#
#   * a feature derived from the outcome scores 1.00, and removing it
#     drops the same model to an honest 0.64;
#   * a scaler fitted on all the data before the split buys essentially
#     nothing (measured over 200 splits), while a group-mean imputer
#     fitted the same way buys eight points;
#   * a target encoding computed before the split beats an out-of-fold
#     one, and restricting it to the training rows removes the gap;
#   * a random split scores 0.88 on time-ordered data where a
#     time-ordered split scores 0.07 -- below the majority baseline;
#   * hour 23 and hour 0 sit 23 apart as integers and 2*sin(pi/24) apart
#     on the circle, exactly, as does every other adjacent pair;
#   * an ordinal code forces monotone predictions where one-hot does not;
#   * a ratio separates classes that neither of its components separates;
#   * a vocabulary chosen on all the documents beats one chosen on the
#     training documents, and unseen test words are dropped not crashed on;
#   * the leakage audit catches both planted leaks and flags none of the
#     four honest columns;
#   * the reference suite (`examples/`) passes in full;
#   * the exercise suite (`starter/`) is all-skip on an untouched
#     checkout, and the harness proves the suite can genuinely FAIL by
#     solving every exercise in a scratch copy, breaking one assertion on
#     purpose, confirming a non-zero exit, then restoring it;
#   * nothing is left behind anywhere.
#
# Everything after the one-time install runs offline. Nothing binds a
# port, nothing needs a key. Deterministic, non-interactive, exits 0 only
# if every check passes.
set -u

export PYTHONDONTWRITEBYTECODE=1

lab_dir="$(cd "$(dirname "${BASH_SOURCE[0]}")/.." && pwd)"

find "${lab_dir}" -name '.venv' -prune -o -type d -name '__pycache__' -exec rm -rf {} + 2>/dev/null || true
find "${lab_dir}" -name '.venv' -prune -o -type d -name '.pytest_cache' -exec rm -rf {} + 2>/dev/null || true

failures=0
checks=0

check() {
  local label="$1" ok="$2"
  checks=$((checks + 1))
  if [ "${ok}" = "yes" ]; then
    echo "  ok: ${label}"
  else
    echo "  FAIL: ${label}"
    failures=$((failures + 1))
  fi
}

resolve_tool() {
  local tool="$1" override="$2"
  if [ -n "${override}" ] && [ -x "${override}" ]; then echo "${override}"; return 0; fi
  if [ -x "${lab_dir}/.venv/bin/${tool}" ]; then echo "${lab_dir}/.venv/bin/${tool}"; return 0; fi
  if command -v "${tool}" >/dev/null 2>&1; then command -v "${tool}"; return 0; fi
  return 1
}

pytest_bin="$(resolve_tool pytest "${PYTEST:-}")" || {
  echo "FAIL: pytest not found." >&2
  echo "  Install the lab's dependencies with:" >&2
  echo "    python3 -m venv .venv" >&2
  echo "    .venv/bin/pip install -r requirements/requirements.txt" >&2
  echo "  Or point this suite at an existing pytest:" >&2
  echo "    PYTEST=/path/to/pytest bash tests/run_tests.sh" >&2
  exit 1
}

python_bin="$(dirname "${pytest_bin}")/python3"
if [ ! -x "${python_bin}" ]; then
  python_bin="$(command -v python3 || true)"
fi
if [ -z "${python_bin}" ]; then
  echo "FAIL: python3 not found on PATH." >&2
  exit 1
fi

if ! "${python_bin}" -c "import pandas, numpy" >/dev/null 2>&1; then
  echo "FAIL: pandas or numpy is not importable from ${python_bin}." >&2
  echo "  Install the lab's dependencies with:" >&2
  echo "    python3 -m venv .venv" >&2
  echo "    .venv/bin/pip install -r requirements/requirements.txt" >&2
  exit 1
fi

echo "Day 137 — Features That Do Not Cheat"
echo

# --------------------------------------------------------------------------
echo "1. The tools and the versions this lab was written against"
# --------------------------------------------------------------------------

versions="$("${python_bin}" - <<'PY'
import platform
from importlib.metadata import version

print(f"python     {platform.python_version()}")
for name in ("pandas", "numpy", "pytest"):
    try:
        print(f"{name:<10} {version(name)}")
    except Exception as exc:  # pragma: no cover
        print(f"{name:<10} NOT INSTALLED ({exc})")
try:
    import sklearn  # noqa: F401

    print("sklearn    INSTALLED — this lab was written without it")
except ImportError:
    print("sklearn    not installed (expected; every model here is written out)")
PY
)"
echo "${versions}"
echo

check "scikit-learn is absent, as the lab's text states" \
  "$( echo "${versions}" | grep -q 'sklearn    not installed' && echo yes || echo no )"

mismatch=0
while IFS= read -r line; do
  [ -z "${line}" ] && continue
  pkg="${line%%==*}"
  pinned="${line#*==}"
  installed="$("${python_bin}" -c "from importlib.metadata import version; print(version('${pkg}'))" 2>/dev/null || echo "MISSING")"
  if [ "${installed}" != "${pinned}" ]; then
    mismatch=1
    echo "  version mismatch: ${pkg} pinned ${pinned}, installed ${installed}"
  fi
done < "${lab_dir}/requirements/requirements.txt"
check "installed packages match requirements.txt exactly" "$( [ ${mismatch} -eq 0 ] && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "2. The nine measurements, run for real"
# --------------------------------------------------------------------------

measured="$(cd "${lab_dir}/examples" && "${python_bin}" - <<'PY'
"""Run every experiment and print one machine-readable line per number."""
import numpy as np

import data
import experiments as E
import features as F

out = {}


def record(key, value):
    if isinstance(value, bool):
        out[key] = "yes" if value else "no"
    elif isinstance(value, float):
        out[key] = f"{value:.4f}"
    else:
        out[key] = str(value)


leak = E.target_leakage()
record("leak_with", leak["with_leak"])
record("leak_without", leak["without_leak"])
record("leak_gap_points", leak["gap_points"])

scaler = E.scaling_contamination()
record("scaler_contaminated", scaler["contaminated"])
record("scaler_correct", scaler["correct"])
record("scaler_optimism_points", scaler["optimism_points"])
record("scaler_trials", int(scaler["trials"]))

imputer = E.imputer_contamination()
record("imputer_contaminated", imputer["contaminated"])
record("imputer_correct", imputer["correct"])
record("imputer_optimism_points", imputer["optimism_points"])

encoding = E.target_encoding()
record("te_naive_all", encoding["naive_all_data"])
record("te_naive_train", encoding["naive_train_only"])
record("te_out_of_fold", encoding["out_of_fold"])
record("te_gap_points", encoding["gap_all_vs_oof_points"])

temporal = E.temporal_leakage()
record("time_random", temporal["random_split"])
record("time_ordered", temporal["time_ordered_split"])
record("time_gap_points", temporal["gap_points"])
record("time_majority_baseline", temporal["majority_rate_in_ordered_test"])
record("time_unseen_batches", int(temporal["test_batches_unseen_by_ordered_train"]))

cyc = E.cyclical_distances()
record("cyc_raw_23_0", cyc["raw_23_to_0"])
record("cyc_circle_23_0", cyc["cyclical_23_to_0"])
record("cyc_circle_3_4", cyc["cyclical_3_to_4"])
record("cyc_spread_below_1e12", cyc["cyclical_adjacent_spread"] < 1e-12)
record("cyc_expected", cyc["expected_adjacent"])

colours = E.ordinal_versus_one_hot()
record("ordinal_monotone", bool(colours["ordinal_is_monotone"]))
record("one_hot_monotone", bool(colours["one_hot_is_monotone"]))
record("ordinal_max_error", colours["ordinal_max_error"])
record("one_hot_max_error", colours["one_hot_max_error"])
record("ordinal_accuracy", colours["ordinal_accuracy"])
record("one_hot_accuracy", colours["one_hot_accuracy"])

inter = E.interaction()
record("income_only", inter["income_only"])
record("spend_only", inter["spend_only"])
record("ratio_only", inter["ratio_only"])
record("income_and_spend", inter["income_and_spend"])

vocab = E.vocabulary_contamination()
record("vocab_all_data", vocab["fitted_on_all_data"])
record("vocab_train_only", vocab["fitted_on_train_only"])
record("vocab_gap_points", vocab["gap_points"])
record("vocab_unseen_words", int(vocab["unseen_test_words"]))
record("vocab_columns", int(vocab["test_matrix_columns"]))

audit = E.audit_result()
record("audit_flagged", ",".join(audit["flagged"]))
record("audit_rules", ",".join(audit["rules"][c] for c in audit["flagged"]))
honest = ["visits", "minutes_on_site", "discount_pct", "channel"]
record("audit_honest_clean", all(c not in audit["flagged"] for c in honest))

frame = E.audit_table()
y = frame["converted"].to_numpy(dtype=float)
leak_column = frame["days_to_first_invoice"].to_numpy(dtype=float)
record("audit_leak_correlation", abs(float(np.corrcoef(leak_column, y)[0, 1])))
strict = F.leakage_audit(frame, "converted", corr_threshold=1.01)
record("audit_strict_flagged", ",".join(f.column for f in strict))

binning = E.binning_decision()
record("bin_width_top_rows", binning["equal_width_top_bin_rows"])
record("bin_count_top_rows", binning["equal_count_top_bin_rows"])
record("bin_width_top_rate", binning["equal_width_top_bin_rate"])
record("bin_count_top_rate", binning["equal_count_top_bin_rate"])

# Two determinism checks: the same generator twice, and the same
# experiment twice, must agree to the last bit.
record("data_is_deterministic", data.signups().equals(data.signups()))
record("experiment_is_deterministic", E.target_leakage() == leak)

for key, value in out.items():
    print(f"{key}={value}")
PY
)"
measured_status=$?
echo "${measured}"
echo

value_of() { echo "${measured}" | grep "^$1=" | cut -d= -f2-; }
at_least() { "${python_bin}" -c "import sys; sys.exit(0 if float('$1') >= float('$2') else 1)" && echo yes || echo no; }
below() { "${python_bin}" -c "import sys; sys.exit(0 if float('$1') < float('$2') else 1)" && echo yes || echo no; }

check "every experiment ran without error" "$( [ ${measured_status} -eq 0 ] && echo yes || echo no )"

echo
echo "  -- 1. target leakage"
check "the leaking feature scores 1.0000" "$( [ "$(value_of leak_with)" = "1.0000" ] && echo yes || echo no )"
check "removing it drops the score into the honest band 0.55-0.80" \
  "$( [ "$(at_least "$(value_of leak_without)" 0.55)" = yes ] && [ "$(below "$(value_of leak_without)" 0.80)" = yes ] && echo yes || echo no )"
check "the gap is at least 25 points" "$(at_least "$(value_of leak_gap_points)" 25)"

echo "  -- 2. a statistic fitted before the split"
check "a contaminated scaler is worth less than 1 point either way" \
  "$( "${python_bin}" -c "import sys; sys.exit(0 if abs(float('$(value_of scaler_optimism_points)')) < 1.0 else 1)" && echo yes || echo no )"
check "the scaler comparison averaged 200 splits" "$( [ "$(value_of scaler_trials)" = "200" ] && echo yes || echo no )"
check "a contaminated group-mean imputer is worth at least 5 points" "$(at_least "$(value_of imputer_optimism_points)" 5)"
check "the contaminated imputer scores above the correct one" \
  "$(at_least "$(value_of imputer_contaminated)" "$(value_of imputer_correct)")"

echo "  -- 3. target encoding"
check "the naive all-data encoding beats out-of-fold by at least 4 points" "$(at_least "$(value_of te_gap_points)" 4)"
check "out-of-fold is at least as good as naive-on-training-rows" \
  "$(at_least "$(value_of te_out_of_fold)" "$(value_of te_naive_train)")"

echo "  -- 4. temporal leakage"
check "the random split scores at least 0.80" "$(at_least "$(value_of time_random)" 0.80)"
check "the time-ordered split scores at most 0.20" "$(below "$(value_of time_ordered)" 0.20)"
check "the time-ordered score is below the majority-class baseline" \
  "$(below "$(value_of time_ordered)" "$(value_of time_majority_baseline)")"
check "exactly one test batch was unseen by the time-ordered training rows" \
  "$( [ "$(value_of time_unseen_batches)" = "1" ] && echo yes || echo no )"

echo "  -- 5. cyclical encoding"
check "raw hours put 23 and 0 exactly 23 apart" "$( [ "$(value_of cyc_raw_23_0)" = "23.0000" ] && echo yes || echo no )"
check "on the circle 23-to-0 equals 3-to-4" \
  "$( [ "$(value_of cyc_circle_23_0)" = "$(value_of cyc_circle_3_4)" ] && echo yes || echo no )"
check "and both equal 2*sin(pi/24)" \
  "$( [ "$(value_of cyc_circle_23_0)" = "$(value_of cyc_expected)" ] && echo yes || echo no )"
check "all 24 adjacent pairs agree to within 1e-12" "$( [ "$(value_of cyc_spread_below_1e12)" = yes ] && echo yes || echo no )"

echo "  -- 6. ordinal versus one-hot"
check "the ordinal model's predictions are monotone in the code" "$( [ "$(value_of ordinal_monotone)" = yes ] && echo yes || echo no )"
check "the one-hot model's predictions are not" "$( [ "$(value_of one_hot_monotone)" = no ] && echo yes || echo no )"
check "one-hot reproduces the observed rates to within 1e-5" "$(below "$(value_of one_hot_max_error)" 0.00001)"
check "the ordinal model is out by more than 0.30 somewhere" "$(at_least "$(value_of ordinal_max_error)" 0.30)"

echo "  -- 7. an interaction"
check "the ratio separates the classes perfectly" "$( [ "$(value_of ratio_only)" = "1.0000" ] && echo yes || echo no )"
check "income alone does not (under 0.60)" "$(below "$(value_of income_only)" 0.60)"
check "spend alone does not (under 0.75)" "$(below "$(value_of spend_only)" 0.75)"

echo "  -- 8. vocabulary"
check "a vocabulary chosen on all documents beats one chosen on training documents" \
  "$(at_least "$(value_of vocab_gap_points)" 1.5)"
check "test documents contain words the training vocabulary never saw" \
  "$(at_least "$(value_of vocab_unseen_words)" 1)"
check "the test matrix keeps the training vocabulary's width" \
  "$( [ "$(value_of vocab_columns)" = "30" ] && echo yes || echo no )"

echo "  -- 9. the audit"
check "both planted leaks are flagged, and only those" \
  "$( [ "$(value_of audit_flagged)" = "days_to_first_invoice,email_template" ] && echo yes || echo no )"
check "each is flagged by the expected rule" \
  "$( [ "$(value_of audit_rules)" = "separable,pure_category" ] && echo yes || echo no )"
check "no honest column is flagged" "$( [ "$(value_of audit_honest_clean)" = yes ] && echo yes || echo no )"
check "the numeric leak's correlation is UNDER the 0.90 threshold" "$(below "$(value_of audit_leak_correlation)" 0.90)"
check "disabling the correlation rule entirely still catches both leaks" \
  "$( [ "$(value_of audit_strict_flagged)" = "days_to_first_invoice,email_template" ] && echo yes || echo no )"

echo "  -- extra: a bin boundary is a decision"
check "equal-width and equal-count binning disagree about the top bin's size" \
  "$( [ "$(value_of bin_width_top_rows)" != "$(value_of bin_count_top_rows)" ] && echo yes || echo no )"
check "and about the rate you would quote for it" \
  "$( [ "$(value_of bin_width_top_rate)" != "$(value_of bin_count_top_rate)" ] && echo yes || echo no )"

echo "  -- determinism"
check "the generators return identical frames on a second call" "$( [ "$(value_of data_is_deterministic)" = yes ] && echo yes || echo no )"
check "the experiments return identical numbers on a second call" "$( [ "$(value_of experiment_is_deterministic)" = yes ] && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "3. Reference suite -- examples/ must pass in full"
# --------------------------------------------------------------------------

examples_output="$(cd "${lab_dir}" && "${pytest_bin}" examples -q 2>&1)"
examples_status=$?
echo "${examples_output}" | tail -3
check "examples/ exits 0" "$( [ ${examples_status} -eq 0 ] && echo yes || echo no )"
check "examples/ reports 9 passed, 0 failed" \
  "$( echo "${examples_output}" | grep -qE '^9 passed' && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "4. Exercise suite -- starter/ is all-skip on an untouched checkout"
# --------------------------------------------------------------------------

starter_output="$(cd "${lab_dir}" && "${pytest_bin}" starter -q 2>&1)"
starter_status=$?
echo "${starter_output}" | tail -3
check "starter/ (untouched) exits 0" "$( [ ${starter_status} -eq 0 ] && echo yes || echo no )"
check "starter/ (untouched) reports 9 skipped, 0 failed" \
  "$( echo "${starter_output}" | grep -qE '^9 skipped' && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "5. Never run 'pytest examples starter' in one invocation -- both"
echo "   directories define a module named test_features.py, and pytest"
echo "   collects by dotted module name. Documented, and run as two commands."
# --------------------------------------------------------------------------

combined_output="$(cd "${lab_dir}" && "${pytest_bin}" examples starter -q 2>&1)"
combined_status=$?
check "'pytest examples starter' aborts rather than silently passing" \
  "$( [ ${combined_status} -ne 0 ] && echo yes || echo no )"
check "the collision is reported as an import file mismatch" \
  "$( echo "${combined_output}" | grep -qi 'import file mismatch' && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "6. Prove the suite can genuinely FAIL: solve every exercise in a"
echo "   scratch copy, confirm green, break one assertion on purpose,"
echo "   confirm a non-zero exit and a printed failure, then restore."
# --------------------------------------------------------------------------

scratch_dir="$(mktemp -d "${TMPDIR:-/tmp}/d137-scratch.XXXXXX")"
cleanup_scratch() { rm -rf "${scratch_dir}"; }
trap cleanup_scratch EXIT

for module in test_features.py experiments.py features.py models.py data.py conftest.py; do
  cp "${lab_dir}/examples/${module}" "${scratch_dir}/${module}"
done

solved_output="$("${pytest_bin}" "${scratch_dir}" -q 2>&1)"
solved_status=$?
check "scratch copy of the solved suite exits 0" "$( [ ${solved_status} -eq 0 ] && echo yes || echo no )"
check "scratch copy reports 9 passed" "$( echo "${solved_output}" | grep -qE '^9 passed' && echo yes || echo no )"

# Break exercise 1's leaking-score assertion on purpose.
sed -i.bak 's/assert leakage\["with_leak"\] == 1.0/assert leakage["with_leak"] == 0.5/' \
  "${scratch_dir}/test_features.py"

broken_output="$("${pytest_bin}" "${scratch_dir}" -q 2>&1)"
broken_status=$?
check "broken scratch copy exits non-zero" "$( [ ${broken_status} -ne 0 ] && echo yes || echo no )"
check "broken scratch copy prints a failure" \
  "$( echo "${broken_output}" | grep -qiE 'failed|assert' && echo yes || echo no )"

mv "${scratch_dir}/test_features.py.bak" "${scratch_dir}/test_features.py"
restored_output="$("${pytest_bin}" "${scratch_dir}" -q 2>&1)"
restored_status=$?
check "restored scratch copy exits 0 again" "$( [ ${restored_status} -eq 0 ] && echo yes || echo no )"
check "restored scratch copy reports 9 passed again" \
  "$( echo "${restored_output}" | grep -qE '^9 passed' && echo yes || echo no )"

cleanup_scratch
trap - EXIT
echo

# --------------------------------------------------------------------------
echo "7. Offline, and nothing left behind"
# --------------------------------------------------------------------------

url_hits="$(grep -rEl 'https?://|ftp://' "${lab_dir}/examples" "${lab_dir}/starter" 2>/dev/null || true)"
check "no URLs inside examples/ or starter/" "$( [ -z "${url_hits}" ] && echo yes || echo no )"

network_hits="$(grep -rElE '\b(requests|urllib|socket|httpx)\b' "${lab_dir}/examples" "${lab_dir}/starter" 2>/dev/null || true)"
check "no networking module is imported anywhere in the lab" "$( [ -z "${network_hits}" ] && echo yes || echo no )"

find "${lab_dir}" -name '.venv' -prune -o -type d -name '__pycache__' -exec rm -rf {} + 2>/dev/null || true
find "${lab_dir}" -name '.venv' -prune -o -type d -name '.pytest_cache' -exec rm -rf {} + 2>/dev/null || true

stray="$(find "${lab_dir}" -name '.venv' -prune -o \( -type d -name '__pycache__' -print -o -type d -name '.pytest_cache' -print \) 2>/dev/null || true)"
check "no __pycache__ or .pytest_cache left behind" "$( [ -z "${stray}" ] && echo yes || echo no )"

leftover_tmp="$(find "${TMPDIR:-/tmp}" -maxdepth 1 -name 'd137-*' -print 2>/dev/null || true)"
check "no d137 temporary directory left in the system temp directory" "$( [ -z "${leftover_tmp}" ] && echo yes || echo no )"
echo

echo "-------------------------------------------------------------"
echo "${checks} checks, ${failures} failure(s)"
if [ "${failures}" -gt 0 ]; then
  exit 1
fi
exit 0