Compute Continuous Ranked Probability Score for forecast

Given a forecast and set of observations, compute the

Usage

crps(fcst, obs, at = NULL, after = NULL, summarize = TRUE)

Arguments

fcst: A forecast object (see output of create_forecast()).
obs: An observations data frame.
at: (Optional) A time (must be compatible with fcst and obs). If specified, the score for this time point will be returned. Mutually exclusive with after.
after: (Optional) A number. If specified, the score at time fcst$forecast_time + after will be returned. Mutually exclusive with at.
summarize: A boolean, defaults to TRUE. If TRUE, a single number will be returned as the score for the forecast. If FALSE, a data frame with columns named time, val_obs, and score will be returned, containing the scores for each individual time point. This can be used by plotting functions to colour-code observations, for example.

Value

If summarize is FALSE, a data frame containing times, observations, and scores. Otherwise, the score at the time speficied by either at or after.

Details

Continuous Ranked Probability Score (CRPS) for every time point. The CRPS for a given distribution f and observation y, let F be the CDF of f. Then the CRPS is the integral of the square of F(x) - H(x - y), where H is the Heaviside function.

The CDF is obtained using an emperical distribution function on the forecast data.

Examples

fc <- create_forecast(list(
  time=c(1,2,3),
  vals=list(c(1,1,1), c(2,2,2), c(3,3,3), c(4,4,4), c(5,5,5)), forecast_time=1)
)

crps(fc, data.frame(time=1:3, val_obs=c(3,4,5)), summarize=FALSE)
#> # A tibble: 3 × 3
#>    time val_obs score
#>   <dbl>   <dbl> <dbl>
#> 1     1       3   0.4
#> 2     2       4   0.6
#> 3     3       5   1.2

crps(fc, data.frame(time=1:3, val_obs=c(3,4,5)), at=2)
#> [1] 0.6