CRPS reliability decomposition — crps_decomposition • modelskill

Decomposes mean ensemble CRPS into a reliability component and potential CRPS following Hersbach (2000). Predictive-distribution columns are treated as equally likely ensemble members. This function accepts predictive samples only, not a predictive mean and standard deviation. Its decomposition is defined for finite ensembles.

Usage

crps_decomposition(obs, distribution, na.rm = TRUE)

Arguments

obs: Numeric observation vector.
distribution: Numeric matrix/data frame of equally weighted predictive samples, one row per observation.
na.rm: Logical; remove incomplete observation/distribution rows?

Value

One-row data frame with crps, reliability, and potential_crps.

Details

$$\mathrm{CRPS} = \mathrm{RELI} + \mathrm{potential\ CRPS}$$

reliability (RELI) is non-negative and zero is ideal; smaller values indicate better distributional calibration. potential_crps is the remainder after removing reliability error. The decomposition is applicable here only to equally weighted predictive samples.

For each retained case $i=1,\ldots,n$, let $x_{i,1}\leq\cdots\leq x_{i,m}$ be the sorted ensemble members and $p_j=j/m$. For interior bins $j=1,\ldots,m-1$, define $$\alpha_{i,j}=\max\{0,\min(obs_i,x_{i,j+1})-x_{i,j}\}, \qquad \beta_{i,j}=\max\{0,x_{i,j+1}-\max(obs_i,x_{i,j})\}.$$ These are the portions of a bin below and above the observation. Equality at a bin edge assigns the full width to the appropriate portion; tied ensemble members have zero width. With bars denoting means over cases, $$g_j=\bar{\alpha}_j+\bar{\beta}_j,\qquad o_j=\frac{\bar{\beta}_j}{g_j}.$$ If $g_j=0$, set $o_j=p_j$; this bin contributes zero. The two exterior bins use $$o_0=\frac{1}{n}\sum_{i=1}^n I(obs_i\leq x_{i,1}),\qquad o_m=\frac{1}{n}\sum_{i=1}^n I(obs_i\leq x_{i,m}),$$ $$g_0=\frac{\frac{1}{n}\sum_{i=1}^n\max(x_{i,1}-obs_i,0)}{o_0}, \qquad g_m=\frac{\frac{1}{n}\sum_{i=1}^n\max(obs_i-x_{i,m},0)}{1-o_m}.$$ Set $g_0=0$ when $o_0=0$ and $g_m=0$ when $o_m=1$. The components, in response units, are then $$\mathrm{RELI}=\sum_{j=0}^m g_j(o_j-p_j)^2,\qquad \mathrm{potential\ CRPS}=\sum_{j=0}^m g_j o_j(1-o_j).$$

References

Hersbach, H. (2000). Decomposition of the continuous ranked probability score for ensemble prediction systems. Weather and Forecasting, 15, 559-570. doi:10.1175/1520-0434(2000)015%3C0559:DOTCRP%3E2.0.CO;2

Examples

crps_decomposition(c(0, 1), matrix(c(-1, 1, 0, 2), nrow = 2))
#>   crps reliability potential_crps
#> 1 0.25           0           0.25