What is Continuous Ranked Probability Score (CRPS)?

Question

I came across some evolution metric at Kaggle:

Continuous Ranked Probability Score (CRPS):

Mathematically,

$C = \frac{1}{199N} \sum_{m=1}^{N} \sum_{n=-99}^{99} (P(y \le n) -H(n - Y_m))^2,$

where P is the predicted distribution, N is the number of plays in the test set, Y is the actual yardage and H(x) is the Heaviside step function $(H(x)=1$ for $x≥0$ and zero otherwise).

Can anyone be kind enough to explain this evolution metric and implementation with Python?

Answer 1

CRPS is in a sense just the mean square error (MSE) of your predicted cumulative density function (CDF) and the true CDF.

The CRPS generalizes the MAE (Mean Absolute Error) to the case of probabilistic forecasts. The CPRS is one of the most widely used accuracy metrics where probabilistic forecasts are involved.

The CRPS is frequently used in order to assess the respective accuracy of two probabilistic forecasting models. In particular, this metric can be combined with a backtesting process in order to stabilize the accuracy assessment by leveraging multiple measurements over the same dataset.

This metric notably differs from simpler metrics such as MAE because of its asymmetric expression: while the forecasts are probabilistic, the observations are deterministic. Unlike the pinball loss function, the CPRS does not focus on any specific point of the probability distribution but considers the distribution of the forecasts as a whole.

Find More Info

Paper

Five ways to calculate CRPS

Using python library properscoring 0.1

pip install properscoring

This library focuses on the closely related Continuous Ranked Probability Score (CRPS) and Brier Score. We like these scores because they are both interpretable (e.g., CRPS is a generalization of mean absolute error) and easily calculated from a finite number of samples of a probability distribution.

We will illustrate how to calculate CRPS against a forecast given by a Gaussian random variable. To begin, import properscoring:

import numpy as np
import properscoring as ps
from scipy.stats import norm

• Exact calculation using crps_gaussian (this is the fastest method):

ps.crps_gaussian(0, mu=0, sig=1)

• Numerical integration with crps_quadrature:

ps.crps_quadrature(0, norm)

• From a finite sample with crps_ensemble:

ensemble = np.random.RandomState(0).randn(1000)
ps.crps_ensemble(0, ensemble)

• Weighted by PDF values with crps_ensemble:

x = np.linspace(-5, 5, num=1000)
ps.crps_ensemble(0, x, weights=norm.pdf(x))

• Based on the threshold decomposition of CRPS with threshold_brier_score:

threshold_scores = ps.threshold_brier_score(0, ensemble, threshold=x)
(x[1] - x[0]) * threshold_scores.sum(axis=-1)

More at