1

Public Documentation

Documentation for ScoringRules's public interface.

ScoringRules.ScoringRules Module

julia

ScoringRules

Proper scoring rules for probabilistic forecasts in Julia.

A port of the R package scoringRules by Alexander I. Jordan, Fabian Krüger, Sebastian Lerch and Sam Allen. See the package README for attribution and provenance. Distributed under GPL-2.0-or-later.

Scores follow the lower-is-better convention (negative orientation) and are broadcast-friendly. The three univariate entry points

crps — continuous ranked probability score
logs — logarithmic score
dss — Dawid–Sebastiani score

dispatch on Distributions.jl types for parametric forecasts and on AbstractVectors for simulated (ensemble) forecasts. Multivariate ensemble forecasts are scored with es (energy score), vs (variogram score) and mmds (maximum-mean-discrepancy score).

Example

julia

using Distributions, ScoringRules
crps(Normal(0, 1), 0.5)

Contents

public api

Index

ScoringRules.ScoringRules
ScoringRules.LogLaplace
ScoringRules.TwoPieceExponential
ScoringRules.TwoPieceNormal
ScoringRules.clogs
ScoringRules.crps
ScoringRules.dss
ScoringRules.dss_moments
ScoringRules.es
ScoringRules.ess_moments
ScoringRules.interval_score
ScoringRules.logs
ScoringRules.mmds
ScoringRules.owcrps
ScoringRules.owes
ScoringRules.owmmds
ScoringRules.owvs
ScoringRules.quantile_score
ScoringRules.rps
ScoringRules.twcrps
ScoringRules.twes
ScoringRules.twmmds
ScoringRules.twvs
ScoringRules.vs

Public API

ScoringRules.LogLaplace Type

julia

struct LogLaplace{T<:Real} <: Distributions.Distribution{Distributions.Univariate, Distributions.Continuous}

μ::Real
σ::Real

julia

LogLaplace(μ, σ)

The log-Laplace distribution: X ~ LogLaplace(μ, σ) iff log(X) ~ Laplace(μ, σ), with μ the location and σ > 0 the scale, both on the log scale. Supported on the positive reals. The mean exists for σ < 1 and the variance for σ < 1/2. Not part of Distributions.jl; provided here so it flows through crps/logs/dss dispatch.

ScoringRules.TwoPieceExponential Type

julia

struct TwoPieceExponential{T<:Real} <: Distributions.Distribution{Distributions.Univariate, Distributions.Continuous}

location::Real
scale1::Real
scale2::Real

julia

TwoPieceExponential(location, scale1, scale2)

The two-piece (double / asymmetric) exponential distribution: back-to-back exponential tails with scale scale1 below location and scale2 above it. Not part of Distributions.jl; provided here so it flows through crps/logs/dss dispatch.

ScoringRules.TwoPieceNormal Type

julia

struct TwoPieceNormal{T<:Real} <: Distributions.Distribution{Distributions.Univariate, Distributions.Continuous}

location::Real
scale1::Real
scale2::Real

julia

TwoPieceNormal(location, scale1, scale2)

The two-piece (split) normal distribution: a normal density with scale scale1 below location and scale2 above it, renormalised to a proper density. With scale1 == scale2 it reduces to Normal(location, scale1). Not part of Distributions.jl; provided here so it flows through crps/logs/dss dispatch.

ScoringRules.clogs Function

julia

clogs(
    dat::AbstractVector{<:Real},
    y::Real;
    a,
    b,
    bw,
    cens
) -> Any

julia

clogs(dat::AbstractVector{<:Real}, y::Real; a=-Inf, b=Inf, bw=nothing, cens=true)

Censored likelihood score (cens = true, the default) or conditional likelihood score (cens = false) of an ensemble forecast dat at observation y, emphasising the window (a, b) (Diks, Panchenko & van Dijk 2011).

Both scores are based on the Gaussian kernel density estimate of the ensemble, with its probability mass in the window and the indicator (strict inequalities):

censored: when y is inside the window and otherwise, so outcomes outside the window enter only through their total probability.
conditional: when y is inside the window and 0 otherwise.

With the default unbounded window both variants reduce to logs(dat, y; bw). If bw is nothing, Silverman's rule-of-thumb bandwidth is used (matching R's bw.nrd). Lower is better.

Arguments

dat: ensemble of simulation draws.
y: scalar observation.

Keyword Arguments

a: lower end of the window (default -Inf).
b: upper end of the window (default Inf).
bw: optional bandwidth; defaults to Silverman's rule-of-thumb.
cens: true for the censored score, false for the conditional score.

Provenance

Ported from clogs_sample in R scoringRules (scores_sample_univ_weighted.R; mixn.cpp; Jordan, Krüger, Lerch, Allen).

Example

julia

using ScoringRules
dat = randn(100)
clogs(dat, 0.5; a = 0.0, b = 1.0)

ScoringRules.crps Function

julia

crps(
    d::Distributions.Distribution{Distributions.Univariate, Distributions.Continuous},
    y::Real
) -> Any

julia

crps(d::UnivariateDistribution, y)

Continuous ranked probability score of the forecast distribution d at the observation y,

Distribution-specific methods provide closed forms; this generic method is a numerical fallback (adaptive quadrature for continuous d). Lower is better.

Arguments

d: forecast distribution (any UnivariateDistribution).
y: scalar observation.

Example

julia

using Distributions, ScoringRules
crps(Normal(0, 1), 0.5)

julia

crps(
    d::Distributions.MixtureModel{Distributions.Univariate, Distributions.Continuous, <:Distributions.Normal},
    y::Real
) -> Any

CRPS of a MixtureModel of normal components at observation y.

Dispatches on MixtureModel{Univariate,Continuous,<:Normal}, extracting component parameters via components(d) and probs(d).

julia

crps(
    dat::AbstractVector{<:Real},
    y::Real;
    method,
    w,
    bw
) -> Any

julia

crps(dat::AbstractVector{<:Real}, y::Real; method=:edf, w=nothing, bw=nothing)

CRPS of an ensemble forecast dat (a vector of m simulation draws) at observation y.

Two approximation methods are available via method:

:edf (default) — empirical distribution function approximation using the quantile decomposition of Laio & Tamea (2007). Optional finite, non-negative weights w (length m) are normalised to sum to one internally.
:kde — Gaussian kernel density estimate with bandwidth bw. If bw is nothing, Silverman's rule-of-thumb is applied (matching R's bw.nrd). The w argument is ignored for this method.

Lower is better.

Arguments

dat: ensemble of m simulation draws.
y: scalar observation.

Keyword Arguments

method: :edf (default) or :kde.
w: optional finite, non-negative weight vector (length m) with a positive sum; only used for :edf.
bw: optional bandwidth; only used for :kde.

Provenance

Ported from crps_sample / crps_edf / crps_kdens in R scoringRules (scores_sample_univ.R; mixn.cpp; Jordan, Krüger, Lerch, Allen).

Example

julia

using ScoringRules
dat = randn(100)
crps(dat, 0.5)

ScoringRules.dss Function

julia

dss(d::Distributions.UnivariateDistribution, y::Real) -> Any

julia

dss(d::UnivariateDistribution, y)

Dawid–Sebastiani score of the forecast distribution d at the observation y,

where and are the mean and variance of d. Only the first two moments of the forecast enter. Lower is better.

Returns NaN when the forecast variance is not finite and positive (for example a Student-t with df ≤ 2, a GEV/GPD with shape ≥ 1/2, a LogLaplace with σ ≥ 1/2 or a log-logistic with β ≤ 2), matching the dss_* functions in R scoringRules.

Arguments

d: forecast distribution (any UnivariateDistribution).
y: scalar observation.

Example

julia

using Distributions, ScoringRules
dss(Normal(0, 1), 0.5)

julia

dss(dat::AbstractVector{<:Real}, y::Real; w) -> Any

julia

dss(dat::AbstractVector{<:Real}, y::Real; w=nothing)

Dawid–Sebastiani score of an ensemble forecast dat (a vector of m simulation draws) at observation y:

where is the sample mean and is the population variance (R uses mean(dat^2) - mean(dat)^2, i.e. the biased estimator). Optional finite, non-negative member weights w (length m) are normalised to sum to one internally and replace the mean and variance with their weighted versions. Lower is better.

Arguments

dat: ensemble of m simulation draws.
y: scalar observation.

Keyword Arguments

w: optional finite, non-negative member weight vector (length m) with a positive sum.

Provenance

Ported from dss_sample / dss_edf in R scoringRules (scores_sample_univ.R; Jordan, Krüger, Lerch, Allen), including the member-weight handling of dss_edf.

Example

julia

using ScoringRules
dat = randn(100)
dss(dat, 0.5)

ScoringRules.dss_moments Function

julia

dss_moments(y::Real, mean::Real, var::Real) -> Any

julia

dss_moments(y, mean, var)

Dawid–Sebastiani score from a moment forecast given directly as its mean and variance, without constructing a distribution. Since the DSS depends on the forecast only through its first two moments, this equals dss(d, y) for any d with that mean and variance. This is the moment-based forecast input mode. Lower is better.

Arguments

y: scalar observation.
mean: forecast mean.
var: forecast variance.

Example

julia

using ScoringRules
dss_moments(0.5, 0.0, 1.0)

ScoringRules.es Function

julia

es(X::AbstractMatrix, y::AbstractVector; w) -> Any

julia

es(X, y; w=nothing)

Energy score of the ensemble forecast X (a d × m matrix, each column one member) at the d-dimensional observation y:

Optional per-member weights w (length m); they are normalised to sum to one internally. If w is nothing, uniform weights are used. Lower is better.

Provenance

Ported from es_sample in R scoringRules (scores_sample_multiv.R; procs_es.cpp; Jordan, Krüger, Lerch, Allen). Gneiting et al. (2008), TEST 17, 211–235.

Example

julia

using ScoringRules
X = randn(2, 50)
y = [0.0, 0.0]
es(X, y)

ScoringRules.ess_moments Function

julia

ess_moments(
    y::Real,
    mean::Real,
    var::Real,
    skew::Real
) -> Any

julia

ess_moments(y, mean, var, skew)

Error-spread score (Christensen, Moroz and Palmer 2015) from a moment forecast given directly as its mean, variance and skewness,

where , and are the forecast mean, variance and skewness. The score assesses whether the spread and skewness of an ensemble forecast are consistent with its error. This is the moment-based forecast input mode. Lower is better.

Arguments

y: scalar observation.
mean: forecast mean.
var: forecast variance; must be non-negative (a negative value throws a DomainError, where R returns NaN with a warning).
skew: forecast skewness.

Example

julia

using ScoringRules
ess_moments(0.5, 0.0, 1.0, 0.5)

ScoringRules.interval_score Function

julia

interval_score(
    lower::Real,
    upper::Real,
    y::Real,
    level::Real
) -> Any

julia

interval_score(lower, upper, y, level)

Interval score (Gneiting & Raftery 2007) for a central prediction interval at nominal coverage level (e.g. level = 0.8 for an 80 % interval), with lower endpoint lower, upper endpoint upper, and observation y.

Arguments

lower — lower quantile forecast (the α₁ = (1-level)/2 quantile).
upper — upper quantile forecast (the α₂ = 1 - (1-level)/2 quantile).
y — scalar observation.
level — nominal coverage, 0 < level < 1.

Returns

A scalar Float64 score.

The score is

julia

IS = (upper - lower)
       + (2 / (1-level)) * (lower - y) * 1{y < lower}
       + (2 / (1-level)) * (y - upper) * 1{y > upper}

This convention uses coverage (not alpha = 1 - coverage), matching ints_quantiles in scores_quantiles.R. Lower values indicate better forecasts.

Example

julia

using ScoringRules
interval_score(-1.64, 1.64, 0.5, 0.9)

julia

interval_score(dat::AbstractVector, y::Real; level, type)

julia

interval_score(dat, y; level, type=7)

Interval score estimated from an ensemble dat, at nominal coverage level. The α₁ and α₂ quantiles are computed from dat with interpolation type type (default 7, matching R), then interval_score(lower, upper, y, level) is called. Mirrors ints_sample in scores_quantiles.R.

Arguments

dat — AbstractVector of ensemble draws.
y — scalar observation.
level — nominal coverage, 0 < level < 1.
type — quantile interpolation type (1–9); default 7.

Returns

A scalar Float64 score.

Example

julia

using ScoringRules
dat = randn(100)
interval_score(dat, 0.5; level = 0.9)

ScoringRules.logs Function

julia

logs(
    d::Distributions.UnivariateDistribution,
    y::Real
) -> Any

julia

logs(d::UnivariateDistribution, y)

Logarithmic score of the forecast distribution d at the observation y, equal to the negative log-likelihood -logpdf(d, y) (or -logpmf for discrete d). Lower is better.

Arguments

d: forecast distribution (any UnivariateDistribution).
y: scalar observation.

Example

julia

using Distributions, ScoringRules
logs(Normal(0, 1), 0.5)

julia

logs(dat::AbstractVector{<:Real}, y::Real; bw) -> Any

julia

logs(dat::AbstractVector{<:Real}, y::Real; bw=nothing)

Logarithmic score of an ensemble forecast dat (a vector of m simulation draws) at observation y using Gaussian kernel density estimation.

If bw is nothing, Silverman's rule-of-thumb bandwidth is used (matching R's bw.nrd). Lower is better.

Arguments

dat: ensemble of m simulation draws.
y: scalar observation.

Keyword Arguments

bw: optional bandwidth; defaults to Silverman's rule-of-thumb.

Provenance

Ported from logs_sample / lsmixnC in R scoringRules (scores_sample_univ.R; mixn.cpp; Jordan, Krüger, Lerch, Allen).

Example

julia

using ScoringRules
dat = randn(100)
logs(dat, 0.5)

ScoringRules.mmds Function

julia

mmds(X::AbstractMatrix, y::AbstractVector; w) -> Any

julia

mmds(X, y; w=nothing)

Maximum-mean-discrepancy score (Gaussian kernel with σ = 1) of the ensemble forecast X (a d × m matrix) at the d-dimensional observation y:

Optional per-member weights w (length m); they are normalised to sum to one internally. If w is nothing, uniform weights are used. Lower is better.

Provenance

Ported from mmds_sample in R scoringRules (scores_sample_multiv.R; procs_es.cpp; Jordan, Krüger, Lerch, Allen). Gneiting & Raftery (2007), JASA 102, 359–378; Gretton et al. (2012), JMLR 13, 723–773.

Example

julia

using ScoringRules
X = randn(2, 50)
y = [0.0, 0.0]
mmds(X, y)

ScoringRules.owcrps Function

julia

owcrps(
    dat::AbstractVector,
    y::Real;
    a,
    b,
    weight_func,
    w
) -> Any

julia

owcrps(dat, y; a=-Inf, b=Inf, weight_func=nothing, w=nothing)

Outcome-weighted CRPS of the ensemble dat at observation y, emphasising outcomes in the interval (a, b).

Each ensemble member xᵢ receives weight w(xᵢ) = 1{a < xᵢ < b}; the observation gets weight w(y) = 1{a < y < b}. The weighted EDF-CRPS (normalised by the sum of member weights) is then multiplied by w(y). Optional member weights w multiply the outcome weights, as in R.

Alternatively, supply a custom vectorised weight_func (takes a Real, returns a non-negative Real); supplying weight_func ignores a and b.

Returns NaN when all combined weights are zero (no member in the region). Lower is better.

Arguments

dat: ensemble of simulation draws.
y: scalar observation.

Keyword Arguments

a: lower threshold (default -Inf).
b: upper threshold (default Inf).
weight_func: custom weight function; overrides a and b when supplied.
w: optional non-negative member weights (length m); multiplied with the outcome weights.

Provenance

Ported from owcrps_sample in R scoringRules (scores_sample_univ_weighted.R; Allen 2024, JSS 110(8)).

Example

julia

using ScoringRules
dat = randn(100)
owcrps(dat, 0.5; a = 0.0, b = 1.0)

ScoringRules.owes Function

julia

owes(
    X::AbstractMatrix,
    y::AbstractVector;
    a,
    b,
    weight_func,
    w
) -> Any

julia

owes(X, y; a=-Inf, b=Inf, weight_func=nothing, w=nothing)

Outcome-weighted energy score of the ensemble X (a d × m matrix) at the d-dimensional observation y.

Each column Xᵢ of X receives weight w(Xᵢ) where the default weight function is w(z) = 1{∀k: a[k] < z[k] < b[k]}. The observation weight is w(y). The energy score is computed with the normalised column weights and then multiplied by w(y). Optional member weights w (length m) multiply the outcome weights, as in R.

Returns NaN when all combined weights are zero. A custom weight_func (takes a length-d vector, returns a non-negative scalar) overrides a and b. Lower is better.

Provenance

Ported from owes_sample in R scoringRules (scores_sample_multiv_weighted.R; Allen 2024, JSS 110(8)).

Example

julia

using ScoringRules
X = randn(2, 50)
y = [0.0, 0.0]
owes(X, y; a = -1.0, b = 1.0)

ScoringRules.owmmds Function

julia

owmmds(
    X::AbstractMatrix,
    y::AbstractVector;
    a,
    b,
    weight_func,
    w
) -> Any

julia

owmmds(X, y; a=-Inf, b=Inf, weight_func=nothing, w=nothing)

Outcome-weighted MMD score (Gaussian kernel, σ = 1) of the ensemble X (a d × m matrix) at the d-dimensional observation y.

The MMD score is computed using the normalised per-column weights w(Xᵢ), then multiplied by w(y). Optional member weights w (length m) multiply the outcome weights, as in R. Returns NaN when all combined weights are zero. See owes for conventions on a, b, and weight_func. Lower is better.

Provenance

Ported from owmmds_sample in R scoringRules (scores_sample_multiv_weighted.R; Allen 2024, JSS 110(8)).

Example

julia

using ScoringRules
X = randn(2, 50)
y = [0.0, 0.0]
owmmds(X, y; a = -1.0, b = 1.0)

ScoringRules.owvs Function

julia

owvs(
    X::AbstractMatrix,
    y::AbstractVector;
    p,
    a,
    b,
    weight_func,
    w
) -> Any

julia

owvs(X, y; p=0.5, a=-Inf, b=Inf, weight_func=nothing, w=nothing)

Outcome-weighted variogram score of order p of the ensemble X (a d × m matrix) at the d-dimensional observation y.

The variogram score is computed using the normalised per-column weights w(Xᵢ), then multiplied by w(y). Optional member weights w (length m) multiply the outcome weights, as in R. Returns NaN when all combined weights are zero. See owes for conventions on a, b, and weight_func. Lower is better.

Provenance

Ported from owvs_sample in R scoringRules (scores_sample_multiv_weighted.R; Allen 2024, JSS 110(8)).

Example

julia

using ScoringRules
X = randn(2, 50)
y = [0.0, 0.0]
owvs(X, y; a = -1.0, b = 1.0)

ScoringRules.quantile_score Function

julia

quantile_score(
    q_levels::AbstractVector,
    q_forecasts::AbstractVector,
    y::Real
) -> Any

julia

quantile_score(q_levels, q_forecasts, y)

Quantile score (pinball loss) for a vector of quantile levels and the corresponding quantile forecasts, evaluated at observation y.

Arguments

q_levels — AbstractVector of quantile levels α ∈ (0, 1).
q_forecasts — AbstractVector of the corresponding quantile forecast values; must have the same length as q_levels.
y — scalar observation.

Returns

A Vector{Float64} of per-level scores, one entry per element of q_levels. To obtain a single summary value use mean(quantile_score(...)).

The score at level α with forecast q and observation y is

julia

score_α = ((y < q) - α) * (q - y)

which equals (1 - α)(q - y) when y < q and α(y - q) when y ≥ q (Koenker & Bassett 1978). Lower values indicate better forecasts.

Example

julia

using ScoringRules
levels = [0.1, 0.5, 0.9]
q = [-1.28, 0.0, 1.28]
quantile_score(levels, q, 1.0)

julia

quantile_score(dat::AbstractVector, y::Real; alpha, type)

julia

quantile_score(dat, y; alpha, type=7)

Quantile score for a single quantile level alpha, where the quantile is estimated from an ensemble dat. Mirrors qs_sample in scores_quantiles.R.

Arguments

dat — AbstractVector of ensemble draws.
y — scalar observation.
alpha — quantile level α ∈ (0, 1).
type — quantile interpolation type (1–9), passed to quantile; default 7 matches R's default.

Returns

A scalar Float64 score.

The empirical α-quantile q̂ is computed from dat, and the score is

julia

score = ((y < q̂) - α) * (q̂ - y)

Example

julia

using ScoringRules
dat = randn(100)
quantile_score(dat, 0.5; alpha = 0.9)

ScoringRules.rps Function

julia

rps(p::AbstractVector, y::Integer) -> Any

julia

rps(p, y)

Ranked Probability Score (RPS; Epstein 1969) for a categorical / ordinal forecast.

Arguments

p — AbstractVector of forecast probabilities over K ≥ 2 ordered outcome categories. Must be non-negative and sum to 1.
y — observed outcome category as an Integer in {1, …, K}, where y = 1 is the smallest and y = K is the largest category.

Returns

A scalar Float64 score. Lower values indicate better forecasts.

The RPS is defined as

julia

RPS = Σ_{k=1}^{K} (P_k - 1{y ≤ k})²

where P_k = Σ_{j=1}^{k} p_j is the cumulative forecast probability through category k. This sums K Brier scores across all category boundaries, reflecting the ordinal structure of the outcome.

The convention matches rps_probs in rps.R (Jordan, Krüger, Lerch, Allen): outcome y is 1-indexed and the sum runs over all K categories (the last term is always zero when probabilities sum to 1, but is retained for consistency with the R source).

Example

julia

using ScoringRules
p = [0.3, 0.2, 0.5]
rps(p, 2)

ScoringRules.twcrps Function

julia

twcrps(
    dat::AbstractVector,
    y::Real;
    a,
    b,
    chain_func,
    w
) -> Any

julia

twcrps(dat, y; a=-Inf, b=Inf, chain_func=nothing, w=nothing)

Threshold-weighted CRPS of the ensemble dat (a vector of m simulation draws) at observation y, emphasising outcomes in the interval (a, b).

The chaining function v(z) = clamp(z, a, b) is applied to both y and every ensemble member; the standard EDF-based CRPS of the transformed forecast is returned. Alternatively, supply a custom vectorised chain_func (takes a Real, returns a Real); supplying chain_func ignores a and b.

Lower is better.

Arguments

dat: ensemble of simulation draws.
y: scalar observation.

Keyword Arguments

a: lower threshold (default -Inf).
b: upper threshold (default Inf).
chain_func: custom chaining function; overrides a and b when supplied.
w: optional non-negative member weights (length m); normalised to sum to one internally.

Provenance

Ported from twcrps_sample in R scoringRules (scores_sample_univ_weighted.R; Allen 2024, JSS 110(8)).

Example

julia

using ScoringRules
dat = randn(100)
twcrps(dat, 0.5; a = 0.0, b = 1.0)

ScoringRules.twes Function

julia

twes(
    X::AbstractMatrix,
    y::AbstractVector;
    a,
    b,
    chain_func,
    w
) -> Any

julia

twes(X, y; a=-Inf, b=Inf, chain_func=nothing, w=nothing)

Threshold-weighted energy score of the ensemble X (a d × m matrix, each column one member) at the d-dimensional observation y.

The chaining function v(z) = clamp.(z, a, b) (element-wise) is applied to y and every column of X; the standard energy score of the transformed forecast is returned. a and b may be scalars (broadcast to all dimensions) or length-d vectors. A custom chain_func (takes a length-d vector, returns a length-d vector) overrides a and b. Optional member weights w (length m) are normalised to sum to one internally.

Lower is better.

Provenance

Ported from twes_sample in R scoringRules (scores_sample_multiv_weighted.R; Allen 2024, JSS 110(8)).

Example

julia

using ScoringRules
X = randn(2, 50)
y = [0.0, 0.0]
twes(X, y; a = -1.0, b = 1.0)

ScoringRules.twmmds Function

julia

twmmds(
    X::AbstractMatrix,
    y::AbstractVector;
    a,
    b,
    chain_func,
    w
) -> Any

julia

twmmds(X, y; a=-Inf, b=Inf, chain_func=nothing, w=nothing)

Threshold-weighted MMD score (Gaussian kernel, σ = 1) of the ensemble X (a d × m matrix) at the d-dimensional observation y.

The chaining function is applied to y and each column of X; the standard MMD score is evaluated on the transformed forecast. Optional member weights w (length m) are normalised to sum to one internally. See twes for conventions on a, b, and chain_func. Lower is better.

Provenance

Ported from twmmds_sample in R scoringRules (scores_sample_multiv_weighted.R; Allen 2024, JSS 110(8)).

Example

julia

using ScoringRules
X = randn(2, 50)
y = [0.0, 0.0]
twmmds(X, y; a = -1.0, b = 1.0)

ScoringRules.twvs Function

julia

twvs(
    X::AbstractMatrix,
    y::AbstractVector;
    p,
    a,
    b,
    chain_func,
    w
) -> Any

julia

twvs(X, y; p=0.5, a=-Inf, b=Inf, chain_func=nothing, w=nothing)

Threshold-weighted variogram score of order p of the ensemble X (a d × m matrix) at the d-dimensional observation y.

The chaining function is applied to y and each column of X; the standard variogram score is then evaluated on the transformed forecast. Optional member weights w (length m) are normalised to sum to one internally. See twes for the conventions on a, b, and chain_func. Lower is better.

Provenance

Ported from twvs_sample in R scoringRules (scores_sample_multiv_weighted.R; Allen 2024, JSS 110(8)).

Example

julia

using ScoringRules
X = randn(2, 50)
y = [0.0, 0.0]
twvs(X, y; a = -1.0, b = 1.0)

ScoringRules.vs Function

julia

vs(X::AbstractMatrix, y::AbstractVector; p, w, w_vs) -> Any

julia

vs(X, y; p=0.5, w=nothing, w_vs=nothing)

Variogram score of order p of the ensemble forecast X (a d × m matrix) at the d-dimensional observation y:

w is an optional per-member weight vector (length m, as in es); it is normalised to sum to one internally. w_vs is an optional d × d non-negative symmetric pairwise weight matrix (R's w_vs; defaults to all ones). When w is supplied the score uses the weighted-ensemble form of the variogram kernel. R's vs_sample silently ignores w_vs when w is given; here both are honoured together. Lower is better.

Provenance

Ported from vs_sample in R scoringRules (scores_sample_multiv.R; procs_es.cpp; Jordan, Krüger, Lerch, Allen). Scheuerer & Hamill (2015), MWR 143, 1321–1334.

Example

julia

using ScoringRules
X = randn(2, 50)
y = [0.0, 0.0]
vs(X, y)