es: Energy Score

In jobstdavid/eppverification: Verification Tools for the Statistical Postprocessing of Ensemble Forecasts

Energy Score

Description

This function calculates the Energy Score (ES) given observations of a multivariate variable and samples of a predictive distribution.

Usage

es(y, x, method = "ens", na.action = na.omit, aggregate = FALSE, ...)

Arguments

y	matrix of observations (see details)
x	3-dimensional array of samples of a multivariate predictive distribution (depending on y; see details)
method	character; "ens" and "mc"; default: "ens" (see details)
na.action	function to handle the NA's. Default: na.omit.
aggregate	logical or function for aggregating the single scores, e.g. sum, mean, weighted.mean, .... Default: FALSE, i.e. no aggregation function.
...	further arguments passed to the aggregate function.

Details

The observations are given in the matrix y with n rows, where each column belongs to an univariate observation variable. The i-th row of matrix y belongs to the i-th third dimension entry of the array x. The i-th third dimension entry must be a matrix with n rows, having the same structure as y, filled with the samples of a multivariate predictive distribution.

If method "ens" is specified, the ES values are calculated for given samples x of a multivariate predictive distribution (Gneiting et al., 2008). If method "mc" is specified, the ES values are calculated by a Monte-Carlo approximation using samples x of a predictive distribution (Gneiting et al., 2008). In the latter case, the number of samples should be "large", e.g. 10.000.

A lower ES indicates a better forecast.

Value

Vector of the score value(s).

Author(s)

David Jobst

References

Gneiting, T., Stanberry, L., Grimit, E., Held, L. and Johnson, N. (2008). Assessing probabilistic forecasts of multivariate quantities, with an application to ensemble predictions of surface winds. Test, 17, 211-264.

Examples

# simulated data
n <- 30
m1 <- 50
m2 <- 10000
y <- cbind(rnorm(n), rgamma(n, shape = 1))
x1 <- array(NA, dim = c(m1, 2, n))
x2 <- array(NA, dim = c(m2, 2, n))
x1[, 1, ] <- rnorm(n*m1)
x1[, 2, ] <- rgamma(n*m1, shape = 1)
x2[, 1, ] <- rnorm(n*m2)
x2[, 2, ] <- rgamma(n*m2, shape = 1)

# es calculation
es(y = y, x = x1, method = "ens")
es(y = y, x = x1, method = "ens", aggregate = mean)

es(y = y, x = x2, method = "mc")
es(y = y, x = x2, method = "mc", aggregate = mean)

jobstdavid/eppverification documentation built on May 13, 2024, 5:20 p.m.