Weighted interval score (WIS), a well-known quantile-based approximation of the commonly-used continuous ranked probability score (CRPS). WIS is a proper score, and can be thought of as a distributional generalization of absolute error. For example, see Bracher et al. (2020) for discussion in the context of COVID-19 forecasting.
Usage
weighted_interval_score(x, actual, quantile_levels = NULL, ...)
# S3 method for class 'dist_quantiles'
weighted_interval_score(
x,
actual,
quantile_levels = NULL,
na_handling = c("impute", "drop", "propagate", "fail"),
...
)
Arguments
- x
distribution. A vector of class distribution. Ideally, this vector contains
dist_quantiles()
, though other distributions are supported whenquantile_levels
is specified. See below.- actual
double. Actual value(s)
- quantile_levels
probabilities. If specified, the score will be computed at this set of levels.
- ...
not used
- na_handling
character. Determines how
quantile_levels
without a correspondingvalue
are handled. For"impute"
, missing values will be calculated if possible using the available quantiles. For"drop"
, explicitly missing values are ignored in the calculation of the score, but implicitly missing values are imputed if possible. For"propogate"
, the resulting score will beNA
if any missing values exist in the originalquantile_levels
. Finally, ifquantile_levels
is specified,"fail"
will result in the score beingNA
when any required quantile levels (implicit or explicit) are do not have corresponding values.
Methods (by class)
weighted_interval_score(dist_quantiles)
: Weighted interval score withdist_quantiles
allows for differentNA
behaviours.
Examples
quantile_levels <- c(.2, .4, .6, .8)
predq_1 <- 1:4 #
predq_2 <- 8:11
dstn <- dist_quantiles(list(predq_1, predq_2), quantile_levels)
actual <- c(3.3, 7.1)
weighted_interval_score(dstn, actual)
#> [1] 0.65 1.90
weighted_interval_score(dstn, actual, c(.25, .5, .75))
#> [1] 0.6833333 1.9833333
library(distributional)
dstn <- dist_normal(c(.75, 2))
weighted_interval_score(dstn, 1, c(.25, .5, .75))
#> [1] 0.3081633 0.7751701
# Missing value behaviours
dstn <- dist_quantiles(c(1, 2, NA, 4), 1:4 / 5)
weighted_interval_score(dstn, 2.5)
#> [1] 0.5
weighted_interval_score(dstn, 2.5, 1:9 / 10)
#> [1] 0.455656
weighted_interval_score(dstn, 2.5, 1:9 / 10, na_handling = "drop")
#> [1] 0.462613
weighted_interval_score(dstn, 2.5, na_handling = "propagate")
#> [1] NA
weighted_interval_score(dist_quantiles(1:4, 1:4 / 5), 2.5, 1:9 / 10,
na_handling = "fail"
)
#> [1] NA
# Using some actual forecasts --------
library(dplyr)
jhu <- covid_case_death_rates %>%
filter(time_value >= "2021-10-01", time_value <= "2021-12-01")
preds <- flatline_forecaster(
jhu, "death_rate",
flatline_args_list(quantile_levels = c(.01, .025, 1:19 / 20, .975, .99))
)$predictions
actuals <- covid_case_death_rates %>%
filter(time_value == as.Date("2021-12-01") + 7) %>%
select(geo_value, time_value, actual = death_rate)
preds <- left_join(preds, actuals,
by = c("target_date" = "time_value", "geo_value")
) %>%
mutate(wis = weighted_interval_score(.pred_distn, actual))
preds
#> # A tibble: 56 × 7
#> geo_value .pred .pred_distn forecast_date target_date actual wis
#> <chr> <dbl> <dist> <date> <date> <dbl> <dbl>
#> 1 ak 0.217 quantiles(0.22)[23] 2021-12-01 2021-12-08 0.0988 0.0673
#> 2 al 0.119 quantiles(0.12)[23] 2021-12-01 2021-12-08 0.174 0.0364
#> 3 ar 0.207 quantiles(0.21)[23] 2021-12-01 2021-12-08 0.514 0.196
#> 4 as 0 quantiles(0)[23] 2021-12-01 2021-12-08 0 0.0146
#> 5 az 0.485 quantiles(0.49)[23] 2021-12-01 2021-12-08 0.826 0.223
#> 6 ca 0.173 quantiles(0.17)[23] 2021-12-01 2021-12-08 0.183 0.0272
#> 7 co 0.521 quantiles(0.52)[23] 2021-12-01 2021-12-08 0.539 0.0302
#> 8 ct 0.177 quantiles(0.18)[23] 2021-12-01 2021-12-08 0.149 0.0302
#> 9 dc 0 quantiles(0)[23] 2021-12-01 2021-12-08 0.0200 0.0166
#> 10 de 0.217 quantiles(0.22)[23] 2021-12-01 2021-12-08 0.391 0.101
#> # ℹ 46 more rows