test-probabilist-bootstrapping.R

context("test bootstrapping for probabilistic evaluation")
library(testthat)
library(R.utils)
library(evalhyd)

all_metrics = c(
  # threshold-based
  'BS', 'BSS', 'BS_CRD', 'BS_LBD', 'REL_DIAG', 'CRPS_FROM_BS', 
  # CDF-based
  'CRPS_FROM_ECDF',
  # quantile-based
  'QS', 'CRPS_FROM_QS',
  # contingency table-based
  'POD', 'POFD', 'FAR', 'CSI', 'ROCSS',
  # ranks-based
# /!\ EXCLUDE BECAUSE OF RANDOM ELEMENT TO DEAL WITH TIES ----------------------
# 'RANK_HIST', 'DS', 'AS',
# ------------------------------------------------------------------------------
  # intervals
  'CR', 'AW', 'AWN', 'WS',
  # multivariate
  'ES'
)

# load some predicted and observed streamflow
dts_1yr = array(
  data=unlist(
    read.csv("./data/q_obs_1yr.csv", header=FALSE, nrows=1)
  ), 
  dim=c(366)
)
obs_1yr = array(
  data=unlist(
    read.csv("./data/q_obs_1yr.csv", header=TRUE, colClasses="numeric")
  ), 
  dim=c(1, 366)
)
prd_1yr = array(
  data=unlist(
    read.csv("./data/q_prd_1yr.csv", header=TRUE, colClasses="numeric")
  ),
  dim=c(50, 366)
)

thr = array(
  c(690, 534, 445), 
  dim=c(1, 3)
)

lvl = c(30., 80.)

# replicate year of data three times
obs_3yrs = cbind(obs_1yr, obs_1yr, obs_1yr)
prd_3yrs = cbind(prd_1yr, prd_1yr, prd_1yr)

# add dimensions for sites and lead times
dim(prd_1yr) = c(1,1,50,366)
dim(prd_3yrs) = c(1,1,50,366*3)


# ------------------------------------------------------------------------------
# compare bootstrap of three block years with three identical years
# ------------------------------------------------------------------------------

for (metric in all_metrics)
{
  testthat::test_that(
    metric,
    {
      # bootstrap with only one year of data