mixle.inference.scoring module

Proper scoring rules for probabilistic forecasts.

A scoring rule S(F, y) measures how well a predictive distribution F matched the value y that actually occurred; it is proper when the forecaster minimises the expected score by reporting their true belief (and strictly proper when that minimiser is unique). These give a single, fair currency for comparing probabilistic forecasts and interval methods at matched coverage – a 90% interval that is wide everywhere and a 90% interval that is tight where it can be both cover 90% of the time, but the second wins on the interval score.

Lower is better for every score here (they are penalties / losses):

  • log_score() – the log loss -log p(y) (local, strictly proper).

  • brier_score() / brier_decomposition() – squared-error score for categorical forecasts, with the reliability/resolution/uncertainty decomposition.

  • crps_ensemble() / crps_gaussian() – the Continuous Ranked Probability Score, the CDF-space analogue of absolute error; reduces to absolute error for a point forecast.

  • interval_score() (a.k.a. winkler_score()) – the proper score for a central interval forecast: width plus an out-of-bounds penalty scaled by the miscoverage level.

  • pinball_loss() – the check loss whose expectation is minimised by the true quantile.

  • energy_score() – the multivariate generalisation of CRPS.

  • skill_score() – a score expressed as fractional improvement over a reference.

All functions are pure NumPy/SciPy and operate on plain arrays; by default they return the mean score over observations, with mean=False returning the per-observation vector for paired comparisons (see mixle.inference.resampling() for CIs on score differences).

log_score(prob, *, mean=True)[source]

Logarithmic score (log loss) -log p(y) of the probability assigned to the outcome.

The log score is local (it depends only on the density at the realised value) and strictly proper. Pass the predictive probability/density evaluated at the observed outcome.

Parameters:

  • prob (ndarray) – array of predictive probabilities (or densities) at the realised outcomes. Values are clipped away from zero before the log so a single zero-probability event yields a large but finite penalty rather than inf.

  • mean (bool) – if True (default) return the mean log loss; otherwise the per-observation vector.

Returns:

Mean log loss (float) or the per-observation array.

Return type:

ndarray | float

brier_score(prob, outcome, *, mean=True)[source]

Brier score (mean squared error of probabilistic classification).

For binary forecasts pass 1-D prob (predicted probability of the positive class) and 0/1 outcome. For K-class forecasts pass prob shaped (n, K) and outcome either as integer class labels (length n) or as a one-hot (n, K) matrix; the multiclass score sums the squared error across classes, so it ranges in [0, 2].

Parameters:

  • prob (ndarray) – (n,) positive-class probabilities, or (n, K) class probabilities.

  • outcome (ndarray) – (n,) 0/1 or integer labels, or (n, K) one-hot.

  • mean (bool) – if True return the mean; otherwise the per-observation vector.

Returns:

Mean Brier score (float) or the per-observation array.

Return type:

ndarray | float

brier_decomposition(prob, outcome, *, bins=10)[source]

Murphy’s reliability–resolution–uncertainty decomposition of the (binary) Brier score.

Bins the forecast probabilities into bins equal-width bins on [0, 1] and computes

Brier = reliability - resolution + uncertainty,

where reliability (lower is better) is the mean squared gap between forecast probability and observed frequency within a bin, resolution (higher is better) is how far bin frequencies move from the base rate, and uncertainty is the base-rate variance p_bar (1 - p_bar) (a property of the data, not the forecaster).

Parameters:

  • prob (ndarray) – (n,) predicted probabilities of the positive class.

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

  • bins (int) – number of equal-width probability bins.

Returns:

{'reliability', 'resolution', 'uncertainty', 'brier'}. The identity reliability - resolution + uncertainty == brier holds up to binning of the score’s in-bin variance term.

Return type:

dict[str, float]

crps_ensemble(forecasts, y, *, fair=False, mean=True)[source]

Continuous Ranked Probability Score from a finite predictive ensemble (sample).

Uses the energy form CRPS = E|X - y| - 1/2 E|X - X'| estimated from the ensemble draws. The CRPS generalises absolute error to distributional forecasts (a point forecast recovers |x-y|) and is reported in the same units as y.

Parameters:

  • forecasts (ndarray) – predictive draws. Either (n, m) (m draws for each of n observations) or (m,) for a single observation. Ragged ensembles are not supported – pad or call per observation.

  • y (ndarray) – (n,) realised values (or a scalar for the (m,) case).

  • fair (bool) – if True use the unbiased 1/(m(m-1)) spread estimator (the “fair”/almost-unbiased CRPS); if False (default) the standard 1/m^2 estimator.

  • mean (bool) – if True return the mean CRPS; otherwise the per-observation vector.

Returns:

Mean CRPS (float) or the per-observation array.

Return type:

ndarray | float

crps_gaussian(mu, sigma, y, *, mean=True)[source]

Closed-form CRPS for a Gaussian predictive distribution N(mu, sigma^2).

CRPS = sigma * [ z (2 Phi(z) - 1) + 2 phi(z) - 1/sqrt(pi) ] with z = (y - mu) / sigma (Gneiting & Raftery 2007). Exact, so it is the right reference when the forecast is Gaussian.

Parameters:

  • mu (ndarray) – predictive means, broadcastable with y.

  • sigma (ndarray) – predictive standard deviations (> 0), broadcastable with y.

  • y (ndarray) – realised values.

  • mean (bool) – if True return the mean; otherwise the per-observation vector.

Returns:

Mean CRPS (float) or the per-observation array.

Return type:

ndarray | float

interval_score(lower, upper, y, alpha, *, mean=True)[source]

Winkler interval score for a central (1 - alpha) prediction interval.

IS = (u - l) + (2/alpha)(l - y) 1[y < l] + (2/alpha)(y - u) 1[y > u] (Gneiting & Raftery 2007). It rewards narrow intervals but penalises misses by 2/alpha times the shortfall, so at matched nominal coverage the tighter, better-placed interval scores lower. This is the proper way to rank interval methods (conformal, quantile-regression, GP credible bands) head to head.

Parameters:

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

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

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

  • alpha (float) – miscoverage level (e.g. 0.1 for a 90% interval), in (0, 1).

  • mean (bool) – if True return the mean; otherwise the per-observation vector.

Returns:

Mean interval score (float) or the per-observation array.

Return type:

ndarray | float

winkler_score(lower, upper, y, alpha, *, mean=True)[source]

Alias for interval_score() (the score is due to Winkler 1972).

Parameters:
Return type:

ndarray | float

pinball_loss(pred, y, tau, *, mean=True)[source]

Pinball (quantile / check) loss for a quantile forecast.

rho_tau(u) = max(tau * u, (tau - 1) * u) with u = y - pred; its expectation is minimised by the true tau-quantile, so it is the proper score for quantile forecasts and the loss that quantile regression and conformalised quantiles optimise.

Parameters:

  • pred (ndarray) – predicted tau-quantiles, shape (n,) (single level) or (n, q) (several levels, one per column of tau).

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

  • tau (float | ndarray) – quantile level(s) in (0, 1) – a scalar for (n,) pred, or a length-q array matching the columns of a (n, q) pred.

  • mean (bool) – if True return the mean over all entries; otherwise the per-observation (or per-observation-per-level) array.

Returns:

Mean pinball loss (float) or the per-observation array.

Return type:

ndarray | float

energy_score(forecasts, y, *, fair=False, mean=True)[source]

Energy score: the multivariate generalisation of CRPS for vector-valued forecasts.

ES = E||X - y|| - 1/2 E||X - X'|| with Euclidean norms; for scalar y it equals crps_ensemble().

Parameters:

  • forecasts (ndarray) – predictive draws shaped (n, m, d) (m d-vectors per observation) or (m, d) for a single observation.

  • y (ndarray) – (n, d) realised vectors (or (d,) for the single-observation case).

  • fair (bool) – use the unbiased 1/(m(m-1)) spread estimator if True.

  • mean (bool) – if True return the mean; otherwise the per-observation vector.

Returns:

Mean energy score (float) or the per-observation array.

Return type:

ndarray | float

skill_score(score, reference, *, perfect=0.0)[source]

Skill score: fractional improvement of score over a reference score.

skill = (reference - score) / (reference - perfect) – 1 means the forecast hit the perfect score, 0 means it only matched the reference (e.g. climatology), and negative means it did worse than the reference. Works for any negatively-oriented score (lower is better) such as those above.

Parameters:

  • score (float) – the forecast’s score (lower is better).

  • reference (float) – the baseline/reference forecast’s score.

  • perfect (float) – the score of a perfect forecast (0.0 for the rules in this module).

Returns:

The skill score (float); nan if the reference already equals the perfect score.

Return type:

float