c4h_verify

Description

Computes quality assessment metrics of climate forecasts.

c4h_verify(
  exp,
  obs,
  metrics,
  ref = NULL,
  brier_thresholds = c(0.1, 0.9),
  prob_thresholds = c(1/3, 2/3),
  alpha = 0.05,
  ncores = NULL,
  N.eff = NA,
  indices_for_clim = NULL,
  cross.val = FALSE,
  clim.cross.val = TRUE,
  weights_exp = NULL,
  weights_ref = NULL,
  comp_dim = NULL,
  limits = NULL,
  conf = TRUE,
  na.rm = FALSE,
  sign = TRUE,
  pval = FALSE,
  sig_method = NULL,
  sig_method.type = "two.sided.approx",
  rocss_cat = 1
)

Arguments

  • exp: An s2dv_cube object of forecast data with at least the ‘sdate’ and ‘member’ dimensions.
  • obs: An s2dv_cube object of observation data with at least the ‘sdate’ and ‘member’ dimensions.
  • metrics: A character vector of the verification metrics to be calculated. Valid options are: “BSS”, “RPSS”, “CRPSS”, “AbsBiasSS”, “MSE”, “MSSS”, “RMSE”, “RMSSS”, “ROCSS”.
  • ref: An optional s2dv_cube object of reference data with at least the ‘sdate’ and ‘ensemble’ dimensions. The dimensions must be the same as ‘exp’ except ‘ensemble’. If ‘ref’ is NULL, the climatological forecast is used as reference forecast. The default value is NULL.
  • brier_thresholds: Optional numeric vector of thresholds for the Brier Skill Score (BSS) calculation. The default is c(0.1, 0.9). Only applicable if “BSS” is in metrics.
  • prob_thresholds: Optional numeric vector of thresholds for the Ranked Probability Skill Score (RPSS) calculation, and for the Relative Operating Characteristic Skill Score (ROCSS). The default is c(1/3, 2/3). Only applicable if “RPSS” or “ROCSS” is in metrics.
  • alpha: A numeric of the significance level to be used in the statistical significance test. The default value is 0.05. Relevant for all metrics except “ROCSS”. In some cases, a warning will be returned if alpha and sig_method.type are not compatible (e.g. if sig_method = “RandomWalk”, sig_method.type is “two.sided.approx” and alpha is not equal to 0.05), and the significance level of the results will be indicated.
  • ncores: An integer indicating the number of cores to use for parallel computation. The default value is NULL.
  • N.eff: Effective sample size to be used in the statistical significance test. It can be NA (and it will be computed with the s2dv:::.Eno), FALSE (and it will use the length of “obs” along “time”, so the autocorrelation is not taken into account), a numeric (which is used for all cases), or an array with the same dimensions as “obs” except “time” (for a particular N.eff to be used for each case). The default value is NA. Only relevant if “BSS”, “RPSS”, “CRPSS”, or “RMSSS” is in metrics.
  • indices_for_clim: A vector of the indices to be taken along ‘sdate’ for computing the thresholds between the probabilistic categories. If NULL, the whole period is used. The default value is NULL. Only relevant if “BSS”, “RPSS”, or “ROCSS” is in metrics.
  • cross.val: A logical indicating whether to compute the thresholds between probabilistics categories in cross-validation. The default value is FALSE. Only relevant if “BSS”, “RPSS”, or “ROCSS” is in metrics.
  • clim.cross.val: A logical indicating whether to build the climatological forecast in cross-validation (i.e. excluding the observed value of the time step when building the probabilistic distribution function for that particular time step). The default value is TRUE. Only used if ‘ref’ is NULL and “CRPSS” is in metrics.
  • weights_exp: A named numerical array of the forecast ensemble weights for probability calculation. The dimensions should include ‘ensemble’, ‘sdate’, and ‘dataset’ if there are multiple datasets. All dimension lengths must be equal to ‘exp’ dimension lengths. The default value is NULL, which means no weighting is applied. The ensemble should have at least 70 members or span at least 10 time steps and have more than 45 members if consistency between the weighted and unweighted methodologies is desired. Only relevant if “BSS” or “RPSS” is in metrics.
  • weights_ref: Same as ‘weights_exp’ but for the reference forecast.
  • comp_dim: A character string indicating the name of dimension along which obs is taken into account only if it is complete. The default value is NULL. Only relevant if “MSE” or “RMSE” is in metrics.
  • limits: A vector of two integers indicating the range along comp_dim to be completed. The default value is c(1, length(comp_dim dimension)). Only relevant if “MSE” or “RMSE” is in metrics.
  • conf: A logical value indicating whether to retrieve the confidence intervals or not. The default value is TRUE. Only relevant if “MSE” or “RMSE” is in metrics.
  • na.rm: A logical. FALSE means that NA values will not be removed, and the function will return ‘NA’ if they are present. TRUE means that NA values will be removed before computation. The default value is FALSE. Only relevant if “BSS”, “RPSS”, or “AbsBiasSS” is in metrics.
  • sign: A logical value indicating whether to compute or not the statistical significance of the test Ho: (R)MSSS = 0. The default value is TRUE. Only relevant if “MSSS” or “RMSSS” is in metrics.
  • pval: A logical value indicating whether to compute or not the p-value of the test Ho: (R)MSSS = 0. The default value is FALSE. Only relevant if “MSSS” or “RMSSS” is in metrics.
  • sig_method: A character string indicating the significance method. The options are NULL, “one-sided Fisher”, “RandomWalk”, and “DieboldMariano”. Relevant for metrics “BSS”, “RPSS”, “CRPSS”, “AbsBiasSS”, “MSSS” and “RMSSS”. The default is NULL. In this case, “RandomWalk” is applied for “BSS”, “RPSS”, “CRPSS”, and “AbsBiasSS”, and “one-sided Fisher” is applied for “MSSS” and “RMSSS”.
  • sig_method.type: A character string indicating the test type of the significance method if sig_method is “RandomWalk” or “DieboldMariano”. Valid options are “two.sided.approx”, “two.sided”, “greater”, “less”. The default is “two.sided.approx”.
  • rocss_cat: An integer indicating the category to be returned when “ROCSS” is in metrics. The default value is 1, which means that the ROCSS for the first probability category (i.e. the lowest category) will be returned.

Returns

a list of lists with the calculated metrics and their significances as s2dv_cube objects.

Note

See clim4health_verification.Rmd in the package vignettes for a detailed tutorial on function arguments and usage.

See Also

c4h_plotskill for visualizing the forecast skill.

Examples

# DISCLAIMER: The purpose of these examples is to show the functionality
# of the different functions. To perform a meaningful verification, a much
# larger time series is required.
# Load data
hcst <- c4h_load(
  system.file("extdata", "hindcast", package = "clim4health"),
  variable = "t2m",
  year = 2010:2012,
  month = 1,
  leadtime_month = 1:3,
  ext = "nc"
)
lats <- hcst$coords$latitude
lons <- hcst$coords$longitude
obs <- c4h_load(
  system.file("extdata", "reanalysis", package = "clim4health"),
  variable = "t2m",
  year = 2010:2012,
  month = 1,
  leadtime_month = 1:3,
  bbox = c(max(lats), min(lons), min(lats), max(lons)),
  ext = "nc"
)

# Example with gridded data
dwn <- c4h_downscale(
  "Interpolation",
  exp = hcst,
  obs = obs,
  method_remap = "bilinear",
  target_grid = system.file(
    "extdata",
    "reanalysis/t2m_201001.nc",
    package = "clim4health"
  )
)
hcst_dwn <- dwn$exp
skill_grid <- c4h_verify(exp = hcst_dwn, obs = obs, metrics = "BSS")
c4h_plotskill(skill_grid$BSS10$bss, skill_grid$BSS10$sign)

# Example with polygon data
areas <- sf::read_sf(system.file("extdata", "areas", "munip_vallecauca.gpkg",
                                 package = "clim4health"))
obs_area <- c4h_space(obs,
                      areas = areas,
                      fun = "mean",
                      areas_id = "munip_code")
hcst_area <- c4h_space(hcst,
                       areas = areas,
                       fun = "mean",
                       areas_id = "munip_code")
skill_area <- c4h_verify(exp = hcst_area,
                         obs = obs_area,
                         metrics = "CRPSS",
                         brier_thresholds = c(0.1, 0.9))
c4h_plotskill(skill_area$CRPSS$crpss, skill_area$CRPSS$sign)

# Example with station data
point_locs <- clim4health::c4h_load(system.file("extdata",
                                               "stations",
                                               package = "clim4health"),
                                    ext = "csv",
                                    year = 2010:2012,
                                    month = 1,
                                    leadtime_month = 1:3,
                                    variable = "temp_mean")
points <- list(latitude = point_locs$attrs$location$latitude,
               longitude = point_locs$attrs$location$longitude)
point_locs <- c4h_time(point_locs,
                       time_aggregation = "monthly",
                       dim_aggregation = "sdate",
                       fun = "mean")
dwn <- c4h_downscale("Interpolation", exp = hcst,
                     points = points, method_point_interp = "bilinear")
hcst_dwn <- dwn$exp
skill_point <- c4h_verify(exp = hcst_dwn,
                          obs = point_locs,
                          metrics = "BSS")
c4h_plotskill(skill = skill_point$BSS10$bss, boundaries = areas)