Add paper-ready benchmark evals + plots as `mhcflurry` subcommands (clean-room, no Keras/TF)

[iskandr]

Summary

The published-paper analysis layer — the multi-tool accuracy benchmark, the significance testing, the headline figures, and the supplemental tables (Tables S1/S2/S6/S7, Data S1–S6) — currently exists only as a stack of one-off Jupyter notebooks (the 2023-retraining set) that never landed in the repo. Those notebooks are TF1/Keras-era, depend on hardcoded cluster paths (/home/odonnt02/..., /data/mhcflurry-evaluation/), pull inputs from Springer/Dropbox URLs, use unseeded RNG, and reference three different download versions. They will not run against the current PyTorch 2.3 tree and are not something we should port line-for-line.

Meanwhile the in-repo reproducible pipeline stops at scored prediction tables (downloads-generation/data_evaluation) and a 2-way new-vs-public release-regression evaluator (mhcflurry compare-models / plot-model-comparison). There is no reproducible in-tree path that answers the paper's actual question — "does MHCflurry beat the field (netMHCpan, MixMHCpred), and is the difference significant?" — nor one that regenerates the figures/tables.

This issue proposes rebuilding that whole layer clean-room as first-class mhcflurry subcommands plus one orchestration script, reusing the existing metric core and CLI conventions, and carrying none of the historical Keras/TF baggage.

Design principles (the "no baggage" contract)

1. Clean-room, current-era inputs only. Consume only current artifacts: the data_evaluation scored benchmark, data_predictions (external-tool scores via mhctools, already produced upstream — never re-invoke netMHCpan/MixMHCpred inline the way the notebooks do), analysis_predictor_info, data_mass_spec_annotated, data_curated, data_references (FM-index for peptide→protein mapping instead of pyensembl), and current PyTorch model bundles. No keras/tensorflow imports anywhere.
2. No hidden state, no hardcoded paths. Every input is a CLI flag with a mhcflurry-downloads-resolvable default. Datasets the notebooks hardcoded as sibling files (expression_groups.xlsx, the immunogenicity screening set, Wolf-Levy S2) become explicit --dataset/--input args or documented download items.
3. Deterministic. All decoy resampling and bootstrap steps take --seed; default seeded. Golden-CSV tests become possible.
4. Follow the established metric/plot split. Metric subcommands write CSV/JSON only (torch/pandas/sklearn, no matplotlib) exactly like compare-models; plot subcommands read those artifacts and render, with matplotlib/logomaker lazy-imported, exactly like plot-model-comparison.
5. One metric implementation. Refactor the shared scoring core (ROC-AUC / PR-AUC / PPV@N, per-allele/per-sample/per-length aggregation) out of mhcflurry/cli/compare_models.py into a library (mhcflurry/benchmarks/) and have both compare-models and the new paper commands call it, so the release-regression numbers and the paper numbers can never silently diverge.
6. Register through the existing lazy dispatcher. Add each new subcommand to _SUBCOMMANDS + _HELP_GROUPS in mhcflurry/cli/main.py (new group: "Paper benchmarks"), each exposing run_argv(argv); add matching mhcflurry-* console-script shims in setup.py for parity. Builds directly on Unify mhcflurry CLI under a single mhcflurry parent command #291.

Proposed shared library: mhcflurry/benchmarks/

Pure-python, matplotlib-free, importable by both the regression and paper commands:

• metrics.py — roc_auc, pr_auc, ppv_at_n, per-{allele,sample,length} aggregation. (Lifted from compare_models.py; compare-models switches to importing it.)
• stats.py — seeded bootstrap CIs, paired ttest_rel, binomial win-counts, percent-change-vs-baseline normalization.
• assembly.py — benchmark table assembly: predictor-column merge, train-leak assertion, seeded per-sample decoy resampling, sample-group tagging, external-tool masking.
• predictor_registry.py + predictor_info.csv — replaces notebook 0 aesthetics: a versioned data file mapping predictor id → display name, plot color, sort order, primary flag. Single source of truth for every plot command (kills the runtime color-hacking in the old aesthetics cell).

Proposed subcommands

Metric / eval (matplotlib-free → CSV/JSON)

Subcommand	Replaces notebook(s)	Reads	Writes
benchmark-assemble	1 prepare benchmark dataset - {mono,multi}allelic	data_evaluation, data_mass_spec_annotated	benchmark.{mono,multi}allelic.parquet, manifest.json
benchmark-score	2 multiallelic accuracy, 2 monoallelic accuracy plots	assembled benchmark	accuracy_scores.{mono,multi}.csv, significance.csv, summary.json
benchmark-novel-alleles	2 monoallelic analysis of novel alleles	assembled benchmark, models_class1_pan_variants	novel_alleles.csv
benchmark-immunogenicity	immunogenicity experiments	--dataset (screening set), data_predictions	immunogenicity_scores.csv
benchmark-processing	1 proteasome mass spec, 2 processing predictor motifs	--dataset (Wolf-Levy S2), data_references FM-index, presentation predictor	processing_proteasome.csv, processing_correlation.parquet, processing_motifs.parquet
model-selection-accuracy	1 model selection predictor accuracy (retires the …UNUSED, NOW WE PRECOMPUTE THIS notebook)	analysis_predictor_info	model_selection_accuracy.csv
sample-table	2 sample table, 4 supplemental sample table with accuracies	data_mass_spec_annotated, assembled benchmark, --expression-groups	sample_table.csv, sample_table_with_accuracies.csv

benchmark-score detail (the scientific core): per-sample × per-length × per-predictor ROC-AUC / PR-AUC / PPV@N over all predictors including external tools; percent-change vs configurable baseline predictors; seeded bootstrap CIs; paired significance (ttest_rel + binomial win-count) for each predictor against each baseline; optional cysteine-removed variant (--exclude-residues C). This is the N-way, statistics-bearing generalization of what compare-models' affinity/presentation components do 2-way — hence the shared metrics.py.

Plot (matplotlib/logomaker isolated → figures)

Subcommand	Replaces	Renders
plot-benchmark	3 scores plots - {mono,multi}allelic	multi-tool AUC/PPV scatter grids, per-length bars, presentation-score panels, mean-PPV summary; --panels selects subsets
plot-processing	4 processing predictor plots (+ its "redo post resubmission" dupe)	AP sequence logos, processing-vs-affinity correlation panels w/ significance annotations

(Novel-allele and immunogenicity panels fold into plot-benchmark --panels rather than getting their own commands.)

Bespoke release packaging → sub-scripts (not core subcommands)

Supplemental assembly is release-specific and file-shuffling-heavy, so it stays a script rather than a general tool:

• scripts/release/make_supplemental.sh — replaces 1 training datasets copy, 3 make additional files, 5 make supplemental data: assembles Tables S1/S2/S6/S7 + Data S1–S6 with publication column renaming from the eval outputs + training-data copies. No ln -s symlink deliverables (the notebooks' approach); real file copies.

Tie-it-together orchestration

scripts/release/paper_eval.sh, mirroring the existing scripts/release/retrain_evaluate_deploy.sh conventions (set -euo pipefail, run_cmd command-echo, --dry-run, phase --skip-* flags):

scripts/release/paper_eval.sh \
    --run-dir /path/to/release-run \
    --out /path/to/paper-eval \
    [--data-dir <data_evaluation path>] \
    [--skip-assemble] [--skip-score] [--skip-specialized] \
    [--skip-plots] [--skip-supplemental] [--dry-run]

Flow: mhcflurry-downloads fetch … → benchmark-assemble → benchmark-score → (novel-alleles, immunogenicity, processing, model-selection-accuracy, sample-table) → plot-benchmark + plot-processing → make_supplemental.sh.

Then wire a --paper-eval phase into retrain_evaluate_deploy.sh that simply invokes this script after the existing compare-models/plot-model-comparison step, so a release run can optionally produce the full paper bundle from the same entry point.

Explicitly out of scope (the baggage we are not carrying)

• No line-for-line port of the Keras/TF notebooks; the notebooks are reference only.
• No inline re-invocation of external predictor binaries — external scores come from data_predictions.
• No hardcoded machine paths, sibling-dir inputs, or Dropbox/Springer URL fetches.
• pyensembl + Ensembl-release dependency dropped in favor of the data_references FM-index already in the pipeline. logomaker/statannot become optional plotting extras, not hard deps of the metric path.

Suggested phasing

1. Extract mhcflurry/benchmarks/ (metrics/stats/assembly/registry); refactor compare-models onto it (no behavior change; golden-CSV regression test).
2. benchmark-assemble + benchmark-score + plot-benchmark — the headline mono/multi benchmark. This alone reproduces the main figures.
3. sample-table + model-selection-accuracy (Table S1/S2 substrate).
4. benchmark-processing + plot-processing; then benchmark-novel-alleles, benchmark-immunogenicity.
5. make_supplemental.sh + paper_eval.sh; wire --paper-eval into retrain_evaluate_deploy.sh.

Testing

Tiny fixtures under test/data/, golden CSVs for each metric command, and a --limit-samples/--limit-files fast path (as compare-models already has --limit-files) so the whole chain runs in CI on a handful of samples.

Related

• Unify mhcflurry CLI under a single mhcflurry parent command #291 (unify CLI under a single mhcflurry parent) — this builds on that dispatcher.
• Roadmap: 2.3.0, 2.4.0, and 3.x model releases #287 (release roadmap) — needed to hold the 2.3.0+ releases to the paper's evidentiary bar.
• Fix unused --proteome-peptides argument in analysis predictor info decoy generation #289 (analysis_predictor_info decoy arg) — model-selection-accuracy consumes that pipeline's output.

---

Source for the notebook→subcommand mapping: audit of the untracked notebooks/2023-retraining/ set against downloads-generation/, scripts/, and mhcflurry/cli/.

[iskandr]

Additional output requirement from the retraining experiments:

The end-to-end retraining / paper-eval workflow should emit compact release-gate tables comparing the newly trained weights against public weights. These should be generated as machine-readable CSV/JSON and optionally rendered as Markdown in the run summary, so the release decision does not depend on manually assembled notebook output.

Concretely, paper_eval.sh / benchmark-score should be able to generate tables like the following from the retraining artifacts plus public-weight baselines. These are also good inputs for paper-ready plots because they provide the same metric/family/flank-mode breakdown the figures need.

Affinity

Eval	Metric	Average	New	Public	Diff	% change
affinity	AUROC	Macro	0.9818	0.9738	+0.0080	+0.82%
affinity	AUPRC	Macro	0.6837	0.6571	+0.0266	+4.05%
affinity	PPV@N	Macro	0.6691	0.6453	+0.0238	+3.68%
affinity	AUROC	Micro	0.9868	0.9812	+0.0056	+0.57%
affinity	AUPRC	Micro	0.6654	0.6100	+0.0554	+9.09%
affinity	PPV@N	Micro	0.6824	0.6461	+0.0363	+5.62%

Presentation

Flank mode	Metric	Average	New	Public	Diff	% change
with_flanks	AUROC	Macro	0.9252	0.9178	+0.0074	+0.80%
with_flanks	AUPRC	Macro	0.3965	0.3657	+0.0308	+8.42%
with_flanks	PPV@N	Macro	0.4318	0.4058	+0.0260	+6.42%
with_flanks	AUROC	Micro	0.9425	0.9376	+0.0049	+0.53%
with_flanks	AUPRC	Micro	0.4407	0.4092	+0.0315	+7.70%
with_flanks	PPV@N	Micro	0.4733	0.4451	+0.0282	+6.33%
without_flanks	AUROC	Macro	0.9216	0.9142	+0.0074	+0.81%
without_flanks	AUPRC	Macro	0.3725	0.3463	+0.0262	+7.57%
without_flanks	PPV@N	Macro	0.4150	0.3920	+0.0229	+5.85%
without_flanks	AUROC	Micro	0.9397	0.9346	+0.0051	+0.55%
without_flanks	AUPRC	Micro	0.4165	0.3878	+0.0287	+7.41%
without_flanks	PPV@N	Micro	0.4589	0.4342	+0.0247	+5.68%

Processing

Flank mode	Metric	Average	New	Public	Diff	% change
with_flanks	AUROC	Macro	0.8491	0.8491	-0.0000	-0.00%
with_flanks	AUPRC	Macro	0.1150	0.1320	-0.0170	-12.88%
with_flanks	PPV@N	Macro	0.1844	0.2015	-0.0170	-8.46%
with_flanks	AUROC	Micro	0.8748	0.8747	+0.0002	+0.02%
with_flanks	AUPRC	Micro	0.1200	0.1420	-0.0219	-15.46%
with_flanks	PPV@N	Micro	0.2036	0.2246	-0.0209	-9.32%
no_flank	AUROC	Macro	0.8560	0.8525	+0.0035	+0.41%
no_flank	AUPRC	Macro	0.1489	0.1508	-0.0019	-1.27%
no_flank	PPV@N	Macro	0.2224	0.2209	+0.0015	+0.67%
no_flank	AUROC	Micro	0.8827	0.8788	+0.0039	+0.44%
no_flank	AUPRC	Micro	0.1668	0.1645	+0.0023	+1.40%
no_flank	PPV@N	Micro	0.2534	0.2503	+0.0031	+1.25%
short_flanks	AUROC	Macro	0.8554	0.8574	-0.0021	-0.24%
short_flanks	AUPRC	Macro	0.1277	0.1389	-0.0112	-8.04%
short_flanks	PPV@N	Macro	0.1962	0.2080	-0.0118	-5.69%
short_flanks	AUROC	Micro	0.8810	0.8825	-0.0015	-0.17%
short_flanks	AUPRC	Micro	0.1358	0.1503	-0.0145	-9.65%
short_flanks	PPV@N	Micro	0.2185	0.2320	-0.0134	-5.80%

Acceptance addendum: after a full retraining run, the workflow should produce these summary tables automatically and then feed the same aggregate outputs into the paper-ready plot commands. That means some of the missing figures should become direct products of retraining rather than separate notebook-only work.