nested-cv-analysis

Analysing nested cross-validation, with transformations that never leak between folds.

Nested CV gives you one winner per outer fold. When those winners disagree, "so which one do I use?" is not answered by the mean score. This package keeps every number the selection already computed and turns it into readings that answer it — does the win hold up, is the model predictable, what does it cost — while making leakage structurally impossible along the way.

import nested_cv_analysis as nca
from sklearn.preprocessing import StandardScaler
from sklearn.svm import SVC

nca.cross_val_score(SVC(), X, y, cv=5, transform=StandardScaler())

Installation

pip install nested-cv-analysis              # core
pip install nested-cv-analysis[plot]        # + figures
pip install nested-cv-analysis[imbalanced]  # + SMOTE and other resamplers

The import name is nested_cv_analysis; the examples abbreviate it as nca.

Why transformations have to go inside the fold

Anything that learns from data — scaling, feature selection, PCA, imputation, SMOTE — has to be refitted inside each fold. Fitting it over the whole set before validating contaminates the test split and inflates the score. The most treacherous case is feature selection: on purely random data, a leaked selection produces 0.85 accuracy where the truth is 0.50.

The two categories of transformation

The package distinguishes two roles, because they follow different rules:

parameter	when it applies	examples
transform	fitted on the fold's training split, applied to both train and test	StandardScaler, SelectKBest, PCA, SimpleImputer
resampler	applied to the training split only, never to test	SMOTE, RandomUnderSampler, ADASYN

The order is transform → resampler → estimator: the transformation is fitted on real data, not on synthetic points.

from imblearn.over_sampling import SMOTE

nca.cross_val_score(
    SVC(), X, y, cv=5, scoring="roc_auc",
    transform=StandardScaler(),
    resampler=SMOTE(random_state=0),
)

Missing values

Imputation is a transformation like any other, and a leak vector like any other: a median computed over the whole dataset carries information from the test rows into training. So it belongs in transform=, refitted per fold.

from sklearn.impute import SimpleImputer

nca.nested_cross_validate(models, X, y, transform=SimpleImputer(strategy="median"))

There is no up-front check that X is free of NaN, and deliberately so: several scikit-learn estimators accept missing values natively, so NaN with no imputer is a legitimate configuration rather than an error.

accepts NaN	does not
DecisionTree*, ExtraTree*, RandomForest*, ExtraTrees*, HistGradientBoosting*	GradientBoostingClassifier/Regressor

GradientBoosting is built on trees but validates its input first, so "the tree family" is not a safe rule of thumb — check the estimator, not the family.

Trees do not impute. At each split they learn which side the missing values fall on, choosing the (threshold, direction) pair that minimises impurity. When missingness is the stronger signal the winning pair has an infinite threshold, so the split is on absence alone. That means missingness itself becomes a predictor: SimpleImputer erases that information by filling the hole, which is a modelling choice rather than mere hygiene.

What is not legitimate is a candidate that cannot cope. If one fails on every inner fold, the run stops and says which candidate it was, adding the NaN diagnosis when X has holes and no transform was given:

RuntimeError: candidate 'logit' failed on every inner fold of outer fold 1/3.

X contains NaN and no transform was given. Either impute inside the fold,
with transform=SimpleImputer(...), or use estimators that accept NaN
natively, such as HistGradientBoostingClassifier/Regressor.

The underlying scikit-learn error is chained, not swallowed.

API

Plain validation

cross_val_score, cross_validate and cross_val_predict mirror the scikit-learn functions of the same name and accept transform= and resampler=. Every other argument — groups, multi-metric scoring, n_jobs, error_score, return_train_score, params — is forwarded unchanged.

Nested cross-validation

The part scikit-learn does not cover directly: choosing among heterogeneous models and still getting an honest estimate.

result = nca.nested_cross_validate(
    {"rf": RandomForestClassifier(), "svc": SVC(), "logit": LogisticRegression()},
    X, y,
    scoring=["accuracy", "roc_auc"],   # several metrics at once
    select_on="accuracy",              # which one drives the selection
    cv_outer=4, cv_inner=5,
    transform=StandardScaler(),
    n_jobs=-1,
)

result.generalization_score()   # the unbiased estimate
result.inner_scores             # matrix (outer folds x models)
result.selection_frequency()    # was the choice stable?
print(result.report())          # the full report

generalization_score() measures the selection procedure, not a specific model. Choosing the best of 30 candidates and reporting the winner's score is optimistic: the maximum of 30 noisy estimates is biased upward. To produce the final model, refit the winner on all the data with nca.build_pipeline(...).

Use compute_holdout=False while exploring — every look at the holdout while you tune the search space reintroduces selection bias.

What a nested CV actually hands you

Not a model. n_outer candidates — one winner per outer fold — plus an honest estimate of the procedure that produced them.

Often fewer than n_outer distinct ones; selection_frequency() says how many. If a single model wins every fold, the decision is already made and the rest of this section is unnecessary. When the folds disagree, what you are holding is a shortlist.

And this data cannot rank that shortlist. Any number computed here to crown the best of the winners is a number the selection already used, and choosing by it reintroduces precisely the bias nested CV was set up to remove. Settling it outright takes a genuinely fresh dataset.

That is a real limit. It is not a dead end. Between "pick one arbitrarily" and "go collect more data" there is a great deal of legitimate evidence — legitimate because none of it reads the holdout for anyone but each fold's own winner.

The evidence you are allowed to use

question	where
How does a winner behave in the folds it did not win?	winners_across_folds(), rank_bumps, score_heatmap
Does it hold its position under the other metrics?	ranking_by_metric(), rank_across_metrics
Is it predictable, or high on average and erratic?	inner_dispersion(), score_vs_dispersion
Is it memorising? What does it cost?	overfit_gap(), complexity_bars

All of it comes from the inner CV — the partitions on which the candidates were actually compared, identical for every candidate within a fold. None of it needs a per-candidate holdout, which is why the package refuses to compute one.

The moves this opens up

Eliminate a fold-lucky winner. A model that takes one fold and collapses in the others won by the split, not on merit. winners_across_folds() puts the two side by side: its score in the fold it won, and its mean and worst position everywhere else. A winner that never leaves the top three elsewhere is a different proposition from one that drops to last — and the second can be dropped in favour of the consistent one, even though both "won a fold".

Ensemble the winners instead of choosing. If several are defensible, you do not have to pick. Combining them is often the better answer to a genuine tie, and the shortlist is already the natural membership list:

from sklearn.ensemble import VotingClassifier

frequency = result.selection_frequency()
winners = list(dict.fromkeys(result.best_model_names))   # distinct, in win order

ensemble = VotingClassifier(
    [(name, nca.build_pipeline(models[name], transform=StandardScaler()))
     for name in winners],
    weights=[frequency[name] for name in winners],       # folds won as the weight
)
ensemble.fit(X, y)

Note what is being ensembled: the winning model types, each refit on all the data. The fitted objects in best_estimators each saw only one outer fold's training split — they are what the estimate was computed from, not what you deploy.

The default voting="hard" is used deliberately: "soft" needs predict_proba on every member, which an SVC() without probability=True does not have. And if the shortlist has one entry, this builds a one-member ensemble — check len(winners) first, or just ship the model.

Or both. Use the evidence to decide which winners are strong enough to stand alone and which only earn their place inside an ensemble. A model that is second or third under every metric and in every fold may never win outright and still be the one you would rather ship.

The one thing to keep straight

generalization_score() estimates the procedure that actually ran: fit the candidates, select by the inner CV under this rule. Judgement applied afterwards — dropping a fold-lucky winner, ensembling the rest — is a step the estimate does not cover, so it is no longer a description of what you did. The bias is far milder than selecting by the holdout, since the evidence is the same inner-CV evidence the selection already used, but it is not zero.

If you want the estimate to cover your rule, encode the rule in select= so it runs inside every fold and the nested CV measures it. one_standard_error is exactly that: a preference for parsimony, expressed as a selector rather than applied by hand afterwards.

Deciding what to do with the best models

With profile_all_models=True (the default), each candidate is fitted on each outer fold's training split to collect training score, complexity and timings — the holdout is not touched.

Does the win hold up?

The outer fold is the unit of analysis. Each model has one value per fold — the mean of its inner CV there, on the same partitions for every candidate — and the question "how does one fold's winner behave in the folds it did not win?" reads straight off it.

result.winners_across_folds()        # each winner in the other folds
result.rank_table()                  # positions fold by fold, best and worst
nca.plots.score_heatmap(result)      # performance per fold, ● on the fold it won
nca.plots.rank_bumps(result)         # the same reading, in positions

  model          fold1   fold2   fold3   fold4   fold5  mean elsewhere  worst pos
  svc-rbf       0.805*   0.795  0.790*  0.805*   0.805           0.800   2 of 11
  forest-100     0.805  0.800*   0.790   0.795   0.765           0.789   5 of 11
  knn-5          0.770   0.780   0.765   0.790  0.820*           0.776   7 of 11

knn-5 won fold 5 with 0.820 — and in the other four it never got past 0.790, dropping to 7th of 11. The split won, not the model. svc-rbf, even in the folds it lost, never left second place.

No model is re-evaluated on someone else's holdout: everything comes from the scores the selection already computed, which is why the reading exists even with profile_all_models=False.

How much did the choosing flatter itself?

The winner's inner score is the maximum of n_models noisy estimates, and the maximum of noise is biased upward. The holdout is the honest read of that same model. The distance between them, fold by fold, is the optimism of the selection — and the picture behind generalization_score().

result.selection_gap()             # per fold: winner, inner, holdout, gap
nca.plots.selection_gap(result)    # the same, as dumbbells

  fold  winner          inner CV   holdout       gap
     1  forest-100        0.8525    0.8250   +0.0275
     2  svc-rbf           0.8562    0.8300   +0.0262
     3  svc-rbf           0.8350    0.8750   -0.0400
     4  gradient-boost    0.8475    0.8250   +0.0225
     5  svc-rbf           0.8562    0.8300   +0.0262
  mean gap +0.0125

Fold 3 came out negative — a single fold's gap is noisy, because the holdout is one split and often a small one. The mean is the reading that carries weight; it is exactly the winners' mean inner score minus generalization_score(). Were the pairs to coincide, reporting the winner's CV score would be safe.

Would another metric have chosen differently?

The ranking is not a property of the models alone — it is a property of the models and the metric. Since every candidate was already scored on every metric in every fold, replaying the selection over a different column costs no refitting.

result.selection_by_metric()             # per metric: who would win each fold
result.ranking_by_metric()               # per model: its position under each metric
nca.plots.rank_across_metrics(result)    # all of it, consolidated

  model           accuracy pos  roc_auc pos  f1 pos  swing
  forest-100                1*           1*       3      2
  gradient-boost            2*            3      1*      2
  svc-rbf                   3*            2       5      3
  knn-15                     4           4*       6      2
  naive-bayes                8            6      2*      6

  accuracy   gradient-boost, forest-100, svc-rbf, forest-100, forest-100  <- drove the selection
  roc_auc    forest-100, forest-100, knn-15, forest-100, forest-100   (3/5 same as accuracy)
  f1         gradient-boost, naive-bayes, gradient-boost, ...          (1/5 same as accuracy)

naive-bayes is 8th on accuracy and 2nd on f1, and would take two folds there. Three metrics, three different selections: on this data the metric is part of the choice, so it has to be argued for rather than defaulted into.

The position is defined within each metric, on its own: take the model's inner-CV score in each outer fold under that metric, average over the outer folds, and count how many models beat that mean — position = 1 + that count. So 1 is best, ties share a position (1, 2, 2, 4, never 1, 2, 2, 3), and a model never scored under a metric has no position there. Scores follow the sklearn convention that higher is better, error metrics already negated, which is why one comparison serves every metric. Rows are ordered by the position under select_on, and swing is worst position minus best.

The figure is rank_bumps with metrics on the x axis instead of folds: a flat line is a model the metric does not move, a diving line is one you would only pick by having picked its metric first.

The * above — a ring on the marker in the figure — is not "position 1". It marks a model that would take at least one outer fold under that metric, and the replay uses the run's own selection rule, not argmax: a one_standard_error run replays as one_standard_error. Position 1 without a star, or a star away from position 1, is expected — the position ranks a mean over folds, a win is decided fold by fold, and the rule may prefer a simpler model within noise.

These re-selections carry no holdout. In each outer fold only the model the real selection chose was ever scored on the held-out data — scoring the others there would create the number it is tempting to select by. So this answers would my choice have changed? but not how would that choice have generalised. For the second question, re-run with select_on= set to that metric.

Dispersion, overfitting and cost

result.inner_dispersion()    # how much each model swings between outer folds
result.overfit_gap()         # training score minus inner-CV score
result.decision_table()      # score, dispersion, positions, gap, cost
result.complexity_ranking()  # most to least complex
result.to_frame()            # the same, as a DataFrame

Many candidates

With dozens, hundreds or thousands of models, every figure and table accepts top and never discards a winner — the winner is exactly what you came looking for. The subtitle declares what was left out; nothing is truncated silently.

nca.plots.rank_bumps(result, top=25)              # default
nca.plots.rank_bumps(result, only_winners=True)   # recommended above ~100
result.report(top=20)

The scatter thins its marks and transparency as the cloud grows, caps the direct labels (max_labels=8) and pushes apart the ones that would overlap.

The holdout is seen by one model only. In each outer fold, only the inner CV's winner is scored on the held-out fold. There is no per-candidate holdout matrix — and that is deliberate: such a number invites you to select by it, and selecting by the holdout voids the guarantee it offers. To compare candidates there are the inner-CV readings, which are the partitions on which they were actually compared.

Parsimony: the one-standard-error rule

Prefer the simplest, or the least dispersed, among the statistically tied:

nca.nested_cross_validate(
    models, X, y,
    select=nca.one_standard_error(prefer="complexity"),
)

It accepts any candidate within one standard error of the best and, among those, picks the one with the lowest complexity (or the lowest dispersion, the default). Complexity is measured as two quantities that are deliberately not merged: effective parameters (coefficients, tree nodes, support vectors — NaN for unrecognised families) and fit and predict time, which is universal.

Figures

pip install nested-cv-analysis[plot]

nca.plots.rank_bumps(result)            # does the win hold up?
nca.plots.winner_facets(result)         # one panel per winner
nca.plots.score_vs_dispersion(result)   # high AND predictable?
nca.plots.metric_panels(result)         # the same, per metric
nca.plots.score_heatmap(result)         # per fold (or source="train")
nca.plots.selection_gap(result)         # inner CV vs holdout, per fold
nca.plots.rank_across_metrics(result)   # does the metric decide the model?
nca.plots.complexity_bars(result)       # parameters and time
nca.plots.decision_dashboard(result)    # the four readings together

All return a Figure and never display it. They accept dark=True. No value exists in colour alone: every figure has a textual equivalent in report() and decision_table().

The same call in a script and in a notebook

The figures are built straight from matplotlib.figure.Figure, never through pyplot. That one decision is what makes both contexts work without flags:

	what you get
Notebook / Colab	returning a figure from a cell displays it once. A pyplot-managed figure is flushed by the inline backend and rendered again by the return value — the same chart twice. They also never accumulate, so no plt.close() and no "More than 20 figures" warning
Script / CLI	no backend is selected on import, so matplotlib.use("Agg") is unnecessary — savefig works headless, over SSH, in CI or in a container

Batch work, in either context:

figures = nca.plots.build_all(result)          # {name: Figure}, nothing written
figures["selection_optimism"]                  # displays inline

nca.plots.save_all(result, "output")           # PNGs, both themes
nca.plots.save_all(result, "output", themes=("light",), names=["dashboard"])

build_all skips what a run cannot support — cost needs profile_all_models=True, selection_optimism needs compute_holdout=True — so list(figures) tells you what that run actually had. Pass skip_unavailable=False to get the error instead.

With several metrics

Every score-based figure takes metric=, and the batch can repeat them:

nca.plots.selection_gap(result, "roc_auc")        # one figure, one metric

figures = nca.plots.build_all(result, metrics="all")
figures["selection_optimism_roc_auc"]             # keyed {figure}_{metric}
nca.plots.save_all(result, "output", metrics=["roc_auc", "f1"])

python examples/generate_output.py --metrics all

Two figures are not repeated (nca.plots.PER_METRIC lists the ones that are): cost reads complexity and fit time, which no metric touches, and metric_panels already puts every metric side by side. A per-metric copy of either would be a byte-identical file under a name claiming otherwise.

The winners never change with metrics — selection is driven by select_on alone. Asking for roc_auc shows the models the accuracy run chose, scored on roc_auc, which is the useful question: the model I picked, how does it look on my other metric?

Mind the count: 3 metrics across both themes is 46 files rather than 18.

The one thing detaching from pyplot costs is plt.show(): it does not know about these figures. nca.plots.show(fig) covers that, for a script with an interactive backend.

examples/05_figures_in_a_notebook.ipynb walks through all of it.

Colours. The identity palette has four slices (nca.plots.MAX_SERIES), found by enumerating combinations of the reference palette against a colour-blindness validator, in the all-pairs pattern: three slices give ΔE 13.0 in both themes, four give 13.0/6.9, and with five no subset passes. Above that, colour is not stretched — winner_facets gives each winner a panel and scales without limit.

For your own colours:

nca.plots.rank_bumps(
    result,
    colors={"svc-rbf": "#7b2d8e", "knn-5": "#0a7d7d"},   # per model
    palette=["#123456", "#654321"],                      # replaces the slices
)

An explicit colour is never cut by the ceiling: whoever you name, you see. Custom colours do not go through the colour-blindness validation — the responsibility passes to whoever chooses them.

Objects outside the scikit-learn API

If your transformer uses different method names, wrap it:

from nested_cv_analysis import MethodAdapter, ResamplerAdapter

nca.cross_val_score(
    SVC(), X, y, cv=5,
    transform=MethodAdapter(my_obj, fit_transform_call="fit_it", transform_call="apply_it"),
    resampler=ResamplerAdapter(other_obj, fit_resample_call="rebalance"),
)

They are two distinct classes on purpose: imblearn identifies samplers by hasattr(step, "fit_resample") and rejects a step that has both transform and fit_resample.

How it works

The package does not reimplement the fold loop. build_pipeline assembles a Pipeline (imblearn's when there is a resampler, sklearn's otherwise) with everything cloned, and delegates the rest to scikit-learn.

That is a design decision, not saved effort: because the same constructor is used for the inner selection and for the holdout refit, there is no execution path where the transformations can diverge. Leakage is not prevented by a check — it is inexpressible.

If you only need plain validation, it is worth knowing that imblearn.pipeline.Pipeline already solves that case on its own (it restricts fit_resample to the training split). This package saves the boilerplate and adds the nested CV and the analysis on top of it.

Examples

python examples/01_transforms_per_fold.py    # the leak, demonstrated
python examples/02_nested_cv.py              # selecting among models
python examples/03_regression.py             # regression and negated metrics
python examples/04_analysing_the_best.py     # reports and figures
python examples/generate_output.py [dir]     # every figure, in both themes

Plus 05_figures_in_a_notebook.ipynb, for Jupyter and Colab.

generate_output.py produces the complete set — each figure in light and dark, plus the textual reports — into a directory (output/ by default).

Development

python -m venv .venv
.venv/bin/pip install -e ".[dev]"     # Windows: .venv\Scripts\pip
.venv/bin/python -m pytest

CI runs the suite on Python 3.9–3.13.

Migrating from the old API

before	now
mycross_val_score(...)	nca.cross_val_score(...)
my_nestedcross_val(...)	nca.nested_cross_validate(...)
train_transform=	resampler=
mlmodel(model, "name")	{"name": model}
fit_transform_call= / transform_call=	MethodAdapter(obj, ...)
train_transform_call=	ResamplerAdapter(obj, ...)
show_all_scores=True	verbose=2
hide_holdout_scores=True	compute_holdout=False
score_strategy_to_sort=	select= (default nca.highest_score, correct for every sklearn scorer)
score_strategy_to_sort='nearest_to_zero_is_better'	select=nca.nearest_to_zero

The old default score_strategy_to_sort="nearest_to_zero_is_better" selected the model with the worst accuracy. Every scikit-learn scorer is "higher is better" by convention — including error metrics, which come negated (neg_mean_squared_error) — so taking the highest is always right. The old strategy is still available as select=nca.nearest_to_zero, for metrics that are genuinely centred on zero.

my_nestedcross_val returned [(name, model), ...]; the return is now a NestedCVResult. The equivalent is list(zip(r.best_model_names, r.best_estimators)).

generalization_score and generalization_std went from properties to methods, so they can take the metric argument: use r.generalization_score().

See CHANGELOG.md for the full history, including the rename from foldsafe.