Compute weighted interval score — weighted_interval

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 when quantile_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 corresponding value 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 be NA if any missing values exist in the original quantile_levels. Finally, if quantile_levels is specified, "fail" will result in the score being NA when any required quantile levels (implicit or explicit) are do not have corresponding values.

Value

a vector of nonnegative scores.

Methods (by class)

weighted_interval_score(dist_quantiles): Weighted interval score with dist_quantiles allows for different NA 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 <- case_death_rate_subset %>%
  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 <- case_death_rate_subset %>%
  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