Merge pull request #525 from epiforecasts/drop-interval-functions

Issue #443 Drop interval functions

NAMESPACE

@@ -23,7 +23,6 @@ export(avail_forecasts)
 export(available_forecasts)
 export(available_metrics)
 export(bias_quantile)
-export(bias_range)
 export(bias_sample)
 export(brier_score)
 export(correlation)
@@ -36,7 +35,6 @@ export(get_forecast_unit)
 export(interval_coverage_dev_quantile)
 export(interval_coverage_quantile)
 export(interval_coverage_sample)
-export(interval_score)
 export(log_shift)
 export(logs_binary)
 export(logs_sample)

NEWS.md

@@ -23,6 +23,9 @@ The update introduces breaking changes. If you want to keep using the older vers
 - `check_forecasts()` was replaced by a new function `validate()`. `validate()` validates the input and in that sense fulfills the purpose of `check_forecasts()`. It has different methods: `validate.default()` assigns the input a class based on their forecast type. Other methods validate the input specifically for the various forecast types.
 - The functionality for computing pairwise comparisons was now split from `summarise_scores()`. Instead of doing pairwise comparisons as part of summarising scores, a new function, `add_pairwise_comparison()`, was introduced that takes summarised scores as an input and adds pairwise comparisons to it. 
 - `add_coverage()` was reworked completely. It's new purpose is now to add coverage information to the raw forecast data (essentially fulfilling some of the functionality that was previously covered by `score_quantile()`)
+- Support for the interval format was mostly dropped (see PR #525 by @nikosbosse and reviewed by @seabbs)
+    - The function `bias_range()` was removed (users should now use `bias_quantile()` instead)
+    - The function `interval_score()` was made an internal function rather than being exported to users. We recommend using `wis()` instead. 
 - The function `find_duplicates()` was renamed to `get_duplicate_forecasts()`
 - Changes to `avail_forecasts()` and `plot_avail_forecasts()`:
   - The function `avail_forecasts()` was renamed to `available_forecasts()` for consistency with `available_metrics()`. The old function, `avail_forecasts()` is still available as an alias, but will be removed in the future.

R/metrics-quantile.R

@@ -97,6 +97,15 @@
 #' or a list with separate entries if `separate_results` is `TRUE`.
 #' @export
 #' @keywords metric
+#' @examples
+#' observed <- c(1, -15, 22)
+#' predicted <- rbind(
+#'   c(-1, 0, 1, 2, 3),
+#'   c(-2, 1, 2, 2, 4),
+#'   c(-2, 0, 3, 3, 4)
+#' )
+#' quantile <- c(0.1, 0.25, 0.5, 0.75, 0.9)
+#' wis(observed, predicted, quantile)
 wis <- function(observed,
                 predicted,
                 quantile,
@@ -362,12 +371,12 @@ interval_coverage_dev_quantile <- function(observed, predicted, quantile) {
 #'
 #' Bias can assume values between
 #' -1 and 1 and is 0 ideally (i.e. unbiased).
+#' @param observed a single number representing the observed value
 #' @param predicted vector of length corresponding to the number of quantiles
 #' that holds predictions
 #' @param quantile vector of corresponding size with the quantile levels for
 #' which predictions were made. If this does not contain the median (0.5) then
 #' the median is imputed as being the mean of the two innermost quantiles.
-#' @inheritParams bias_range
 #' @param na.rm logical. Should missing values be removed?
 #' @return scalar with the quantile bias for a single quantile prediction
 #' @export
@@ -522,7 +531,7 @@ ae_median_quantile <- function(observed, predicted, quantile) {
 #' @description
 #' Proper Scoring Rule to score quantile predictions. Smaller values are better.
 #' The quantile score is
-#' closely related to the Interval score (see [interval_score()]) and is
+#' closely related to the Interval score (see [wis()]) and is
 #' the quantile equivalent that works with single quantiles instead of
 #' central prediction intervals.
 #'

R/metrics-range.R

@@ -65,25 +65,26 @@
 #' lower <- qnorm(alpha / 2, rnorm(30, mean = 1:30))
 #' upper <- qnorm((1 - alpha / 2), rnorm(30, mean = 11:40))
 #'
-#' interval_score(
+#' scoringutils:::interval_score(
 #'   observed = observed,
 #'   lower = lower,
 #'   upper = upper,
 #'   interval_range = interval_range
 #' )
 #'
 #' # gives a warning, as the interval_range should likely be 50 instead of 0.5
-#' interval_score(observed = 4, upper = 8, lower = 2, interval_range = 0.5)
+#' scoringutils:::interval_score(
+#'   observed = 4, upper = 8, lower = 2, interval_range = 0.5
+#' )
 #'
 #' # example with missing values and separate results
-#' interval_score(
+#' scoringutils:::interval_score(
 #'   observed = c(observed, NA),
 #'   lower = c(lower, NA),
 #'   upper = c(NA, upper),
 #'   separate_results = TRUE,
 #'   interval_range = 90
 #' )
-#' @export
 #' @keywords metric
 #' @references Strictly Proper Scoring Rules, Prediction,and Estimation,
 #' Tilmann Gneiting and Adrian E. Raftery, 2007, Journal of the American
@@ -132,126 +133,3 @@ interval_score <- function(observed,
     return(score)
   }
 }
-
-
-################################################################################
-# Metrics with a many-to-one relationship between input and score
-################################################################################
-
-#' @title Determines Bias of Quantile Forecasts based on the range of the
-#' prediction intervals
-#'
-#' @description
-#' Determines bias from quantile forecasts based on the range of the
-#' prediction intervals. For an increasing number of quantiles this measure
-#' converges against the sample based bias version for integer and continuous
-#' forecasts.
-#'
-#' @details
-#' For quantile forecasts, bias is measured as
-#'
-#' \deqn{
-#' B_t = (1 - 2 \cdot \max \{i | q_{t,i} \in Q_t \land q_{t,i} \leq x_t\})
-#' \mathbf{1}( x_t \leq q_{t, 0.5}) \\
-#' + (1 - 2 \cdot \min \{i | q_{t,i} \in Q_t \land q_{t,i} \geq x_t\})
-#'  \mathbf{1}( x_t \geq q_{t, 0.5}),
-#' }{
-#' B_t = (1 - 2 * max(i | q_{t,i} in Q_t and q_{t,i} <= x_t\))
-#' 1( x_t <= q_{t, 0.5}) + (1 - 2 * min(i | q_{t,i} in Q_t and q_{t,i} >= x_t))
-#'  1( x_t >= q_{t, 0.5}),
-#' }
-#'
-#' where \eqn{Q_t} is the set of quantiles that form the predictive
-#' distribution at time \eqn{t}. They represent our
-#' belief about what the later observed value \eqn{x_t} will be. For
-#' consistency, we define
-#' \eqn{Q_t} such that it always includes the element
-#' \eqn{q_{t, 0} = - \infty} and \eqn{q_{t,1} = \infty}.
-#' \eqn{\mathbf{1}()}{1()} is the indicator function that is \eqn{1} if the
-#' condition is satisfied and $0$ otherwise. In clearer terms, \eqn{B_t} is
-#' defined as the maximum percentile rank for which the corresponding quantile
-#' is still below the observed value, if the observed value is smaller than the
-#' median of the predictive distribution. If the observed value is above the
-#' median of the predictive distribution, then $B_t$ is the minimum percentile
-#' rank for which the corresponding quantile is still larger than the true
-#' value. If the observed value is exactly the median, both terms cancel out and
-#' \eqn{B_t} is zero. For a large enough number of quantiles, the
-#' percentile rank will equal the proportion of predictive samples below the
-#' observed value, and this metric coincides with the one for
-#' continuous forecasts.
-#'
-#' Bias can assume values between
-#' -1 and 1 and is 0 ideally.
-#' @param lower vector of length corresponding to the number of central
-#' prediction intervals that holds predictions for the lower bounds of a
-#' prediction interval
-#' @param upper vector of length corresponding to the number of central
-#' prediction intervals that holds predictions for the upper bounds of a
-#' prediction interval
-#' @param range vector of corresponding size with information about the width
-#' of the central prediction interval
-#' @param observed a single observed value
-#' @return scalar with the quantile bias for a single quantile prediction
-#' @seealso bias_quantile bias_sample
-#' @author Nikos Bosse \email{nikosbosse@@gmail.com}
-#' @examples
-#'
-#' lower <- c(
-#'   6341.000, 6329.500, 6087.014, 5703.500,
-#'   5451.000, 5340.500, 4821.996, 4709.000,
-#'   4341.500, 4006.250, 1127.000, 705.500
-#' )
-#'
-#' upper <- c(
-#'   6341.000, 6352.500, 6594.986, 6978.500,
-#'   7231.000, 7341.500, 7860.004, 7973.000,
-#'   8340.500, 8675.750, 11555.000, 11976.500
-#' )
-#'
-#' range <- c(0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 95, 98)
-#'
-#' observed <- 8062
-#'
-#' bias_range(
-#'   lower = lower, upper = upper,
-#'   range = range, observed = observed
-#' )
-#' @export
-#' @keywords metric
-
-bias_range <- function(lower, upper, range, observed) {
-
-  if (anyNA(range)) {
-    if (is.na(range[1]) && !any(range[-1] == 0)) {
-      range[1] <- 0
-    }
-    range <- range[!is.na(range)]
-    lower <- lower[!is.na(range)]
-    upper <- upper[!is.na(range)]
-  }
-
-  if (length(range) > 1 && !all(diff(range) > 0)) {
-    stop("Range must be increasing")
-  }
-
-  if (length(lower) != length(upper) || length(range) != length(lower)) {
-    stop("Inputs must have same length")
-  }
-
-  check_quantiles(range, name = "range", range = c(0, 100))
-
-  # Convert range to quantiles
-  quantile <- c(
-    rev(abs(100 - range) / (2 * 100)),
-    abs(100 + range[range != 0]) / (2 * 100)
-  )
-
-  # Combine predictions
-  upper_without_median <- upper[range != 0]
-  predicted <- c(rev(lower), upper_without_median)
-
-  # Call bias_quantile
-  bias <- bias_quantile(observed, predicted, quantile)
-
-  return(bias)
-}

tests/testthat/test-metrics-quantile.R

@@ -798,3 +798,10 @@ test_that("bias_quantile only produces one message", {
     "Median not available, computing bias as mean of the two innermost quantiles in order to compute bias."
   )
 })
+
+test_that("bias_quantile() works with point forecasts", {
+  predicted <- 1
+  observed <- 1
+  quantile <- 0.5
+  expect_equal(bias_quantile(observed, predicted, quantile), 0)
+})