"""Does retraining actually recover performance under drift? Measure it. The rest of this sub-project establishes that SECOM drifts hard — the weekly fail rate falls from ~23% to 1-5%, and 55-66% of sensors clear the PSI action threshold in every test segment — and then *asserts* that a monitor-and-retrain loop matters more than model choice. That assertion is untested: the main pipeline fits one model on one split and stops. This script tests it. It walks forward through the wafer stream the way a deployed system would, and compares four retraining policies at identical evaluation windows: frozen train once at the start, never again — the usual portfolio setup expanding retrain on all history before each window sliding retrain on the most recent WINDOW_TRAIN wafers only triggered retrain on all history, but only when the PSI monitor fires — the realistic production policy, and what common/drift.py is for `triggered` deliberately uses the expanding window rather than the sliding one. An earlier version had it retrain on recent wafers, which made it identical to `sliding` by construction whenever the monitor fired — and since the monitor fires in 7 of 8 windows, it tested nothing of its own. Each retrain also records how many failures were in its training block, because with a rare positive class that count, not recency, turns out to be what decides the outcome. Nothing at window w may touch data from window w or later — that includes column selection, imputation medians, the scaler, the model, and the decision threshold. Each retrain refits all of them. Run: python 01_secom/rolling_retrain.py """ from __future__ import annotations import sys from pathlib import Path import matplotlib matplotlib.use("Agg") import matplotlib.pyplot as plt import numpy as np import pandas as pd from xgboost import XGBClassifier HERE = Path(__file__).parent sys.path.insert(0, str(HERE.parent)) sys.path.insert(0, str(HERE)) from common import drift, metrics # noqa: E402 from run_pipeline import (COST_FN, COST_FP, load, preprocess_apply, # noqa: E402 preprocess_fit) FIGS = HERE / "figures" OUT = HERE / "outputs" INITIAL_TRAIN = 600 # wafers used before the simulation starts N_WINDOWS = 8 # evaluation windows carved from the remainder WINDOW_TRAIN = 600 # sliding-window memory, in wafers VAL_FRAC = 0.30 # tail of each training block reserved for the threshold SEEDS = [0, 1, 2, 3, 4] # the spread across seeds is the noise floor PSI_TRIGGER_FRAC = 0.10 # retrain when >10% of sensors reach the action band def fit_model(train: pd.DataFrame, sensor_cols: list[str], seed: int): """Refit everything — preprocessing, model, threshold — on `train` only. Returns a callable scoring new wafers, plus the chosen threshold. The threshold comes from the tail of the training block, standing in for "the most recent labelled data I had when I deployed this". """ n_val = max(40, int(VAL_FRAC * len(train))) fit_part, val_part = train.iloc[:-n_val], train.iloc[-n_val:] prep = preprocess_fit(fit_part, sensor_cols) X_fit = preprocess_apply(fit_part, prep) X_val = preprocess_apply(val_part, prep) y_fit = fit_part["y"].to_numpy() y_val = val_part["y"].to_numpy() spw = (y_fit == 0).sum() / max(1, (y_fit == 1).sum()) model = XGBClassifier( n_estimators=300, max_depth=3, learning_rate=0.05, subsample=0.8, colsample_bytree=0.5, scale_pos_weight=spw, eval_metric="aucpr", reg_lambda=2.0, n_jobs=-1, random_state=seed) model.fit(X_fit, y_fit) if y_val.sum() > 0: thr, _ = metrics.best_threshold( y_val, model.predict_proba(X_val)[:, 1], COST_FN, COST_FP) else: thr = 0.5 # refit on the full block so the deployed model sees the freshest wafers prep_full = preprocess_fit(train, sensor_cols) X_all = preprocess_apply(train, prep_full) model.fit(X_all, train["y"].to_numpy()) def score(df: pd.DataFrame) -> np.ndarray: return model.predict_proba(preprocess_apply(df, prep_full))[:, 1] info = {"n_train": len(train), "n_fail_train": int(train["y"].sum()), "train_prevalence": float(train["y"].mean())} return score, float(thr), prep_full, info def simulate(df: pd.DataFrame, sensor_cols: list[str], seed: int) -> pd.DataFrame: """One walk-forward pass. Returns per-window, per-strategy results.""" bounds = np.linspace(INITIAL_TRAIN, len(df), N_WINDOWS + 1).astype(int) strategies = ["frozen", "expanding", "sliding", "triggered"] # frozen: fitted once, before the simulation starts, and never touched again frozen_score, frozen_thr, frozen_prep, frozen_info = fit_model( df.iloc[:INITIAL_TRAIN], sensor_cols, seed) # triggered: carries its own state — reference distribution and last model trig_score, trig_thr, trig_prep, trig_info = ( frozen_score, frozen_thr, frozen_prep, frozen_info) trig_ref = df.iloc[:INITIAL_TRAIN] trig_retrains = 0 rows = [] for w in range(N_WINDOWS): lo, hi = bounds[w], bounds[w + 1] window = df.iloc[lo:hi] y_true = window["y"].to_numpy() if y_true.sum() == 0: # no failures to detect; metric undefined continue # --- decide whether the triggered policy fires, using only past data recent = df.iloc[max(0, lo - WINDOW_TRAIN // 2):lo] psi_df = drift.psi_report(trig_ref[trig_prep["cols"]], recent[trig_prep["cols"]]) fired = drift.drift_trigger(psi_df, PSI_TRIGGER_FRAC)["retrain"] if fired: # all history, not just recent — see the module docstring trig_score, trig_thr, trig_prep, trig_info = fit_model( df.iloc[:lo], sensor_cols, seed) trig_ref = df.iloc[max(0, lo - WINDOW_TRAIN):lo] trig_retrains += 1 for strat in strategies: if strat == "frozen": score, thr, info, n_fit = frozen_score, frozen_thr, frozen_info, 0 elif strat == "expanding": score, thr, _, info = fit_model(df.iloc[:lo], sensor_cols, seed) n_fit = 1 elif strat == "sliding": score, thr, _, info = fit_model( df.iloc[max(0, lo - WINDOW_TRAIN):lo], sensor_cols, seed) n_fit = 1 else: score, thr, info, n_fit = trig_score, trig_thr, trig_info, int(fired) s = score(window) pred = (s >= thr).astype(int) rows.append({ "seed": seed, "window": w, "strategy": strat, "n": len(window), "n_fail": int(y_true.sum()), "prevalence": float(y_true.mean()), "start": window["timestamp"].min(), "end": window["timestamp"].max(), "pr_auc": metrics.pr_auc(y_true, s), "capture@20%": metrics.decile_capture(y_true, s, 0.2), "expected_cost": metrics.expected_cost(y_true, pred, COST_FN, COST_FP), "retrained": n_fit, "psi_fired": bool(fired), "threshold": thr, **info, "_score": s, "_y": y_true, }) for r in rows: if r["strategy"] == "triggered": r["total_retrains"] = trig_retrains return pd.DataFrame(rows) def pooled(res: pd.DataFrame) -> pd.DataFrame: """Pool every window's predictions per strategy per seed. Per-window PR-AUC on ~120 wafers with ~7 failures is far too noisy to read directly; pooling the whole deployment period is the fair comparison and answers the operational question — over these three months, which policy caught more failures? """ out = [] for (seed, strat), g in res.groupby(["seed", "strategy"]): y = np.concatenate(g["_y"].to_list()) s = np.concatenate(g["_score"].to_list()) out.append({ "seed": seed, "strategy": strat, "pr_auc": metrics.pr_auc(y, s), "capture@20%": metrics.decile_capture(y, s, 0.2), "capture@5%": metrics.decile_capture(y, s, 0.05), "expected_cost": float(g["expected_cost"].mean()), "retrains": int(g["retrained"].sum()), # the mechanism: how many failures did this policy get to learn from? "mean_fails_in_train": float(g["n_fail_train"].mean()), "mean_train_size": float(g["n_train"].mean()), }) return pd.DataFrame(out) def paired_test(res: pd.DataFrame) -> pd.DataFrame: """Is any retraining policy actually better than frozen, window by window? The pooled score is fragile here: a single 121-wafer window can hold 5 failures and swing pooled PR-AUC on its own. The honest comparison is paired — same window, same seed, different policy — and the independent unit is the *window*, not the seed, since seeds re-score identical wafers. So seeds are averaged first, leaving 8 paired differences and, with them, very little statistical power. That is the real state of the evidence. """ from scipy import stats per_window = (res.groupby(["strategy", "window"])["pr_auc"].mean().unstack(0)) rows = [] for strat in per_window.columns: if strat == "frozen": continue diff = (per_window[strat] - per_window["frozen"]).dropna() diff = diff[diff.abs() > 1e-12] # drop exact ties (window 0) if len(diff) < 3: continue stat, p = stats.wilcoxon(diff) rows.append({ "strategy": strat, "n_windows": len(diff), "wins_vs_frozen": int((diff > 0).sum()), "losses_vs_frozen": int((diff < 0).sum()), "median_diff": float(diff.median()), "mean_diff": float(diff.mean()), "wilcoxon_p": float(p), }) return pd.DataFrame(rows) def main(): FIGS.mkdir(exist_ok=True) OUT.mkdir(exist_ok=True) df = load() sensor_cols = [c for c in df.columns if c.startswith("s")] print(f"{len(df)} wafers, overall fail rate {df['y'].mean():.3%}") print(f"initial training block {INITIAL_TRAIN}, " f"{N_WINDOWS} evaluation windows over the remaining {len(df) - INITIAL_TRAIN}") frames = [] for seed in SEEDS: print(f" seed {seed} ...", flush=True) frames.append(simulate(df, sensor_cols, seed)) res = pd.concat(frames, ignore_index=True) pool = pooled(res) summary = (pool.groupby("strategy") .agg(pr_auc=("pr_auc", "mean"), pr_auc_sd=("pr_auc", "std"), capture20=("capture@20%", "mean"), capture20_sd=("capture@20%", "std"), capture5=("capture@5%", "mean"), cost=("expected_cost", "mean"), retrains=("retrains", "mean"), fails_in_train=("mean_fails_in_train", "mean"), train_size=("mean_train_size", "mean")) .reindex(["frozen", "triggered", "sliding", "expanding"])) base = summary.loc["frozen", "pr_auc"] summary["vs_frozen_%"] = (summary["pr_auc"] / base - 1) * 100 res.drop(columns=["_score", "_y"]).to_csv(OUT / "rolling_retrain_windows.csv", index=False) summary.to_csv(OUT / "rolling_retrain_summary.csv") prevalence = res.groupby("window")["prevalence"].first() print("\n=== pooled over the whole deployment period " f"(mean of {len(SEEDS)} seeds) ===") print(summary.round(4).to_string()) print(f"\nnoise floor for reference: seed-to-seed sd of PR-AUC is " f"{summary['pr_auc_sd'].mean():.4f} " f"({summary['pr_auc_sd'].mean() / base * 100:.0f}% of the frozen baseline)") print("\n=== window prevalence (the drift the policies face) ===") print(prevalence.round(4).to_string()) # how concentrated is the pooled difference? one window can carry it pw = res.groupby(["strategy", "window"])["pr_auc"].mean().unstack(0) gap = (pw["frozen"] - pw["sliding"]) top = gap.abs().idxmax() print(f"\nconcentration check: window {top} alone accounts for " f"{gap[top] / gap.abs().sum() * 100:.0f}% of the frozen-vs-sliding gap, " f"and holds only {int(res[res.window == top]['n_fail'].iloc[0])} failures") paired = paired_test(res) paired.to_csv(OUT / "rolling_retrain_paired_test.csv", index=False) print("\n=== paired per-window test against frozen " "(seeds averaged; independent unit is the window) ===") print(paired.round(4).to_string(index=False)) make_figure(res, summary, paired, prevalence) def make_figure(res: pd.DataFrame, summary: pd.DataFrame, paired: pd.DataFrame, prevalence: pd.Series) -> None: """Draw from tables only, so `--plot-only` can redraw without refitting.""" colors = {"frozen": "#c0392b", "triggered": "#1f5c8b", "sliding": "#1e7d54", "expanding": "#d98c00"} base = summary.loc["frozen", "pr_auc"] fig, axes = plt.subplots(1, 3, figsize=(19, 5.4)) ax = axes[0] per_win = res.groupby(["strategy", "window"])["pr_auc"].mean().unstack(0) for strat in ["frozen", "triggered", "sliding", "expanding"]: if strat in per_win: ax.plot(per_win.index, per_win[strat], marker="o", ms=4, label=strat, color=colors[strat], lw=2.2 if strat in ("frozen", "triggered") else 1.4) ax2 = ax.twinx() ax2.bar(prevalence.index, prevalence.to_numpy() * 100, alpha=0.13, color="gray", width=0.65) ax2.set_ylabel("window fail rate (%) — grey bars", color="gray", fontsize=9) ax2.tick_params(axis="y", labelcolor="gray", labelsize=8) ax.set(xlabel="deployment window (time →)", ylabel="PR-AUC in that window") ax.set_title("Per-window performance as the process drifts\n" "~5 failures per window; window 2 carries the pooled gap", fontsize=11) ax.legend(fontsize=8, loc="upper right") ax = axes[1] order = ["frozen", "triggered", "sliding", "expanding"] x = np.arange(len(order)) vals = summary.loc[order, "pr_auc"].to_numpy() errs = summary.loc[order, "pr_auc_sd"].to_numpy() bars = ax.bar(x, vals, yerr=errs, capsize=5, color=[colors[s] for s in order], edgecolor="black") for xi, v, r in zip(x, vals, summary.loc[order, "retrains"]): ax.text(xi, v + 0.002, f"{v:.3f}\n({r:.0f} retrains)", ha="center", va="bottom", fontsize=8.5) ax.axhline(base, color="#c0392b", ls="--", lw=1, alpha=0.6) ax.set_xticks(x, order) # the pooled bars invite a conclusion the paired test does not support, # so the p-values go on the chart itself pvals = paired.set_index("strategy")["wilcoxon_p"].to_dict() for xi, strat in zip(x, order): if strat in pvals: ax.text(xi, 0.004, f"p={pvals[strat]:.2f}", ha="center", va="bottom", fontsize=8.5, color="white", fontweight="bold") ax.set_ylabel("PR-AUC pooled over the whole deployment period") ax.set_title("No policy beats frozen at this sample size\n" f"bars = sd over {len(SEEDS)} seeds; p = paired Wilcoxon", fontsize=11) ax.set_ylim(0, max(vals + errs) * 1.30) # --- the mechanism: how many failures did each policy get to learn from? ax = axes[2] fails = summary.loc[order, "fails_in_train"].to_numpy() b = ax.bar(x, fails, color=[colors[s] for s in order], edgecolor="black") for xi, v, sz in zip(x, fails, summary.loc[order, "train_size"]): ax.text(xi, v + 0.6, f"{v:.0f} fails\nof {sz:.0f} wafers", ha="center", va="bottom", fontsize=8.5) ax.set_xticks(x, order) ax.set_ylabel("failures in the training block (mean over windows)") ax.set_title("Untested hypothesis for the trend: sliding\n" "discards the failures the model needs", fontsize=11) ax.set_ylim(0, fails.max() * 1.32) fig.tight_layout() fig.savefig(FIGS / "rolling_retrain.png", dpi=150) plt.close(fig) print(f"\nwrote {FIGS / 'rolling_retrain.png'}") def replot() -> None: """Redraw the figure from saved tables — the simulation takes ~12 minutes, the figure should not have to.""" res = pd.read_csv(OUT / "rolling_retrain_windows.csv") summary = pd.read_csv(OUT / "rolling_retrain_summary.csv", index_col=0) paired = pd.read_csv(OUT / "rolling_retrain_paired_test.csv") prevalence = res.groupby("window")["prevalence"].first() make_figure(res, summary, paired, prevalence) if __name__ == "__main__": if "--plot-only" in sys.argv: replot() else: main()