I would like to share with you some doubts that I have for the Hands-on 2. Below there's only the code and then there are some questions:
Code
# Step (0): Source the packages
library(CSTools) # CST_Calibration
library(startR) # Start
library(s2dv)
# Step (0): Load the data
# Start call [...]
# Step (1): Create 's2dv_cube'
# hcst <- as.s2dv_cube(hcst)
# fcst <- as.s2dv_cube(fcst)
# obs <- as.s2dv_cube(obs)
load("/esarchive/scratch/erifarov/rpackages/PATC2023/data-ho2/data_ho2.Rdata")
# Step (2): Calibrate the data
source("https://earth.bsc.es/gitlab/external/cstools/-/raw/master/R/CST_Calibration.R")
library(multiApply)
hcst_cal <- CST_Calibration(exp = hcst, obs = obs,
cal.method = "evmos",
eval.method = "leave-one-out",
multi.model = FALSE,
na.fill = TRUE,
na.rm = TRUE,
apply_to = NULL,
alpha = NULL,
memb_dim = "ensemble",
sdate_dim = "syear",
ncores = 10)
fcst_cal <- CST_Calibration(exp = hcst, obs = obs, exp_cor = fcst,
cal.method = "evmos",
eval.method = "leave-one-out",
multi.model = FALSE,
na.fill = TRUE,
na.rm = TRUE,
apply_to = NULL,
alpha = NULL,
memb_dim = "ensemble",
sdate_dim = "syear",
ncores = 10)
# Step (3.1): Compute anomalies from hindcast and obs
source("https://earth.bsc.es/gitlab/external/cstools/-/raw/master/R/CST_Anomaly.R")
source("https://earth.bsc.es/gitlab/es/s2dv/-/raw/master/R/Ano_CrossValid.R")
hcst_anom <- CST_Anomaly(exp = hcst_cal, obs = obs,
cross = TRUE,
memb = TRUE,
memb_dim = 'ensemble',
dim_anom = 'syear',
dat_dim = c('dat', 'ensemble'),
ftime_dim = 'time',
ncores = 10)
# Step (3.2): Compute forecast anomalies
clim <- s2dv::Clim(exp = hcst_cal$data, obs = obs$data,
time_dim = "syear",
dat_dim = c("dat", "ensemble"),
memb = FALSE,
memb_dim = "ensemble",
ftime_dim = "time",
ncores = 10)
clim_hcst <- InsertDim(clim$clim_exp, posdim = 1, lendim = 1, name = "syear")
dims <- dim(clim_hcst)
clim_hcst <- rep(clim_hcst, fcst$dim[['ensemble']])
dim(clim_hcst) <- c(dims, ensemble = fcst$dim[['ensemble']])
clim_hcst <- Reorder(clim_hcst, order = names(fcst$dim))
fcst_anom <- fcst_cal$data - clim_hcst
# Step (4): Compute RPSS
skill_data <- RPSS(exp = hcst_cal$data, obs = obs$data, memb_dim = 'ensemble',
time_dim = 'syear', Fair = FALSE, cross.val = TRUE,
ncores = 10)
skill_anom <- RPSS(exp = hcst_anom$exp$data, obs = hcst_anom$obs$data,
time_dim = 'syear',
memb_dim = 'ensemble',
Fair = FALSE,
cross.val = TRUE,
ncores = 10)
Questions
- I've decided to leave the loading part with Start (it will be changed to CST_Start with CSTools new release) if someone wants to run the exercise alone. Do you agree with this @aho? Is it mentioned that if the data is already loaded, go directly to the next step, that it is how it will be in the PATC 2023.
- For the Anomalies section (Step 3.2), I've copied from SUNSET the part of the forecast anomalies. Do you think it's necessary to compute it?
- For the Skill section (Step 4), now there is the RPSS computation for anomalies and for the data. Is it good to have both? Also, do you know if there is a way to include the RPSS using also the forecast? Now, it's not possible since the dimensions of the data are:
> dim(hcst$data)
dat var syear ensemble time latitude longitude
1 1 24 25 2 60 60
> dim(fcst$data)
dat var syear ensemble time latitude longitude
1 1 1 51 2 60 60
- Something missing?
Thanks in advance,
Eva