Develop into main

.Rbuildignore

@@ -12,6 +12,7 @@
 ^Meta$
 ^_pkgdown\.yml$
 ^inst/manuscript/manuscript_cache$
+^inst/manuscript/.trackdown$
 ^\.lintr$
 ^docs$
 ^\.devcontainer$

.github/dependabot.yml

@@ -0,0 +1,11 @@
+# To get started with Dependabot version updates, you'll need to specify which
+# package ecosystems to update and where the package manifests are located.
+# Please see the documentation for all configuration options:
+# https://docs.github.com/github/administering-a-repository/configuration-options-for-dependency-updates
+
+version: 2
+updates:
+  - package-ecosystem: "github-actions" # See documentation for possible values
+    directory: "/" # Location of package manifests
+    schedule:
+      interval: "weekly"

.github/workflows/render_readme.yaml

@@ -9,6 +9,7 @@ on:
   push:
     paths:
       - 'README.Rmd'
+      - DESCRIPTION
 
 jobs:
   render-readme:

.gitignore

@@ -8,9 +8,12 @@ doc
 Meta
 inst/manuscript/manuscript_cache/
 inst/manuscript/manuscript.log
+inst/manuscript/manuscript.aux
+inst/manuscript/manuscript.blg
 inst/manuscript/manuscript.pdf
 inst/manuscript/manuscript.tex
 inst/manuscript/manuscript_files/
+inst/manuscript/.trackdown
 docs
 ..bfg-report/
 .DS_Store

.lintr

@@ -11,6 +11,6 @@ linters: linters_with_tags(
 exclusions: c(
     list.files("tests", recursive = TRUE, full.names = TRUE),
     list.files("inst", recursive = TRUE, full.names = TRUE),
-    "vignettes/metric-details.Rmd"
+    list.files("vignettes", pattern = ".R$", full.names = TRUE)
   )
 exclude: "# nolint"

DESCRIPTION

@@ -1,6 +1,6 @@
 Package: scoringutils
 Title: Utilities for Scoring and Assessing Predictions
-Version: 1.2.2
+Version: 1.2.2.9000
 Language: en-GB
 Authors@R: c(
     person(given = "Nikos",
@@ -47,11 +47,13 @@ License: MIT + file LICENSE
 Encoding: UTF-8
 LazyData: true
 Imports: 
+    checkmate, 
     data.table,
     ggdist (>= 3.2.0),
     ggplot2 (>= 3.4.0),
     lifecycle,
     methods,
+    Metrics,
     rlang,
     scoringRules,
     stats

NAMESPACE

@@ -1,34 +1,51 @@
 # Generated by roxygen2: do not edit by hand
 
+S3method(as_forecast,default)
 S3method(print,scoringutils_check)
-export(abs_error)
+S3method(quantile_to_interval,data.frame)
+S3method(quantile_to_interval,numeric)
+S3method(score,default)
+S3method(score,forecast_binary)
+S3method(score,forecast_point)
+S3method(score,forecast_quantile)
+S3method(score,forecast_sample)
+S3method(validate_forecast,forecast_binary)
+S3method(validate_forecast,forecast_point)
+S3method(validate_forecast,forecast_quantile)
+S3method(validate_forecast,forecast_sample)
 export(add_coverage)
+export(add_pairwise_comparison)
 export(ae_median_quantile)
 export(ae_median_sample)
-export(avail_forecasts)
+export(as_forecast)
 export(available_metrics)
 export(bias_quantile)
-export(bias_range)
 export(bias_sample)
 export(brier_score)
-export(check_forecasts)
 export(correlation)
 export(crps_sample)
+export(dispersion)
 export(dss_sample)
-export(find_duplicates)
-export(interval_score)
+export(get_duplicate_forecasts)
+export(get_forecast_counts)
+export(get_forecast_type)
+export(get_forecast_unit)
+export(interval_coverage)
+export(interval_coverage_deviation)
 export(log_shift)
 export(logs_binary)
 export(logs_sample)
 export(mad_sample)
 export(make_NA)
 export(make_na)
 export(merge_pred_and_obs)
+export(new_forecast)
+export(overprediction)
 export(pairwise_comparison)
 export(pit)
 export(pit_sample)
-export(plot_avail_forecasts)
 export(plot_correlation)
+export(plot_forecast_counts)
 export(plot_heatmap)
 export(plot_interval_coverage)
 export(plot_pairwise_comparison)
@@ -39,15 +56,49 @@ export(plot_ranges)
 export(plot_score_table)
 export(plot_wis)
 export(quantile_score)
+export(quantile_to_interval)
+export(rules_binary)
+export(rules_point)
+export(rules_quantile)
+export(rules_sample)
+export(run_safely)
 export(sample_to_quantile)
 export(score)
 export(se_mean_sample)
+export(select_rules)
 export(set_forecast_unit)
-export(squared_error)
 export(summarise_scores)
 export(summarize_scores)
 export(theme_scoringutils)
 export(transform_forecasts)
+export(underprediction)
+export(validate_forecast)
+export(validate_general)
+export(wis)
+importFrom(Metrics,ae)
+importFrom(Metrics,ape)
+importFrom(Metrics,se)
+importFrom(checkmate,assert)
+importFrom(checkmate,assert_character)
+importFrom(checkmate,assert_data_frame)
+importFrom(checkmate,assert_data_table)
+importFrom(checkmate,assert_factor)
+importFrom(checkmate,assert_list)
+importFrom(checkmate,assert_logical)
+importFrom(checkmate,assert_number)
+importFrom(checkmate,assert_numeric)
+importFrom(checkmate,assert_string)
+importFrom(checkmate,assert_subset)
+importFrom(checkmate,assert_vector)
+importFrom(checkmate,check_atomic_vector)
+importFrom(checkmate,check_data_frame)
+importFrom(checkmate,check_function)
+importFrom(checkmate,check_matrix)
+importFrom(checkmate,check_numeric)
+importFrom(checkmate,check_vector)
+importFrom(checkmate,test_factor)
+importFrom(checkmate,test_list)
+importFrom(checkmate,test_numeric)
 importFrom(data.table,"%like%")
 importFrom(data.table,':=')
 importFrom(data.table,.I)
@@ -59,8 +110,11 @@ importFrom(data.table,data.table)
 importFrom(data.table,dcast)
 importFrom(data.table,is.data.table)
 importFrom(data.table,melt)
+importFrom(data.table,nafill)
 importFrom(data.table,rbindlist)
 importFrom(data.table,setDT)
+importFrom(data.table,setattr)
+importFrom(data.table,setcolorder)
 importFrom(data.table,setnames)
 importFrom(ggdist,geom_lineribbon)
 importFrom(ggplot2,.data)
@@ -103,7 +157,6 @@ importFrom(ggplot2,xlab)
 importFrom(ggplot2,ylab)
 importFrom(lifecycle,deprecated)
 importFrom(methods,hasArg)
-importFrom(methods,is)
 importFrom(rlang,enexprs)
 importFrom(rlang,warn)
 importFrom(scoringRules,crps_sample)
@@ -121,5 +174,6 @@ importFrom(stats,rbinom)
 importFrom(stats,reorder)
 importFrom(stats,runif)
 importFrom(stats,sd)
+importFrom(stats,weighted.mean)
 importFrom(stats,wilcox.test)
 importFrom(utils,combn)

NEWS.md

@@ -1,3 +1,48 @@
+# scoringutils 1.2.2.9000
+
+This major update and addresses a variety of comments made by reviewers from the Journal of Statistical Software (see preprint of the manuscript [here](https://arxiv.org/abs/2205.07090)).
+
+The update introduces breaking changes. If you want to keep using the older version, you can download it using `remotes::install_github("epiforecasts/[email protected]")`. 
+
+## Package updates
+- In `score()`, required columns "true_value" and "prediction" were renamed and replaced by required columns "observed" and "predicted". Scoring functions now also use the function arguments "observed" and "predicted" everywhere consistently. 
+- The overall scoring workflow was updated. `score()` is now a generic function that dispatches the correct method based on the forecast type. forecast types currently supported are "binary", "point", "sample" and "quantile" with corresponding classes "forecast_binary", "forecast_point", "forecast_sample" and "forecast_quantile". An object of class `forecast_*` can be created using the function `as_forecast()`, which also replaces the previous function `check_forecasts()` (see more information below). 
+- Scoring rules (functions used for scoring) received a consistent interface and input checks:
+  - Scoring rules for binary forecasts:
+    - `observed`: factor with exactly 2 levels
+    - `predicted`: numeric, vector with probabilities
+  - Scoring rules for point forecasts:
+    - `observed`: numeric vector
+    - `predicted`: numeric vector
+  - Scoring rules for sample-based forecasts:
+    - `observed`: numeric, either a scalar or a vector
+    - `predicted`: numeric, a vector (if `observed` is a scalar) or a matrix (if `observed` is a vector)
+  - Scoring rules for quantile-based forecasts:
+    - `observed`: numeric, either a scalar or a vector
+    - `predicted`: numeric, a vector (if `observed` is a scalar) or a matrix (if `observed` is a vector)
+    - `quantile`: numeric, a vector with quantile-levels. Can alternatively be a matrix of the same shape as `predicted`.
+- Users can now supply their own scoring rules to `score()` as a list of functions. Default scoring rules can be accessed using the functions `rules_point()`, `rules_sample()`, `rules_quantile()` and `rules_binary()`, which return a list of scoring rules suitable for the respective forecast type.
+- `check_forecasts()` was replaced by a different workflow. There now is a function, `as_forecast()`, that determines forecast type of the data, constructs a forecasting object and validates it using the function `validate_forecast()` (a generic that dispatches the correct method based on the forecast type). Objects of class `forecast_binary`, `forecast_point`, `forecast_sample` and `forecast_quantile` have print methods that fulfill the functionality of `check_forecasts()`.
+- 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 `get_forecast_counts()`. This represents a change in the naming convention where we aim to name functions that provide the user with additional useful information about the data with a prefix "get_". Sees Issue #403 and #521 and PR #511 by @nikosbosse and reviewed by @seabbs for details. 
+  - For clarity, the output column in `get_forecast_counts()` was renamed from "Number forecasts" to "count".
+  - `get_forecast_counts()` now also displays combinations where there are 0 forecasts, instead of silently dropping corresponding rows.
+  - `plot_avail_forecasts()` was renamed `plot_forecast_counts()` in line with the change in the function name. The `x` argument no longer has a default value, as the value will depend on the data provided by the user.
+- The deprecated `..density..` was replaced with `after_stat(density)` in ggplot calls.
+- Files ending in ".Rda" were renamed to ".rds" where appropriate when used together with `saveRDS()` or `readRDS()`.
+- `score()` now calls `na.omit()` on the data, instead of only removing rows with missing values in the columns `observed` and `predicted`. This is because `NA` values in other columns can also mess up e.g. grouping of forecasts according to the unit of a single forecast.
+- added documentation for the return value of `summarise_scores()`.
+- Removed abs_error and squared_error from the package in favour of `Metrics::ae` and `Metrics::se`.
+- Renamed `interval_coverage_quantile()` and `interval_coverage_dev_quantile()` to `interval_coverage()` and `interval_coverage_deviation()`, respectively. Removed `interval_coverage_sample()` as users are now expected to convert to a quantile format first before scoring.
+- Added unit tests for `interval_coverage_quantile()` and `interval_coverage_dev_quantile()` in order to make sure that the functions provide the correct warnings when insufficient quantiles are provided.
+- Documentation pkgdown pages are now created both for the stable and dev versions.
+
 # scoringutils 1.2.2
 
 ## Package updates
@@ -6,7 +51,7 @@
 
 ## Bug fixes
 - Fixes a bug with `set_forecast_unit()` where the function only workded with a data.table, but not a data.frame as an input. 
-- The metrics table in the vignette [Details on the metrics implemented in `scoringutils`](https://epiforecasts.io/scoringutils/articles/metric-details.html) had duplicated entries. This was fixed by removing the duplicated rows. 
+- The metrics table in the vignette [Details on the metrics implemented in `scoringutils`](https://epiforecasts.io/scoringutils/articles/metric-details.html) had duplicated entries. This was fixed by removing the duplicated rows.
 
 # scoringutils 1.2.1
 

R/add_coverage.R

@@ -0,0 +1,83 @@
+#' @title Add Coverage Values to Quantile-Based Forecasts
+#'
+#' @description Adds interval coverage of central prediction intervals,
+#' quantile coverage for predictive quantiles, as well as the deviation between
+#' desired and actual coverage to a data.table. Forecasts should be in a
+#' quantile format (following the input requirements of `score()`).
+#'
+#' **Interval coverage**
+#'
+#' Interval coverage for a given interval range is defined as the proportion of
+#' observations that fall within the corresponding central prediction intervals.
+#' Central prediction intervals are symmetric around the median and and formed
+#' by two quantiles that denote the lower and upper bound. For example, the 50%
+#' central prediction interval is the interval between the 0.25 and 0.75
+#' quantiles of the predictive distribution.
+#'
+#' The function `add_coverage()` computes the coverage per central prediction
+#' interval, so the interval coverage will always be either `TRUE`
+#' (observed value falls within the interval) or `FALSE`  (observed value falls
+#' outside the interval). You can summarise the interval coverage values to get
+#' the proportion of observations that fall within the central prediction
+#' intervals.
+#'
+#' **Quantile coverage**
+#'
+#' Quantile coverage for a given quantile is defined as the proportion of
+#' observed values that are smaller than the corresponding predictive quantile.
+#' For example, the 0.5 quantile coverage is the proportion of observed values
+#' that are smaller than the 0.5 quantile of the predictive distribution.
+#' Just as above, for a single observation and the quantile of a single
+#' predictive distribution, the value will either be `TRUE` or `FALSE`.
+#'
+#' **Coverage deviation**
+#'
+#' The coverage deviation is the difference between the desired coverage
+#' (can be either interval or quantile coverage) and the
+#' actual coverage. For example, if the desired coverage is 90% and the actual
+#' coverage is 80%, the coverage deviation is -0.1.
+#'
+#' @inheritParams score
+#' @return a data.table with the input and columns "interval_coverage",
+#' "interval_coverage_deviation", "quantile_coverage",
+#' "quantile_coverage_deviation" added.
+#' @importFrom data.table setcolorder
+#' @examples
+#' library(magrittr) # pipe operator
+#' example_quantile %>%
+#'   add_coverage()
+#' @export
+#' @keywords scoring
+#' @export
+add_coverage <- function(data) {
+  stored_attributes <- get_scoringutils_attributes(data)
+  data <- as_forecast(data)
+  forecast_unit <- get_forecast_unit(data)
+  data_cols <- colnames(data) # store so we can reset column order later
+
+  interval_data <- quantile_to_interval(data, format = "wide")
+  interval_data[,
+    interval_coverage := (observed <= upper) & (observed >= lower)
+  ][, c("lower", "upper", "observed") := NULL]
+
+  data[, range := get_range_from_quantile(quantile)]
+
+  data <- merge(interval_data, data, by = unique(c(forecast_unit, "range")))
+  data[, interval_coverage_deviation := interval_coverage - range / 100]
+  data[, quantile_coverage := observed <= predicted]
+  data[, quantile_coverage_deviation := quantile_coverage - quantile]
+
+  # reset column order
+  new_metrics <- c("interval_coverage", "interval_coverage_deviation",
+                   "quantile_coverage", "quantile_coverage_deviation")
+  setcolorder(data, unique(c(data_cols, "range", new_metrics)))
+
+  # add coverage "metrics" to list of stored metrics
+  # this makes it possible to use `summarise_scores()` later on
+  stored_attributes[["metric_names"]] <- c(
+    stored_attributes[["metric_names"]],
+    new_metrics
+  )
+  data <- assign_attributes(data, stored_attributes)
+  return(data[])
+}