mixle.inference.calibration module

Calibration diagnostics: “is my probability / interval actually calibrated?”

A forecast is calibrated when its stated probabilities match observed frequencies: events given 70% should happen ~70% of the time, and 90% intervals should contain the truth ~90% of the time. This is distinct from sharpness (how concentrated the forecast is) and from accuracy – a forecaster can be perfectly calibrated while carrying little decision value (for example, always predicting the base rate), so calibration is a necessary, not sufficient, condition that you check separately. These diagnostics are model-free: they look only at predicted probabilities/intervals and what happened.

Three families, by forecast type:

  • Probability classifiers – reliability_curve() (the reliability diagram), and the expected_calibration_error() / maximum_calibration_error() summaries. For multiclass use top_label_confidence() to reduce to the (confidence, correct) calibration problem.

  • Full predictive distributions – the Probability Integral Transform: pit_values() / pit_ensemble(), the pit_histogram(), and pit_calibration_error(). Under a calibrated forecast the PIT values are Uniform(0, 1); a U-shaped histogram means under-dispersion, a hump means over-dispersion, a slope means bias.

  • Intervals / quantiles – interval_coverage() (coverage and mean width at one level) and coverage_curve() (empirical-vs-nominal coverage across a grid of levels).

Several functions take ci=True to attach nonparametric bootstrap intervals. These are PER-BIN POINTWISE intervals, not a simultaneous band: each bin’s interval covers that bin’s own frequency at the stated level, but the chance that at least one of bins intervals excludes its truth is far higher (measured 30-45% familywise exclusion at a nominal 5% per bin on reliability diagrams) – so one bin poking outside its interval is expected under perfect calibration, not evidence against it. The binned summaries also have a NULL FLOOR: with finitely many outcomes per bin, ECE and MCE are strictly positive even for a perfectly calibrated forecaster; compare them to calibration_null_expectation(), not to zero.

reliability_curve(prob, outcome, *, bins=10, strategy='uniform', ci=False, n_boot=1000, ci_level=0.95, seed=0)[source]

Reliability diagram: observed frequency vs mean forecast probability, per bin.

Bins the forecasts, then within each bin compares the mean predicted probability to the observed event frequency. A perfectly calibrated forecaster lies on the diagonal observed == predicted.

Parameters:

  • prob (ndarray) – (n,) predicted probabilities of the positive class (or top-label confidences).

  • outcome (ndarray) – (n,) 0/1 outcomes (or correctness indicators).

  • bins (int) – number of bins.

  • strategy (str) – "uniform" (equal-width) or "quantile" (equal-count) bins.

  • ci (bool) – if True attach a percentile bootstrap interval on the observed frequency in each bin.

  • n_boot (int) – bootstrap resamples when ci is True.

  • ci_level (float) – central probability of each PER-BIN interval (e.g. 0.95).

  • seed (int | RandomState | None) – RNG seed for the bootstrap.

Returns:

{'mean_pred', 'obs_freq', 'count', 'bin_edges'} (one entry per non-empty bin), plus 'obs_lo' / 'obs_hi' / 'obs_ci_effective_boot' when ci is True.

Return type:

dict[str, ndarray]

The intervals are POINTWISE, one bin at a time – with 10 bins at 95% each, at least one bin escaping its interval under perfect calibration is a ~30-45% event, so judge single-bin excursions accordingly. A resample in which a bin comes up EMPTY contributes nothing to that bin’s interval; the quantiles are taken over the non-empty resamples only, which conditions the interval on the bin being observed (a real bias for sparse bins). obs_ci_effective_boot reports, per bin, how many of the n_boot resamples actually informed the interval – treat an interval resting on far fewer than n_boot replicates as unstable rather than narrow.

expected_calibration_error(prob, outcome, *, bins=10, strategy='uniform', norm='l1', ci=False, n_boot=1000, ci_level=0.95, seed=0)[source]

Expected Calibration Error: count-weighted average gap between confidence and accuracy.

ECE = sum_b (n_b / n) |obs_b - pred_b| over bins (norm='l2' uses the squared gap, square- rooted). For multiclass classifiers reduce with top_label_confidence() first.

Zero is NOT the null value. The plug-in ECE takes an absolute value of binomial noise in every bin, so a PERFECTLY calibrated forecaster still scores roughly sqrt(bins / n) in expectation – and a percentile bootstrap interval around the plug-in estimate inherits that upward shift, so under perfect calibration it excludes the true value 0 essentially always. Judge the point estimate against calibration_null_expectation() (the measured perfect-calibration distribution of this same statistic at your prob profile and binning), and read the bootstrap interval as sampling uncertainty around the BIASED plug-in functional, not as a test of calibration.

Parameters:

  • prob (ndarray) – (n,) predicted probabilities / confidences.

  • outcome (ndarray) – (n,) 0/1 outcomes / correctness indicators.

  • bins (int) – number of bins.

  • strategy (str) – "uniform" or "quantile" binning.

  • norm (str) – "l1" (mean absolute gap) or "l2" (root mean squared gap).

  • ci (bool) – if True also return a percentile bootstrap interval.

  • n_boot (int) – bootstrap controls.

  • ci_level (float) – bootstrap controls.

  • seed (int | RandomState | None) – bootstrap controls.

Returns:

The ECE (float), or (ece, lo, hi) when ci is True.

Return type:

float | tuple[float, float, float]

maximum_calibration_error(prob, outcome, *, bins=10, strategy='uniform')[source]

Maximum Calibration Error: the worst per-bin gap max_b |obs_b - pred_b|.

Unlike expected_calibration_error() this is not count-weighted, so it surfaces a small but badly-miscalibrated region that the average would hide – and for exactly that reason its null value GROWS with bins: a max over more (and therefore sparser) bins of pure binomial noise rises from ~0.05 to ~0.6 as bins goes 10 -> 500 under perfect calibration, with the argmax typically landing on the sparsest bin. A large MCE from a thin bin is what noise looks like at that binning; compare to calibration_null_expectation() before reading it as a defect, and prefer more data per bin (fewer bins, or quantile binning) when the max matters.

Parameters:
Return type:

float

calibration_null_expectation(prob, *, bins=10, strategy='uniform', norm='l1', n_sim=500, seed=0)[source]

What ECE and MCE look like for a PERFECTLY calibrated forecaster at this prob profile.

Simulates outcomes y*_i ~ Bernoulli(prob_i) – the definition of perfect calibration at the observed forecasts – and recomputes the binned statistics each time, with the same bins the real statistics would use. The plug-in ECE/MCE are strictly positive under this null (absolute values of binomial noise), so THESE numbers, not zero, are the calibrated baseline: an observed ECE near ece here is consistent with perfect calibration, and evidence of miscalibration starts around the *_q95 quantiles, which give a one-sided 5% test by construction.

Parameters:

  • prob (ndarray) – (n,) the forecaster’s predicted probabilities (outcomes are not needed – the null generates them).

  • bins (int) – the same binning controls you pass to the real statistics.

  • strategy (str) – the same binning controls you pass to the real statistics.

  • norm (str) – the same binning controls you pass to the real statistics.

  • n_sim (int) – null simulations; at least 20. The *_q95 fields are the CONSERVATIVE order statistic X_(ceil(0.95 (n_sim+1))), calibrated for the STRICT comparison observed > q95: for a continuous statistic its fresh-null strict exceedance is exactly 1 - k/(n_sim + 1) <= 5% at every accepted n_sim (4.76% at the floor of 20). With DISCRETE ties the two operators split – a single forecast at p = 0.5 makes every null ECE equal the threshold, so strict exceedance is 0 while >= reads 1 (STAT-RR22-10) – so compare with > and read the bound as <= (ties only make it more conservative). The interpolated quantile realized 9.83% at the floor, and n_sim=1 once labeled a single draw a “q95” with 49.8% exceedance (STAT-RR19-10/STAT-RR21-09). The default 500 keeps the threshold’s own MC noise small.

  • seed (int | RandomState | None) – RNG seed.

Returns:

{'ece', 'ece_q95', 'mce', 'mce_q95', 'n_sim'} – null means and 95th percentiles of the plug-in ECE (under norm) and MCE at this forecast profile and binning.

Return type:

dict[str, float]

top_label_confidence(prob, labels)[source]

Reduce a multiclass classifier to the (confidence, correct) top-label calibration problem.

This checks calibration OF THE ARGMAX CLASS ONLY – a deliberately narrow reduction, not THE multiclass calibration check. Probabilities the model puts on classes it never predicts are invisible to it: a forecaster that always assigns 0.4 to a class that never occurs (while its argmax confidences are accurate) passes top-label calibration with ECE 0. Classwise or full- vector calibration is a strictly stronger property and is not implemented here; when the non-argmax probabilities carry decisions (thresholded alerts on a minority class, expected-cost rankings), check those probabilities against their own outcomes directly.

Parameters:

  • prob (ndarray) – (n, K) class-probability matrix.

  • labels (ndarray) – (n,) integer true labels.

Returns:

the max predicted probability per row and a 0/1 indicator of whether that argmax class was the true label. Feed these to reliability_curve() / expected_calibration_error().

Return type:

(confidence, correct)

pit_values(y, cdf)[source]

Probability Integral Transform values u_i = F_i(y_i).

Under a calibrated continuous predictive distribution the PIT values are Uniform(0, 1). Pass either the precomputed CDF values F_i(y_i) or a callable cdf(y) -> F(y); a scalar-only callable such as a fitted univariate distribution’s .cdf is evaluated element-wise (mixtures and other structured families define no cdf; pass their CDF as your own callable).

Parameters:
Returns:

(n,) PIT values in [0, 1]. Invalid CDF values raise rather than being clipped.

Return type:

ndarray

pit_ensemble(y, forecasts, *, randomize=True, seed=0)[source]

Rank-based PIT from a finite predictive ensemble.

u_i is the fraction of ensemble members <= y_i. With randomize=True ties are broken by a uniform jitter within the rank gap, which makes the PIT exactly Uniform(0, 1) under calibration even for discrete ensembles (the randomized PIT of Czado et al. 2009).

Parameters:

  • y (ndarray) – (n,) realised values.

  • forecasts (ndarray) – (n, m) ensemble (m draws per observation).

  • randomize (bool) – jitter ties for an exactly-uniform PIT.

  • seed (int | RandomState | None) – RNG seed when randomize is True.

Returns:

(n,) PIT values in [0, 1].

Return type:

ndarray

pit_histogram(pit, *, bins=10)[source]

Histogram of PIT values with the uniform reference level.

Parameters:

  • pit (ndarray) – (n,) PIT values.

  • bins (int) – number of equal-width bins on [0, 1].

Returns:

{'counts', 'density', 'edges', 'uniform'} where density integrates to 1 and uniform is the flat reference density (1.0) a calibrated forecast would match.

Return type:

dict[str, ndarray]

pit_calibration_error(pit, *, bins=10)[source]

Calibration error of a PIT histogram: TOTAL absolute deviation from uniform mass.

sum_b |count_b/n - 1/bins| – the SUM over bins (twice the total-variation distance to uniform, range [0, 2(1 - 1/bins)]), not a per-bin mean: a threshold sized for a mean would be off by a factor of bins. 0 when the PIT histogram is perfectly flat (calibrated), larger when it is U-shaped (under-dispersed) or humped (over-dispersed).

Parameters:
Return type:

float

interval_coverage(lower, upper, y)[source]

Empirical coverage and mean width of a set of prediction intervals.

Parameters:

  • lower (ndarray) – (n,) lower endpoints.

  • upper (ndarray) – (n,) upper endpoints.

  • y (ndarray) – (n,) realised values.

Returns:

{'coverage', 'mean_width'} – the fraction of y inside [lower, upper] and the mean interval width. Compare coverage to the nominal level the interval was built for.

Return type:

dict[str, float]

coverage_curve(forecasts, y, *, levels=None)[source]

Empirical-vs-nominal coverage of central intervals across a grid of nominal levels.

For each nominal central level c the per-observation central interval [quantile((1-c)/2), quantile((1+c)/2)] is read off the predictive ensemble and its empirical coverage is measured.

Judge empirical against null_expectation, NOT the diagonal. Plug-in quantiles of a finite ensemble cover BELOW nominal even when the ensemble is perfect: with m = 5 members, a perfectly calibrated forecast’s 95% interval covers only ~0.64 – the interval simply cannot reach past the sample extremes. null_expectation is the expected coverage of a perfectly calibrated m-member ensemble at each level, computed for a uniform reference with the same quantile rule. It is almost – not exactly – distribution-free: linear quantile interpolation leaves a small dependence on the ensemble’s shape at very small m (measured 0.632 / 0.644 / 0.636 for uniform / Gaussian / exponential at m = 5 and nominal 0.95, converging by m ~ 50), which is an order of magnitude smaller than the diagonal gap it corrects. empirical below the DIAGONAL at small m is expected; empirical meaningfully below null_expectation indicates over-confidence.

Parameters:

  • forecasts (ndarray) – (n, m) predictive ensemble.

  • y (ndarray) – (n,) realised values.

  • levels (ndarray | None) – nominal central coverage levels in (0, 1); defaults to 0.05 .. 0.95 by 0.05.

Returns:

{'nominal', 'empirical', 'null_expectation'} arrays of equal length – null_expectation is the finite-m perfect-calibration benchmark per level.

Return type:

dict[str, ndarray]

class ProbabilityCalibrator(method='isotonic')[source]

Bases: object

Map raw scores to calibrated probabilities – fit against binary outcomes.

A raw score (a model’s confidence, a self-consistency fraction, a token likelihood) need not be a probability of anything: it can be monotone-but-miscalibrated, or have no relationship to the outcome at all. This learns the transform score -> P(outcome = 1 | score) from labeled data, so the output is a probability of the event you calibrated against.

  • method="isotonic" – monotone, non-parametric (pool-adjacent-violators). Assumes higher score => not-lower probability; flexible, needs enough calibration points.

  • method="platt" – logistic sigmoid(a * score + b). Two parameters, usable on little data, but assumes a sigmoidal relationship. Fit against Platt’s smoothed targets (N+ + 1)/(N+ + 2) and 1/(N- + 2) rather than raw 0/1: with raw targets, perfectly separable calibration data drives a to infinity and the “calibrated” outputs to exact 0/1 – maximal overconfidence from the method that was supposed to remove it. The smoothing keeps the map finite and the outputs strictly inside (0, 1).

The guarantee is conditional on HELD-OUT evaluation. The isotonic fit interpolates its own calibration outcomes, so re-checking calibration on the SAME data it was fit on passes by construction (a reliability curve of the fitted values against the fitting outcomes sits on the diagonal no matter how little signal the score carried); any honest check of a fitted calibrator – isotonic especially – uses data it never saw. Both methods assume the score/outcome relationship is stable between fitting and serving.

A near-flat fitted curve is itself the finding: it means the raw score carried little information about the outcome (its “likelihood” was unrelated to the event).

Parameters:

method (str)

fit(scores, outcomes)[source]

Fit the score->probability map on scores with binary outcomes (0/1).

Parameters:
Return type:

ProbabilityCalibrator

predict(scores)[source]

Calibrated probabilities for scores (clamped to [0, 1]).

Parameters:

scores (Any)

Return type:

ndarray

calibrate_probabilities(scores, outcomes, *, method='isotonic')[source]

Fit a ProbabilityCalibrator mapping scores to P(outcome=1 | score).

Parameters:
Return type:

ProbabilityCalibrator