diff --git a/.gitignore b/.gitignore
index cfb074d..fb80953 100644
--- a/.gitignore
+++ b/.gitignore
@@ -28,3 +28,9 @@ examples/gusto/out/
 validation/predrupdate/gusto_y_p.csv
 validation/predrupdate/synpm_estimates.csv
 validation/predrupdate/gusto_estimates.csv
+
+# Same reasoning for the IMPROVE neoepitope fixture: build_fixture.py regenerates both parquets
+# deterministically from the released IMPROVE_paper results archive.
+examples/improve/cohort.parquet
+examples/improve/predictions.parquet
+examples/improve/out/
diff --git a/examples/improve/README.md b/examples/improve/README.md
new file mode 100644
index 0000000..b21f5b8
--- /dev/null
+++ b/examples/improve/README.md
@@ -0,0 +1,120 @@
+# IMPROVE neoepitope immunogenicity — predval example cohort
+
+This directory applies predval to a problem outside the clinical-risk-score setting it was
+built in: deciding which candidate neoepitopes from a patient's tumour are worth putting into a
+personalised vaccine or a TCR-T product.
+
+## Why this cohort
+
+Every published neoepitope immunogenicity predictor reports **ranking** metrics — AUROC, AUPRC,
+top-N recall, fraction-ranked. A literature check (2026-10-04) found no published work reporting
+**calibration** for peptide-level immunogenicity predictors: whether a peptide scored 0.8
+validates 80% of the time. The nearest work is adjacent rather than overlapping and is cited
+below.
+
+Calibration is not a cosmetic concern here, because selection is budget-constrained. A
+personalised vaccine encodes on the order of 20 neoepitopes; a TCR-T product carries 3-5. Those
+slots are the scarce resource, and a probability that over-promises by a large factor misprices
+every one of them.
+
+## Provenance
+
+Source: Borch et al. 2024, *IMPROVE: a feature model to predict neoepitope immunogenicity
+through broad-scale validation of T cell recognition*, Frontiers in Immunology. Data release:
+[github.com/SRHgroup/IMPROVE_paper](https://github.com/SRHgroup/IMPROVE_paper), `results.zip`.
+
+- 17,520 candidate neoepitopes, 70 patients, three accrual cohorts (bladder 6,237 peptides / 24
+  patients; melanoma 5,921 / 26; a mixed "Basket" group 5,362 / 20).
+- 467 peptides are experimentally immunogenic — **2.666% prevalence**.
+- Predictions are the paper's own **out-of-fold** scores from its 5-fold cross-validation. That
+  CV was partitioned **by patient**: every patient sits in exactly one partition, verified in
+  `build_fixture.py`. No score here comes from a model that had seen another peptide from the
+  same patient, which is what makes the file admissible as a predictions-only validation input
+  rather than an apparent-fit artefact.
+- Seven models: the IMPROVE random forest in its three released feature configurations
+  (`Simple`; `TME_excluded`, adding cellular prevalence and a prioritisation score;
+  `TME_included`, further adding cytolytic activity, HLA expression and MCP-counter estimates),
+  each in full and in the `wo_prime` ablation that drops the PRIME feature, plus the paper's
+  `nnalign` comparator.
+
+Nothing in this directory re-trains, re-tunes or re-ranks anything. predval reads probabilities
+and labels and nothing else.
+
+## Unit of analysis
+
+The cluster is the **patient**, not the peptide. Peptides from one patient share that patient's
+HLA genotype, mutational process, tumour transcriptome and microenvironment — the dominant
+sources of variance in every feature these models consume. Mean cluster size is 250 peptides,
+an order of magnitude larger than the PCam fixture's slides, so the cost of getting this wrong
+is correspondingly larger. Every interval in the report resamples whole patients.
+
+## Prespecified subgroups
+
+Both were declared before any metric was computed on this cohort.
+
+- `tumour_cohort` — the release's own three accrual groups. A model that holds in melanoma, the
+  indication these predictors were overwhelmingly trained and tuned on, and fails elsewhere is
+  precisely the failure mode that matters for carrying a verdict to another disease. This is the
+  closest available proxy for that transfer question, and its result bounds how far anything
+  here may be carried. MSS colorectal cancer is not in this cohort at all.
+- `peptide_length` — 9-mers versus everything else. The binding predictors that supply the
+  dominant features are trained overwhelmingly on 9-mers and are known to degrade off that
+  length, so the off-length stratum is where a calibration fault would be expected *a priori*.
+
+## Reproducing
+
+```sh
+# 1. fetch the release and extract results.zip
+git clone https://github.com/SRHgroup/IMPROVE_paper.git
+unzip -q IMPROVE_paper/results.zip -x '__MACOSX/*' -d improve_results
+
+# 2. build the fixture (writes cohort.parquet and predictions.parquet here)
+uv run python examples/improve/build_fixture.py --src 'improve_results/results copy'
+
+# 3. the main evaluation -- B=2000 patient-clustered bootstrap, several minutes
+uv run python examples/improve/run_evaluation.py
+
+# 4. the budgeted-yield sidecar
+uv run python examples/improve/budget_yield.py
+```
+
+`run_evaluation.py` writes `out/report.html`, `out/findings.json` and the parquet artefacts.
+`budget_yield.py` writes `out/budget_yield.parquet`.
+
+## The budgeted-yield sidecar
+
+`budget_yield.py` is deliberately *not* part of predval's contract. predval reports calibration;
+the mapping from calibration to a decision is cohort-specific, and this is that mapping for this
+cohort. For each model and each budget k ∈ {3, 5, 10, 20} it takes every patient's own top-k
+candidates by score and reports, per patient and averaged with a patient-level bootstrap:
+
+- `observed` — immunogenic peptides actually in that top-k;
+- `promised` — the sum of the model's published probabilities over the *same* top-k set, which is
+  what those probabilities assert the yield will be;
+- `recalibrated` — the same sum after predval's cross-fitted rung1 (intercept-only) correction,
+  fitted on other patients.
+
+Because all three columns score the identical top-k set, `promised − observed` isolates the
+miscalibration cost in units of antigens. It is not a ranking error.
+
+## Prior art this sits next to
+
+- **ImmUQBench** (Qayyum et al. 2026, *Oxford Open Immunology* iqag003, PMID 41859687) reports
+  ECE, Brier and NLL for uncertainty-quantification methods on **whole-protein** antigen
+  immunogenicity. Complementary scope: no peptide-HLA level, no patient clustering, no budget.
+- **arXiv 2604.13254** — calibrated conformal abstention for TCR-pMHC **binding**.
+- **CaliPPer** (arXiv 2606.07258) — distance-aware Platt recalibration, largely for binding and
+  presentation rather than immunogenicity.
+- **TESLA** (Wells et al. 2020, *Cell*) and the shortcut-bias analysis in Zhang et al. 2026
+  (*Cell Genomics*, PMC13347938) report ranking metrics only.
+
+## Limits worth stating plainly
+
+- Three tumour types, none of them colorectal. Nothing here licenses a claim about MSS CRC.
+- 70 patients is a small number of independent units for interval estimation, which is exactly
+  why the patient-clustered intervals are reported and the naive per-peptide ones are shown
+  beside them labelled incorrect.
+- The experimental labels are themselves a screen: "not immunogenic" means "not detected under
+  the assays applied to that patient", which is not the same as "not immunogenic". A calibration
+  result is a statement about agreement with these labels.
+- A calibration test can fail to reject "calibrated". It can never establish it.
diff --git a/examples/improve/budget_yield.py b/examples/improve/budget_yield.py
new file mode 100644
index 0000000..7348aba
--- /dev/null
+++ b/examples/improve/budget_yield.py
@@ -0,0 +1,167 @@
+#!/usr/bin/env python
+"""Budgeted-yield sidecar for the IMPROVE fixture: what a slot budget actually buys.
+
+    uv run python examples/improve/budget_yield.py [--out DIR] [--n-boot N]
+
+Calibration is a property of probabilities; antigen selection is a property of a *budget*. A
+personalised vaccine encodes on the order of 20 neoepitopes and a TCR-T product carries 3-5, so
+the operative question is not "is peptide X immunogenic" but "of the k slots I am going to fill
+for this patient, how many carry a real target, and how many did the score promise?"
+
+For each model and each budget k this script takes every patient's own top-k candidates by
+score and reports three numbers per patient, averaged over patients with a patient-level
+bootstrap interval:
+
+  observed   -- immunogenic peptides actually in the top-k (the truth)
+  promised   -- the sum of the model's own published probabilities over that same top-k, which
+                is what those probabilities assert the yield will be
+  recalibrated -- the same sum after a cross-fitted intercept correction (predval's rung1),
+                fitted on other patients and applied here
+
+`promised` minus `observed` is the slot budget's miscalibration cost, in units of antigens. It
+is not a ranking error: the same top-k set is used for all three columns.
+
+This is a sidecar, not part of predval's contract: predval reports calibration, and the mapping
+from calibration to a decision is cohort-specific. It writes budget_yield.parquet and prints a
+summary. No network calls.
+"""
+
+from __future__ import annotations
+
+import argparse
+from pathlib import Path
+
+import numpy as np
+import pandas as pd
+
+from predval import load_cohort, load_predictions
+
+# predval's own cross-fitting, reused rather than reimplemented so the correction applied here is
+# byte-for-byte the one the main report's rung1 column is scored on. Private because it is not a
+# contract surface; this example is allowed to reach for it, library code is not.
+from predval.recalibrate import _crossfit_predictions
+
+HERE = Path(__file__).resolve().parent
+BUDGETS = (3, 5, 10, 20)
+W = 26
+
+
+def _per_patient_topk(
+    frame: pd.DataFrame, k: int, score: str, extra: tuple[str, ...]
+) -> pd.DataFrame:
+    """Each patient's own top-k rows by `score`, ties broken deterministically by subject_id."""
+    ordered = frame.sort_values([score, "subject_id"], ascending=[False, True])
+    top = ordered.groupby("patient", sort=True).head(k)
+    cols = ["patient", "response", score, *extra]
+    return top[cols]
+
+
+def _cluster_bootstrap(
+    per_patient: pd.DataFrame, columns: list[str], n_boot: int, seed: int, ci: float
+) -> dict[str, tuple[float, float]]:
+    """Resample whole patients with replacement; return a percentile interval per column."""
+    rng = np.random.default_rng(seed)
+    values = per_patient[columns].to_numpy(dtype=float)
+    n = len(per_patient)
+    draws = np.empty((n_boot, len(columns)), dtype=float)
+    for b in range(n_boot):
+        idx = rng.integers(0, n, size=n)
+        draws[b] = values[idx].mean(axis=0)
+    lo = (1.0 - ci) / 2.0 * 100.0
+    return {
+        col: (float(np.percentile(draws[:, j], lo)), float(np.percentile(draws[:, j], 100.0 - lo)))
+        for j, col in enumerate(columns)
+    }
+
+
+def main(argv: list[str] | None = None) -> int:
+    ap = argparse.ArgumentParser(description=__doc__)
+    ap.add_argument("--out", type=Path, default=HERE / "out")
+    ap.add_argument("--n-boot", type=int, default=None)
+    args = ap.parse_args(argv)
+
+    cohort = load_cohort(HERE / "cohort.yaml")
+    predictions = load_predictions(HERE / "predictions.parquet")
+    spec = cohort.spec
+    n_boot = args.n_boot if args.n_boot is not None else spec.uncertainty.n_boot
+    seed = spec.uncertainty.seed
+    ci = spec.uncertainty.ci_level
+
+    table = cohort.table[[spec.subject_key, "patient", spec.outcome.field]].rename(
+        columns={spec.subject_key: "subject_id", spec.outcome.field: "response"}
+    )
+    prevalence = float(table["response"].mean())
+
+    rows = []
+    for model_id, grp in predictions.frame.groupby("model_id"):
+        merged = table.merge(
+            grp[["subject_id", "predicted"]], on="subject_id", how="inner", validate="1:1"
+        )
+        y = merged["response"].to_numpy(dtype=float)
+        p = merged["predicted"].to_numpy(dtype=float)
+        groups = merged["patient"].to_numpy()
+        k_folds = min(5, np.unique(groups).size)
+        cf, reason = _crossfit_predictions("rung1", y, p, groups, k_folds)
+        if cf is None:
+            raise SystemExit(f"{model_id}: rung1 cross-fit unavailable ({reason})")
+        merged["recalibrated"] = cf
+
+        for k in BUDGETS:
+            top = _per_patient_topk(merged, k, "predicted", ("recalibrated",))
+            per_patient = top.groupby("patient", sort=True).agg(
+                observed=("response", "sum"),
+                promised=("predicted", "sum"),
+                recalibrated=("recalibrated", "sum"),
+                slots=("response", "size"),
+            )
+            cols = ["observed", "promised", "recalibrated"]
+            intervals = _cluster_bootstrap(per_patient, cols, n_boot, seed, ci)
+            means = {c: float(per_patient[c].mean()) for c in cols}
+            rows.append(
+                {
+                    "model_id": model_id,
+                    "budget_k": k,
+                    "n_patients": int(len(per_patient)),
+                    "slots_median": float(per_patient["slots"].median()),
+                    **means,
+                    **{f"{c}_ci_low": intervals[c][0] for c in cols},
+                    **{f"{c}_ci_high": intervals[c][1] for c in cols},
+                    "patients_with_zero_hits": int((per_patient["observed"] == 0).sum()),
+                    "random_baseline": k * prevalence,
+                    "n_boot": n_boot,
+                    "seed": seed,
+                    "ci_level": ci,
+                }
+            )
+
+    out = pd.DataFrame(rows)
+    args.out.mkdir(parents=True, exist_ok=True)
+    path = args.out / "budget_yield.parquet"
+    out.to_parquet(path, index=False)
+
+    print(f"budget yield -> {path}")
+    print(
+        f"cohort prevalence {prevalence:.4%}; {out['n_patients'].iloc[0]} patients; "
+        f"B={n_boot} patient-level resamples, seed {seed}\n"
+    )
+    for k in BUDGETS:
+        sub = out[out["budget_k"] == k].sort_values("observed", ascending=False)
+        print(f"budget k={k} slots per patient  (random pick would yield {k * prevalence:.3f})")
+        print(
+            f"  {'model':{W}} {'observed':>24} {'promised':>10} {'recalib':>9} "
+            f"{'over by':>8} {'0-hit pts':>10}"
+        )
+        for _, r in sub.iterrows():
+            print(
+                f"  {r['model_id']:{W}} {r['observed']:8.3f} "
+                f"[{r['observed_ci_low']:.3f}, {r['observed_ci_high']:.3f}]"
+                f" {r['promised']:10.3f} {r['recalibrated']:9.3f}"
+                f" {r['promised'] / max(r['observed'], 1e-9):7.1f}x"
+                f" {r['patients_with_zero_hits']:6d}/{r['n_patients']}"
+            )
+        print()
+    return 0
+
+
+if __name__ == "__main__":
+    raise SystemExit(main())
diff --git a/examples/improve/build_fixture.py b/examples/improve/build_fixture.py
new file mode 100644
index 0000000..e3f2f1b
--- /dev/null
+++ b/examples/improve/build_fixture.py
@@ -0,0 +1,174 @@
+#!/usr/bin/env python
+"""Build the IMPROVE neoantigen-immunogenicity fixture for predval.
+
+Source: github.com/SRHgroup/IMPROVE_paper (Borch et al. 2024, Front Immunol), results.zip,
+the 5-fold cross-validation out-of-fold prediction tables plus the NNAlign comparator.
+
+    uv run python examples/improve/build_fixture.py --src <extracted results dir>
+
+Writes cohort.parquet and predictions.parquet next to cohort.yaml. No network calls.
+"""
+
+from __future__ import annotations
+
+import argparse
+from pathlib import Path
+
+import pandas as pd
+
+HERE = Path(__file__).resolve().parent
+
+# model_id -> (subdirectory, filename) under the extracted `results copy/5_fold_CV/` tree.
+RF_MODELS = {
+    "improve_simple": ("Simple", "pred_df_Simple.txt"),
+    "improve_simple_wo_prime": ("Simple", "pred_df_wo_primeSimple.txt"),
+    "improve_tme_excl": ("TME_excluded", "pred_df_TME_excluded.txt"),
+    "improve_tme_excl_wo_prime": ("TME_excluded", "pred_df_wo_primeTME_excluded.txt"),
+    "improve_tme_incl": ("TME_included", "pred_df_TME_included.txt"),
+    "improve_tme_incl_wo_prime": ("TME_included", "pred_df_wo_primeTME_included.txt"),
+}
+
+
+def _sid(frame: pd.DataFrame) -> pd.Series:
+    return (
+        frame["Patient"].astype(str)
+        + "|"
+        + frame["Mut_peptide"].astype(str)
+        + "|"
+        + frame["HLA_allele"].astype(str)
+    )
+
+
+def load_rf(src: Path) -> dict[str, pd.DataFrame]:
+    out = {}
+    for model_id, (sub, name) in RF_MODELS.items():
+        path = src / "5_fold_CV" / sub / name
+        frame = pd.read_csv(path, sep=r"\s+")
+        frame["subject_id"] = _sid(frame)
+        out[model_id] = frame
+    return out
+
+
+def load_nnalign(src: Path) -> pd.DataFrame:
+    """Parse the NNAlign comparator table.
+
+    The released file's header line crams its first five field names into one tab-delimited cell
+    ("Peptide Target Split HLA donor") while the data rows are space-delimited there, so pandas
+    reads peptide/target/split/HLA as a four-level index and shifts every named column one to
+    the left: `donor` lands in the first named column and `Measure`/`Prediction` are the target
+    and the score. Read the identifiers off the index levels rather than trusting the header,
+    and drop the four trailing malformed rows (they carry column *names* as values).
+    """
+    path = src / "NNalign" / "all.test_pred_cmb_ny"
+    raw = pd.read_csv(path, sep="\t")
+    donor_col = raw.columns[0]
+    frame = pd.DataFrame(
+        {
+            "Mut_peptide": pd.Index(raw.index.get_level_values(0)).astype(str).str.strip(),
+            "target": pd.to_numeric(
+                pd.Index(raw.index.get_level_values(1)).astype(str), errors="coerce"
+            ),
+            "HLA_allele": pd.Index(raw.index.get_level_values(3)).astype(str).str.strip(),
+            "Patient": raw[donor_col].astype(str).str.strip().to_numpy(),
+            "predicted": pd.to_numeric(raw["Prediction"], errors="coerce").to_numpy(),
+        }
+    )
+    frame = frame.dropna(subset=["target", "predicted"]).reset_index(drop=True)
+    frame["subject_id"] = _sid(frame)
+    return frame
+
+
+def main(argv: list[str] | None = None) -> int:
+    ap = argparse.ArgumentParser(description=__doc__)
+    ap.add_argument(
+        "--src",
+        type=Path,
+        required=True,
+        help="extracted IMPROVE_paper results.zip directory (the 'results copy' folder)",
+    )
+    ap.add_argument("--out", type=Path, default=HERE)
+    args = ap.parse_args(argv)
+
+    rf = load_rf(args.src)
+    spine = rf["improve_simple"]
+
+    # One row per (patient, peptide, HLA). The released tables already satisfy this; assert it
+    # rather than trusting it, because a silent duplicate would be averaged away downstream.
+    if spine["subject_id"].duplicated().any():
+        raise SystemExit("spine has duplicate subject_id values")
+
+    pep_len = spine["Mut_peptide"].str.len()
+    cohort = pd.DataFrame(
+        {
+            "subject_id": spine["subject_id"],
+            "patient": spine["Patient"].astype(str),
+            "tumour_cohort": spine["cohort"].astype(str),
+            # Predictors are trained overwhelmingly on 9-mers; off-length behaviour is the
+            # prespecified question, so the split is 9 versus everything else.
+            "peptide_length": pep_len.map(lambda n: "9mer" if n == 9 else "non_9mer"),
+            "cv_partition": spine["Partition"].astype(int),
+            "mut_peptide": spine["Mut_peptide"].astype(str),
+            "hla_allele": spine["HLA_allele"].astype(str),
+            "response": spine["response"].astype(int),
+        }
+    ).reset_index(drop=True)
+
+    rows = []
+    for model_id, frame in rf.items():
+        # Outcome must agree across tables; they are six models over one labelled row set.
+        merged = frame[["subject_id", "response", "prediction_rf"]].merge(
+            cohort[["subject_id", "response"]], on="subject_id", suffixes=("_m", "_c")
+        )
+        if len(merged) != len(cohort) or (merged["response_m"] != merged["response_c"]).any():
+            raise SystemExit(f"{model_id}: row set or labels disagree with the spine")
+        rows.append(
+            pd.DataFrame(
+                {
+                    "subject_id": frame["subject_id"],
+                    "model_id": model_id,
+                    "predicted": frame["prediction_rf"].astype(float),
+                }
+            )
+        )
+
+    nn = load_nnalign(args.src)
+    nn = nn[nn["subject_id"].isin(set(cohort["subject_id"]))].copy()
+    nn = nn.drop_duplicates(subset="subject_id", keep=False)
+    labels = cohort.set_index("subject_id")["response"]
+    nn_label = (nn["target"] > 0.5).astype(int).to_numpy()
+    if (labels.loc[nn["subject_id"]].to_numpy() != nn_label).any():
+        raise SystemExit("nnalign: labels disagree with the spine")
+    rows.append(
+        pd.DataFrame(
+            {
+                "subject_id": nn["subject_id"],
+                "model_id": "nnalign",
+                "predicted": nn["predicted"].clip(0.0, 1.0).astype(float),
+            }
+        )
+    )
+
+    predictions = pd.concat(rows, ignore_index=True)
+
+    args.out.mkdir(parents=True, exist_ok=True)
+    cohort.to_parquet(args.out / "cohort.parquet", index=False)
+    predictions.to_parquet(args.out / "predictions.parquet", index=False)
+
+    print(f"cohort:      {len(cohort):,} peptides, {cohort['patient'].nunique()} patients")
+    print(
+        f"outcome:     {int(cohort['response'].sum())} immunogenic "
+        f"({cohort['response'].mean():.3%} prevalence)"
+    )
+    print(f"predictions: {len(predictions):,} rows, {predictions['model_id'].nunique()} models")
+    for model_id, grp in predictions.groupby("model_id"):
+        print(
+            f"  {model_id:28} n={len(grp):6,}  "
+            f"coverage={len(grp) / len(cohort):.4f}  "
+            f"mean={grp['predicted'].mean():.4f}  "
+            f"range=[{grp['predicted'].min():.4f}, {grp['predicted'].max():.4f}]"
+        )
+    return 0
+
+
+if __name__ == "__main__":
+    raise SystemExit(main())
diff --git a/examples/improve/cohort.yaml b/examples/improve/cohort.yaml
new file mode 100644
index 0000000..371595c
--- /dev/null
+++ b/examples/improve/cohort.yaml
@@ -0,0 +1,108 @@
+# IMPROVE neoepitope immunogenicity cohort -- predval's first non-imaging, non-clinical fixture.
+#
+# 17,520 candidate neoepitopes from 70 patients (bladder, melanoma and a mixed "Basket" cohort),
+# each with an experimentally determined T-cell recognition outcome; 467 are immunogenic
+# (2.666% prevalence). Source: Borch et al. 2024, "IMPROVE: a feature model to predict
+# neoepitope immunogenicity through broad-scale validation of T-cell recognition",
+# Frontiers in Immunology, data release github.com/SRHgroup/IMPROVE_paper (results.zip).
+#
+# The predictions are the out-of-fold scores from that paper's own 5-fold cross-validation,
+# which was partitioned BY PATIENT: every patient sits in exactly one partition, so no
+# prediction here was produced by a model that had seen another peptide from the same patient.
+# That is what makes the file admissible as a predictions-only validation input rather than
+# an apparent-fit artefact.
+#
+# Regenerate the parquet files with:
+#     uv run python examples/improve/build_fixture.py --src <extracted 'results copy' dir>
+
+cohort_id: improve-neoepitope-cv-patient-grouped
+version: 0
+subject_key: subject_id
+data: cohort.parquet
+
+outcome:
+  type: binary
+  field: response
+  positive_label: 1
+
+# A candidate neoepitope is not a patient. These 17,520 peptides come from 70 patients, and
+# peptides from one patient share that patient's HLA genotype, mutational process, tumour
+# transcriptome and TME -- the dominant sources of variance in every feature these models use.
+# Uncertainty therefore has to be computed over patients. This is the direct analogue of PCam's
+# 22 slides, with a much larger mean cluster size (250 peptides per patient).
+#
+# unit_noun is the unit of OBSERVATION (one row = one candidate peptide-HLA pair);
+# clustering.name is the unit of INDEPENDENCE (one patient).
+unit_noun: peptide
+
+clustering:
+  field: patient
+  name: patient
+
+coverage:
+  min_fraction: 0.95
+  compare_on: both
+
+completeness:
+  require_outcome: true
+  on_violation: drop_and_report
+
+# Prespecified, both chosen before any metric was computed on this cohort.
+#
+# tumour_cohort: the release's own three accrual groups (bladder n=6,237/24 patients;
+#   melanoma n=5,921/26; Basket n=5,362/20). A selection model that holds in melanoma -- the
+#   indication these predictors were overwhelmingly trained and tuned on -- and fails elsewhere
+#   is precisely the failure mode that matters for MSS colorectal disease, which is not in this
+#   cohort at all. This subgroup is the closest available proxy for that transfer question and
+#   its result bounds how far any verdict here may be carried.
+#
+# peptide_length: 9-mers versus everything else. Peptide-HLA binding predictors, which supply
+#   the dominant features (RankEL, RankBA, Stability), are trained overwhelmingly on 9-mers and
+#   are known to degrade off that length. Declared because the off-length stratum is where a
+#   calibration fault would be expected a priori, not because it was observed.
+subgroups:
+  - name: tumour_cohort
+    field: tumour_cohort
+  - name: peptide_length
+    field: peptide_length
+
+# 0.5 is the nominal probability-scale decision threshold and is reported because it is what a
+# naive consumer of these scores would use. It is the wrong threshold for this problem -- see
+# the budgeted-yield sidecar (examples/improve/budget_yield.py), which treats selection as
+# top-k under a 3-20 antigen slot budget rather than as thresholding -- but a threshold the
+# field would actually reach for is worth measuring before replacing it.
+thresholds: [0.5]
+
+uncertainty:
+  n_boot: 2000
+  seed: 1337
+  ci_level: 0.95
+  show_naive_ci: true
+
+# Roster: 7 declared, 7 expected present.
+#
+# Six are the IMPROVE random forest itself, in the three feature configurations the paper
+# releases (Simple; TME_excluded, which adds cellular prevalence and a prioritisation score;
+# TME_included, which further adds cytolytic activity, HLA expression and MCP-counter TME
+# estimates), each in its full form and in the `wo_prime` ablation that drops the PRIME feature.
+# The seventh, `nnalign`, is the paper's NNAlign comparator -- a different model class trained
+# against a 0.001/1.0 target encoding, which is why its scores live on a different numeric scale
+# and why including it is informative rather than redundant.
+expected_models:
+  - improve_simple
+  - improve_simple_wo_prime
+  - improve_tme_excl
+  - improve_tme_excl_wo_prime
+  - improve_tme_incl
+  - improve_tme_incl_wo_prime
+  - nnalign
+on_missing_model: fail
+
+# Empty, deliberately, and the reasoning is worth recording because a random forest invites the
+# opposite answer. predval's `ensemble_members` flags a model_id whose predictions are a blend
+# of SEVERAL SEPARATELY TRAINED MODELS, because such a blend's apparent performance is inflated
+# by the selection of its own weights. An IMPROVE random forest is one trained model whose
+# internal bagging is part of its fitting procedure, not a post-hoc blend over a roster -- the
+# same reason a gradient-boosted model is not declared here either. No member of this roster is
+# a blend, so the ensemble-construction-bias caution correctly does not appear.
+ensemble_members: []
diff --git a/examples/improve/run_evaluation.py b/examples/improve/run_evaluation.py
new file mode 100644
index 0000000..8f6023d
--- /dev/null
+++ b/examples/improve/run_evaluation.py
@@ -0,0 +1,198 @@
+#!/usr/bin/env python
+"""Evaluate the IMPROVE neoepitope fixture end to end and write the artefacts.
+
+    uv run python examples/improve/run_evaluation.py [--out DIR] [--n-boot N]
+
+Writes metrics.parquet, fragility.parquet, coverage.parquet, manifest.json, report.html and
+findings.json, then prints the summary: prevalence anchoring, discrimination, the
+unit-of-analysis exhibit over patients, calibration slope and intercept, the recalibration
+ladder (as-published versus cross-fitted) and fragility to dropping any single patient.
+At the declared B=2000 over 7 models this takes a few minutes; pass --n-boot for a smoke run.
+"""
+
+from __future__ import annotations
+
+import argparse
+from pathlib import Path
+
+import pandas as pd
+
+from predval import evaluate, load_cohort, load_predictions, write_findings, write_report
+from predval.metrics import optimism
+
+HERE = Path(__file__).resolve().parent
+W = 26  # model ids here are long; widen every column that prints one
+
+
+def _overall(metrics: pd.DataFrame, metric: str, ci_method: str | None = None) -> pd.DataFrame:
+    sel = (
+        (metrics["stratum_kind"] == "overall")
+        & (metrics["subset"] == "common")
+        & (metrics["metric"] == metric)
+    )
+    if ci_method is not None:
+        sel &= metrics["ci_method"] == ci_method
+    return metrics[sel]
+
+
+def main(argv: list[str] | None = None) -> int:
+    ap = argparse.ArgumentParser(description=__doc__)
+    ap.add_argument("--out", type=Path, default=HERE / "out")
+    ap.add_argument("--n-boot", type=int, default=None, help="override uncertainty.n_boot")
+    args = ap.parse_args(argv)
+
+    cohort = load_cohort(HERE / "cohort.yaml")
+    predictions = load_predictions(HERE / "predictions.parquet")
+
+    if args.n_boot is not None:
+        unc = cohort.spec.uncertainty.model_copy(update={"n_boot": args.n_boot})
+        spec = cohort.spec.model_copy(update={"uncertainty": unc})
+        cohort = type(cohort)(**{**cohort.__dict__, "spec": spec})
+
+    result = evaluate(cohort, predictions)
+    paths = result.write(args.out)
+    report_path = write_report(result, args.out / "report.html")
+    findings_path = write_findings(result, args.out / "findings.json")
+
+    print(f"cohort:      {cohort.n_subjects:,} peptides")
+    print(f"models:      {len(predictions.model_ids)}")
+    print(f"metrics:     {len(result.metrics):,} rows -> {paths['metrics']}")
+    print(f"report:      {report_path}")
+    print(f"findings:    {findings_path}")
+    print()
+
+    print("FLAGS")
+    for flag in result.flags:
+        print(f"  [{flag.severity:7}] {flag.code}: {flag.message}")
+    print()
+
+    # ---- prevalence anchoring: the first thing to look at under 2.7% positives -------------
+    frame = predictions.frame
+    prev = float(cohort.table[cohort.spec.outcome.field].mean())
+    print(f"Prevalence anchoring -- observed immunogenic rate {prev:.4%}")
+    print(f"  {'model':{W}} {'mean pred':>10} {'ratio':>7} {'min':>8} {'max':>8}")
+    for model, grp in frame.groupby("model_id"):
+        mean_p = float(grp["predicted"].mean())
+        print(
+            f"  {model:{W}} {mean_p:10.4f} {mean_p / prev:6.1f}x "
+            f"{grp['predicted'].min():8.4f} {grp['predicted'].max():8.4f}"
+        )
+    print()
+
+    auroc = _overall(result.metrics, "auroc", "cluster_bootstrap").sort_values(
+        "value", ascending=False
+    )
+    print("AUROC by model (common subset, patient-clustered bootstrap 95% CI)")
+    for _, r in auroc.iterrows():
+        print(f"  {r['model_id']:{W}} {r['value']:.4f}  [{r['ci_low']:.4f}, {r['ci_high']:.4f}]")
+    print()
+
+    ap_rows = _overall(result.metrics, "average_precision", "cluster_bootstrap").sort_values(
+        "value", ascending=False
+    )
+    if len(ap_rows):
+        print(f"Average precision (no-skill baseline = prevalence = {prev:.4f})")
+        for _, r in ap_rows.iterrows():
+            print(
+                f"  {r['model_id']:{W}} {r['value']:.4f}  [{r['ci_low']:.4f}, {r['ci_high']:.4f}]"
+            )
+        print()
+
+    naive = _overall(result.metrics, "auroc")
+    print("Unit-of-analysis exhibit -- AUROC interval width, patient vs naive per-peptide")
+    for model, grp in naive.groupby("model_id"):
+        widths = {r["ci_method"]: r["ci_high"] - r["ci_low"] for _, r in grp.iterrows()}
+        if "naive_row_bootstrap" in widths and widths["naive_row_bootstrap"] > 0:
+            ratio = widths["cluster_bootstrap"] / widths["naive_row_bootstrap"]
+            print(
+                f"  {model:{W}} patient {widths['cluster_bootstrap']:.4f}  "
+                f"naive {widths['naive_row_bootstrap']:.4f}  ratio {ratio:5.1f}x"
+            )
+    print()
+
+    for metric, blurb in (
+        ("calibration_slope", "1.0 = correct spread; <1 = over-spread/over-confident"),
+        ("calibration_intercept", "0.0 = correctly anchored; <0 = systematically too high"),
+    ):
+        rows = _overall(result.metrics, metric, "cluster_robust_t").sort_values("value")
+        if not len(rows):
+            continue
+        print(f"{metric.replace('_', ' ').capitalize()} ({blurb})")
+        for _, r in rows.iterrows():
+            interval = f"[{r['ci_low']:8.3f}, {r['ci_high']:8.3f}]"
+            print(f"  {r['model_id']:{W}} {r['value']:8.3f}  {interval}")
+        print()
+
+    lad = _overall(result.metrics, "brier").pivot_table(
+        index="model_id", columns=["rung", "fit_mode"], values="value"
+    )
+
+    def cell(model: str, rung: str, mode: str) -> float:
+        col = (rung, mode)
+        return float(lad.loc[model, col]) if col in lad.columns else float("nan")
+
+    def cf_gain(model: str) -> float:
+        r0 = cell(model, "rung0", "apparent")
+        cf = min(cell(model, r, "crossfit") for r in ("rung1", "rung2", "rung3"))
+        return r0 - cf
+
+    order = sorted(lad.index, key=cf_gain, reverse=True)
+    print("Recalibration ladder -- Brier, as-published vs cross-fitted by held-out patients")
+    print(f"  reference: a constant predictor at the prevalence scores {prev * (1 - prev):.4f}")
+    print("  optim_r2 > 0 == the apparent rung2 fit flattered itself versus held-out patients")
+    print(
+        f"  {'model':{W}} {'rung0':>8} {'r1_cf':>8} {'r2_cf':>8} {'r3_cf':>8} "
+        f"{'gain_cf':>8} {'optim_r2':>9}"
+    )
+    for model in order:
+        r0 = cell(model, "rung0", "apparent")
+        r1_cf = cell(model, "rung1", "crossfit")
+        r2_cf = cell(model, "rung2", "crossfit")
+        r3_cf = cell(model, "rung3", "crossfit")
+        gain = r0 - min(r1_cf, r2_cf, r3_cf)
+        optim = optimism("brier", cell(model, "rung2", "apparent"), r2_cf)
+        print(
+            f"  {model:{W}} {r0:8.4f} {r1_cf:8.4f} {r2_cf:8.4f} {r3_cf:8.4f} "
+            f"{gain:8.4f} {optim:+9.4f}"
+        )
+    print()
+
+    gains = _overall(result.metrics, "paired_gain_brier")
+    if len(gains):
+        print("Paired Brier gain, cross-fitted (CI excluding 0 == a demonstrated repair)")
+        for model, grp in gains.groupby("model_id"):
+            parts = []
+            for _, r in grp.sort_values("rung").iterrows():
+                mark = " "
+                if pd.notna(r["ci_low"]) and pd.notna(r["ci_high"]):
+                    if r["ci_low"] > 0:
+                        mark = "+"
+                    elif r["ci_high"] < 0:
+                        mark = "-"
+                parts.append(f"{r['rung']}{mark}{r['value']:+.5f}")
+            print(f"  {model:{W}} " + "  ".join(parts))
+        print("  + = repair demonstrated   - = recalibration reliably HARMS   blank = unproven")
+        print()
+
+    frag = result.fragility[
+        (result.fragility["stratum_kind"] == "overall")
+        & (result.fragility["subset"] == "common")
+        & (result.fragility["metric"] == "auroc")
+    ].sort_values("max_abs_delta", ascending=False)
+    print("Fragility -- largest AUROC change from dropping any single patient")
+    for _, r in frag.head(7).iterrows():
+        print(f"  {r['model_id']:{W}} {r['max_abs_delta']:.4f}  (patient {r['culprit_cluster']})")
+
+    boundary = _overall(result.metrics, "boundary_count")
+    boundary = boundary[boundary["value"] > 0]
+    if len(boundary):
+        print()
+        print("Data quality -- predictions exactly at 0 or 1 (eps-clip load)")
+        for _, r in boundary.iterrows():
+            print(f"  {r['model_id']:{W}} {int(r['value']):5d}")
+
+    return 0
+
+
+if __name__ == "__main__":
+    raise SystemExit(main())