From 208ed9743b1a5ceb4501ee6081d7654c25e46b52 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Wed, 29 Jun 2022 16:12:54 +0200 Subject: [PATCH 01/47] Allow dat_dim = NULL in Corr() function --- R/Corr.R | 39 +++++++++++++++++++++++++++------------ 1 file changed, 27 insertions(+), 12 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index 0382a39..3cf640d 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -23,7 +23,8 @@ #'@param time_dim A character string indicating the name of dimension along #' which the correlations are computed. The default value is 'sdate'. #'@param dat_dim A character string indicating the name of dataset (nobs/nexp) -#' dimension. The default value is 'dataset'. +#' dimension. The default value is 'dataset'. If there is no dataset +#' dimension, set NULL. #'@param 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. @@ -51,7 +52,8 @@ #' c(nexp, nobs, exp_memb, obs_memb, all other dimensions of exp except #' time_dim and memb_dim).\cr #'nexp is the number of experiment (i.e., 'dat_dim' in exp), and nobs is the -#'number of observation (i.e., 'dat_dim' in obs). exp_memb is the number of +#'number of observation (i.e., 'dat_dim' in obs). If +#' dat_dim is NULL, nexp and nobs are omitted. exp_memb is the number of #'member in experiment (i.e., 'memb_dim' in exp) and obs_memb is the number of #'member in observation (i.e., 'memb_dim' in obs).\cr\cr #'\item{$corr}{ @@ -129,11 +131,14 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + "Set it as NULL if there is no dataset dimension.") + } } ## comp_dim if (!is.null(comp_dim)) { @@ -194,8 +199,10 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] + if (!is.null(dat_dim)) { + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] + } if (!is.null(memb_dim)) { name_exp <- name_exp[-which(name_exp == memb_dim)] name_obs <- name_obs[-which(name_obs == memb_dim)] @@ -327,15 +334,19 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', return(res) } -.Corr <- function(exp, obs, time_dim = 'sdate', method = 'pearson', +.Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', method = 'pearson', conf = TRUE, pval = TRUE, conf.lev = 0.95, ncores_input = NULL) { if (length(dim(exp)) == 2) { # exp: [sdate, dat_exp] # obs: [sdate, dat_obs] + if (is.null(dat_dim)) { + nexp <- 1 + nobs <- 1 + } else { nexp <- as.numeric(dim(exp)[2]) nobs <- as.numeric(dim(obs)[2]) - + } # NOTE: Use sapply to replace the for loop CORR <- sapply(1:nobs, function(i) { sapply(1:nexp, function (x) { @@ -356,11 +367,15 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # exp: [sdate, dat_exp, memb_exp] # obs: [sdate, dat_obs, memb_obs] + if (is.null(dat_dim)) { + nexp <- 1 + nobs <- 1 + } else { nexp <- as.numeric(dim(exp)[2]) nobs <- as.numeric(dim(obs)[2]) exp_memb <- as.numeric(dim(exp)[3]) obs_memb <- as.numeric(dim(obs)[3]) - + } CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) for (j in 1:obs_memb) { -- GitLab From 5743f5818bb9aa785472830ab1b54c97a64aff01 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Wed, 29 Jun 2022 17:18:10 +0200 Subject: [PATCH 02/47] Allow dat_dim = NULL in Consist_Trend() function --- R/Consist_Trend.R | 21 ++++++++++++++------- 1 file changed, 14 insertions(+), 7 deletions(-) diff --git a/R/Consist_Trend.R b/R/Consist_Trend.R index b02aa5f..5ac797e 100644 --- a/R/Consist_Trend.R +++ b/R/Consist_Trend.R @@ -14,7 +14,8 @@ #'@param dat_dim A character string indicating the name of the dataset #' dimensions. If data at some point of 'time_dim' are not complete along #' 'dat_dim' in both 'exp' and 'obs', this point in all 'dat_dim' will be -#' discarded. The default value is 'dataset'. +#' discarded. The default value is 'dataset'. If there is no dataset +#' dimension, set NULL. #'@param time_dim A character string indicating the name of dimension along #' which the trend is computed. The default value is 'sdate'. #'@param interval A positive numeric indicating the unit length between two @@ -29,7 +30,8 @@ #' with dimensions c(stats = 2, nexp + nobs, the rest dimensions of 'exp' and #' 'obs' except time_dim), where 'nexp' is the length of 'dat_dim' in 'exp' #' and 'nobs' is the length of 'dat_dim' in 'obs. The 'stats' dimension -#' contains the intercept and the slope. +#' contains the intercept and the slope. If dat_dim is NULL, nexp and nobs are +#' omitted. #'} #'\item{$conf.lower}{ #' A numeric array of the lower limit of 95\% confidence interval with @@ -109,19 +111,24 @@ Consist_Trend <- function(exp, obs, dat_dim = 'dataset', time_dim = 'sdate', int stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.character(dat_dim)) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!all(dat_dim %in% names(dim(exp))) | !all(dat_dim %in% names(dim(obs)))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + "Set it as NULL if there is no dataset dimension.") + } } ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) + if (!is.null(dat_dim)) { for (i in 1:length(dat_dim)) { name_exp <- name_exp[-which(name_exp == dat_dim[i])] name_obs <- name_obs[-which(name_obs == dat_dim[i])] } + } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", "all dimension expect 'dat_dim'.")) -- GitLab From c9b969fd231bb1bfc61023d31e51d550685ba9e8 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Wed, 29 Jun 2022 17:26:18 +0200 Subject: [PATCH 03/47] Allow dat_dim = NULL in RMS() function --- R/RMS.R | 18 ++++++++++++------ 1 file changed, 12 insertions(+), 6 deletions(-) diff --git a/R/RMS.R b/R/RMS.R index b3c8ad4..194a3f6 100644 --- a/R/RMS.R +++ b/R/RMS.R @@ -21,7 +21,8 @@ #'@param time_dim A character string indicating the name of dimension along #' which the correlations are computed. The default value is 'sdate'. #'@param dat_dim A character string indicating the name of member (nobs/nexp) -#' dimension. The default value is 'dataset'. +#' dimension. The default value is 'dataset'. If there is no dataset +#' dimension, set NULL. #'@param 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. @@ -107,11 +108,14 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + "Set it as NULL if there is no dataset dimension.") + } } ## comp_dim if (!is.null(comp_dim)) { @@ -151,8 +155,10 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) + if (!is.null(dat_dim)) { name_exp <- name_exp[-which(name_exp == dat_dim)] name_obs <- name_obs[-which(name_obs == dat_dim)] + } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", "all dimension expect 'dat_dim'.")) -- GitLab From 29f39928b9e03057faefa51446834c3cffe5ae2f Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Wed, 29 Jun 2022 17:29:21 +0200 Subject: [PATCH 04/47] Allow dat_dim = NULL in RMSSS() function --- R/RMSSS.R | 18 ++++++++++++------ 1 file changed, 12 insertions(+), 6 deletions(-) diff --git a/R/RMSSS.R b/R/RMSSS.R index 5fa9659..3c3e006 100644 --- a/R/RMSSS.R +++ b/R/RMSSS.R @@ -23,7 +23,8 @@ #' 'exp', then the vector will automatically be 'time_dim' and 'dat_dim' will #' be 1. #'@param dat_dim A character string indicating the name of dataset (nobs/nexp) -#' dimension. The default value is 'dataset'. +#' dimension. The default value is 'dataset'. If there is no dataset +#' dimension, set NULL. #'@param time_dim A character string indicating the name of dimension along #' which the RMSSS are computed. The default value is 'sdate'. #'@param pval A logical value indicating whether to compute or not the p-value @@ -96,11 +97,14 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + "Set it as NULL if there is no dataset dimension.") + } } ## pval if (!is.logical(pval) | length(pval) > 1) { @@ -116,8 +120,10 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) + if (!is.null(dat_dim)) { name_exp <- name_exp[-which(name_exp == dat_dim)] name_obs <- name_obs[-which(name_obs == dat_dim)] + } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", "all dimension expect 'dat_dim'.")) -- GitLab From 846c8d6cf55381d0c318e411d2fb22c12a61185b Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Wed, 29 Jun 2022 17:33:38 +0200 Subject: [PATCH 05/47] Allow dat_dim = NULL in UlitmateBrier() function --- R/UltimateBrier.R | 15 ++++++++++----- 1 file changed, 10 insertions(+), 5 deletions(-) diff --git a/R/UltimateBrier.R b/R/UltimateBrier.R index aeaddcd..9e7698f 100644 --- a/R/UltimateBrier.R +++ b/R/UltimateBrier.R @@ -102,11 +102,14 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti stop("Parameter 'exp' and 'obs' must have dimension names.") } ## dat_dim - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + "Set it as NULL if there is no dataset dimension.") + } } ## memb_dim if (!is.character(memb_dim) | length(memb_dim) > 1) { @@ -133,8 +136,10 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti name_obs <- sort(names(dim(obs))) name_exp <- name_exp[-which(name_exp == memb_dim)] name_obs <- name_obs[-which(name_obs == memb_dim)] + if (!is.null(dat_dim)) { name_exp <- name_exp[-which(name_exp == dat_dim)] name_obs <- name_obs[-which(name_obs == dat_dim)] + } if (any(name_exp != name_obs)) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", "of all the dimensions expect 'dat_dim' and 'memb_dim'.")) -- GitLab From 37a85620434cd3fa7f8ae97df952d5290a41a5bd Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 30 Jun 2022 13:19:17 +0200 Subject: [PATCH 06/47] Corrected typo of output text - expect for except - in some functions --- R/ACC.R | 2 +- R/BrierScore.R | 10 +++++----- R/Clim.R | 2 +- R/Consist_Trend.R | 2 +- R/Corr.R | 2 +- R/NAO.R | 2 +- R/RMS.R | 2 +- R/RMSSS.R | 2 +- R/RPS.R | 2 +- R/RPSS.R | 4 ++-- R/RatioRMS.R | 2 +- R/RatioSDRMS.R | 2 +- R/ResidualCorr.R | 2 +- R/UltimateBrier.R | 4 ++-- 14 files changed, 20 insertions(+), 20 deletions(-) diff --git a/R/ACC.R b/R/ACC.R index 44c724c..ceb1e94 100644 --- a/R/ACC.R +++ b/R/ACC.R @@ -276,7 +276,7 @@ ACC <- function(exp, obs, dat_dim = 'dataset', lat_dim = 'lat', lon_dim = 'lon', } if (!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all the dimensions expect 'dat_dim' and 'memb_dim'.")) + "all the dimensions except 'dat_dim' and 'memb_dim'.")) } #----------------------------------------------------------------- diff --git a/R/BrierScore.R b/R/BrierScore.R index c3673ad..22f497d 100644 --- a/R/BrierScore.R +++ b/R/BrierScore.R @@ -57,13 +57,13 @@ #''bs_check_gres', 'bss_gres', 'rel_bias_corrected', 'gres_bias_corrected', #''unc_bias_corrected', and 'bss_bias_corrected' are (a) a number (b) an array #'with dimensions c(nexp, nobs, all the rest dimensions in 'exp' and 'obs' -#'expect 'time_dim' and 'memb_dim') (c) an array with dimensions of +#'except 'time_dim' and 'memb_dim') (c) an array with dimensions of #''exp' and 'obs' except 'time_dim' and 'memb_dim'\cr #'Items 'nk', 'fkbar', and 'okbar' are (a) a vector of length of bin number #'determined by 'threshold' (b) an array with dimensions c(nexp, nobs, -#'no. of bins, all the rest dimensions in 'exp' and 'obs' expect 'time_dim' and +#'no. of bins, all the rest dimensions in 'exp' and 'obs' except 'time_dim' and #''memb_dim') (c) an array with dimensions c(no. of bin, all the rest dimensions -#'in 'exp' and 'obs' expect 'time_dim' and 'memb_dim') +#'in 'exp' and 'obs' except 'time_dim' and 'memb_dim') #' #'@references #'Wilks (2006) Statistical Methods in the Atmospheric Sciences.\cr @@ -176,11 +176,11 @@ BrierScore <- function(exp, obs, thresholds = seq(0.1, 0.9, 0.1), time_dim = 'sd } if (any(!name_exp %in% name_obs) | any(!name_obs %in% name_exp)) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.")) + "of all the dimensions except 'dat_dim' and 'memb_dim'.")) } if (!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.")) + "of all the dimensions except 'dat_dim' and 'memb_dim'.")) } ## ncores if (!is.null(ncores)) { diff --git a/R/Clim.R b/R/Clim.R index ce9ba2e..c144025 100644 --- a/R/Clim.R +++ b/R/Clim.R @@ -173,7 +173,7 @@ Clim <- function(exp, obs, time_dim = 'sdate', dat_dim = c('dataset', 'member'), } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have the same dimensions ", - "expect 'dat_dim'.")) + "except 'dat_dim'.")) } ############################### diff --git a/R/Consist_Trend.R b/R/Consist_Trend.R index 5ac797e..84b7699 100644 --- a/R/Consist_Trend.R +++ b/R/Consist_Trend.R @@ -131,7 +131,7 @@ Consist_Trend <- function(exp, obs, dat_dim = 'dataset', time_dim = 'sdate', int } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension expect 'dat_dim'.")) + "all dimension except 'dat_dim'.")) } ## interval if (!is.numeric(interval) | interval <= 0 | length(interval) > 1) { diff --git a/R/Corr.R b/R/Corr.R index 3cf640d..628511d 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -209,7 +209,7 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension expect 'dat_dim' and 'memb_dim'.")) + "all dimension except 'dat_dim' and 'memb_dim'.")) } if (dim(exp)[time_dim] < 3) { stop("The length of time_dim must be at least 3 to compute correlation.") diff --git a/R/NAO.R b/R/NAO.R index af4893a..6d48dba 100644 --- a/R/NAO.R +++ b/R/NAO.R @@ -194,7 +194,7 @@ NAO <- function(exp = NULL, obs = NULL, lat, lon, time_dim = 'sdate', } if (throw_error) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'memb_dim'.")) + "of all the dimensions except 'memb_dim'.")) } } ## ftime_avg diff --git a/R/RMS.R b/R/RMS.R index 194a3f6..20da981 100644 --- a/R/RMS.R +++ b/R/RMS.R @@ -161,7 +161,7 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension expect 'dat_dim'.")) + "all dimension except 'dat_dim'.")) } if (dim(exp)[time_dim] < 2) { stop("The length of time_dim must be at least 2 to compute RMS.") diff --git a/R/RMSSS.R b/R/RMSSS.R index 3c3e006..9526517 100644 --- a/R/RMSSS.R +++ b/R/RMSSS.R @@ -126,7 +126,7 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension expect 'dat_dim'.")) + "all dimension except 'dat_dim'.")) } if (dim(exp)[time_dim] <= 2) { stop("The length of time_dim must be more than 2 to compute RMSSS.") diff --git a/R/RPS.R b/R/RPS.R index 8ded53e..200c59c 100644 --- a/R/RPS.R +++ b/R/RPS.R @@ -88,7 +88,7 @@ RPS <- function(exp, obs, time_dim = 'sdate', memb_dim = 'member', if (!identical(length(name_exp), length(name_obs)) | !identical(dim(exp)[name_exp], dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimensions expect 'memb_dim'.")) + "all dimensions except 'memb_dim'.")) } ## prob_thresholds if (!is.numeric(prob_thresholds) | !is.vector(prob_thresholds) | diff --git a/R/RPSS.R b/R/RPSS.R index 3b24777..cbd87bd 100644 --- a/R/RPSS.R +++ b/R/RPSS.R @@ -108,7 +108,7 @@ RPSS <- function(exp, obs, ref = NULL, time_dim = 'sdate', memb_dim = 'member', if (!identical(length(name_exp), length(name_obs)) | !identical(dim(exp)[name_exp], dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimensions expect 'memb_dim'.")) + "all dimensions except 'memb_dim'.")) } if (!is.null(ref)) { name_ref <- sort(names(dim(ref))) @@ -116,7 +116,7 @@ RPSS <- function(exp, obs, ref = NULL, time_dim = 'sdate', memb_dim = 'member', if (!identical(length(name_exp), length(name_ref)) | !identical(dim(exp)[name_exp], dim(ref)[name_ref])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimensions expect 'memb_dim'.")) + "all dimensions except 'memb_dim'.")) } } ## prob_thresholds diff --git a/R/RatioRMS.R b/R/RatioRMS.R index f7e34b4..51f3984 100644 --- a/R/RatioRMS.R +++ b/R/RatioRMS.R @@ -21,7 +21,7 @@ #' computation. The default value is NULL. #' #'@return A list containing the numeric arrays with dimensions identical with -#' 'exp1', 'exp2', and 'obs', expect 'time_dim': +#' 'exp1', 'exp2', and 'obs', except 'time_dim': #'\item{$ratiorms}{ #' The ratio between the RMSE (i.e., RMSE1/RMSE2). #'} diff --git a/R/RatioSDRMS.R b/R/RatioSDRMS.R index 544ca6f..d527625 100644 --- a/R/RatioSDRMS.R +++ b/R/RatioSDRMS.R @@ -109,7 +109,7 @@ RatioSDRMS <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', name_obs <- name_obs[-which(name_obs == memb_dim)] if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all the dimensions expect 'dat_dim' and 'memb_dim'.")) + "all the dimensions except 'dat_dim' and 'memb_dim'.")) } ## pval if (!is.logical(pval) | length(pval) > 1) { diff --git a/R/ResidualCorr.R b/R/ResidualCorr.R index 7428d1b..8d9f040 100644 --- a/R/ResidualCorr.R +++ b/R/ResidualCorr.R @@ -125,7 +125,7 @@ ResidualCorr <- function(exp, obs, ref, N.eff = NA, time_dim = 'sdate', if (length(name_exp) != length(name_obs) | length(name_exp) != length(name_ref) | any(dim(exp)[name_exp] != dim(obs)[name_obs]) | any(dim(exp)[name_exp] != dim(ref)[name_ref])) { stop(paste0("Parameter 'exp', 'obs', and 'ref' must have same length of ", - "all dimensions expect 'memb_dim'.")) + "all dimensions except 'memb_dim'.")) } ## method if (!method %in% c("pearson", "kendall", "spearman")) { diff --git a/R/UltimateBrier.R b/R/UltimateBrier.R index 9e7698f..a2c594e 100644 --- a/R/UltimateBrier.R +++ b/R/UltimateBrier.R @@ -142,11 +142,11 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } if (any(name_exp != name_obs)) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.")) + "of all the dimensions except 'dat_dim' and 'memb_dim'.")) } if (!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.")) + "of all the dimensions except 'dat_dim' and 'memb_dim'.")) } ## quantile if (!is.logical(quantile) | length(quantile) > 1) { -- GitLab From cad3511398b10702b1d9867c4881631d1c0f7bc8 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 30 Jun 2022 16:59:15 +0200 Subject: [PATCH 07/47] Correct typo - except - also in tests --- tests/testthat/test-ACC.R | 2 +- tests/testthat/test-BrierScore.R | 2 +- tests/testthat/test-Clim.R | 2 +- tests/testthat/test-Consist_Trend.R | 2 +- tests/testthat/test-Corr.R | 2 +- tests/testthat/test-NAO.R | 4 ++-- tests/testthat/test-RMS.R | 2 +- tests/testthat/test-RMSSS.R | 2 +- tests/testthat/test-RPS.R | 2 +- tests/testthat/test-RPSS.R | 4 ++-- tests/testthat/test-RatioSDRMS.R | 2 +- tests/testthat/test-ResidualCorr.R | 2 +- tests/testthat/test-UltimateBrier.R | 4 ++-- 13 files changed, 16 insertions(+), 16 deletions(-) diff --git a/tests/testthat/test-ACC.R b/tests/testthat/test-ACC.R index 5c36725..04e60ad 100644 --- a/tests/testthat/test-ACC.R +++ b/tests/testthat/test-ACC.R @@ -152,7 +152,7 @@ test_that("1. Input checks", { obs = array(1:4, dim = c(dataset = 2, member = 2, lat = 1, lon = 1)), lat = c(1, 2), lon = c(1), avg_dim = NULL), - "Parameter 'exp' and 'obs' must have same length of all the dimensions expect 'dat_dim' and 'memb_dim'." + "Parameter 'exp' and 'obs' must have same length of all the dimensions except 'dat_dim' and 'memb_dim'." ) }) diff --git a/tests/testthat/test-BrierScore.R b/tests/testthat/test-BrierScore.R index 5668a08..e2c34f9 100644 --- a/tests/testthat/test-BrierScore.R +++ b/tests/testthat/test-BrierScore.R @@ -51,7 +51,7 @@ test_that("1. Input checks", { expect_error( BrierScore(exp3, obs3), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.") + "of all the dimensions except 'dat_dim' and 'memb_dim'.") ) # thresholds expect_error( diff --git a/tests/testthat/test-Clim.R b/tests/testthat/test-Clim.R index 7751afb..f5e288e 100644 --- a/tests/testthat/test-Clim.R +++ b/tests/testthat/test-Clim.R @@ -134,7 +134,7 @@ test_that("1. Input checks", { expect_error( Clim(array(1:10, dim = c(dataset = 2, member = 5, sdate = 4, ftime = 3)), array(1:4, dim = c(dataset = 2, member = 2, sdate = 5, ftime = 3))), - "Parameter 'exp' and 'obs' must have the same dimensions expect 'dat_dim'." + "Parameter 'exp' and 'obs' must have the same dimensions except 'dat_dim'." ) }) diff --git a/tests/testthat/test-Consist_Trend.R b/tests/testthat/test-Consist_Trend.R index aa66f45..91dacf7 100644 --- a/tests/testthat/test-Consist_Trend.R +++ b/tests/testthat/test-Consist_Trend.R @@ -62,7 +62,7 @@ test_that("1. Input checks", { Consist_Trend(array(1:10, dim = c(dataset = 2, member = 5, sdate = 4, ftime = 3)), array(1:4, dim = c(dataset = 2, member = 2, sdate = 5, ftime = 3))), paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension expect 'dat_dim'.") + "all dimension except 'dat_dim'.") ) # interval expect_error( diff --git a/tests/testthat/test-Corr.R b/tests/testthat/test-Corr.R index 9d5d4a3..164e61e 100644 --- a/tests/testthat/test-Corr.R +++ b/tests/testthat/test-Corr.R @@ -132,7 +132,7 @@ test_that("1. Input checks", { expect_error( Corr(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension expect 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." ) expect_error( Corr(exp = array(1:10, dim = c(sdate = 2, dataset = 5, a = 1)), diff --git a/tests/testthat/test-NAO.R b/tests/testthat/test-NAO.R index f3c6d21..2da0d40 100644 --- a/tests/testthat/test-NAO.R +++ b/tests/testthat/test-NAO.R @@ -95,12 +95,12 @@ test_that("1. Input checks", { expect_error( NAO(exp1, array(rnorm(10), dim = c(member = 1, sdate = 3, ftime = 4, lat = 2, lon = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'memb_dim'.") + "of all the dimensions except 'memb_dim'.") ) expect_error( NAO(exp1, array(rnorm(10), dim = c(dataset = 1, member = 1, sdate = 3, ftime = 4, lat = 1, lon = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'memb_dim'.") + "of all the dimensions except 'memb_dim'.") ) # ftime_avg expect_error( diff --git a/tests/testthat/test-RMS.R b/tests/testthat/test-RMS.R index 6e18bee..f1083dd 100644 --- a/tests/testthat/test-RMS.R +++ b/tests/testthat/test-RMS.R @@ -89,7 +89,7 @@ test_that("1. Input checks", { expect_error( RMS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension expect 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." ) expect_error( RMS(exp = array(1:5, dim = c(sdate = 1, dataset = 5, a = 1)), diff --git a/tests/testthat/test-RMSSS.R b/tests/testthat/test-RMSSS.R index 7d2a16d..91ef930 100644 --- a/tests/testthat/test-RMSSS.R +++ b/tests/testthat/test-RMSSS.R @@ -77,7 +77,7 @@ test_that("1. Input checks", { expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension expect 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." ) expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), diff --git a/tests/testthat/test-RPS.R b/tests/testthat/test-RPS.R index 04df1bb..2029fd9 100644 --- a/tests/testthat/test-RPS.R +++ b/tests/testthat/test-RPS.R @@ -56,7 +56,7 @@ test_that("1. Input checks", { # exp, ref, and obs (2) expect_error( RPS(exp1, array(1:9, dim = c(sdate = 9))), - "Parameter 'exp' and 'obs' must have same length of all dimensions expect 'memb_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimensions except 'memb_dim'." ) # prob_thresholds expect_error( diff --git a/tests/testthat/test-RPSS.R b/tests/testthat/test-RPSS.R index c63add0..b8623b6 100644 --- a/tests/testthat/test-RPSS.R +++ b/tests/testthat/test-RPSS.R @@ -63,11 +63,11 @@ test_that("1. Input checks", { # exp, ref, and obs (2) expect_error( RPSS(exp1, array(1:9, dim = c(sdate = 9))), - "Parameter 'exp' and 'obs' must have same length of all dimensions expect 'memb_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimensions except 'memb_dim'." ) expect_error( RPSS(exp1, obs1, ref2), - "Parameter 'exp' and 'obs' must have same length of all dimensions expect 'memb_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimensions except 'memb_dim'." ) # prob_thresholds expect_error( diff --git a/tests/testthat/test-RatioSDRMS.R b/tests/testthat/test-RatioSDRMS.R index 4b93faf..5dbc171 100644 --- a/tests/testthat/test-RatioSDRMS.R +++ b/tests/testthat/test-RatioSDRMS.R @@ -77,7 +77,7 @@ test_that("1. Input checks", { expect_error( RatioSDRMS(exp = exp3, obs = array(1:2, dim = c(member = 1, sdate = 2)), dat_dim = NULL), - "Parameter 'exp' and 'obs' must have same length of all the dimensions expect 'dat_dim' and 'memb_dim'." + "Parameter 'exp' and 'obs' must have same length of all the dimensions except 'dat_dim' and 'memb_dim'." ) }) diff --git a/tests/testthat/test-ResidualCorr.R b/tests/testthat/test-ResidualCorr.R index c15276c..be71b47 100644 --- a/tests/testthat/test-ResidualCorr.R +++ b/tests/testthat/test-ResidualCorr.R @@ -66,7 +66,7 @@ test_that("1. Input checks", { # exp, ref, and obs (2) expect_error( ResidualCorr(exp1, obs1, array(1:10, dim = c(sdate = 10, memb = 1)), memb_dim = 'memb'), - "Parameter 'exp', 'obs', and 'ref' must have same length of all dimensions expect 'memb_dim'." + "Parameter 'exp', 'obs', and 'ref' must have same length of all dimensions except 'memb_dim'." ) # method expect_error( diff --git a/tests/testthat/test-UltimateBrier.R b/tests/testthat/test-UltimateBrier.R index 654aad0..0eb0118 100644 --- a/tests/testthat/test-UltimateBrier.R +++ b/tests/testthat/test-UltimateBrier.R @@ -58,12 +58,12 @@ test_that("1. Input checks", { expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 6, ftime = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.") + "of all the dimensions except 'dat_dim' and 'memb_dim'.") ) expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 5, time = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.") + "of all the dimensions except 'dat_dim' and 'memb_dim'.") ) # quantile expect_error( -- GitLab From 2f9f2206e154714918a05fcef69a2b672f64094f Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Fri, 1 Jul 2022 16:43:25 +0200 Subject: [PATCH 08/47] Upload files in roxygen format --- man/BrierScore.Rd | 6 +++--- man/Consist_Trend.Rd | 6 ++++-- man/Corr.Rd | 6 ++++-- man/RMS.Rd | 3 ++- man/RMSSS.Rd | 3 ++- man/RatioRMS.Rd | 2 +- man/s2dv-package.Rd | 10 +++++++++- man/sampleDepthData.Rd | 6 ++---- man/sampleMap.Rd | 6 ++---- man/sampleTimeSeries.Rd | 6 ++---- 10 files changed, 31 insertions(+), 23 deletions(-) diff --git a/man/BrierScore.Rd b/man/BrierScore.Rd index 9271a2a..76fb27e 100644 --- a/man/BrierScore.Rd +++ b/man/BrierScore.Rd @@ -73,13 +73,13 @@ Items 'rel', 'res', 'unc', 'bs', 'bs_check_res', 'bss_res', 'gres', 'bs_check_gres', 'bss_gres', 'rel_bias_corrected', 'gres_bias_corrected', 'unc_bias_corrected', and 'bss_bias_corrected' are (a) a number (b) an array with dimensions c(nexp, nobs, all the rest dimensions in 'exp' and 'obs' -expect 'time_dim' and 'memb_dim') (c) an array with dimensions of +except 'time_dim' and 'memb_dim') (c) an array with dimensions of 'exp' and 'obs' except 'time_dim' and 'memb_dim'\cr Items 'nk', 'fkbar', and 'okbar' are (a) a vector of length of bin number determined by 'threshold' (b) an array with dimensions c(nexp, nobs, -no. of bins, all the rest dimensions in 'exp' and 'obs' expect 'time_dim' and +no. of bins, all the rest dimensions in 'exp' and 'obs' except 'time_dim' and 'memb_dim') (c) an array with dimensions c(no. of bin, all the rest dimensions -in 'exp' and 'obs' expect 'time_dim' and 'memb_dim') +in 'exp' and 'obs' except 'time_dim' and 'memb_dim') } \description{ Compute the Brier score (BS) and the components of its standard decompostion diff --git a/man/Consist_Trend.Rd b/man/Consist_Trend.Rd index 2ac7d42..b93dad5 100644 --- a/man/Consist_Trend.Rd +++ b/man/Consist_Trend.Rd @@ -23,7 +23,8 @@ parameter 'exp' except along 'dat_dim'.} \item{dat_dim}{A character string indicating the name of the dataset dimensions. If data at some point of 'time_dim' are not complete along 'dat_dim' in both 'exp' and 'obs', this point in all 'dat_dim' will be -discarded. The default value is 'dataset'.} +discarded. The default value is 'dataset'. If there is no dataset +dimension, set NULL.} \item{time_dim}{A character string indicating the name of dimension along which the trend is computed. The default value is 'sdate'.} @@ -41,7 +42,8 @@ A list containing: with dimensions c(stats = 2, nexp + nobs, the rest dimensions of 'exp' and 'obs' except time_dim), where 'nexp' is the length of 'dat_dim' in 'exp' and 'nobs' is the length of 'dat_dim' in 'obs. The 'stats' dimension - contains the intercept and the slope. + contains the intercept and the slope. If dat_dim is NULL, nexp and nobs are + omitted. } \item{$conf.lower}{ A numeric array of the lower limit of 95\% confidence interval with diff --git a/man/Corr.Rd b/man/Corr.Rd index 077791f..54f101b 100644 --- a/man/Corr.Rd +++ b/man/Corr.Rd @@ -31,7 +31,8 @@ parameter 'exp' except along 'dat_dim' and 'memb_dim'.} which the correlations are computed. The default value is 'sdate'.} \item{dat_dim}{A character string indicating the name of dataset (nobs/nexp) -dimension. The default value is 'dataset'.} +dimension. The default value is 'dataset'. If there is no dataset +dimension, set NULL.} \item{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 @@ -68,7 +69,8 @@ A list containing the numeric arrays with dimension:\cr c(nexp, nobs, exp_memb, obs_memb, all other dimensions of exp except time_dim and memb_dim).\cr nexp is the number of experiment (i.e., 'dat_dim' in exp), and nobs is the -number of observation (i.e., 'dat_dim' in obs). exp_memb is the number of +number of observation (i.e., 'dat_dim' in obs). If + dat_dim is NULL, nexp and nobs are omitted. exp_memb is the number of member in experiment (i.e., 'memb_dim' in exp) and obs_memb is the number of member in observation (i.e., 'memb_dim' in obs).\cr\cr \item{$corr}{ diff --git a/man/RMS.Rd b/man/RMS.Rd index 4391df4..a84f965 100644 --- a/man/RMS.Rd +++ b/man/RMS.Rd @@ -31,7 +31,8 @@ length as 'exp', then the vector will automatically be 'time_dim' and which the correlations are computed. The default value is 'sdate'.} \item{dat_dim}{A character string indicating the name of member (nobs/nexp) -dimension. The default value is 'dataset'.} +dimension. The default value is 'dataset'. If there is no dataset +dimension, set NULL.} \item{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 diff --git a/man/RMSSS.Rd b/man/RMSSS.Rd index 9ebcf65..2c23ee4 100644 --- a/man/RMSSS.Rd +++ b/man/RMSSS.Rd @@ -30,7 +30,8 @@ be 1.} which the RMSSS are computed. The default value is 'sdate'.} \item{dat_dim}{A character string indicating the name of dataset (nobs/nexp) -dimension. The default value is 'dataset'.} +dimension. The default value is 'dataset'. If there is no dataset +dimension, set NULL.} \item{pval}{A logical value indicating whether to compute or not the p-value of the test Ho: RMSSS = 0. If pval = TRUE, the insignificant RMSSS will diff --git a/man/RatioRMS.Rd b/man/RatioRMS.Rd index 194c6b9..3330eb5 100644 --- a/man/RatioRMS.Rd +++ b/man/RatioRMS.Rd @@ -30,7 +30,7 @@ computation. The default value is NULL.} } \value{ A list containing the numeric arrays with dimensions identical with - 'exp1', 'exp2', and 'obs', expect 'time_dim': + 'exp1', 'exp2', and 'obs', except 'time_dim': \item{$ratiorms}{ The ratio between the RMSE (i.e., RMSE1/RMSE2). } diff --git a/man/s2dv-package.Rd b/man/s2dv-package.Rd index 3c98a95..ffb6783 100644 --- a/man/s2dv-package.Rd +++ b/man/s2dv-package.Rd @@ -6,7 +6,15 @@ \alias{s2dv-package} \title{s2dv: A Set of Common Tools for Seasonal to Decadal Verification} \description{ -The advanced version of package 's2dverification'. It is intended for 'seasonal to decadal' (s2d) climate forecast verification, but it can also be used in other kinds of forecasts or general climate analysis. This package is specially designed for the comparison between the experimental and observational datasets. The functionality of the included functions covers from data retrieval, data post-processing, skill scores against observation, to visualization. Compared to 's2dverification', 's2dv' is more compatible with the package 'startR', able to use multiple cores for computation and handle multi-dimensional arrays with a higher flexibility. +The advanced version of package 's2dverification'. It is + intended for 'seasonal to decadal' (s2d) climate forecast verification, but + it can also be used in other kinds of forecasts or general climate analysis. + This package is specially designed for the comparison between the experimental + and observational datasets. The functionality of the included functions covers + from data retrieval, data post-processing, skill scores against observation, + to visualization. Compared to 's2dverification', 's2dv' is more compatible + with the package 'startR', able to use multiple cores for computation and + handle multi-dimensional arrays with a higher flexibility. } \references{ \url{https://earth.bsc.es/gitlab/es/s2dv/} diff --git a/man/sampleDepthData.Rd b/man/sampleDepthData.Rd index 47e2f1b..77e4a7a 100644 --- a/man/sampleDepthData.Rd +++ b/man/sampleDepthData.Rd @@ -5,8 +5,7 @@ \alias{sampleDepthData} \title{Sample of Experimental Data for Forecast Verification In Function Of Latitudes And Depths} -\format{ -The data set provides with a variable named 'sampleDepthData'.\cr\cr +\format{The data set provides with a variable named 'sampleDepthData'.\cr\cr sampleDepthData$exp is an array that contains the experimental data and the dimension meanings and values are:\cr @@ -19,8 +18,7 @@ but in this sample is not defined (NULL).\cr\cr sampleDepthData$depths is an array with the 7 longitudes covered by the data.\cr\cr -sampleDepthData$lat is an array with the 21 latitudes covered by the data.\cr\cr -} +sampleDepthData$lat is an array with the 21 latitudes covered by the data.\cr\cr} \usage{ data(sampleDepthData) } diff --git a/man/sampleMap.Rd b/man/sampleMap.Rd index e4ec5a5..eaf8aa5 100644 --- a/man/sampleMap.Rd +++ b/man/sampleMap.Rd @@ -4,8 +4,7 @@ \name{sampleMap} \alias{sampleMap} \title{Sample Of Observational And Experimental Data For Forecast Verification In Function Of Longitudes And Latitudes} -\format{ -The data set provides with a variable named 'sampleMap'.\cr\cr +\format{The data set provides with a variable named 'sampleMap'.\cr\cr sampleMap$mod is an array that contains the experimental data and the dimension meanings and values are:\cr c(# of experimental datasets, # of members, # of starting dates, # of lead-times, # of latitudes, # of longitudes)\cr @@ -17,8 +16,7 @@ sampleMap$obs is an array that contains the observational data and the dimension sampleMap$lat is an array with the 2 latitudes covered by the data (see examples on Load() for details on why such low resolution).\cr\cr - sampleMap$lon is an array with the 3 longitudes covered by the data (see examples on Load() for details on why such low resolution). -} + sampleMap$lon is an array with the 3 longitudes covered by the data (see examples on Load() for details on why such low resolution).} \usage{ data(sampleMap) } diff --git a/man/sampleTimeSeries.Rd b/man/sampleTimeSeries.Rd index 7f058e2..05a8e79 100644 --- a/man/sampleTimeSeries.Rd +++ b/man/sampleTimeSeries.Rd @@ -4,8 +4,7 @@ \name{sampleTimeSeries} \alias{sampleTimeSeries} \title{Sample Of Observational And Experimental Data For Forecast Verification As Area Averages} -\format{ -The data set provides with a variable named 'sampleTimeSeries'.\cr\cr +\format{The data set provides with a variable named 'sampleTimeSeries'.\cr\cr sampleTimeSeries$mod is an array that contains the experimental data and the dimension meanings and values are:\cr c(# of experimental datasets, # of members, # of starting dates, # of lead-times)\cr @@ -17,8 +16,7 @@ sampleTimeSeries$obs is an array that contains the observational data and the di sampleTimeSeries$lat is an array with the 2 latitudes covered by the data that was area averaged to calculate the time series (see examples on Load() for details on why such low resolution).\cr\cr -sampleTimeSeries$lon is an array with the 3 longitudes covered by the data that was area averaged to calculate the time series (see examples on Load() for details on why such low resolution). -} +sampleTimeSeries$lon is an array with the 3 longitudes covered by the data that was area averaged to calculate the time series (see examples on Load() for details on why such low resolution).} \usage{ data(sampleTimeSeries) } -- GitLab From ea68861e43836af438f429a6c6ce5b6db367b66d Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Tue, 5 Jul 2022 12:50:31 +0200 Subject: [PATCH 09/47] Corrected development of dat_dim = NULL --- R/Corr.R | 160 +++++++++++++++++++++++++++---------- tests/testthat/test-Corr.R | 59 +++++++++++++- 2 files changed, 177 insertions(+), 42 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index 628511d..2e3ea17 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -242,6 +242,16 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } if (is.null(memb_dim)) { + if (is.null(dat_dim)) { + res <- Apply(list(exp, obs), + target_dims = list(c(time_dim), + c(time_dim)), + fun = .Corr, + dat_dim = dat_dim, memb_dim = memb_dim, + time_dim = time_dim, method = method, + pval = pval, conf = conf, conf.lev = conf.lev, + ncores = ncores) + } else { # Define output_dims if (conf & pval) { output_dims <- list(corr = c('nexp', 'nobs'), @@ -264,11 +274,12 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', c(time_dim, dat_dim)), output_dims = output_dims, fun = .Corr, + dat_dim = dat_dim, memb_dim = memb_dim, time_dim = time_dim, method = method, pval = pval, conf = conf, conf.lev = conf.lev, ncores = ncores) - } else { + }} else { if (!memb) { #ensemble mean name_exp <- names(dim(exp)) margin_dims_ind <- c(1:length(name_exp))[-which(name_exp == memb_dim)] @@ -297,11 +308,40 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', c(time_dim, dat_dim)), output_dims = output_dims, fun = .Corr, + dat_dim = dat_dim, memb_dim = NULL, time_dim = time_dim, method = method, pval = pval, conf = conf, conf.lev = conf.lev, ncores_input = ncores, ncores = ncores) } else { # correlation for each member + if (is.null(dat_dim)) { + # Define output_dims + if (conf & pval) { + output_dims <- list(corr = c('exp_memb', 'obs_memb'), + p.val = c('exp_memb', 'obs_memb'), + conf.lower = c('exp_memb', 'obs_memb'), + conf.upper = c('exp_memb', 'obs_memb')) + } else if (conf & !pval) { + output_dims <- list(corr = c('exp_memb', 'obs_memb'), + conf.lower = c('exp_memb', 'obs_memb'), + conf.upper = c('exp_memb', 'obs_memb')) + } else if (!conf & pval) { + output_dims <- list(corr = c('exp_memb', 'obs_memb'), + p.val = c('exp_memb', 'obs_memb')) + } else { + output_dims <- list(corr = c('exp_memb', 'obs_memb')) + } + + res <- Apply(list(exp, obs), + target_dims = list(c(time_dim, memb_dim), + c(time_dim, memb_dim)), + output_dims = output_dims, + fun = .Corr, + dat_dim = dat_dim, memb_dim = memb_dim, + time_dim = time_dim, method = method, + pval = pval, conf = conf, conf.lev = conf.lev, ncores_input = ncores, + ncores = ncores) + } else { # Define output_dims if (conf & pval) { @@ -325,28 +365,33 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', c(time_dim, dat_dim, memb_dim)), output_dims = output_dims, fun = .Corr, + dat_dim = dat_dim, memb_dim = memb_dim, time_dim = time_dim, method = method, pval = pval, conf = conf, conf.lev = conf.lev, ncores_input = ncores, ncores = ncores) } + } } return(res) } -.Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', method = 'pearson', +.Corr <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', time_dim = 'sdate', method = 'pearson', conf = TRUE, pval = TRUE, conf.lev = 0.95, ncores_input = NULL) { - if (length(dim(exp)) == 2) { + if (is.null(memb_dim)) { + if (is.null(dat_dim)) { + # exp: [sdate] + # obs: [sdate] + nexp <- 1 + nobs <- 1 + CORR <- cor(exp, obs,use = "pairwise.complete.obs",method = method) + } else { # exp: [sdate, dat_exp] # obs: [sdate, dat_obs] - if (is.null(dat_dim)) { - nexp <- 1 - nobs <- 1 - } else { - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) - } + nexp <- as.numeric(dim(exp)[dat_dim]) + nobs <- as.numeric(dim(obs)[dat_dim]) + # NOTE: Use sapply to replace the for loop CORR <- sapply(1:nobs, function(i) { sapply(1:nexp, function (x) { @@ -361,40 +406,58 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', }) if (is.null(dim(CORR))) { CORR <- array(CORR, dim = c(1, 1)) - } + }} } else { # member - - # exp: [sdate, dat_exp, memb_exp] - # obs: [sdate, dat_obs, memb_obs] if (is.null(dat_dim)) { - nexp <- 1 - nobs <- 1 - } else { - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) - exp_memb <- as.numeric(dim(exp)[3]) - obs_memb <- as.numeric(dim(obs)[3]) - } - CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) - - for (j in 1:obs_memb) { - for (y in 1:exp_memb) { - CORR[, , y, j] <- sapply(1:nobs, function(i) { - sapply(1:nexp, function (x) { - if (any(!is.na(exp[, x, y])) && sum(!is.na(obs[, i, j])) > 2) { - cor(exp[, x, y], obs[, i, j], - use = "pairwise.complete.obs", - method = method) - } else { - NA #CORR[, i] <- NA - } - }) - }) + # exp: [sdate, memb_exp] + # obs: [sdate, memb_obs] + nexp <- 1 + nobs <- 1 + exp_memb <- as.numeric(dim(exp)[memb_dim]) # memb_dim + obs_memb <- as.numeric(dim(obs)[memb_dim]) + CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) + + for (j in 1:obs_memb) { + for (y in 1:exp_memb) { + + if (any(!is.na(exp[,y])) && sum(!is.na(obs[, j])) > 2) { + CORR[, , y, j] <- cor(exp[, y], obs[, j], + use = "pairwise.complete.obs", + method = method) + } else { + NA #CORR[, i] <- NA + } + + } + } + } else { + # exp: [sdate, dat_exp, memb_exp] + # obs: [sdate, dat_obs, memb_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) + exp_memb <- as.numeric(dim(exp)[3]) + obs_memb <- as.numeric(dim(obs)[3]) + + CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) + + for (j in 1:obs_memb) { + for (y in 1:exp_memb) { + CORR[, , y, j] <- sapply(1:nobs, function(i) { + sapply(1:nexp, function (x) { + if (any(!is.na(exp[, x, y])) && sum(!is.na(obs[, i, j])) > 2) { + cor(exp[, x, y], obs[, i, j], + use = "pairwise.complete.obs", + method = method) + } else { + NA #CORR[, i] <- NA + } + }) + }) + } + } } - } - } @@ -420,7 +483,7 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', eno <- Eno(obs, time_dim, ncores = ncores_input) } - if (length(dim(exp)) == 2) { + if (is.null(memb_dim)) { eno_expand <- array(dim = c(nexp = nexp, nobs = nobs)) for (i in 1:nexp) { eno_expand[i, ] <- eno @@ -450,6 +513,21 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', confhigh <- tanh(atanh(CORR) + qnorm(conf.upper) / sqrt(eno_expand - 3)) } +################################### NEW +# Remove nexp and nobs if dat_dim = NULL + if (is.null(dat_dim) & !is.null(memb_dim)) { + dim(CORR) <- dim(CORR)[3:length(dim(CORR))] + if (pval) { + dim(p.val) <- dim(p.val)[3:length(dim(p.val))] + } + if (conf) { + dim(conflow) <- dim(conflow)[3:length(dim(conflow))] + dim(confhigh) <- dim(confhigh)[3:length(dim(confhigh))] + } + } + +################################### + if (pval & conf) { res <- list(corr = CORR, p.val = p.val, conf.lower = conflow, conf.upper = confhigh) @@ -462,6 +540,6 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', res <- list(corr = CORR) } - return(res) + return(res) } diff --git a/tests/testthat/test-Corr.R b/tests/testthat/test-Corr.R index 164e61e..65407d1 100644 --- a/tests/testthat/test-Corr.R +++ b/tests/testthat/test-Corr.R @@ -40,6 +40,22 @@ context("s2dv::Corr tests") obs4 <- array(rnorm(10), dim = c(member = 1, dataset = 1, sdate = 5, lat = 2)) + # dat5: exp and obs have memb_dim and dataset = NULL + set.seed(1) + exp5 <- array(rnorm(90), dim = c(member = 3, sdate = 5, + lat = 2, lon = 3)) + + set.seed(2) + obs5 <- array(rnorm(30), dim = c(member = 1, sdate = 5, + lat = 2, lon = 3)) + + # dat6: exp and obs have memb_dim = NULL and dataset = NULL + set.seed(1) + exp6 <- array(rnorm(90), dim = c(sdate = 5, lat = 2, lon = 3)) + + set.seed(2) + obs6 <- array(rnorm(30), dim = c(sdate = 5, lat = 2, lon = 3)) + ############################################## test_that("1. Input checks", { @@ -132,7 +148,7 @@ test_that("1. Input checks", { expect_error( Corr(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim' and 'memb_dim'." ) expect_error( Corr(exp = array(1:10, dim = c(sdate = 2, dataset = 5, a = 1)), @@ -392,4 +408,45 @@ test_that("5. Output checks: dat4", { ) }) + +############################################## +test_that("6. Output checks: dat5", { + expect_equal( + dim(Corr(exp5, obs5, dat_dim = NULL, memb_dim = 'member')$corr), + c(exp_memb = 3, obs_memb = 1, lat = 2, lon = 3) + ) + expect_equal( + names(Corr(exp5, obs5, dat_dim = NULL, memb_dim = 'member')), + c("corr", "p.val", "conf.lower", "conf.upper") + ) + expect_equal( + names(Corr(exp5, obs5, dat_dim = NULL, memb_dim = 'member', pval = FALSE, conf = FALSE)), + c("corr") + ) + expect_equal( + mean(Corr(exp5, obs5, dat_dim = NULL, memb_dim = 'member', pval = FALSE, conf = FALSE)$corr), + 0.1880204, + tolerance = 0.0001 + ) +}) ############################################## +test_that("7. Output checks: dat6", { + expect_equal( + dim(Corr(exp6, obs6, dat_dim = NULL, memb_dim = NULL)$corr), + c(lat = 2, lon = 3) + ) + expect_equal( + names(Corr(exp6, obs6, dat_dim = NULL, memb_dim = NULL)), + c("corr", "p.val", "conf.lower", "conf.upper") + ) + expect_equal( + names(Corr(exp6, obs6, dat_dim = NULL, memb_dim = NULL, pval = FALSE, conf = FALSE)), + c("corr") + ) + expect_equal( + mean(Corr(exp6, obs6, dat_dim = NULL, memb_dim = NULL, pval = FALSE, conf = FALSE)$corr), + 0.1084488, + tolerance = 0.0001 + ) +}) +############################################## \ No newline at end of file -- GitLab From 2a693808e7e5309abf888d3e3a450eb4aa67b105 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Tue, 5 Jul 2022 16:25:46 +0200 Subject: [PATCH 10/47] Added test for new development --- tests/testthat/test-Corr.R | 15 +++++++++++++++ 1 file changed, 15 insertions(+) diff --git a/tests/testthat/test-Corr.R b/tests/testthat/test-Corr.R index 65407d1..7d477ad 100644 --- a/tests/testthat/test-Corr.R +++ b/tests/testthat/test-Corr.R @@ -56,6 +56,13 @@ context("s2dv::Corr tests") set.seed(2) obs6 <- array(rnorm(30), dim = c(sdate = 5, lat = 2, lon = 3)) + # dat7: exp and obs have memb_dim = NULL and dataset = 1 + set.seed(1) + exp7 <- array(rnorm(90), dim = c(dataset = 1, sdate = 5, lat = 2, lon = 3)) + + set.seed(2) + obs7 <- array(rnorm(30), dim = c(dataset = 1, sdate = 5, lat = 2, lon = 3)) + ############################################## test_that("1. Input checks", { @@ -449,4 +456,12 @@ test_that("7. Output checks: dat6", { tolerance = 0.0001 ) }) +############################################## +test_that("8. Output checks: dat6 and dat7", { + expect_equal( + mean(Corr(exp6, obs6, dat_dim = NULL, memb_dim = NULL, pval = FALSE, conf = FALSE)$corr), + mean(Corr(exp7, obs7, memb_dim = NULL, pval = FALSE, conf = FALSE)$corr), + tolerance = 0.0001 + ) +}) ############################################## \ No newline at end of file -- GitLab From f4e53ae2f7cd5cb3ce5dccb00e41bc7af70e43c6 Mon Sep 17 00:00:00 2001 From: aho Date: Tue, 5 Jul 2022 18:24:42 +0200 Subject: [PATCH 11/47] Refine documentation --- R/Corr.R | 32 ++++++++++++++++---------------- man/Corr.Rd | 32 ++++++++++++++++---------------- 2 files changed, 32 insertions(+), 32 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index 2e3ea17..d9e8144 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -2,22 +2,22 @@ #' #'Calculate the correlation coefficient (Pearson, Kendall or Spearman) for #'an array of forecast and an array of observation. The correlations are -#'computed along time_dim, the startdate dimension. If comp_dim is given, -#'the correlations are computed only if obs along the comp_dim dimension are -#'complete between limits[1] and limits[2], i.e., there is no NA between -#'limits[1] and limits[2]. This option can be activated if the user wants to -#'account only for the forecasts which the corresponding observations are -#'available at all leadtimes.\cr +#'computed along 'time_dim' that usually refers to the start date dimension. If +#''comp_dim' is given, the correlations are computed only if obs along comp_dim +#'dimension are complete between limits[1] and limits[2], i.e., there is no NA +#'between limits[1] and limits[2]. This option can be activated if the user +#'wants to account only for the forecasts which the corresponding observations +#'are available at all leadtimes.\cr #'The confidence interval is computed by the Fisher transformation and the #'significance level relies on an one-sided student-T distribution.\cr -#'If the dataset has more than one member, ensemble mean is necessary necessary -#'before using this function since it only allows one dimension 'dat_dim' to -#'have inconsistent length between 'exp' and 'obs'. If all the dimensions of -#''exp' and 'obs' are identical, you can simply use apply() and cor() to +#'The function can calculate ensemble mean before correlation by 'memb_dim' +#'specified and 'memb = F'. If ensemble mean is not calculated, correlation will +#'be calculated for each member. +#'If there is only one dataset for exp and obs, you can simply use cor() to #'compute the correlation. #' -#'@param exp A named numeric array of experimental data, with at least two -#' dimensions 'time_dim' and 'dat_dim'. +#'@param exp A named numeric array of experimental data, with at least dimension +#' 'time_dim'. #'@param obs A named numeric array of observational data, same dimensions as #' parameter 'exp' except along 'dat_dim' and 'memb_dim'. #'@param time_dim A character string indicating the name of dimension along @@ -52,10 +52,10 @@ #' c(nexp, nobs, exp_memb, obs_memb, all other dimensions of exp except #' time_dim and memb_dim).\cr #'nexp is the number of experiment (i.e., 'dat_dim' in exp), and nobs is the -#'number of observation (i.e., 'dat_dim' in obs). If -#' dat_dim is NULL, nexp and nobs are omitted. exp_memb is the number of -#'member in experiment (i.e., 'memb_dim' in exp) and obs_memb is the number of -#'member in observation (i.e., 'memb_dim' in obs).\cr\cr +#'number of observation (i.e., 'dat_dim' in obs). If dat_dim is NULL, nexp and +#'nobs are omitted. exp_memb is the number of member in experiment (i.e., +#''memb_dim' in exp) and obs_memb is the number of member in observation (i.e., +#''memb_dim' in obs). If memb = F, exp_memb and obs_memb are omitted.\cr\cr #'\item{$corr}{ #' The correlation coefficient. #'} diff --git a/man/Corr.Rd b/man/Corr.Rd index 54f101b..4cbd098 100644 --- a/man/Corr.Rd +++ b/man/Corr.Rd @@ -21,8 +21,8 @@ Corr( ) } \arguments{ -\item{exp}{A named numeric array of experimental data, with at least two -dimensions 'time_dim' and 'dat_dim'.} +\item{exp}{A named numeric array of experimental data, with at least dimension +'time_dim'.} \item{obs}{A named numeric array of observational data, same dimensions as parameter 'exp' except along 'dat_dim' and 'memb_dim'.} @@ -69,10 +69,10 @@ A list containing the numeric arrays with dimension:\cr c(nexp, nobs, exp_memb, obs_memb, all other dimensions of exp except time_dim and memb_dim).\cr nexp is the number of experiment (i.e., 'dat_dim' in exp), and nobs is the -number of observation (i.e., 'dat_dim' in obs). If - dat_dim is NULL, nexp and nobs are omitted. exp_memb is the number of -member in experiment (i.e., 'memb_dim' in exp) and obs_memb is the number of -member in observation (i.e., 'memb_dim' in obs).\cr\cr +number of observation (i.e., 'dat_dim' in obs). If dat_dim is NULL, nexp and +nobs are omitted. exp_memb is the number of member in experiment (i.e., +'memb_dim' in exp) and obs_memb is the number of member in observation (i.e., +'memb_dim' in obs). If memb = F, exp_memb and obs_memb are omitted.\cr\cr \item{$corr}{ The correlation coefficient. } @@ -89,18 +89,18 @@ member in observation (i.e., 'memb_dim' in obs).\cr\cr \description{ Calculate the correlation coefficient (Pearson, Kendall or Spearman) for an array of forecast and an array of observation. The correlations are -computed along time_dim, the startdate dimension. If comp_dim is given, -the correlations are computed only if obs along the comp_dim dimension are -complete between limits[1] and limits[2], i.e., there is no NA between -limits[1] and limits[2]. This option can be activated if the user wants to -account only for the forecasts which the corresponding observations are -available at all leadtimes.\cr +computed along 'time_dim' that usually refers to the start date dimension. If +'comp_dim' is given, the correlations are computed only if obs along comp_dim +dimension are complete between limits[1] and limits[2], i.e., there is no NA +between limits[1] and limits[2]. This option can be activated if the user +wants to account only for the forecasts which the corresponding observations +are available at all leadtimes.\cr The confidence interval is computed by the Fisher transformation and the significance level relies on an one-sided student-T distribution.\cr -If the dataset has more than one member, ensemble mean is necessary necessary -before using this function since it only allows one dimension 'dat_dim' to -have inconsistent length between 'exp' and 'obs'. If all the dimensions of -'exp' and 'obs' are identical, you can simply use apply() and cor() to +The function can calculate ensemble mean before correlation by 'memb_dim' +specified and 'memb = F'. If ensemble mean is not calculated, correlation will +be calculated for each member. +If there is only one dataset for exp and obs, you can simply use cor() to compute the correlation. } \examples{ -- GitLab From 5bf47ac9a85c39f090101d0241bc98273e9659b5 Mon Sep 17 00:00:00 2001 From: aho Date: Tue, 5 Jul 2022 18:24:52 +0200 Subject: [PATCH 12/47] Formatting --- man/s2dv-package.Rd | 10 +--------- man/sampleDepthData.Rd | 6 ++++-- man/sampleMap.Rd | 6 ++++-- man/sampleTimeSeries.Rd | 6 ++++-- 4 files changed, 13 insertions(+), 15 deletions(-) diff --git a/man/s2dv-package.Rd b/man/s2dv-package.Rd index ffb6783..3c98a95 100644 --- a/man/s2dv-package.Rd +++ b/man/s2dv-package.Rd @@ -6,15 +6,7 @@ \alias{s2dv-package} \title{s2dv: A Set of Common Tools for Seasonal to Decadal Verification} \description{ -The advanced version of package 's2dverification'. It is - intended for 'seasonal to decadal' (s2d) climate forecast verification, but - it can also be used in other kinds of forecasts or general climate analysis. - This package is specially designed for the comparison between the experimental - and observational datasets. The functionality of the included functions covers - from data retrieval, data post-processing, skill scores against observation, - to visualization. Compared to 's2dverification', 's2dv' is more compatible - with the package 'startR', able to use multiple cores for computation and - handle multi-dimensional arrays with a higher flexibility. +The advanced version of package 's2dverification'. It is intended for 'seasonal to decadal' (s2d) climate forecast verification, but it can also be used in other kinds of forecasts or general climate analysis. This package is specially designed for the comparison between the experimental and observational datasets. The functionality of the included functions covers from data retrieval, data post-processing, skill scores against observation, to visualization. Compared to 's2dverification', 's2dv' is more compatible with the package 'startR', able to use multiple cores for computation and handle multi-dimensional arrays with a higher flexibility. } \references{ \url{https://earth.bsc.es/gitlab/es/s2dv/} diff --git a/man/sampleDepthData.Rd b/man/sampleDepthData.Rd index 77e4a7a..47e2f1b 100644 --- a/man/sampleDepthData.Rd +++ b/man/sampleDepthData.Rd @@ -5,7 +5,8 @@ \alias{sampleDepthData} \title{Sample of Experimental Data for Forecast Verification In Function Of Latitudes And Depths} -\format{The data set provides with a variable named 'sampleDepthData'.\cr\cr +\format{ +The data set provides with a variable named 'sampleDepthData'.\cr\cr sampleDepthData$exp is an array that contains the experimental data and the dimension meanings and values are:\cr @@ -18,7 +19,8 @@ but in this sample is not defined (NULL).\cr\cr sampleDepthData$depths is an array with the 7 longitudes covered by the data.\cr\cr -sampleDepthData$lat is an array with the 21 latitudes covered by the data.\cr\cr} +sampleDepthData$lat is an array with the 21 latitudes covered by the data.\cr\cr +} \usage{ data(sampleDepthData) } diff --git a/man/sampleMap.Rd b/man/sampleMap.Rd index eaf8aa5..e4ec5a5 100644 --- a/man/sampleMap.Rd +++ b/man/sampleMap.Rd @@ -4,7 +4,8 @@ \name{sampleMap} \alias{sampleMap} \title{Sample Of Observational And Experimental Data For Forecast Verification In Function Of Longitudes And Latitudes} -\format{The data set provides with a variable named 'sampleMap'.\cr\cr +\format{ +The data set provides with a variable named 'sampleMap'.\cr\cr sampleMap$mod is an array that contains the experimental data and the dimension meanings and values are:\cr c(# of experimental datasets, # of members, # of starting dates, # of lead-times, # of latitudes, # of longitudes)\cr @@ -16,7 +17,8 @@ sampleMap$obs is an array that contains the observational data and the dimension sampleMap$lat is an array with the 2 latitudes covered by the data (see examples on Load() for details on why such low resolution).\cr\cr - sampleMap$lon is an array with the 3 longitudes covered by the data (see examples on Load() for details on why such low resolution).} + sampleMap$lon is an array with the 3 longitudes covered by the data (see examples on Load() for details on why such low resolution). +} \usage{ data(sampleMap) } diff --git a/man/sampleTimeSeries.Rd b/man/sampleTimeSeries.Rd index 05a8e79..7f058e2 100644 --- a/man/sampleTimeSeries.Rd +++ b/man/sampleTimeSeries.Rd @@ -4,7 +4,8 @@ \name{sampleTimeSeries} \alias{sampleTimeSeries} \title{Sample Of Observational And Experimental Data For Forecast Verification As Area Averages} -\format{The data set provides with a variable named 'sampleTimeSeries'.\cr\cr +\format{ +The data set provides with a variable named 'sampleTimeSeries'.\cr\cr sampleTimeSeries$mod is an array that contains the experimental data and the dimension meanings and values are:\cr c(# of experimental datasets, # of members, # of starting dates, # of lead-times)\cr @@ -16,7 +17,8 @@ sampleTimeSeries$obs is an array that contains the observational data and the di sampleTimeSeries$lat is an array with the 2 latitudes covered by the data that was area averaged to calculate the time series (see examples on Load() for details on why such low resolution).\cr\cr -sampleTimeSeries$lon is an array with the 3 longitudes covered by the data that was area averaged to calculate the time series (see examples on Load() for details on why such low resolution).} +sampleTimeSeries$lon is an array with the 3 longitudes covered by the data that was area averaged to calculate the time series (see examples on Load() for details on why such low resolution). +} \usage{ data(sampleTimeSeries) } -- GitLab From 4664a335c692ef327e5b544436589ae14366ee73 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Wed, 6 Jul 2022 16:22:03 +0200 Subject: [PATCH 13/47] Correct case memb parameter and delete repeated lines --- R/Corr.R | 163 +++++++++---------------------------------------------- 1 file changed, 25 insertions(+), 138 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index d9e8144..124e525 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -241,137 +241,23 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', rm(obs_sub, outrows) } - if (is.null(memb_dim)) { - if (is.null(dat_dim)) { - res <- Apply(list(exp, obs), - target_dims = list(c(time_dim), - c(time_dim)), - fun = .Corr, - dat_dim = dat_dim, memb_dim = memb_dim, - time_dim = time_dim, method = method, - pval = pval, conf = conf, conf.lev = conf.lev, - ncores = ncores) - } else { - # Define output_dims - if (conf & pval) { - output_dims <- list(corr = c('nexp', 'nobs'), - p.val = c('nexp', 'nobs'), - conf.lower = c('nexp', 'nobs'), - conf.upper = c('nexp', 'nobs')) - } else if (conf & !pval) { - output_dims <- list(corr = c('nexp', 'nobs'), - conf.lower = c('nexp', 'nobs'), - conf.upper = c('nexp', 'nobs')) - } else if (!conf & pval) { - output_dims <- list(corr = c('nexp', 'nobs'), - p.val = c('nexp', 'nobs')) - } else { - output_dims <- list(corr = c('nexp', 'nobs')) - } - - res <- Apply(list(exp, obs), - target_dims = list(c(time_dim, dat_dim), - c(time_dim, dat_dim)), - output_dims = output_dims, - fun = .Corr, - dat_dim = dat_dim, memb_dim = memb_dim, - time_dim = time_dim, method = method, - pval = pval, conf = conf, conf.lev = conf.lev, - ncores = ncores) - - }} else { + if (!is.null(memb_dim)) { if (!memb) { #ensemble mean name_exp <- names(dim(exp)) margin_dims_ind <- c(1:length(name_exp))[-which(name_exp == memb_dim)] exp <- apply(exp, margin_dims_ind, mean, na.rm = TRUE) #NOTE: remove NAs here obs <- apply(obs, margin_dims_ind, mean, na.rm = TRUE) - - # Define output_dims - if (conf & pval) { - output_dims <- list(corr = c('nexp', 'nobs'), - p.val = c('nexp', 'nobs'), - conf.lower = c('nexp', 'nobs'), - conf.upper = c('nexp', 'nobs')) - } else if (conf & !pval) { - output_dims <- list(corr = c('nexp', 'nobs'), - conf.lower = c('nexp', 'nobs'), - conf.upper = c('nexp', 'nobs')) - } else if (!conf & pval) { - output_dims <- list(corr = c('nexp', 'nobs'), - p.val = c('nexp', 'nobs')) - } else { - output_dims <- list(corr = c('nexp', 'nobs')) - } - - res <- Apply(list(exp, obs), - target_dims = list(c(time_dim, dat_dim), - c(time_dim, dat_dim)), - output_dims = output_dims, - fun = .Corr, - dat_dim = dat_dim, memb_dim = NULL, - time_dim = time_dim, method = method, - pval = pval, conf = conf, conf.lev = conf.lev, ncores_input = ncores, - ncores = ncores) - - } else { # correlation for each member - if (is.null(dat_dim)) { - # Define output_dims - if (conf & pval) { - output_dims <- list(corr = c('exp_memb', 'obs_memb'), - p.val = c('exp_memb', 'obs_memb'), - conf.lower = c('exp_memb', 'obs_memb'), - conf.upper = c('exp_memb', 'obs_memb')) - } else if (conf & !pval) { - output_dims <- list(corr = c('exp_memb', 'obs_memb'), - conf.lower = c('exp_memb', 'obs_memb'), - conf.upper = c('exp_memb', 'obs_memb')) - } else if (!conf & pval) { - output_dims <- list(corr = c('exp_memb', 'obs_memb'), - p.val = c('exp_memb', 'obs_memb')) - } else { - output_dims <- list(corr = c('exp_memb', 'obs_memb')) - } - - res <- Apply(list(exp, obs), - target_dims = list(c(time_dim, memb_dim), - c(time_dim, memb_dim)), - output_dims = output_dims, - fun = .Corr, - dat_dim = dat_dim, memb_dim = memb_dim, - time_dim = time_dim, method = method, - pval = pval, conf = conf, conf.lev = conf.lev, ncores_input = ncores, - ncores = ncores) - } else { - - # Define output_dims - if (conf & pval) { - output_dims <- list(corr = c('nexp', 'nobs', 'exp_memb', 'obs_memb'), - p.val = c('nexp', 'nobs', 'exp_memb', 'obs_memb'), - conf.lower = c('nexp', 'nobs', 'exp_memb', 'obs_memb'), - conf.upper = c('nexp', 'nobs', 'exp_memb', 'obs_memb')) - } else if (conf & !pval) { - output_dims <- list(corr = c('nexp', 'nobs', 'exp_memb', 'obs_memb'), - conf.lower = c('nexp', 'nobs', 'exp_memb', 'obs_memb'), - conf.upper = c('nexp', 'nobs', 'exp_memb', 'obs_memb')) - } else if (!conf & pval) { - output_dims <- list(corr = c('nexp', 'nobs', 'exp_memb', 'obs_memb'), - p.val = c('nexp', 'nobs', 'exp_memb', 'obs_memb')) - } else { - output_dims <- list(corr = c('nexp', 'nobs', 'exp_memb', 'obs_memb')) - } - - res <- Apply(list(exp, obs), + memb_dim <- NULL + } + } + res <- Apply(list(exp, obs), target_dims = list(c(time_dim, dat_dim, memb_dim), c(time_dim, dat_dim, memb_dim)), - output_dims = output_dims, fun = .Corr, dat_dim = dat_dim, memb_dim = memb_dim, time_dim = time_dim, method = method, pval = pval, conf = conf, conf.lev = conf.lev, ncores_input = ncores, - ncores = ncores) - } - } - } + ncores = ncores) return(res) } @@ -385,28 +271,29 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # obs: [sdate] nexp <- 1 nobs <- 1 + CORR <- array(dim = c(nexp = nexp, nobs = nobs)) CORR <- cor(exp, obs,use = "pairwise.complete.obs",method = method) } else { # exp: [sdate, dat_exp] # obs: [sdate, dat_obs] nexp <- as.numeric(dim(exp)[dat_dim]) nobs <- as.numeric(dim(obs)[dat_dim]) + CORR <- array(dim = c(nexp = nexp, nobs = nobs)) + + for (j in 1:nobs) { + for (y in 1:nexp) { + + if (any(!is.na(exp[,y])) && sum(!is.na(obs[, j])) > 2) { + CORR[y, j] <- cor(exp[, y], obs[, j], + use = "pairwise.complete.obs", + method = method) + } else { + NA #CORR[, i] <- NA + } -# NOTE: Use sapply to replace the for loop - CORR <- sapply(1:nobs, function(i) { - sapply(1:nexp, function (x) { - if (any(!is.na(exp[, x])) && sum(!is.na(obs[, i])) > 2) { #NOTE: Is this necessary? - cor(exp[, x], obs[, i], - use = "pairwise.complete.obs", - method = method) - } else { - NA #CORR[, i] <- NA + } } - }) - }) - if (is.null(dim(CORR))) { - CORR <- array(CORR, dim = c(1, 1)) - }} + } } else { # member if (is.null(dat_dim)) { @@ -434,10 +321,10 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } else { # exp: [sdate, dat_exp, memb_exp] # obs: [sdate, dat_obs, memb_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) - exp_memb <- as.numeric(dim(exp)[3]) - obs_memb <- as.numeric(dim(obs)[3]) + nexp <- as.numeric(dim(exp)[dat_dim]) + nobs <- as.numeric(dim(obs)[dat_dim]) + exp_memb <- as.numeric(dim(exp)[memb_dim]) + obs_memb <- as.numeric(dim(obs)[memb_dim]) CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) -- GitLab From e2d0bb064eee6c786cc332db0a004dc2b67b7c99 Mon Sep 17 00:00:00 2001 From: aho Date: Fri, 8 Jul 2022 16:44:33 +0200 Subject: [PATCH 14/47] Formatting and code improvement --- R/Corr.R | 99 ++++++++++++++++++++++++++++++++------------------------ 1 file changed, 56 insertions(+), 43 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index 124e525..efea0c8 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -243,13 +243,16 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', if (!is.null(memb_dim)) { if (!memb) { #ensemble mean - name_exp <- names(dim(exp)) - margin_dims_ind <- c(1:length(name_exp))[-which(name_exp == memb_dim)] - exp <- apply(exp, margin_dims_ind, mean, na.rm = TRUE) #NOTE: remove NAs here - obs <- apply(obs, margin_dims_ind, mean, na.rm = TRUE) + exp <- MeanDims(exp, memb_dim, na.rm = TRUE) + obs <- MeanDims(obs, memb_dim, na.rm = TRUE) +# name_exp <- names(dim(exp)) +# margin_dims_ind <- c(1:length(name_exp))[-which(name_exp == memb_dim)] +# exp <- apply(exp, margin_dims_ind, mean, na.rm = TRUE) #NOTE: remove NAs here +# obs <- apply(obs, margin_dims_ind, mean, na.rm = TRUE) memb_dim <- NULL } } + res <- Apply(list(exp, obs), target_dims = list(c(time_dim, dat_dim, memb_dim), c(time_dim, dat_dim, memb_dim)), @@ -264,7 +267,6 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', .Corr <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', time_dim = 'sdate', method = 'pearson', conf = TRUE, pval = TRUE, conf.lev = 0.95, ncores_input = NULL) { - if (is.null(memb_dim)) { if (is.null(dat_dim)) { # exp: [sdate] @@ -272,59 +274,70 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', nexp <- 1 nobs <- 1 CORR <- array(dim = c(nexp = nexp, nobs = nobs)) - CORR <- cor(exp, obs,use = "pairwise.complete.obs",method = method) + if (any(!is.na(exp)) && sum(!is.na(obs)) > 2) { + CORR <- cor(exp, obs, use = "pairwise.complete.obs", method = method) + } } else { - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[dat_dim]) - nobs <- as.numeric(dim(obs)[dat_dim]) - CORR <- array(dim = c(nexp = nexp, nobs = nobs)) - + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[dat_dim]) + nobs <- as.numeric(dim(obs)[dat_dim]) + CORR <- array(dim = c(nexp = nexp, nobs = nobs)) for (j in 1:nobs) { - for (y in 1:nexp) { - - if (any(!is.na(exp[,y])) && sum(!is.na(obs[, j])) > 2) { - CORR[y, j] <- cor(exp[, y], obs[, j], - use = "pairwise.complete.obs", - method = method) - } else { - NA #CORR[, i] <- NA - } - - } + for (y in 1:nexp) { + if (any(!is.na(exp[, y])) && sum(!is.na(obs[, j])) > 2) { + CORR[y, j] <- cor(exp[, y], obs[, j], + use = "pairwise.complete.obs", + method = method) + } + } } - } +#---------------------------------------- +# Same as above calculation. +#TODO: Compare which is faster. +# CORR <- sapply(1:nobs, function(i) { +# sapply(1:nexp, function (x) { +# if (any(!is.na(exp[, x])) && sum(!is.na(obs[, i])) > 2) { +# cor(exp[, x], obs[, i], +# use = "pairwise.complete.obs", +# method = method) +# } else { +# NA +# } +# }) +# }) +#----------------------------------------- + } - } else { # member - if (is.null(dat_dim)) { + } else { # memb_dim != NULL + exp_memb <- as.numeric(dim(exp)[memb_dim]) # memb_dim + obs_memb <- as.numeric(dim(obs)[memb_dim]) + + if (is.null(dat_dim)) { # exp: [sdate, memb_exp] # obs: [sdate, memb_obs] nexp <- 1 nobs <- 1 - exp_memb <- as.numeric(dim(exp)[memb_dim]) # memb_dim - obs_memb <- as.numeric(dim(obs)[memb_dim]) CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) for (j in 1:obs_memb) { - for (y in 1:exp_memb) { + for (y in 1:exp_memb) { if (any(!is.na(exp[,y])) && sum(!is.na(obs[, j])) > 2) { - CORR[, , y, j] <- cor(exp[, y], obs[, j], - use = "pairwise.complete.obs", - method = method) - } else { - NA #CORR[, i] <- NA - } - + CORR[, , y, j] <- cor(exp[, y], obs[, j], + use = "pairwise.complete.obs", + method = method) + } else { + NA #CORR[, i] <- NA } + + } } } else { # exp: [sdate, dat_exp, memb_exp] # obs: [sdate, dat_obs, memb_obs] nexp <- as.numeric(dim(exp)[dat_dim]) nobs <- as.numeric(dim(obs)[dat_dim]) - exp_memb <- as.numeric(dim(exp)[memb_dim]) - obs_memb <- as.numeric(dim(obs)[memb_dim]) CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) @@ -400,16 +413,16 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', confhigh <- tanh(atanh(CORR) + qnorm(conf.upper) / sqrt(eno_expand - 3)) } -################################### NEW -# Remove nexp and nobs if dat_dim = NULL +################################### + # Remove nexp and nobs if dat_dim = NULL if (is.null(dat_dim) & !is.null(memb_dim)) { dim(CORR) <- dim(CORR)[3:length(dim(CORR))] if (pval) { - dim(p.val) <- dim(p.val)[3:length(dim(p.val))] + dim(p.val) <- dim(p.val)[3:length(dim(p.val))] } if (conf) { - dim(conflow) <- dim(conflow)[3:length(dim(conflow))] - dim(confhigh) <- dim(confhigh)[3:length(dim(confhigh))] + dim(conflow) <- dim(conflow)[3:length(dim(conflow))] + dim(confhigh) <- dim(confhigh)[3:length(dim(confhigh))] } } -- GitLab From 81272e111de70302341d0bc24e503bb8bb499644 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 14 Jul 2022 10:06:51 +0200 Subject: [PATCH 15/47] Fix code corr method 'kendall' and 'spearman' --- R/Corr.R | 21 ++++----------------- 1 file changed, 4 insertions(+), 17 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index efea0c8..3e400b3 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -289,7 +289,7 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', CORR[y, j] <- cor(exp[, y], obs[, j], use = "pairwise.complete.obs", method = method) - } + } } } #---------------------------------------- @@ -319,18 +319,13 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', nexp <- 1 nobs <- 1 CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) - for (j in 1:obs_memb) { for (y in 1:exp_memb) { - if (any(!is.na(exp[,y])) && sum(!is.na(obs[, j])) > 2) { CORR[, , y, j] <- cor(exp[, y], obs[, j], use = "pairwise.complete.obs", method = method) - } else { - NA #CORR[, i] <- NA } - } } } else { @@ -338,9 +333,7 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # obs: [sdate, dat_obs, memb_obs] nexp <- as.numeric(dim(exp)[dat_dim]) nobs <- as.numeric(dim(obs)[dat_dim]) - CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) - for (j in 1:obs_memb) { for (y in 1:exp_memb) { CORR[, , y, j] <- sapply(1:nobs, function(i) { @@ -349,19 +342,14 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', cor(exp[, x, y], obs[, i, j], use = "pairwise.complete.obs", method = method) - } else { - NA #CORR[, i] <- NA } }) }) - } } } - } - # if (pval) { # for (i in 1:nobs) { # p.val[, i] <- try(sapply(1:nexp, @@ -375,14 +363,13 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # } if (pval | conf) { - if (method == "kendall" | method == "spearman") { + if ((method == "kendall" | method == "spearman") && (!is.null(dat_dim) | !is.null(memb_dim))) { tmp <- apply(obs, c(1:length(dim(obs)))[-1], rank) # for memb_dim = NULL, 2; for memb_dim, c(2, 3) names(dim(tmp))[1] <- time_dim eno <- Eno(tmp, time_dim, ncores = ncores_input) - } else if (method == "pearson") { - eno <- Eno(obs, time_dim, ncores = ncores_input) + } else { # (dat_dim = NULL and memb_dim = NULL) or method = 'pearson' + eno <- Eno(obs, time_dim, ncores = ncores_input) } - if (is.null(memb_dim)) { eno_expand <- array(dim = c(nexp = nexp, nobs = nobs)) for (i in 1:nexp) { -- GitLab From e7f4ec1b7cb43572601b7a3d82682d25cdd5f791 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 14 Jul 2022 13:13:34 +0200 Subject: [PATCH 16/47] Revert changes due to failed pipeline --- R/Corr.R | 21 +++++++++++++++++---- 1 file changed, 17 insertions(+), 4 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index 3e400b3..efea0c8 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -289,7 +289,7 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', CORR[y, j] <- cor(exp[, y], obs[, j], use = "pairwise.complete.obs", method = method) - } + } } } #---------------------------------------- @@ -319,13 +319,18 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', nexp <- 1 nobs <- 1 CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) + for (j in 1:obs_memb) { for (y in 1:exp_memb) { + if (any(!is.na(exp[,y])) && sum(!is.na(obs[, j])) > 2) { CORR[, , y, j] <- cor(exp[, y], obs[, j], use = "pairwise.complete.obs", method = method) + } else { + NA #CORR[, i] <- NA } + } } } else { @@ -333,7 +338,9 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # obs: [sdate, dat_obs, memb_obs] nexp <- as.numeric(dim(exp)[dat_dim]) nobs <- as.numeric(dim(obs)[dat_dim]) + CORR <- array(dim = c(nexp = nexp, nobs = nobs, exp_memb = exp_memb, obs_memb = obs_memb)) + for (j in 1:obs_memb) { for (y in 1:exp_memb) { CORR[, , y, j] <- sapply(1:nobs, function(i) { @@ -342,14 +349,19 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', cor(exp[, x, y], obs[, i, j], use = "pairwise.complete.obs", method = method) + } else { + NA #CORR[, i] <- NA } }) }) + } } } + } + # if (pval) { # for (i in 1:nobs) { # p.val[, i] <- try(sapply(1:nexp, @@ -363,13 +375,14 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # } if (pval | conf) { - if ((method == "kendall" | method == "spearman") && (!is.null(dat_dim) | !is.null(memb_dim))) { + if (method == "kendall" | method == "spearman") { tmp <- apply(obs, c(1:length(dim(obs)))[-1], rank) # for memb_dim = NULL, 2; for memb_dim, c(2, 3) names(dim(tmp))[1] <- time_dim eno <- Eno(tmp, time_dim, ncores = ncores_input) - } else { # (dat_dim = NULL and memb_dim = NULL) or method = 'pearson' - eno <- Eno(obs, time_dim, ncores = ncores_input) + } else if (method == "pearson") { + eno <- Eno(obs, time_dim, ncores = ncores_input) } + if (is.null(memb_dim)) { eno_expand <- array(dim = c(nexp = nexp, nobs = nobs)) for (i in 1:nexp) { -- GitLab From 053a121892109017305f94ec7f5843124e1adda4 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 14 Jul 2022 13:25:46 +0200 Subject: [PATCH 17/47] Correct case corr method 'kendall' or 'spearman' --- R/Corr.R | 8 ++------ 1 file changed, 2 insertions(+), 6 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index efea0c8..fb24787 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -327,8 +327,6 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', CORR[, , y, j] <- cor(exp[, y], obs[, j], use = "pairwise.complete.obs", method = method) - } else { - NA #CORR[, i] <- NA } } @@ -349,8 +347,6 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', cor(exp[, x, y], obs[, i, j], use = "pairwise.complete.obs", method = method) - } else { - NA #CORR[, i] <- NA } }) }) @@ -375,11 +371,11 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # } if (pval | conf) { - if (method == "kendall" | method == "spearman") { + if ((method == "kendall" | method == "spearman") && (!is.null(dat_dim) | !is.null(memb_dim))) { tmp <- apply(obs, c(1:length(dim(obs)))[-1], rank) # for memb_dim = NULL, 2; for memb_dim, c(2, 3) names(dim(tmp))[1] <- time_dim eno <- Eno(tmp, time_dim, ncores = ncores_input) - } else if (method == "pearson") { + } else { eno <- Eno(obs, time_dim, ncores = ncores_input) } -- GitLab From ee9b8c85ebb1ed8833633a11095cb80809db06b2 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 14 Jul 2022 13:41:03 +0200 Subject: [PATCH 18/47] Revert commit --- R/Corr.R | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index fb24787..efea0c8 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -327,6 +327,8 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', CORR[, , y, j] <- cor(exp[, y], obs[, j], use = "pairwise.complete.obs", method = method) + } else { + NA #CORR[, i] <- NA } } @@ -347,6 +349,8 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', cor(exp[, x, y], obs[, i, j], use = "pairwise.complete.obs", method = method) + } else { + NA #CORR[, i] <- NA } }) }) @@ -371,11 +375,11 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # } if (pval | conf) { - if ((method == "kendall" | method == "spearman") && (!is.null(dat_dim) | !is.null(memb_dim))) { + if (method == "kendall" | method == "spearman") { tmp <- apply(obs, c(1:length(dim(obs)))[-1], rank) # for memb_dim = NULL, 2; for memb_dim, c(2, 3) names(dim(tmp))[1] <- time_dim eno <- Eno(tmp, time_dim, ncores = ncores_input) - } else { + } else if (method == "pearson") { eno <- Eno(obs, time_dim, ncores = ncores_input) } -- GitLab From 7b3ee5f54ede45e03303fda3a60ea7a1809ee9c8 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 14 Jul 2022 17:06:45 +0200 Subject: [PATCH 19/47] Correct case 'kendall' and 'separman' corr type --- R/Corr.R | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index efea0c8..d99a1bf 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -375,12 +375,12 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # } if (pval | conf) { - if (method == "kendall" | method == "spearman") { + if ((method == "kendall" | method == "spearman") && (!is.null(dat_dim) | !is.null(memb_dim))) { tmp <- apply(obs, c(1:length(dim(obs)))[-1], rank) # for memb_dim = NULL, 2; for memb_dim, c(2, 3) names(dim(tmp))[1] <- time_dim eno <- Eno(tmp, time_dim, ncores = ncores_input) - } else if (method == "pearson") { - eno <- Eno(obs, time_dim, ncores = ncores_input) + } else { + eno <- Eno(obs, time_dim, ncores = ncores_input) } if (is.null(memb_dim)) { -- GitLab From 3b4d67a367f3718a000f4ecab292faee213d7bca Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 14 Jul 2022 17:33:58 +0200 Subject: [PATCH 20/47] Delete comments #CORR[, i] <- NA --- R/Corr.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index d99a1bf..b3be256 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -328,7 +328,7 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', use = "pairwise.complete.obs", method = method) } else { - NA #CORR[, i] <- NA + NA } } @@ -350,7 +350,7 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', use = "pairwise.complete.obs", method = method) } else { - NA #CORR[, i] <- NA + NA } }) }) -- GitLab From 662e8393d5462f050fd4f1b0e8ccaeeed25d8bef Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Fri, 15 Jul 2022 11:34:03 +0200 Subject: [PATCH 21/47] Delete else part test --- R/Corr.R | 4 ---- 1 file changed, 4 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index b3be256..245b93a 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -327,8 +327,6 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', CORR[, , y, j] <- cor(exp[, y], obs[, j], use = "pairwise.complete.obs", method = method) - } else { - NA } } @@ -349,8 +347,6 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', cor(exp[, x, y], obs[, i, j], use = "pairwise.complete.obs", method = method) - } else { - NA } }) }) -- GitLab From 92f741219faeb9d93449d64db4be23cd83052965 Mon Sep 17 00:00:00 2001 From: aho Date: Fri, 15 Jul 2022 11:46:23 +0200 Subject: [PATCH 22/47] Add 'else' back for necessary part --- R/Corr.R | 2 ++ 1 file changed, 2 insertions(+) diff --git a/R/Corr.R b/R/Corr.R index 245b93a..d72df52 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -347,6 +347,8 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', cor(exp[, x, y], obs[, i, j], use = "pairwise.complete.obs", method = method) + } else { + NA } }) }) -- GitLab From 099de1fe3a37de68942a34e1bdd76f20e2c76689 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Fri, 15 Jul 2022 15:18:41 +0200 Subject: [PATCH 23/47] Correct case for rank corr method dat_dim and memb_dim NULL --- R/Corr.R | 19 +++++++++++++------ 1 file changed, 13 insertions(+), 6 deletions(-) diff --git a/R/Corr.R b/R/Corr.R index d72df52..ef710b5 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -373,12 +373,19 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', # } if (pval | conf) { - if ((method == "kendall" | method == "spearman") && (!is.null(dat_dim) | !is.null(memb_dim))) { - tmp <- apply(obs, c(1:length(dim(obs)))[-1], rank) # for memb_dim = NULL, 2; for memb_dim, c(2, 3) - names(dim(tmp))[1] <- time_dim - eno <- Eno(tmp, time_dim, ncores = ncores_input) - } else { - eno <- Eno(obs, time_dim, ncores = ncores_input) + if (method == "kendall" | method == "spearman") { + if (!is.null(dat_dim) | !is.null(memb_dim)) { + tmp <- apply(obs, c(1:length(dim(obs)))[-1], rank) # for memb_dim = NULL, 2; for memb_dim, c(2, 3) + names(dim(tmp))[1] <- time_dim + eno <- Eno(tmp, time_dim, ncores = ncores_input) + } else { + tmp <- rank(obs) + tmp <- array(tmp) + names(dim(tmp)) <- time_dim + eno <- Eno(tmp, time_dim, ncores = ncores_input) + } + } else if (method == "pearson") { + eno <- Eno(obs, time_dim, ncores = ncores_input) } if (is.null(memb_dim)) { -- GitLab From 5eb89fd16643d0a66fb90d36c1f3432ee0403c15 Mon Sep 17 00:00:00 2001 From: aho Date: Fri, 15 Jul 2022 20:48:50 +0200 Subject: [PATCH 24/47] Add test for kendall --- tests/testthat/test-Corr.R | 19 ++++++++++++++++++- 1 file changed, 18 insertions(+), 1 deletion(-) diff --git a/tests/testthat/test-Corr.R b/tests/testthat/test-Corr.R index 7d477ad..6e3c016 100644 --- a/tests/testthat/test-Corr.R +++ b/tests/testthat/test-Corr.R @@ -455,6 +455,23 @@ test_that("7. Output checks: dat6", { 0.1084488, tolerance = 0.0001 ) + # kendall + expect_equal( + as.vector(Corr(exp6, obs6, dat_dim = NULL, memb_dim = NULL, method = 'kendall')$corr[2:4]), + c(0.0, 0.6, 0.2), + tolerance = 0.0001 + ) + expect_equal( + as.vector(Corr(exp6, obs6, dat_dim = NULL, memb_dim = NULL, method = 'kendall')$p[2:4]), + c(0.5000000, 0.1423785, 0.3735300), + tolerance = 0.0001 + ) + expect_equal( + as.vector(Corr(exp6, obs6, dat_dim = NULL, memb_dim = NULL, method = 'kendall')$conf.low[2:4]), + c(-0.8822664, -0.5997500, -0.8284490), + tolerance = 0.0001 + ) + }) ############################################## test_that("8. Output checks: dat6 and dat7", { @@ -464,4 +481,4 @@ test_that("8. Output checks: dat6 and dat7", { tolerance = 0.0001 ) }) -############################################## \ No newline at end of file +############################################## -- GitLab From 2a1349d00d41d2ea3da692d392f202886c50c29e Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 12:43:24 +0200 Subject: [PATCH 25/47] Add space in text output --- R/Corr.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/Corr.R b/R/Corr.R index ef710b5..67efb09 100644 --- a/R/Corr.R +++ b/R/Corr.R @@ -137,7 +137,7 @@ Corr <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", - "Set it as NULL if there is no dataset dimension.") + " Set it as NULL if there is no dataset dimension.") } } ## comp_dim -- GitLab From b16946d180c3bae0f0db7d905cfb8cf6c1387b70 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 12:50:06 +0200 Subject: [PATCH 26/47] Allow dat_dim NULL in RMS, RMSSS and added tests for both --- R/RMS.R | 36 +++-- R/RMSSS.R | 259 +++++++++++++++++++++--------------- tests/testthat/test-RMS.R | 50 +++++++ tests/testthat/test-RMSSS.R | 27 ++++ 4 files changed, 258 insertions(+), 114 deletions(-) diff --git a/R/RMS.R b/R/RMS.R index 20da981..7a38a72 100644 --- a/R/RMS.R +++ b/R/RMS.R @@ -156,8 +156,8 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) if (!is.null(dat_dim)) { - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", @@ -202,12 +202,22 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } .RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { - - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) + conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { + if (is.null(dat_dim)) { + # exp: [sdate] + # obs: [sdate] + nexp <- 1 + nobs <- 1 + ini_dims <- dim(exp) + dim(exp) <- c(ini_dims, dat_dim = 1) + dim(obs) <- c(ini_dims, dat_dim = 1) + } else { + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) + } + nsdate <- as.numeric(dim(exp)[1]) dif <- array(dim = c(sdate = nsdate, nexp = nexp, nobs = nobs)) @@ -242,6 +252,16 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', conf.upper <- (eno * rms ** 2 / chi) ** 0.5 } + ################################### + # Remove nexp and nobs if dat_dim = NULL + if (is.null(dat_dim)) { + dim(rms) <- NULL + dim(conf.lower) <- NULL + dim(conf.upper) <- NULL + } + + ################################### + if (conf) { res <- list(rms = rms, conf.lower = conf.lower, conf.upper = conf.upper) } else { diff --git a/R/RMSSS.R b/R/RMSSS.R index 9526517..7a38a72 100644 --- a/R/RMSSS.R +++ b/R/RMSSS.R @@ -1,64 +1,75 @@ -#'Compute root mean square error skill score +#'Compute root mean square error #' -#'Compute the root mean square error skill score (RMSSS) between an array of -#'forecast 'exp' and an array of observation 'obs'. The two arrays should -#'have the same dimensions except along dat_dim, where the length can be -#'different, with the number of experiments/models (nexp) and the number of -#'observational datasets (nobs).\cr -#'RMSSS computes the root mean square error skill score of each jexp in 1:nexp -#'against each jobs in 1:nobs which gives nexp * nobs RMSSS for each other -#'grid point of the array.\cr -#'The RMSSS are computed along the time_dim dimension which should corresponds -#'to the startdate dimension.\cr -#'The p-value is optionally provided by an one-sided Fisher test.\cr +#'Compute the root mean square error for an array of forecasts and an array of +#'observations. The RMSEs are computed along time_dim, the dimension which +#'corresponds to the startdate dimension. If comp_dim is given, the RMSEs are +#'computed only if obs along the comp_dim dimension are complete between +#'limits[1] and limits[2], i.e. there are no NAs between limits[1] and +#'limits[2]. This option can be activated if the user wishes to account only +#'for the forecasts for which the corresponding observations are available at +#'all leadtimes.\cr +#'The confidence interval is computed by the chi2 distribution.\cr #' -#'@param exp A named numeric array of experimental data which contains at least -#' two dimensions for dat_dim and time_dim. It can also be a vector with the +#'@param exp A named numeric array of experimental data, with at least two +#' dimensions 'time_dim' and 'dat_dim'. It can also be a vector with the #' same length as 'obs', then the vector will automatically be 'time_dim' and #' 'dat_dim' will be 1. -#'@param obs A named numeric array of observational data which contains at least -#' two dimensions for dat_dim and time_dim. The dimensions should be the same -#' as paramter 'exp' except the length of 'dat_dim' dimension. The order of -#' dimension can be different. It can also be a vector with the same length as -#' 'exp', then the vector will automatically be 'time_dim' and 'dat_dim' will -#' be 1. -#'@param dat_dim A character string indicating the name of dataset (nobs/nexp) +#'@param obs A named numeric array of observational data, same dimensions as +#' parameter 'exp' except along dat_dim. It can also be a vector with the same +#' length as 'exp', then the vector will automatically be 'time_dim' and +#' 'dat_dim' will be 1. +#'@param time_dim A character string indicating the name of dimension along +#' which the correlations are computed. The default value is 'sdate'. +#'@param dat_dim A character string indicating the name of member (nobs/nexp) #' dimension. The default value is 'dataset'. If there is no dataset #' dimension, set NULL. -#'@param time_dim A character string indicating the name of dimension along -#' which the RMSSS are computed. The default value is 'sdate'. -#'@param pval A logical value indicating whether to compute or not the p-value -#' of the test Ho: RMSSS = 0. If pval = TRUE, the insignificant RMSSS will -#' return NA. The default value is TRUE. +#'@param 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. +#'@param 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)). +#'@param conf A logical value indicating whether to retrieve the confidence +#' intervals or not. The default value is TRUE. +#'@param conf.lev A numeric indicating the confidence level for the +#' regression computation. The default value is 0.95. #'@param ncores An integer indicating the number of cores to use for parallel #' computation. The default value is NULL. #' -#'@return +#'@return #'A list containing the numeric arrays with dimension:\cr #' c(nexp, nobs, all other dimensions of exp except time_dim).\cr #'nexp is the number of experiment (i.e., dat_dim in exp), and nobs is the #'number of observation (i.e., dat_dim in obs).\cr -#'\item{$rmsss}{ -#' The root mean square error skill score. +#'\item{$rms}{ +#' The root mean square error. +#'} +#'\item{$conf.lower}{ +#' The lower confidence interval. Only present if \code{conf = TRUE}. #'} -#'\item{$p.val}{ -#' The p-value. Only present if \code{pval = TRUE}. +#'\item{$conf.upper}{ +#' The upper confidence interval. Only present if \code{conf = TRUE}. #'} #' #'@examples -#' set.seed(1) -#' exp <- array(rnorm(30), dim = c(dataset = 2, time = 3, memb = 5)) -#' set.seed(2) -#' obs <- array(rnorm(15), dim = c(time = 3, memb = 5, dataset = 1)) -#' res <- RMSSS(exp, obs, time_dim = 'time') +#'# Load sample data as in Load() example: +#' set.seed(1) +#' exp1 <- array(rnorm(120), dim = c(dataset = 3, sdate = 5, ftime = 2, lon = 1, lat = 4)) +#' set.seed(2) +#' obs1 <- array(rnorm(80), dim = c(dataset = 2, sdate = 5, ftime = 2, lon = 1, lat = 4)) +#' set.seed(2) +#' na <- floor(runif(10, min = 1, max = 80)) +#' obs1[na] <- NA +#' res <- RMS(exp1, obs1, comp_dim = 'ftime') +#' # Renew example when Ano and Smoothing are ready #' -#'@rdname RMSSS +#'@rdname RMS #'@import multiApply -#'@importFrom stats pf +#'@importFrom ClimProjDiags Subset +#'@importFrom stats qchisq #'@export -RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - pval = TRUE, ncores = NULL) { - +RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', + comp_dim = NULL, limits = NULL, + conf = TRUE, conf.lev = 0.95, ncores = NULL) { # Check inputs ## exp and obs (1) if (is.null(exp) | is.null(obs)) { @@ -87,7 +98,7 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } if(!all(names(dim(exp)) %in% names(dim(obs))) | !all(names(dim(obs)) %in% names(dim(exp)))) { - stop("Parameter 'exp' and 'obs' must have same dimension name.") + stop("Parameter 'exp' and 'obs' must have same dimension name") } ## time_dim if (!is.character(time_dim) | length(time_dim) > 1) { @@ -106,9 +117,33 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', "Set it as NULL if there is no dataset dimension.") } } - ## pval - if (!is.logical(pval) | length(pval) > 1) { - stop("Parameter 'pval' must be one logical value.") + ## comp_dim + if (!is.null(comp_dim)) { + if (!is.character(comp_dim) | length(comp_dim) > 1) { + stop("Parameter 'comp_dim' must be a character string.") + } + if (!comp_dim %in% names(dim(exp)) | !comp_dim %in% names(dim(obs))) { + stop("Parameter 'comp_dim' is not found in 'exp' or 'obs' dimension.") + } + } + ## limits + if (!is.null(limits)) { + if (is.null(comp_dim)) { + stop("Paramter 'comp_dim' cannot be NULL if 'limits' is assigned.") + } + if (!is.numeric(limits) | any(limits %% 1 != 0) | any(limits < 0) | + length(limits) != 2 | any(limits > dim(exp)[comp_dim])) { + stop(paste0("Parameter 'limits' must be a vector of two positive ", + "integers smaller than the length of paramter 'comp_dim'.")) + } + } + ## conf + if (!is.logical(conf) | length(conf) > 1) { + stop("Parameter 'conf' must be one logical value.") + } + ## conf.lev + if (!is.numeric(conf.lev) | conf.lev < 0 | conf.lev > 1 | length(conf.lev) > 1) { + stop("Parameter 'conf.lev' must be a numeric number between 0 and 1.") } ## ncores if (!is.null(ncores)) { @@ -121,17 +156,17 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) if (!is.null(dat_dim)) { - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", "all dimension except 'dat_dim'.")) } - if (dim(exp)[time_dim] <= 2) { - stop("The length of time_dim must be more than 2 to compute RMSSS.") + if (dim(exp)[time_dim] < 2) { + stop("The length of time_dim must be at least 2 to compute RMS.") } - + ############################### # Sort dimension @@ -142,85 +177,97 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ############################### - # Calculate RMSSS + # Calculate RMS + + # Remove data along comp_dim dim if there is at least one NA between limits + if (!is.null(comp_dim)) { + if (is.null(limits)) { + limits <- c(1, dim(obs)[comp_dim]) + } + pos <- which(names(dim(obs)) == comp_dim) + obs_sub <- Subset(obs, pos, list(limits[1]:limits[2])) + outrows <- is.na(MeanDims(obs_sub, pos, na.rm = FALSE)) + outrows <- InsertDim(outrows, pos, dim(obs)[comp_dim]) + obs[which(outrows)] <- NA + } res <- Apply(list(exp, obs), target_dims = list(c(time_dim, dat_dim), c(time_dim, dat_dim)), - fun = .RMSSS, + fun = .RMS, time_dim = time_dim, dat_dim = dat_dim, - pval = pval, ncores_input = ncores, + conf = conf, conf.lev = conf.lev, ncores_input = ncores, ncores = ncores) - return(res) } -.RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', pval = TRUE, - ncores_input = NULL) { - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) - nsdate <- as.numeric(dim(exp)[1]) +.RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', + conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { + if (is.null(dat_dim)) { + # exp: [sdate] + # obs: [sdate] + nexp <- 1 + nobs <- 1 + ini_dims <- dim(exp) + dim(exp) <- c(ini_dims, dat_dim = 1) + dim(obs) <- c(ini_dims, dat_dim = 1) + } else { + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) + } - p_val <- array(dim = c(nexp = nexp, nobs = nobs)) - dif1 <- array(dim = c(nsdate, nexp, nobs)) - names(dim(dif1)) <- c(time_dim, 'nexp', 'nobs') + nsdate <- as.numeric(dim(exp)[1]) -# if (conf) { -# conflow <- (1 - conf.lev) / 2 -# confhigh <- 1 - conflow -# conf_low <- array(dim = c(nexp = nexp, nobs = nobs)) -# conf_high <- array(dim = c(nexp = nexp, nobs = nobs)) -# } + dif <- array(dim = c(sdate = nsdate, nexp = nexp, nobs = nobs)) + chi <- array(dim = c(nexp = nexp, nobs = nobs)) + if (conf) { + conflow <- (1 - conf.lev) / 2 + confhigh <- 1 - conflow + conf.lower <- array(dim = c(nexp = nexp, nobs = nobs)) + conf.upper <- array(dim = c(nexp = nexp, nobs = nobs)) + } - # dif1 + # dif for (i in 1:nobs) { - dif1[, , i] <- sapply(1:nexp, function(x) {exp[, x] - obs[, i]}) + dif[, , i] <- sapply(1:nexp, function(x) {exp[, x] - obs[, i]}) } + rms <- apply(dif^2, c(2, 3), mean, na.rm = TRUE)^0.5 #array(dim = c(_exp, nobs)) - # rms1 and eno1 - rms1 <- apply(dif1^2, c(2, 3), mean, na.rm = TRUE)^0.5 #array(dim = c(nexp, nobs)) - # rms2 and eno2 - rms2 <- array(colMeans(obs^2, na.rm = TRUE)^0.5, dim = c(nobs = nobs)) - rms2[which(abs(rms2) <= (max(abs(rms2), na.rm = TRUE) / 1000))] <- max(abs( - rms2), na.rm = TRUE) / 1000 - #rms2 above: [nobs] - rms2 <- array(rms2, dim = c(nobs = nobs, nexp = nexp)) - #rms2 above: [nobs, nexp] - rms2 <- Reorder(rms2, c(2, 1)) - #rms2 above: [nexp, nobs] + if (conf) { + #eno <- Eno(dif, 1) #count effective sample along sdate. dim = c(nexp, nobs) + eno <- Eno(dif, time_dim, ncores = ncores_input) #change to this line when Eno() is done - # use rms1 and rms2 to calculate rmsss - rmsss <- 1 - rms1/rms2 + # conf.lower + chi <- sapply(1:nobs, function(i) { + qchisq(confhigh, eno[, i] - 1) + }) + conf.lower <- (eno * rms ** 2 / chi) ** 0.5 - ## pval and conf - if (pval) { - eno1 <- Eno(dif1, time_dim, ncores = ncores_input) - eno2 <- Eno(obs, time_dim, ncores = ncores_input) - eno2 <- array(eno2, dim = c(nobs = nobs, nexp = nexp)) - eno2 <- Reorder(eno2, c(2, 1)) + # conf.upper + chi <- sapply(1:nobs, function(i) { + qchisq(conflow, eno[, i] - 1) + }) + conf.upper <- (eno * rms ** 2 / chi) ** 0.5 } - # pval - if (pval) { - - F.stat <- (eno2 * rms2^2 / (eno2- 1)) / ((eno1 * rms1^2 / (eno1- 1))) - tmp <- !is.na(eno1) & !is.na(eno2) & eno1 > 2 & eno2 > 2 - p_val <- 1 - pf(F.stat, eno1 - 1, eno2 - 1) - p_val[which(!tmp)] <- NA - - # change not significant rmsss to NA - rmsss[which(!tmp)] <- NA + ################################### + # Remove nexp and nobs if dat_dim = NULL + if (is.null(dat_dim)) { + dim(rms) <- NULL + dim(conf.lower) <- NULL + dim(conf.upper) <- NULL } - # output - if (pval) { - res <- list(rmsss = rmsss, p.val = p_val) + ################################### + if (conf) { + res <- list(rms = rms, conf.lower = conf.lower, conf.upper = conf.upper) } else { - res <- list(rmsss = rmsss) - } + res <- list(rms = rms) + } return(res) + } diff --git a/tests/testthat/test-RMS.R b/tests/testthat/test-RMS.R index f1083dd..65edfac 100644 --- a/tests/testthat/test-RMS.R +++ b/tests/testthat/test-RMS.R @@ -17,6 +17,13 @@ context("s2dv::RMS tests") set.seed(6) obs2 <- rnorm(10) + # dat3 + set.seed(1) + exp3 <- array(rnorm(120), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) + + set.seed(2) + obs3 <- array(rnorm(80), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) + ############################################## test_that("1. Input checks", { @@ -186,3 +193,46 @@ test_that("3. Output checks: dat2", { }) ############################################## +test_that("4. Output checks: dat3", { + + expect_equal( + dim(RMS(exp3, obs3, dat_dim = NULL)$rms), + c(ftime = 2, lon = 1, lat = 4) + ) + + suppressWarnings( + expect_equal( + RMS(exp3, obs3, dat_dim = NULL)$rms[1:6], + c(1.6458118, 0.8860392, 0.8261295, 1.1681939, 2.1693538, 1.3064454), + tolerance = 0.001 + ) +) +suppressWarnings( + expect_equal( + length(which(is.na(RMS(exp3, obs3, dat_dim = NULL)$conf.lower))), + 0 + ) +) +suppressWarnings( + expect_equal( + max(RMS(exp3, obs3, dat_dim = NULL)$conf.lower, na.rm = T), + 1.453144, + tolerance = 0.001 + ) +) +suppressWarnings( + expect_equal( + min(RMS(exp3, obs3, dat_dim = NULL, conf.lev = 0.99)$conf.upper, na.rm = TRUE), + 2.646274, + tolerance = 0.0001 + ) +) +suppressWarnings( + expect_equal( + length(RMS(exp3, obs3, dat_dim = NULL, conf = FALSE)), + 1 + ) +) +}) + +############################################## \ No newline at end of file diff --git a/tests/testthat/test-RMSSS.R b/tests/testthat/test-RMSSS.R index 91ef930..60b4672 100644 --- a/tests/testthat/test-RMSSS.R +++ b/tests/testthat/test-RMSSS.R @@ -25,6 +25,12 @@ context("s2dv::RMSSS tests") set.seed(6) obs3 <- rnorm(10) + # case 4 + set.seed(7) + exp4 <- array(rnorm(120), dim = c(sdate = 10, dat = 1, lon = 3, lat = 2)) + set.seed(8) + obs4 <- array(rnorm(60), dim = c(dat = 1, sdate = 10, lat = 2, lon = 3)) + ############################################## test_that("1. Input checks", { @@ -157,5 +163,26 @@ test_that("4. Output checks: case 3", { tolerance = 0.00001 ) +############################################## +test_that("4. Output checks: case 4", { + + expect_equal( + dim(RMSSS(exp4, obs4, dat_dim = NULL)$rmsss), + c(dat = 1, lon = 3, lat = 2) + ) + expect_equal( + mean(RMSSS(exp4, obs4, dat_dim = NULL)$rmsss), + -0.3114424, + tolerance = 0.00001 + ) + expect_equal( + range(RMSSS(exp4, obs4, dat_dim = NULL)$p.val), + c(0.3560534, 0.9192801), + tolerance = 0.00001 + ) + +}) + +############################################## }) -- GitLab From 9015838a870424b2e7d08129fab0fac7db232e52 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 12:52:14 +0200 Subject: [PATCH 27/47] Allow dat_dim NULL in UltimateBrier and added test --- R/UltimateBrier.R | 51 +++++++--- tests/testthat/test-UltimateBrier.R | 142 ++++++++++++++++++++++++++++ 2 files changed, 182 insertions(+), 11 deletions(-) diff --git a/R/UltimateBrier.R b/R/UltimateBrier.R index a2c594e..18c8e5b 100644 --- a/R/UltimateBrier.R +++ b/R/UltimateBrier.R @@ -107,10 +107,10 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti stop("Parameter 'dat_dim' must be a character string or NULL.") } if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", - "Set it as NULL if there is no dataset dimension.") + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") } } + ## memb_dim if (!is.character(memb_dim) | length(memb_dim) > 1) { stop("Parameter 'memb_dim' must be a character string.") @@ -136,10 +136,8 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti name_obs <- sort(names(dim(obs))) name_exp <- name_exp[-which(name_exp == memb_dim)] name_obs <- name_obs[-which(name_obs == memb_dim)] - if (!is.null(dat_dim)) { name_exp <- name_exp[-which(name_exp == dat_dim)] name_obs <- name_obs[-which(name_obs == dat_dim)] - } if (any(name_exp != name_obs)) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", "of all the dimensions except 'dat_dim' and 'memb_dim'.")) @@ -188,6 +186,7 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti target_dims = list(c(time_dim, dat_dim, memb_dim), c(time_dim, dat_dim, memb_dim)), fun = .UltimateBrier, + dat_dim = dat_dim, memb_dim = memb_dim, thr = thr, type = type, decomposition = decomposition, ncores = ncores)$output1 @@ -208,6 +207,7 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti target_dims = list(c(time_dim, dat_dim), c(time_dim, dat_dim)), fun = .UltimateBrier, + dat_dim = dat_dim, memb_dim = memb_dim, thr = thr, type = type, decomposition = decomposition, ncores = ncores) @@ -222,8 +222,8 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti return(res) } -.UltimateBrier <- function(exp, obs, thr = c(5/100, 95/100), type = 'BS', - decomposition = TRUE) { +.UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', thr = c(5/100, 95/100), + type = 'BS', decomposition = TRUE) { # If exp and obs are probablistics # exp: [sdate, nexp] # obs: [sdate, nobs] @@ -234,6 +234,10 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti #NOTE: 'thr' is used in 'FairEnsembleBSS' and 'FairEnsembleBS'. But if quantile = F and # thr is real value, does it work? if (type == 'FairEnsembleBSS') { + if (is.null(dat_dim)) { + obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') + exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') + } size_ens_ref <- prod(dim(obs)[c(1, 3)]) res <- array(dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]), @@ -250,6 +254,9 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } } } + if (is.null(dat_dim)) { + dim(res) <- dim(res)[3:length(dim(res))] + } } else if (type == 'FairEnsembleBS') { #NOTE: The calculation in s2dverification::UltimateBrier is wrong. In the final stage, @@ -257,6 +264,10 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti # but the 3rd dim of result is 'bins' instead of decomposition. 'FairEnsembleBS' does # not have decomposition. # The calculation is fixed here. + if (is.null(dat_dim)) { + obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') + exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') + } res <- array(dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]), bin = length(thr) + 1)) @@ -269,10 +280,17 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } } } + if (is.null(dat_dim)) { + dim(res) <- dim(res)[3:length(dim(res))] + } # tmp <- res[, , 1] - res[, , 2] + res[, , 3] # res <- array(tmp, dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]))) } else if (type == 'BS') { + if (is.null(dat_dim)) { + obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') + exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') + } comp <- array(dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]), comp = 3)) @@ -285,17 +303,28 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } if (decomposition) { rel <- comp[, , 1] - dim(rel) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) res <- comp[, , 2] - dim(res) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) unc <- comp[, , 3] - dim(unc) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) bs <- rel - res + unc - dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + if (is.null(dat_dim)) { + dim(rel) <- NULL + dim(res) <- NULL + dim(unc) <- NULL + dim(bs) <- NULL + } else { + dim(rel) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + dim(res) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + dim(unc) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + } res <- list(bs = bs, rel = rel, res = res, unc = unc) } else { bs <- comp[, , 1] - comp[, , 2] + comp[, , 3] - dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + if (is.null(dat_dim)) { + dim(bs) <- NULL + } else { + dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + } res <- list(bs = bs) } diff --git a/tests/testthat/test-UltimateBrier.R b/tests/testthat/test-UltimateBrier.R index 0eb0118..e31f807 100644 --- a/tests/testthat/test-UltimateBrier.R +++ b/tests/testthat/test-UltimateBrier.R @@ -7,6 +7,12 @@ exp1 <- array(rnorm(30), dim = c(dataset = 1, member = 3, sdate = 5, ftime = 2)) set.seed(2) obs1 <- array(round(rnorm(10)), dim = c(dataset = 1, sdate = 5, ftime = 2)) +# dat2 +set.seed(1) +exp2 <- array(rnorm(30), dim = c(member = 3, sdate = 5, ftime = 2)) +set.seed(2) +obs2 <- array(round(rnorm(10)), dim = c(sdate = 5, ftime = 2)) + ############################################## test_that("1. Input checks", { @@ -238,3 +244,139 @@ tolerance = 0.0001 }) +############################################## +test_that("3. Output checks: dat2", { + +# 'BS' +expect_equal( +is.list(UltimateBrier(exp2, obs2, dat_dim = NULL)), +TRUE +) +expect_equal( +names(UltimateBrier(exp2, obs2, dat_dim = NULL)), +c('bs', 'rel', 'res', 'unc') +) +expect_equal( +is.list(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE)), +FALSE +) +expect_equal( +dim(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE)), +c(bin = 3, ftime = 2) +) +expect_equal( +dim(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, thr = c(0.25, 0.5, 0.75))), +c(bin = 4, ftime = 2) +) +expect_equal( +UltimateBrier(exp2, obs2, dat_dim = NULL)$bs, +UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE) +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$bs), +c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), +tolerance = 0.0001 +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$rel), +c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), +tolerance = 0.0001 +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$res), +c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), +tolerance = 0.0001 +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$unc), +c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), +tolerance = 0.0001 +) + +# 'BSS' +expect_equal( +dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'BSS')), +c(bin = 3, ftime = 2) +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'BSS')), +c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), +tolerance = 0.0001 +) + +# 'FairStartDatesBS' +expect_equal( +is.list(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')), +TRUE +) +expect_equal( +names(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')), +c('bs', 'rel', 'res', 'unc') +) +expect_equal( +is.list(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, type = 'FairStartDatesBS')), +FALSE +) +expect_equal( +dim(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, type = 'FairStartDatesBS')), +c(bin = 3, ftime = 2) +) +expect_equal( +UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$bs, +UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, type = 'FairStartDatesBS') +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$bs), +c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), +tolerance = 0.0001 +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$rel), +c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), +tolerance = 0.0001 +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$res), +c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), +tolerance = 0.0001 +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$unc), +c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), +tolerance = 0.0001 +) + +# 'FairStartDatesBSS' +expect_equal( +dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBSS')), +c(bin = 3, ftime = 2) +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBSS')), +c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), +tolerance = 0.0001 +) +# 'FairEnsembleBS' +expect_equal( +dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBS')), +c(bin = 3, ftime = 2) +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBS')), +c(0.1333333, 0.2000000, 0.2000000, 0.1333333, 0.4000000, 0.2000000), +tolerance = 0.0001 +) +# 'FairEnsembleBSS' +expect_equal( +dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBSS')), +c(bin = 3, ftime = 2) +) +expect_equal( +as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBSS')), +c(-0.1111111, -0.6666667, -0.6666667, 0.2592593, -1.2222222, -0.6666667), +tolerance = 0.0001 +) + +}) + + -- GitLab From 966701a0ae263fe90836b02ab38350186fadfe22 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 13:17:41 +0200 Subject: [PATCH 28/47] Revert commits --- R/RMS.R | 56 ++---- R/RMSSS.R | 275 +++++++++++----------------- R/UltimateBrier.R | 62 ++----- tests/testthat/test-RMS.R | 52 +----- tests/testthat/test-RMSSS.R | 29 +-- tests/testthat/test-UltimateBrier.R | 146 +-------------- 6 files changed, 144 insertions(+), 476 deletions(-) diff --git a/R/RMS.R b/R/RMS.R index 7a38a72..b3c8ad4 100644 --- a/R/RMS.R +++ b/R/RMS.R @@ -21,8 +21,7 @@ #'@param time_dim A character string indicating the name of dimension along #' which the correlations are computed. The default value is 'sdate'. #'@param dat_dim A character string indicating the name of member (nobs/nexp) -#' dimension. The default value is 'dataset'. If there is no dataset -#' dimension, set NULL. +#' dimension. The default value is 'dataset'. #'@param 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. @@ -108,14 +107,11 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.null(dat_dim)) { - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string or NULL.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", - "Set it as NULL if there is no dataset dimension.") - } + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") } ## comp_dim if (!is.null(comp_dim)) { @@ -155,13 +151,11 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - if (!is.null(dat_dim)) { - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] - } + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension except 'dat_dim'.")) + "all dimension expect 'dat_dim'.")) } if (dim(exp)[time_dim] < 2) { stop("The length of time_dim must be at least 2 to compute RMS.") @@ -202,22 +196,12 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } .RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { - if (is.null(dat_dim)) { - # exp: [sdate] - # obs: [sdate] - nexp <- 1 - nobs <- 1 - ini_dims <- dim(exp) - dim(exp) <- c(ini_dims, dat_dim = 1) - dim(obs) <- c(ini_dims, dat_dim = 1) - } else { - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) - } - + conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { + + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) nsdate <- as.numeric(dim(exp)[1]) dif <- array(dim = c(sdate = nsdate, nexp = nexp, nobs = nobs)) @@ -252,16 +236,6 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', conf.upper <- (eno * rms ** 2 / chi) ** 0.5 } - ################################### - # Remove nexp and nobs if dat_dim = NULL - if (is.null(dat_dim)) { - dim(rms) <- NULL - dim(conf.lower) <- NULL - dim(conf.upper) <- NULL - } - - ################################### - if (conf) { res <- list(rms = rms, conf.lower = conf.lower, conf.upper = conf.upper) } else { diff --git a/R/RMSSS.R b/R/RMSSS.R index 7a38a72..5fa9659 100644 --- a/R/RMSSS.R +++ b/R/RMSSS.R @@ -1,75 +1,63 @@ -#'Compute root mean square error +#'Compute root mean square error skill score #' -#'Compute the root mean square error for an array of forecasts and an array of -#'observations. The RMSEs are computed along time_dim, the dimension which -#'corresponds to the startdate dimension. If comp_dim is given, the RMSEs are -#'computed only if obs along the comp_dim dimension are complete between -#'limits[1] and limits[2], i.e. there are no NAs between limits[1] and -#'limits[2]. This option can be activated if the user wishes to account only -#'for the forecasts for which the corresponding observations are available at -#'all leadtimes.\cr -#'The confidence interval is computed by the chi2 distribution.\cr +#'Compute the root mean square error skill score (RMSSS) between an array of +#'forecast 'exp' and an array of observation 'obs'. The two arrays should +#'have the same dimensions except along dat_dim, where the length can be +#'different, with the number of experiments/models (nexp) and the number of +#'observational datasets (nobs).\cr +#'RMSSS computes the root mean square error skill score of each jexp in 1:nexp +#'against each jobs in 1:nobs which gives nexp * nobs RMSSS for each other +#'grid point of the array.\cr +#'The RMSSS are computed along the time_dim dimension which should corresponds +#'to the startdate dimension.\cr +#'The p-value is optionally provided by an one-sided Fisher test.\cr #' -#'@param exp A named numeric array of experimental data, with at least two -#' dimensions 'time_dim' and 'dat_dim'. It can also be a vector with the +#'@param exp A named numeric array of experimental data which contains at least +#' two dimensions for dat_dim and time_dim. It can also be a vector with the #' same length as 'obs', then the vector will automatically be 'time_dim' and #' 'dat_dim' will be 1. -#'@param obs A named numeric array of observational data, same dimensions as -#' parameter 'exp' except along dat_dim. It can also be a vector with the same -#' length as 'exp', then the vector will automatically be 'time_dim' and -#' 'dat_dim' will be 1. +#'@param obs A named numeric array of observational data which contains at least +#' two dimensions for dat_dim and time_dim. The dimensions should be the same +#' as paramter 'exp' except the length of 'dat_dim' dimension. The order of +#' dimension can be different. It can also be a vector with the same length as +#' 'exp', then the vector will automatically be 'time_dim' and 'dat_dim' will +#' be 1. +#'@param dat_dim A character string indicating the name of dataset (nobs/nexp) +#' dimension. The default value is 'dataset'. #'@param time_dim A character string indicating the name of dimension along -#' which the correlations are computed. The default value is 'sdate'. -#'@param dat_dim A character string indicating the name of member (nobs/nexp) -#' dimension. The default value is 'dataset'. If there is no dataset -#' dimension, set NULL. -#'@param 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. -#'@param 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)). -#'@param conf A logical value indicating whether to retrieve the confidence -#' intervals or not. The default value is TRUE. -#'@param conf.lev A numeric indicating the confidence level for the -#' regression computation. The default value is 0.95. +#' which the RMSSS are computed. The default value is 'sdate'. +#'@param pval A logical value indicating whether to compute or not the p-value +#' of the test Ho: RMSSS = 0. If pval = TRUE, the insignificant RMSSS will +#' return NA. The default value is TRUE. #'@param ncores An integer indicating the number of cores to use for parallel #' computation. The default value is NULL. #' -#'@return +#'@return #'A list containing the numeric arrays with dimension:\cr #' c(nexp, nobs, all other dimensions of exp except time_dim).\cr #'nexp is the number of experiment (i.e., dat_dim in exp), and nobs is the #'number of observation (i.e., dat_dim in obs).\cr -#'\item{$rms}{ -#' The root mean square error. -#'} -#'\item{$conf.lower}{ -#' The lower confidence interval. Only present if \code{conf = TRUE}. +#'\item{$rmsss}{ +#' The root mean square error skill score. #'} -#'\item{$conf.upper}{ -#' The upper confidence interval. Only present if \code{conf = TRUE}. +#'\item{$p.val}{ +#' The p-value. Only present if \code{pval = TRUE}. #'} #' #'@examples -#'# Load sample data as in Load() example: -#' set.seed(1) -#' exp1 <- array(rnorm(120), dim = c(dataset = 3, sdate = 5, ftime = 2, lon = 1, lat = 4)) -#' set.seed(2) -#' obs1 <- array(rnorm(80), dim = c(dataset = 2, sdate = 5, ftime = 2, lon = 1, lat = 4)) -#' set.seed(2) -#' na <- floor(runif(10, min = 1, max = 80)) -#' obs1[na] <- NA -#' res <- RMS(exp1, obs1, comp_dim = 'ftime') -#' # Renew example when Ano and Smoothing are ready +#' set.seed(1) +#' exp <- array(rnorm(30), dim = c(dataset = 2, time = 3, memb = 5)) +#' set.seed(2) +#' obs <- array(rnorm(15), dim = c(time = 3, memb = 5, dataset = 1)) +#' res <- RMSSS(exp, obs, time_dim = 'time') #' -#'@rdname RMS +#'@rdname RMSSS #'@import multiApply -#'@importFrom ClimProjDiags Subset -#'@importFrom stats qchisq +#'@importFrom stats pf #'@export -RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - comp_dim = NULL, limits = NULL, - conf = TRUE, conf.lev = 0.95, ncores = NULL) { +RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', + pval = TRUE, ncores = NULL) { + # Check inputs ## exp and obs (1) if (is.null(exp) | is.null(obs)) { @@ -98,7 +86,7 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } if(!all(names(dim(exp)) %in% names(dim(obs))) | !all(names(dim(obs)) %in% names(dim(exp)))) { - stop("Parameter 'exp' and 'obs' must have same dimension name") + stop("Parameter 'exp' and 'obs' must have same dimension name.") } ## time_dim if (!is.character(time_dim) | length(time_dim) > 1) { @@ -108,42 +96,15 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.null(dat_dim)) { - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string or NULL.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", - "Set it as NULL if there is no dataset dimension.") - } + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string.") } - ## comp_dim - if (!is.null(comp_dim)) { - if (!is.character(comp_dim) | length(comp_dim) > 1) { - stop("Parameter 'comp_dim' must be a character string.") - } - if (!comp_dim %in% names(dim(exp)) | !comp_dim %in% names(dim(obs))) { - stop("Parameter 'comp_dim' is not found in 'exp' or 'obs' dimension.") - } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") } - ## limits - if (!is.null(limits)) { - if (is.null(comp_dim)) { - stop("Paramter 'comp_dim' cannot be NULL if 'limits' is assigned.") - } - if (!is.numeric(limits) | any(limits %% 1 != 0) | any(limits < 0) | - length(limits) != 2 | any(limits > dim(exp)[comp_dim])) { - stop(paste0("Parameter 'limits' must be a vector of two positive ", - "integers smaller than the length of paramter 'comp_dim'.")) - } - } - ## conf - if (!is.logical(conf) | length(conf) > 1) { - stop("Parameter 'conf' must be one logical value.") - } - ## conf.lev - if (!is.numeric(conf.lev) | conf.lev < 0 | conf.lev > 1 | length(conf.lev) > 1) { - stop("Parameter 'conf.lev' must be a numeric number between 0 and 1.") + ## pval + if (!is.logical(pval) | length(pval) > 1) { + stop("Parameter 'pval' must be one logical value.") } ## ncores if (!is.null(ncores)) { @@ -155,18 +116,16 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - if (!is.null(dat_dim)) { - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] - } + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension except 'dat_dim'.")) + "all dimension expect 'dat_dim'.")) } - if (dim(exp)[time_dim] < 2) { - stop("The length of time_dim must be at least 2 to compute RMS.") + if (dim(exp)[time_dim] <= 2) { + stop("The length of time_dim must be more than 2 to compute RMSSS.") } - + ############################### # Sort dimension @@ -177,97 +136,85 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ############################### - # Calculate RMS - - # Remove data along comp_dim dim if there is at least one NA between limits - if (!is.null(comp_dim)) { - if (is.null(limits)) { - limits <- c(1, dim(obs)[comp_dim]) - } - pos <- which(names(dim(obs)) == comp_dim) - obs_sub <- Subset(obs, pos, list(limits[1]:limits[2])) - outrows <- is.na(MeanDims(obs_sub, pos, na.rm = FALSE)) - outrows <- InsertDim(outrows, pos, dim(obs)[comp_dim]) - obs[which(outrows)] <- NA - } + # Calculate RMSSS res <- Apply(list(exp, obs), target_dims = list(c(time_dim, dat_dim), c(time_dim, dat_dim)), - fun = .RMS, + fun = .RMSSS, time_dim = time_dim, dat_dim = dat_dim, - conf = conf, conf.lev = conf.lev, ncores_input = ncores, + pval = pval, ncores_input = ncores, ncores = ncores) + return(res) } -.RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { - if (is.null(dat_dim)) { - # exp: [sdate] - # obs: [sdate] - nexp <- 1 - nobs <- 1 - ini_dims <- dim(exp) - dim(exp) <- c(ini_dims, dat_dim = 1) - dim(obs) <- c(ini_dims, dat_dim = 1) - } else { - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) - } - +.RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', pval = TRUE, + ncores_input = NULL) { + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) nsdate <- as.numeric(dim(exp)[1]) + + p_val <- array(dim = c(nexp = nexp, nobs = nobs)) + dif1 <- array(dim = c(nsdate, nexp, nobs)) + names(dim(dif1)) <- c(time_dim, 'nexp', 'nobs') - dif <- array(dim = c(sdate = nsdate, nexp = nexp, nobs = nobs)) - chi <- array(dim = c(nexp = nexp, nobs = nobs)) - if (conf) { - conflow <- (1 - conf.lev) / 2 - confhigh <- 1 - conflow - conf.lower <- array(dim = c(nexp = nexp, nobs = nobs)) - conf.upper <- array(dim = c(nexp = nexp, nobs = nobs)) - } +# if (conf) { +# conflow <- (1 - conf.lev) / 2 +# confhigh <- 1 - conflow +# conf_low <- array(dim = c(nexp = nexp, nobs = nobs)) +# conf_high <- array(dim = c(nexp = nexp, nobs = nobs)) +# } - # dif + # dif1 for (i in 1:nobs) { - dif[, , i] <- sapply(1:nexp, function(x) {exp[, x] - obs[, i]}) + dif1[, , i] <- sapply(1:nexp, function(x) {exp[, x] - obs[, i]}) } - rms <- apply(dif^2, c(2, 3), mean, na.rm = TRUE)^0.5 #array(dim = c(_exp, nobs)) - if (conf) { - #eno <- Eno(dif, 1) #count effective sample along sdate. dim = c(nexp, nobs) - eno <- Eno(dif, time_dim, ncores = ncores_input) #change to this line when Eno() is done + # rms1 and eno1 + rms1 <- apply(dif1^2, c(2, 3), mean, na.rm = TRUE)^0.5 #array(dim = c(nexp, nobs)) + # rms2 and eno2 + rms2 <- array(colMeans(obs^2, na.rm = TRUE)^0.5, dim = c(nobs = nobs)) + rms2[which(abs(rms2) <= (max(abs(rms2), na.rm = TRUE) / 1000))] <- max(abs( + rms2), na.rm = TRUE) / 1000 + #rms2 above: [nobs] + rms2 <- array(rms2, dim = c(nobs = nobs, nexp = nexp)) + #rms2 above: [nobs, nexp] + rms2 <- Reorder(rms2, c(2, 1)) + #rms2 above: [nexp, nobs] - # conf.lower - chi <- sapply(1:nobs, function(i) { - qchisq(confhigh, eno[, i] - 1) - }) - conf.lower <- (eno * rms ** 2 / chi) ** 0.5 + # use rms1 and rms2 to calculate rmsss + rmsss <- 1 - rms1/rms2 - # conf.upper - chi <- sapply(1:nobs, function(i) { - qchisq(conflow, eno[, i] - 1) - }) - conf.upper <- (eno * rms ** 2 / chi) ** 0.5 + ## pval and conf + if (pval) { + eno1 <- Eno(dif1, time_dim, ncores = ncores_input) + eno2 <- Eno(obs, time_dim, ncores = ncores_input) + eno2 <- array(eno2, dim = c(nobs = nobs, nexp = nexp)) + eno2 <- Reorder(eno2, c(2, 1)) } - ################################### - # Remove nexp and nobs if dat_dim = NULL - if (is.null(dat_dim)) { - dim(rms) <- NULL - dim(conf.lower) <- NULL - dim(conf.upper) <- NULL + # pval + if (pval) { + + F.stat <- (eno2 * rms2^2 / (eno2- 1)) / ((eno1 * rms1^2 / (eno1- 1))) + tmp <- !is.na(eno1) & !is.na(eno2) & eno1 > 2 & eno2 > 2 + p_val <- 1 - pf(F.stat, eno1 - 1, eno2 - 1) + p_val[which(!tmp)] <- NA + + # change not significant rmsss to NA + rmsss[which(!tmp)] <- NA } - ################################### + # output + if (pval) { + res <- list(rmsss = rmsss, p.val = p_val) - if (conf) { - res <- list(rms = rms, conf.lower = conf.lower, conf.upper = conf.upper) } else { - res <- list(rms = rms) - } + res <- list(rmsss = rmsss) + } return(res) - } diff --git a/R/UltimateBrier.R b/R/UltimateBrier.R index 18c8e5b..aeaddcd 100644 --- a/R/UltimateBrier.R +++ b/R/UltimateBrier.R @@ -102,15 +102,12 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti stop("Parameter 'exp' and 'obs' must have dimension names.") } ## dat_dim - if (!is.null(dat_dim)) { - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string or NULL.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") - } + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") } - ## memb_dim if (!is.character(memb_dim) | length(memb_dim) > 1) { stop("Parameter 'memb_dim' must be a character string.") @@ -140,11 +137,11 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti name_obs <- name_obs[-which(name_obs == dat_dim)] if (any(name_exp != name_obs)) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions except 'dat_dim' and 'memb_dim'.")) + "of all the dimensions expect 'dat_dim' and 'memb_dim'.")) } if (!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions except 'dat_dim' and 'memb_dim'.")) + "of all the dimensions expect 'dat_dim' and 'memb_dim'.")) } ## quantile if (!is.logical(quantile) | length(quantile) > 1) { @@ -186,7 +183,6 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti target_dims = list(c(time_dim, dat_dim, memb_dim), c(time_dim, dat_dim, memb_dim)), fun = .UltimateBrier, - dat_dim = dat_dim, memb_dim = memb_dim, thr = thr, type = type, decomposition = decomposition, ncores = ncores)$output1 @@ -207,7 +203,6 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti target_dims = list(c(time_dim, dat_dim), c(time_dim, dat_dim)), fun = .UltimateBrier, - dat_dim = dat_dim, memb_dim = memb_dim, thr = thr, type = type, decomposition = decomposition, ncores = ncores) @@ -222,8 +217,8 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti return(res) } -.UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', thr = c(5/100, 95/100), - type = 'BS', decomposition = TRUE) { +.UltimateBrier <- function(exp, obs, thr = c(5/100, 95/100), type = 'BS', + decomposition = TRUE) { # If exp and obs are probablistics # exp: [sdate, nexp] # obs: [sdate, nobs] @@ -234,10 +229,6 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti #NOTE: 'thr' is used in 'FairEnsembleBSS' and 'FairEnsembleBS'. But if quantile = F and # thr is real value, does it work? if (type == 'FairEnsembleBSS') { - if (is.null(dat_dim)) { - obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') - exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') - } size_ens_ref <- prod(dim(obs)[c(1, 3)]) res <- array(dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]), @@ -254,9 +245,6 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } } } - if (is.null(dat_dim)) { - dim(res) <- dim(res)[3:length(dim(res))] - } } else if (type == 'FairEnsembleBS') { #NOTE: The calculation in s2dverification::UltimateBrier is wrong. In the final stage, @@ -264,10 +252,6 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti # but the 3rd dim of result is 'bins' instead of decomposition. 'FairEnsembleBS' does # not have decomposition. # The calculation is fixed here. - if (is.null(dat_dim)) { - obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') - exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') - } res <- array(dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]), bin = length(thr) + 1)) @@ -280,17 +264,10 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } } } - if (is.null(dat_dim)) { - dim(res) <- dim(res)[3:length(dim(res))] - } # tmp <- res[, , 1] - res[, , 2] + res[, , 3] # res <- array(tmp, dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]))) } else if (type == 'BS') { - if (is.null(dat_dim)) { - obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') - exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') - } comp <- array(dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]), comp = 3)) @@ -303,28 +280,17 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } if (decomposition) { rel <- comp[, , 1] + dim(rel) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) res <- comp[, , 2] + dim(res) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) unc <- comp[, , 3] + dim(unc) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) bs <- rel - res + unc - if (is.null(dat_dim)) { - dim(rel) <- NULL - dim(res) <- NULL - dim(unc) <- NULL - dim(bs) <- NULL - } else { - dim(rel) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) - dim(res) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) - dim(unc) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) - dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) - } + dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) res <- list(bs = bs, rel = rel, res = res, unc = unc) } else { bs <- comp[, , 1] - comp[, , 2] + comp[, , 3] - if (is.null(dat_dim)) { - dim(bs) <- NULL - } else { - dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) - } + dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) res <- list(bs = bs) } diff --git a/tests/testthat/test-RMS.R b/tests/testthat/test-RMS.R index 65edfac..6e18bee 100644 --- a/tests/testthat/test-RMS.R +++ b/tests/testthat/test-RMS.R @@ -17,13 +17,6 @@ context("s2dv::RMS tests") set.seed(6) obs2 <- rnorm(10) - # dat3 - set.seed(1) - exp3 <- array(rnorm(120), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) - - set.seed(2) - obs3 <- array(rnorm(80), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) - ############################################## test_that("1. Input checks", { @@ -96,7 +89,7 @@ test_that("1. Input checks", { expect_error( RMS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension expect 'dat_dim'." ) expect_error( RMS(exp = array(1:5, dim = c(sdate = 1, dataset = 5, a = 1)), @@ -193,46 +186,3 @@ test_that("3. Output checks: dat2", { }) ############################################## -test_that("4. Output checks: dat3", { - - expect_equal( - dim(RMS(exp3, obs3, dat_dim = NULL)$rms), - c(ftime = 2, lon = 1, lat = 4) - ) - - suppressWarnings( - expect_equal( - RMS(exp3, obs3, dat_dim = NULL)$rms[1:6], - c(1.6458118, 0.8860392, 0.8261295, 1.1681939, 2.1693538, 1.3064454), - tolerance = 0.001 - ) -) -suppressWarnings( - expect_equal( - length(which(is.na(RMS(exp3, obs3, dat_dim = NULL)$conf.lower))), - 0 - ) -) -suppressWarnings( - expect_equal( - max(RMS(exp3, obs3, dat_dim = NULL)$conf.lower, na.rm = T), - 1.453144, - tolerance = 0.001 - ) -) -suppressWarnings( - expect_equal( - min(RMS(exp3, obs3, dat_dim = NULL, conf.lev = 0.99)$conf.upper, na.rm = TRUE), - 2.646274, - tolerance = 0.0001 - ) -) -suppressWarnings( - expect_equal( - length(RMS(exp3, obs3, dat_dim = NULL, conf = FALSE)), - 1 - ) -) -}) - -############################################## \ No newline at end of file diff --git a/tests/testthat/test-RMSSS.R b/tests/testthat/test-RMSSS.R index 60b4672..7d2a16d 100644 --- a/tests/testthat/test-RMSSS.R +++ b/tests/testthat/test-RMSSS.R @@ -25,12 +25,6 @@ context("s2dv::RMSSS tests") set.seed(6) obs3 <- rnorm(10) - # case 4 - set.seed(7) - exp4 <- array(rnorm(120), dim = c(sdate = 10, dat = 1, lon = 3, lat = 2)) - set.seed(8) - obs4 <- array(rnorm(60), dim = c(dat = 1, sdate = 10, lat = 2, lon = 3)) - ############################################## test_that("1. Input checks", { @@ -83,7 +77,7 @@ test_that("1. Input checks", { expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension expect 'dat_dim'." ) expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), @@ -163,26 +157,5 @@ test_that("4. Output checks: case 3", { tolerance = 0.00001 ) -############################################## -test_that("4. Output checks: case 4", { - - expect_equal( - dim(RMSSS(exp4, obs4, dat_dim = NULL)$rmsss), - c(dat = 1, lon = 3, lat = 2) - ) - expect_equal( - mean(RMSSS(exp4, obs4, dat_dim = NULL)$rmsss), - -0.3114424, - tolerance = 0.00001 - ) - expect_equal( - range(RMSSS(exp4, obs4, dat_dim = NULL)$p.val), - c(0.3560534, 0.9192801), - tolerance = 0.00001 - ) - -}) - -############################################## }) diff --git a/tests/testthat/test-UltimateBrier.R b/tests/testthat/test-UltimateBrier.R index e31f807..654aad0 100644 --- a/tests/testthat/test-UltimateBrier.R +++ b/tests/testthat/test-UltimateBrier.R @@ -7,12 +7,6 @@ exp1 <- array(rnorm(30), dim = c(dataset = 1, member = 3, sdate = 5, ftime = 2)) set.seed(2) obs1 <- array(round(rnorm(10)), dim = c(dataset = 1, sdate = 5, ftime = 2)) -# dat2 -set.seed(1) -exp2 <- array(rnorm(30), dim = c(member = 3, sdate = 5, ftime = 2)) -set.seed(2) -obs2 <- array(round(rnorm(10)), dim = c(sdate = 5, ftime = 2)) - ############################################## test_that("1. Input checks", { @@ -64,12 +58,12 @@ test_that("1. Input checks", { expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 6, ftime = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions except 'dat_dim' and 'memb_dim'.") + "of all the dimensions expect 'dat_dim' and 'memb_dim'.") ) expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 5, time = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions except 'dat_dim' and 'memb_dim'.") + "of all the dimensions expect 'dat_dim' and 'memb_dim'.") ) # quantile expect_error( @@ -244,139 +238,3 @@ tolerance = 0.0001 }) -############################################## -test_that("3. Output checks: dat2", { - -# 'BS' -expect_equal( -is.list(UltimateBrier(exp2, obs2, dat_dim = NULL)), -TRUE -) -expect_equal( -names(UltimateBrier(exp2, obs2, dat_dim = NULL)), -c('bs', 'rel', 'res', 'unc') -) -expect_equal( -is.list(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE)), -FALSE -) -expect_equal( -dim(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE)), -c(bin = 3, ftime = 2) -) -expect_equal( -dim(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, thr = c(0.25, 0.5, 0.75))), -c(bin = 4, ftime = 2) -) -expect_equal( -UltimateBrier(exp2, obs2, dat_dim = NULL)$bs, -UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE) -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$bs), -c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$rel), -c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$res), -c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$unc), -c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), -tolerance = 0.0001 -) - -# 'BSS' -expect_equal( -dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'BSS')), -c(bin = 3, ftime = 2) -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'BSS')), -c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), -tolerance = 0.0001 -) - -# 'FairStartDatesBS' -expect_equal( -is.list(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')), -TRUE -) -expect_equal( -names(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')), -c('bs', 'rel', 'res', 'unc') -) -expect_equal( -is.list(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, type = 'FairStartDatesBS')), -FALSE -) -expect_equal( -dim(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, type = 'FairStartDatesBS')), -c(bin = 3, ftime = 2) -) -expect_equal( -UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$bs, -UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, type = 'FairStartDatesBS') -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$bs), -c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$rel), -c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$res), -c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$unc), -c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), -tolerance = 0.0001 -) - -# 'FairStartDatesBSS' -expect_equal( -dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBSS')), -c(bin = 3, ftime = 2) -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBSS')), -c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), -tolerance = 0.0001 -) -# 'FairEnsembleBS' -expect_equal( -dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBS')), -c(bin = 3, ftime = 2) -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBS')), -c(0.1333333, 0.2000000, 0.2000000, 0.1333333, 0.4000000, 0.2000000), -tolerance = 0.0001 -) -# 'FairEnsembleBSS' -expect_equal( -dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBSS')), -c(bin = 3, ftime = 2) -) -expect_equal( -as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBSS')), -c(-0.1111111, -0.6666667, -0.6666667, 0.2592593, -1.2222222, -0.6666667), -tolerance = 0.0001 -) - -}) - - -- GitLab From 85c975bd58c29fcf45ff7367b7ff03d97a146e87 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 14:47:15 +0200 Subject: [PATCH 29/47] Allow dat_dim NULL in RMS() --- R/RMS.R | 56 +++++++++++++++++++++++++++++++++++++++++--------------- 1 file changed, 41 insertions(+), 15 deletions(-) diff --git a/R/RMS.R b/R/RMS.R index b3c8ad4..7a38a72 100644 --- a/R/RMS.R +++ b/R/RMS.R @@ -21,7 +21,8 @@ #'@param time_dim A character string indicating the name of dimension along #' which the correlations are computed. The default value is 'sdate'. #'@param dat_dim A character string indicating the name of member (nobs/nexp) -#' dimension. The default value is 'dataset'. +#' dimension. The default value is 'dataset'. If there is no dataset +#' dimension, set NULL. #'@param 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. @@ -107,11 +108,14 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + "Set it as NULL if there is no dataset dimension.") + } } ## comp_dim if (!is.null(comp_dim)) { @@ -151,11 +155,13 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] + if (!is.null(dat_dim)) { + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] + } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension expect 'dat_dim'.")) + "all dimension except 'dat_dim'.")) } if (dim(exp)[time_dim] < 2) { stop("The length of time_dim must be at least 2 to compute RMS.") @@ -196,12 +202,22 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } .RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { - - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) + conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { + if (is.null(dat_dim)) { + # exp: [sdate] + # obs: [sdate] + nexp <- 1 + nobs <- 1 + ini_dims <- dim(exp) + dim(exp) <- c(ini_dims, dat_dim = 1) + dim(obs) <- c(ini_dims, dat_dim = 1) + } else { + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) + } + nsdate <- as.numeric(dim(exp)[1]) dif <- array(dim = c(sdate = nsdate, nexp = nexp, nobs = nobs)) @@ -236,6 +252,16 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', conf.upper <- (eno * rms ** 2 / chi) ** 0.5 } + ################################### + # Remove nexp and nobs if dat_dim = NULL + if (is.null(dat_dim)) { + dim(rms) <- NULL + dim(conf.lower) <- NULL + dim(conf.upper) <- NULL + } + + ################################### + if (conf) { res <- list(rms = rms, conf.lower = conf.lower, conf.upper = conf.upper) } else { -- GitLab From 4e39b3383a715464cc2e680fdd5bb6adec5587c5 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 15:02:02 +0200 Subject: [PATCH 30/47] revert commit --- R/RMS.R | 56 +++++++++++++++----------------------------------------- 1 file changed, 15 insertions(+), 41 deletions(-) diff --git a/R/RMS.R b/R/RMS.R index 7a38a72..b3c8ad4 100644 --- a/R/RMS.R +++ b/R/RMS.R @@ -21,8 +21,7 @@ #'@param time_dim A character string indicating the name of dimension along #' which the correlations are computed. The default value is 'sdate'. #'@param dat_dim A character string indicating the name of member (nobs/nexp) -#' dimension. The default value is 'dataset'. If there is no dataset -#' dimension, set NULL. +#' dimension. The default value is 'dataset'. #'@param 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. @@ -108,14 +107,11 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.null(dat_dim)) { - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string or NULL.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", - "Set it as NULL if there is no dataset dimension.") - } + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") } ## comp_dim if (!is.null(comp_dim)) { @@ -155,13 +151,11 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - if (!is.null(dat_dim)) { - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] - } + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension except 'dat_dim'.")) + "all dimension expect 'dat_dim'.")) } if (dim(exp)[time_dim] < 2) { stop("The length of time_dim must be at least 2 to compute RMS.") @@ -202,22 +196,12 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } .RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { - if (is.null(dat_dim)) { - # exp: [sdate] - # obs: [sdate] - nexp <- 1 - nobs <- 1 - ini_dims <- dim(exp) - dim(exp) <- c(ini_dims, dat_dim = 1) - dim(obs) <- c(ini_dims, dat_dim = 1) - } else { - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) - } - + conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { + + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) nsdate <- as.numeric(dim(exp)[1]) dif <- array(dim = c(sdate = nsdate, nexp = nexp, nobs = nobs)) @@ -252,16 +236,6 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', conf.upper <- (eno * rms ** 2 / chi) ** 0.5 } - ################################### - # Remove nexp and nobs if dat_dim = NULL - if (is.null(dat_dim)) { - dim(rms) <- NULL - dim(conf.lower) <- NULL - dim(conf.upper) <- NULL - } - - ################################### - if (conf) { res <- list(rms = rms, conf.lower = conf.lower, conf.upper = conf.upper) } else { -- GitLab From c3c1b92524356c44f4e6f47fbce12d1e823f6388 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 15:14:46 +0200 Subject: [PATCH 31/47] allow dat_dim NULL in RMS() --- R/RMS.R | 53 +++++++++++++++++++++++++++++++++++++++-------------- 1 file changed, 39 insertions(+), 14 deletions(-) diff --git a/R/RMS.R b/R/RMS.R index b3c8ad4..cc013e6 100644 --- a/R/RMS.R +++ b/R/RMS.R @@ -107,11 +107,14 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension. ", + "Set it as NULL if there is no dataset dimension.") + } } ## comp_dim if (!is.null(comp_dim)) { @@ -151,11 +154,13 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] + if (!is.null(dat_dim)) { + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] + } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension expect 'dat_dim'.")) + "all dimension except 'dat_dim'.")) } if (dim(exp)[time_dim] < 2) { stop("The length of time_dim must be at least 2 to compute RMS.") @@ -196,12 +201,22 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } .RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { - - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) + conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { + if (is.null(dat_dim)) { + # exp: [sdate] + # obs: [sdate] + nexp <- 1 + nobs <- 1 + ini_dims <- dim(exp) + dim(exp) <- c(ini_dims, dat_dim = 1) + dim(obs) <- c(ini_dims, dat_dim = 1) + } else { + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) + } + nsdate <- as.numeric(dim(exp)[1]) dif <- array(dim = c(sdate = nsdate, nexp = nexp, nobs = nobs)) @@ -236,6 +251,16 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', conf.upper <- (eno * rms ** 2 / chi) ** 0.5 } + ################################### + # Remove nexp and nobs if dat_dim = NULL + if (is.null(dat_dim)) { + dim(rms) <- NULL + dim(conf.lower) <- NULL + dim(conf.upper) <- NULL + } + + ################################### + if (conf) { res <- list(rms = rms, conf.lower = conf.lower, conf.upper = conf.upper) } else { -- GitLab From 138e3be98fb8aef41e6ebf5e06799f64c45e3b24 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 15:31:05 +0200 Subject: [PATCH 32/47] Corrected typo in test file --- tests/testthat/test-RMS.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/testthat/test-RMS.R b/tests/testthat/test-RMS.R index 6e18bee..f1083dd 100644 --- a/tests/testthat/test-RMS.R +++ b/tests/testthat/test-RMS.R @@ -89,7 +89,7 @@ test_that("1. Input checks", { expect_error( RMS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension expect 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." ) expect_error( RMS(exp = array(1:5, dim = c(sdate = 1, dataset = 5, a = 1)), -- GitLab From c3f9d73c63be6a8e6ee73eff031f28973223d02b Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 15:53:12 +0200 Subject: [PATCH 33/47] Revert commit --- R/RMS.R | 53 ++++++++++++++--------------------------------------- 1 file changed, 14 insertions(+), 39 deletions(-) diff --git a/R/RMS.R b/R/RMS.R index cc013e6..b3c8ad4 100644 --- a/R/RMS.R +++ b/R/RMS.R @@ -107,14 +107,11 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.null(dat_dim)) { - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string or NULL.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension. ", - "Set it as NULL if there is no dataset dimension.") - } + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") } ## comp_dim if (!is.null(comp_dim)) { @@ -154,13 +151,11 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - if (!is.null(dat_dim)) { - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] - } + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension except 'dat_dim'.")) + "all dimension expect 'dat_dim'.")) } if (dim(exp)[time_dim] < 2) { stop("The length of time_dim must be at least 2 to compute RMS.") @@ -201,22 +196,12 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } .RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { - if (is.null(dat_dim)) { - # exp: [sdate] - # obs: [sdate] - nexp <- 1 - nobs <- 1 - ini_dims <- dim(exp) - dim(exp) <- c(ini_dims, dat_dim = 1) - dim(obs) <- c(ini_dims, dat_dim = 1) - } else { - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) - } - + conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { + + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) nsdate <- as.numeric(dim(exp)[1]) dif <- array(dim = c(sdate = nsdate, nexp = nexp, nobs = nobs)) @@ -251,16 +236,6 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', conf.upper <- (eno * rms ** 2 / chi) ** 0.5 } - ################################### - # Remove nexp and nobs if dat_dim = NULL - if (is.null(dat_dim)) { - dim(rms) <- NULL - dim(conf.lower) <- NULL - dim(conf.upper) <- NULL - } - - ################################### - if (conf) { res <- list(rms = rms, conf.lower = conf.lower, conf.upper = conf.upper) } else { -- GitLab From beffc1fcbc4a2774ae6758257cb4a5e4f88a381f Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 16:06:03 +0200 Subject: [PATCH 34/47] Allow dat_dim NULL in RMS --- R/RMS.R | 53 +++++++++++++++++++++++++++++++++++++++-------------- 1 file changed, 39 insertions(+), 14 deletions(-) diff --git a/R/RMS.R b/R/RMS.R index b3c8ad4..b3188fc 100644 --- a/R/RMS.R +++ b/R/RMS.R @@ -107,11 +107,14 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + "Set it as NULL if there is no dataset dimension.") + } } ## comp_dim if (!is.null(comp_dim)) { @@ -151,11 +154,13 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - name_exp <- name_exp[-which(name_exp == dat_dim)] - name_obs <- name_obs[-which(name_obs == dat_dim)] + if (!is.null(dat_dim)) { + name_exp <- name_exp[-which(name_exp == dat_dim)] + name_obs <- name_obs[-which(name_obs == dat_dim)] + } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension expect 'dat_dim'.")) + "all dimension except 'dat_dim'.")) } if (dim(exp)[time_dim] < 2) { stop("The length of time_dim must be at least 2 to compute RMS.") @@ -196,12 +201,22 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', } .RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', - conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { - - # exp: [sdate, dat_exp] - # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) + conf = TRUE, conf.lev = 0.95, ncores_input = NULL) { + if (is.null(dat_dim)) { + # exp: [sdate] + # obs: [sdate] + nexp <- 1 + nobs <- 1 + ini_dims <- dim(exp) + dim(exp) <- c(ini_dims, dat_dim = 1) + dim(obs) <- c(ini_dims, dat_dim = 1) + } else { + # exp: [sdate, dat_exp] + # obs: [sdate, dat_obs] + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) + } + nsdate <- as.numeric(dim(exp)[1]) dif <- array(dim = c(sdate = nsdate, nexp = nexp, nobs = nobs)) @@ -236,6 +251,16 @@ RMS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', conf.upper <- (eno * rms ** 2 / chi) ** 0.5 } + ################################### + # Remove nexp and nobs if dat_dim = NULL + if (is.null(dat_dim)) { + dim(rms) <- NULL + dim(conf.lower) <- NULL + dim(conf.upper) <- NULL + } + + ################################### + if (conf) { res <- list(rms = rms, conf.lower = conf.lower, conf.upper = conf.upper) } else { -- GitLab From 7f870ec234b4ca6271365a891b50c5cb4a6afa9b Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 16:18:02 +0200 Subject: [PATCH 35/47] Added tests for RMS and corrected typo in tests ultimatebrier and rmsss --- tests/testthat/test-RMS.R | 50 +++++++++++++++++++++++++++++ tests/testthat/test-RMSSS.R | 2 +- tests/testthat/test-UltimateBrier.R | 4 +-- 3 files changed, 53 insertions(+), 3 deletions(-) diff --git a/tests/testthat/test-RMS.R b/tests/testthat/test-RMS.R index f1083dd..65edfac 100644 --- a/tests/testthat/test-RMS.R +++ b/tests/testthat/test-RMS.R @@ -17,6 +17,13 @@ context("s2dv::RMS tests") set.seed(6) obs2 <- rnorm(10) + # dat3 + set.seed(1) + exp3 <- array(rnorm(120), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) + + set.seed(2) + obs3 <- array(rnorm(80), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) + ############################################## test_that("1. Input checks", { @@ -186,3 +193,46 @@ test_that("3. Output checks: dat2", { }) ############################################## +test_that("4. Output checks: dat3", { + + expect_equal( + dim(RMS(exp3, obs3, dat_dim = NULL)$rms), + c(ftime = 2, lon = 1, lat = 4) + ) + + suppressWarnings( + expect_equal( + RMS(exp3, obs3, dat_dim = NULL)$rms[1:6], + c(1.6458118, 0.8860392, 0.8261295, 1.1681939, 2.1693538, 1.3064454), + tolerance = 0.001 + ) +) +suppressWarnings( + expect_equal( + length(which(is.na(RMS(exp3, obs3, dat_dim = NULL)$conf.lower))), + 0 + ) +) +suppressWarnings( + expect_equal( + max(RMS(exp3, obs3, dat_dim = NULL)$conf.lower, na.rm = T), + 1.453144, + tolerance = 0.001 + ) +) +suppressWarnings( + expect_equal( + min(RMS(exp3, obs3, dat_dim = NULL, conf.lev = 0.99)$conf.upper, na.rm = TRUE), + 2.646274, + tolerance = 0.0001 + ) +) +suppressWarnings( + expect_equal( + length(RMS(exp3, obs3, dat_dim = NULL, conf = FALSE)), + 1 + ) +) +}) + +############################################## \ No newline at end of file diff --git a/tests/testthat/test-RMSSS.R b/tests/testthat/test-RMSSS.R index 7d2a16d..91ef930 100644 --- a/tests/testthat/test-RMSSS.R +++ b/tests/testthat/test-RMSSS.R @@ -77,7 +77,7 @@ test_that("1. Input checks", { expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension expect 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." ) expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), diff --git a/tests/testthat/test-UltimateBrier.R b/tests/testthat/test-UltimateBrier.R index 654aad0..0eb0118 100644 --- a/tests/testthat/test-UltimateBrier.R +++ b/tests/testthat/test-UltimateBrier.R @@ -58,12 +58,12 @@ test_that("1. Input checks", { expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 6, ftime = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.") + "of all the dimensions except 'dat_dim' and 'memb_dim'.") ) expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 5, time = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.") + "of all the dimensions except 'dat_dim' and 'memb_dim'.") ) # quantile expect_error( -- GitLab From 71a0ced64c0391eb2d4ea716fd3876eed6dd2bc0 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 16:31:04 +0200 Subject: [PATCH 36/47] Revert commit --- tests/testthat/test-RMSSS.R | 2 +- tests/testthat/test-UltimateBrier.R | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/tests/testthat/test-RMSSS.R b/tests/testthat/test-RMSSS.R index 91ef930..7d2a16d 100644 --- a/tests/testthat/test-RMSSS.R +++ b/tests/testthat/test-RMSSS.R @@ -77,7 +77,7 @@ test_that("1. Input checks", { expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension expect 'dat_dim'." ) expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), diff --git a/tests/testthat/test-UltimateBrier.R b/tests/testthat/test-UltimateBrier.R index 0eb0118..654aad0 100644 --- a/tests/testthat/test-UltimateBrier.R +++ b/tests/testthat/test-UltimateBrier.R @@ -58,12 +58,12 @@ test_that("1. Input checks", { expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 6, ftime = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions except 'dat_dim' and 'memb_dim'.") + "of all the dimensions expect 'dat_dim' and 'memb_dim'.") ) expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 5, time = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions except 'dat_dim' and 'memb_dim'.") + "of all the dimensions expect 'dat_dim' and 'memb_dim'.") ) # quantile expect_error( -- GitLab From 92a1b16040e7a870b9b73c9fae031cf4498029c6 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 16:47:51 +0200 Subject: [PATCH 37/47] Correct typo --- tests/testthat/test-RMS.R | 50 --------------------------------------- 1 file changed, 50 deletions(-) diff --git a/tests/testthat/test-RMS.R b/tests/testthat/test-RMS.R index 65edfac..f1083dd 100644 --- a/tests/testthat/test-RMS.R +++ b/tests/testthat/test-RMS.R @@ -17,13 +17,6 @@ context("s2dv::RMS tests") set.seed(6) obs2 <- rnorm(10) - # dat3 - set.seed(1) - exp3 <- array(rnorm(120), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) - - set.seed(2) - obs3 <- array(rnorm(80), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) - ############################################## test_that("1. Input checks", { @@ -193,46 +186,3 @@ test_that("3. Output checks: dat2", { }) ############################################## -test_that("4. Output checks: dat3", { - - expect_equal( - dim(RMS(exp3, obs3, dat_dim = NULL)$rms), - c(ftime = 2, lon = 1, lat = 4) - ) - - suppressWarnings( - expect_equal( - RMS(exp3, obs3, dat_dim = NULL)$rms[1:6], - c(1.6458118, 0.8860392, 0.8261295, 1.1681939, 2.1693538, 1.3064454), - tolerance = 0.001 - ) -) -suppressWarnings( - expect_equal( - length(which(is.na(RMS(exp3, obs3, dat_dim = NULL)$conf.lower))), - 0 - ) -) -suppressWarnings( - expect_equal( - max(RMS(exp3, obs3, dat_dim = NULL)$conf.lower, na.rm = T), - 1.453144, - tolerance = 0.001 - ) -) -suppressWarnings( - expect_equal( - min(RMS(exp3, obs3, dat_dim = NULL, conf.lev = 0.99)$conf.upper, na.rm = TRUE), - 2.646274, - tolerance = 0.0001 - ) -) -suppressWarnings( - expect_equal( - length(RMS(exp3, obs3, dat_dim = NULL, conf = FALSE)), - 1 - ) -) -}) - -############################################## \ No newline at end of file -- GitLab From fe9572df91fb98dbe98d58de82e070b607512275 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 17:03:15 +0200 Subject: [PATCH 38/47] Allow dat_dim NULL in RMSSS and changed test --- R/RMSSS.R | 39 +++++++++++++++++++++++++++++-------- tests/testthat/test-RMSSS.R | 5 +++-- 2 files changed, 34 insertions(+), 10 deletions(-) diff --git a/R/RMSSS.R b/R/RMSSS.R index 5fa9659..e44105c 100644 --- a/R/RMSSS.R +++ b/R/RMSSS.R @@ -96,11 +96,14 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + " Set it as NULL if there is no dataset dimension.") + } } ## pval if (!is.logical(pval) | length(pval) > 1) { @@ -116,11 +119,13 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) + if (!is.null(dat_dim)) { name_exp <- name_exp[-which(name_exp == dat_dim)] name_obs <- name_obs[-which(name_obs == dat_dim)] + } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", - "all dimension expect 'dat_dim'.")) + "all dimension except 'dat_dim'.")) } if (dim(exp)[time_dim] <= 2) { stop("The length of time_dim must be more than 2 to compute RMSSS.") @@ -151,10 +156,20 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', .RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', pval = TRUE, ncores_input = NULL) { + if (is.null(dat_dim)) { + # exp: [sdate] + # obs: [sdate] + nexp <- 1 + nobs <- 1 + dim(exp) <- c(dim(exp), nexp = nexp) + dim(obs) <- c(dim(obs), nobs = nobs) + } else { # exp: [sdate, dat_exp] # obs: [sdate, dat_obs] - nexp <- as.numeric(dim(exp)[2]) - nobs <- as.numeric(dim(obs)[2]) + nexp <- as.numeric(dim(exp)[2]) + nobs <- as.numeric(dim(obs)[2]) + } + nsdate <- as.numeric(dim(exp)[1]) p_val <- array(dim = c(nexp = nexp, nobs = nobs)) @@ -208,6 +223,14 @@ RMSSS <- function(exp, obs, time_dim = 'sdate', dat_dim = 'dataset', rmsss[which(!tmp)] <- NA } + ################################### + # Remove extra dimensions if dat_dim = NULL + if (is.null(dat_dim)) { + dim(rmsss) <- NULL + dim(p_val) <- NULL + } + ################################### + # output if (pval) { res <- list(rmsss = rmsss, p.val = p_val) diff --git a/tests/testthat/test-RMSSS.R b/tests/testthat/test-RMSSS.R index 7d2a16d..6064508 100644 --- a/tests/testthat/test-RMSSS.R +++ b/tests/testthat/test-RMSSS.R @@ -64,7 +64,8 @@ test_that("1. Input checks", { ) expect_error( RMSSS(exp0, obs0, dat_dim = 'memb'), - "Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension." + paste0("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + " Set it as NULL if there is no dataset dimension.") ) expect_error( RMSSS(exp0, obs0, pval = c(T, T)), @@ -77,7 +78,7 @@ test_that("1. Input checks", { expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), obs = array(1:4, dim = c(a = 1, sdate = 2, dataset = 2))), - "Parameter 'exp' and 'obs' must have same length of all dimension expect 'dat_dim'." + "Parameter 'exp' and 'obs' must have same length of all dimension except 'dat_dim'." ) expect_error( RMSSS(exp = array(1:10, dim = c(sdate = 1, dataset = 5, a = 1)), -- GitLab From f66380ad20c623d51c011e5bdd835046346c89cf Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 17:21:54 +0200 Subject: [PATCH 39/47] Allow dat_dim NULL and add test --- R/UltimateBrier.R | 62 ++++- tests/testthat/test-UltimateBrier.R | 403 +++++++++++++++++++--------- 2 files changed, 318 insertions(+), 147 deletions(-) diff --git a/R/UltimateBrier.R b/R/UltimateBrier.R index aeaddcd..18c8e5b 100644 --- a/R/UltimateBrier.R +++ b/R/UltimateBrier.R @@ -102,12 +102,15 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti stop("Parameter 'exp' and 'obs' must have dimension names.") } ## dat_dim - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim) | length(dat_dim) > 1) { + stop("Parameter 'dat_dim' must be a character string or NULL.") + } + if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + } } + ## memb_dim if (!is.character(memb_dim) | length(memb_dim) > 1) { stop("Parameter 'memb_dim' must be a character string.") @@ -137,11 +140,11 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti name_obs <- name_obs[-which(name_obs == dat_dim)] if (any(name_exp != name_obs)) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.")) + "of all the dimensions except 'dat_dim' and 'memb_dim'.")) } if (!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.")) + "of all the dimensions except 'dat_dim' and 'memb_dim'.")) } ## quantile if (!is.logical(quantile) | length(quantile) > 1) { @@ -183,6 +186,7 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti target_dims = list(c(time_dim, dat_dim, memb_dim), c(time_dim, dat_dim, memb_dim)), fun = .UltimateBrier, + dat_dim = dat_dim, memb_dim = memb_dim, thr = thr, type = type, decomposition = decomposition, ncores = ncores)$output1 @@ -203,6 +207,7 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti target_dims = list(c(time_dim, dat_dim), c(time_dim, dat_dim)), fun = .UltimateBrier, + dat_dim = dat_dim, memb_dim = memb_dim, thr = thr, type = type, decomposition = decomposition, ncores = ncores) @@ -217,8 +222,8 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti return(res) } -.UltimateBrier <- function(exp, obs, thr = c(5/100, 95/100), type = 'BS', - decomposition = TRUE) { +.UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', thr = c(5/100, 95/100), + type = 'BS', decomposition = TRUE) { # If exp and obs are probablistics # exp: [sdate, nexp] # obs: [sdate, nobs] @@ -229,6 +234,10 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti #NOTE: 'thr' is used in 'FairEnsembleBSS' and 'FairEnsembleBS'. But if quantile = F and # thr is real value, does it work? if (type == 'FairEnsembleBSS') { + if (is.null(dat_dim)) { + obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') + exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') + } size_ens_ref <- prod(dim(obs)[c(1, 3)]) res <- array(dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]), @@ -245,6 +254,9 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } } } + if (is.null(dat_dim)) { + dim(res) <- dim(res)[3:length(dim(res))] + } } else if (type == 'FairEnsembleBS') { #NOTE: The calculation in s2dverification::UltimateBrier is wrong. In the final stage, @@ -252,6 +264,10 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti # but the 3rd dim of result is 'bins' instead of decomposition. 'FairEnsembleBS' does # not have decomposition. # The calculation is fixed here. + if (is.null(dat_dim)) { + obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') + exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') + } res <- array(dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]), bin = length(thr) + 1)) @@ -264,10 +280,17 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } } } + if (is.null(dat_dim)) { + dim(res) <- dim(res)[3:length(dim(res))] + } # tmp <- res[, , 1] - res[, , 2] + res[, , 3] # res <- array(tmp, dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]))) } else if (type == 'BS') { + if (is.null(dat_dim)) { + obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') + exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') + } comp <- array(dim = c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2]), comp = 3)) @@ -280,17 +303,28 @@ UltimateBrier <- function(exp, obs, dat_dim = 'dataset', memb_dim = 'member', ti } if (decomposition) { rel <- comp[, , 1] - dim(rel) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) res <- comp[, , 2] - dim(res) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) unc <- comp[, , 3] - dim(unc) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) bs <- rel - res + unc - dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + if (is.null(dat_dim)) { + dim(rel) <- NULL + dim(res) <- NULL + dim(unc) <- NULL + dim(bs) <- NULL + } else { + dim(rel) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + dim(res) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + dim(unc) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + } res <- list(bs = bs, rel = rel, res = res, unc = unc) } else { bs <- comp[, , 1] - comp[, , 2] + comp[, , 3] - dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + if (is.null(dat_dim)) { + dim(bs) <- NULL + } else { + dim(bs) <- c(nexp = as.numeric(dim(exp)[2]), nobs = as.numeric(dim(obs)[2])) + } res <- list(bs = bs) } diff --git a/tests/testthat/test-UltimateBrier.R b/tests/testthat/test-UltimateBrier.R index 654aad0..09412f0 100644 --- a/tests/testthat/test-UltimateBrier.R +++ b/tests/testthat/test-UltimateBrier.R @@ -7,7 +7,11 @@ exp1 <- array(rnorm(30), dim = c(dataset = 1, member = 3, sdate = 5, ftime = 2)) set.seed(2) obs1 <- array(round(rnorm(10)), dim = c(dataset = 1, sdate = 5, ftime = 2)) - +# dat2 +set.seed(1) +exp2 <- array(rnorm(30), dim = c(member = 3, sdate = 5, ftime = 2)) +set.seed(2) +obs2 <- array(round(rnorm(10)), dim = c(sdate = 5, ftime = 2)) ############################################## test_that("1. Input checks", { # exp and obs @@ -26,7 +30,7 @@ test_that("1. Input checks", { # dat_dim expect_error( UltimateBrier(exp1, obs1, dat_dim = 2), - "Parameter 'dat_dim' must be a character string." + "Parameter 'dat_dim' must be a character string or NULL." ) expect_error( UltimateBrier(exp1, obs1, dat_dim = 'dat'), @@ -58,12 +62,12 @@ test_that("1. Input checks", { expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 6, ftime = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.") + "of all the dimensions except 'dat_dim' and 'memb_dim'.") ) expect_error( UltimateBrier(exp1, array(1:10, dim = c(dataset = 1, sdate = 5, time = 2))), paste0("Parameter 'exp' and 'obs' must have the same names and lengths ", - "of all the dimensions expect 'dat_dim' and 'memb_dim'.") + "of all the dimensions except 'dat_dim' and 'memb_dim'.") ) # quantile expect_error( @@ -105,136 +109,269 @@ test_that("1. Input checks", { ############################################## test_that("2. Output checks: dat1", { -# 'BS' -expect_equal( -is.list(UltimateBrier(exp1, obs1)), -TRUE -) -expect_equal( -names(UltimateBrier(exp1, obs1)), -c('bs', 'rel', 'res', 'unc') -) -expect_equal( -is.list(UltimateBrier(exp1, obs1, decomposition = FALSE)), -FALSE -) -expect_equal( -dim(UltimateBrier(exp1, obs1, decomposition = FALSE)), -c(nexp = 1, nobs = 1, bin = 3, ftime = 2) -) -expect_equal( -dim(UltimateBrier(exp1, obs1, decomposition = FALSE, thr = c(0.25, 0.5, 0.75))), -c(nexp = 1, nobs = 1, bin = 4, ftime = 2) -) -expect_equal( -UltimateBrier(exp1, obs1)$bs, -UltimateBrier(exp1, obs1, decomposition = FALSE) -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1)$bs), -c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1)$rel), -c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1)$res), -c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1)$unc), -c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), -tolerance = 0.0001 -) - -# 'BSS' -expect_equal( -dim(UltimateBrier(exp1, obs1, type = 'BSS')), -c(nexp = 1, nobs = 1, bin = 3, ftime = 2) -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1, type = 'BSS')), -c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), -tolerance = 0.0001 -) - -# 'FairStartDatesBS' -expect_equal( -is.list(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')), -TRUE -) -expect_equal( -names(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')), -c('bs', 'rel', 'res', 'unc') -) -expect_equal( -is.list(UltimateBrier(exp1, obs1, decomposition = FALSE, type = 'FairStartDatesBS')), -FALSE -) -expect_equal( -dim(UltimateBrier(exp1, obs1, decomposition = FALSE, type = 'FairStartDatesBS')), -c(nexp = 1, nobs = 1, bin = 3, ftime = 2) -) -expect_equal( -UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$bs, -UltimateBrier(exp1, obs1, decomposition = FALSE, type = 'FairStartDatesBS') -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$bs), -c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$rel), -c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$res), -c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), -tolerance = 0.0001 -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$unc), -c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), -tolerance = 0.0001 -) - -# 'FairStartDatesBSS' -expect_equal( -dim(UltimateBrier(exp1, obs1, type = 'FairStartDatesBSS')), -c(nexp = 1, nobs = 1, bin = 3, ftime = 2) -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBSS')), -c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), -tolerance = 0.0001 -) -# 'FairEnsembleBS' -expect_equal( -dim(UltimateBrier(exp1, obs1, type = 'FairEnsembleBS')), -c(nexp = 1, nobs = 1, bin = 3, ftime = 2) -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1, type = 'FairEnsembleBS')), -c(0.1333333, 0.2000000, 0.2000000, 0.1333333, 0.4000000, 0.2000000), -tolerance = 0.0001 -) -# 'FairEnsembleBSS' -expect_equal( -dim(UltimateBrier(exp1, obs1, type = 'FairEnsembleBSS')), -c(nexp = 1, nobs = 1, bin = 3, ftime = 2) -) -expect_equal( -as.vector(UltimateBrier(exp1, obs1, type = 'FairEnsembleBSS')), -c(-0.1111111, -0.6666667, -0.6666667, 0.2592593, -1.2222222, -0.6666667), -tolerance = 0.0001 -) + # 'BS' + expect_equal( + is.list(UltimateBrier(exp1, obs1)), + TRUE + ) + expect_equal( + names(UltimateBrier(exp1, obs1)), + c('bs', 'rel', 'res', 'unc') + ) + expect_equal( + is.list(UltimateBrier(exp1, obs1, decomposition = FALSE)), + FALSE + ) + expect_equal( + dim(UltimateBrier(exp1, obs1, decomposition = FALSE)), + c(nexp = 1, nobs = 1, bin = 3, ftime = 2) + ) + expect_equal( + dim(UltimateBrier(exp1, obs1, decomposition = FALSE, thr = c(0.25, 0.5, 0.75))), + c(nexp = 1, nobs = 1, bin = 4, ftime = 2) + ) + expect_equal( + UltimateBrier(exp1, obs1)$bs, + UltimateBrier(exp1, obs1, decomposition = FALSE) + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1)$bs), + c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1)$rel), + c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1)$res), + c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1)$unc), + c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), + tolerance = 0.0001 + ) + + # 'BSS' + expect_equal( + dim(UltimateBrier(exp1, obs1, type = 'BSS')), + c(nexp = 1, nobs = 1, bin = 3, ftime = 2) + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1, type = 'BSS')), + c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), + tolerance = 0.0001 + ) + + # 'FairStartDatesBS' + expect_equal( + is.list(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')), + TRUE + ) + expect_equal( + names(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')), + c('bs', 'rel', 'res', 'unc') + ) + expect_equal( + is.list(UltimateBrier(exp1, obs1, decomposition = FALSE, type = 'FairStartDatesBS')), + FALSE + ) + expect_equal( + dim(UltimateBrier(exp1, obs1, decomposition = FALSE, type = 'FairStartDatesBS')), + c(nexp = 1, nobs = 1, bin = 3, ftime = 2) + ) + expect_equal( + UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$bs, + UltimateBrier(exp1, obs1, decomposition = FALSE, type = 'FairStartDatesBS') + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$bs), + c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$rel), + c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$res), + c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBS')$unc), + c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), + tolerance = 0.0001 + ) + + # 'FairStartDatesBSS' + expect_equal( + dim(UltimateBrier(exp1, obs1, type = 'FairStartDatesBSS')), + c(nexp = 1, nobs = 1, bin = 3, ftime = 2) + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1, type = 'FairStartDatesBSS')), + c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), + tolerance = 0.0001 + ) + # 'FairEnsembleBS' + expect_equal( + dim(UltimateBrier(exp1, obs1, type = 'FairEnsembleBS')), + c(nexp = 1, nobs = 1, bin = 3, ftime = 2) + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1, type = 'FairEnsembleBS')), + c(0.1333333, 0.2000000, 0.2000000, 0.1333333, 0.4000000, 0.2000000), + tolerance = 0.0001 + ) + # 'FairEnsembleBSS' + expect_equal( + dim(UltimateBrier(exp1, obs1, type = 'FairEnsembleBSS')), + c(nexp = 1, nobs = 1, bin = 3, ftime = 2) + ) + expect_equal( + as.vector(UltimateBrier(exp1, obs1, type = 'FairEnsembleBSS')), + c(-0.1111111, -0.6666667, -0.6666667, 0.2592593, -1.2222222, -0.6666667), + tolerance = 0.0001 + ) }) +############################################## +test_that("3. Output checks: dat2", { + # 'BS' + expect_equal( + is.list(UltimateBrier(exp2, obs2, dat_dim = NULL)), + TRUE + ) + expect_equal( + names(UltimateBrier(exp2, obs2, dat_dim = NULL)), + c('bs', 'rel', 'res', 'unc') + ) + expect_equal( + is.list(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE)), + FALSE + ) + expect_equal( + dim(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE)), + c(bin = 3, ftime = 2) + ) + expect_equal( + dim(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, thr = c(0.25, 0.5, 0.75))), + c(bin = 4, ftime = 2) + ) + expect_equal( + UltimateBrier(exp2, obs2, dat_dim = NULL)$bs, + UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE) + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$bs), + c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$rel), + c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$res), + c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL)$unc), + c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), + tolerance = 0.0001 + ) + + # 'BSS' + expect_equal( + dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'BSS')), + c(bin = 3, ftime = 2) + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'BSS')), + c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), + tolerance = 0.0001 + ) + + # 'FairStartDatesBS' + expect_equal( + is.list(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')), + TRUE + ) + expect_equal( + names(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')), + c('bs', 'rel', 'res', 'unc') + ) + expect_equal( + is.list(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, type = 'FairStartDatesBS')), + FALSE + ) + expect_equal( + dim(UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, type = 'FairStartDatesBS')), + c(bin = 3, ftime = 2) + ) + expect_equal( + UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$bs, + UltimateBrier(exp2, obs2, dat_dim = NULL, decomposition = FALSE, type = 'FairStartDatesBS') + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$bs), + c(0.42222222, 0.44444444, 0.02222222, 0.48888889, 0.37777778, 0.02222222), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$rel), + c(0.22222222, 0.31111111, 0.02222222, 0.28888889, 0.24444444, 0.02222222), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$res), + c(0.0400000, 0.1066667, 0.0000000, 0.0400000, 0.1066667, 0.0000000), + tolerance = 0.0001 + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBS')$unc), + c(0.24, 0.24, 0.00, 0.24, 0.24, 0.00), + tolerance = 0.0001 + ) + + # 'FairStartDatesBSS' + expect_equal( + dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBSS')), + c(bin = 3, ftime = 2) + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairStartDatesBSS')), + c(-0.7592593, -0.8518519, -Inf, -1.0370370, -0.5740741, -Inf), + tolerance = 0.0001 + ) + # 'FairEnsembleBS' + expect_equal( + dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBS')), + c(bin = 3, ftime = 2) + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBS')), + c(0.1333333, 0.2000000, 0.2000000, 0.1333333, 0.4000000, 0.2000000), + tolerance = 0.0001 + ) + # 'FairEnsembleBSS' + expect_equal( + dim(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBSS')), + c(bin = 3, ftime = 2) + ) + expect_equal( + as.vector(UltimateBrier(exp2, obs2, dat_dim = NULL, type = 'FairEnsembleBSS')), + c(-0.1111111, -0.6666667, -0.6666667, 0.2592593, -1.2222222, -0.6666667), + tolerance = 0.0001 + ) + +}) -- GitLab From e40e8d8c73b2957ead96c3e49224672a1344bd31 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 17:34:20 +0200 Subject: [PATCH 40/47] Add test for dat_dim NULL --- tests/testthat/test-RMSSS.R | 25 +++++++++++++++++++++++++ 1 file changed, 25 insertions(+) diff --git a/tests/testthat/test-RMSSS.R b/tests/testthat/test-RMSSS.R index 6064508..d592130 100644 --- a/tests/testthat/test-RMSSS.R +++ b/tests/testthat/test-RMSSS.R @@ -25,6 +25,12 @@ context("s2dv::RMSSS tests") set.seed(6) obs3 <- rnorm(10) + # case 4 + set.seed(7) + exp4 <- array(rnorm(120), dim = c(sdate = 10, dat = 1, lon = 3, lat = 2)) + set.seed(8) + obs4 <- array(rnorm(60), dim = c(dat = 1, sdate = 10, lat = 2, lon = 3)) + ############################################## test_that("1. Input checks", { @@ -158,5 +164,24 @@ test_that("4. Output checks: case 3", { tolerance = 0.00001 ) +}) + +############################################## +test_that("5. Output checks: case 4", { + + expect_equal( + dim(RMSSS(exp4, obs4, dat_dim = NULL)$rmsss), + c(dat = 1, lon = 3, lat = 2) + ) + expect_equal( + mean(RMSSS(exp4, obs4, dat_dim = NULL)$rmsss), + -0.3114424, + tolerance = 0.00001 + ) + expect_equal( + range(RMSSS(exp4, obs4, dat_dim = NULL)$p.val), + c(0.3560534, 0.9192801), + tolerance = 0.00001 + ) }) -- GitLab From 10573fe6c2d48b41c948190911f308aee8a72fc9 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 17:44:26 +0200 Subject: [PATCH 41/47] Add test for dat_dim NULL --- tests/testthat/test-RMS.R | 50 +++++++++++++++++++++++++++++++++++++++ 1 file changed, 50 insertions(+) diff --git a/tests/testthat/test-RMS.R b/tests/testthat/test-RMS.R index f1083dd..65edfac 100644 --- a/tests/testthat/test-RMS.R +++ b/tests/testthat/test-RMS.R @@ -17,6 +17,13 @@ context("s2dv::RMS tests") set.seed(6) obs2 <- rnorm(10) + # dat3 + set.seed(1) + exp3 <- array(rnorm(120), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) + + set.seed(2) + obs3 <- array(rnorm(80), dim = c(sdate = 5, ftime = 2, lon = 1, lat = 4)) + ############################################## test_that("1. Input checks", { @@ -186,3 +193,46 @@ test_that("3. Output checks: dat2", { }) ############################################## +test_that("4. Output checks: dat3", { + + expect_equal( + dim(RMS(exp3, obs3, dat_dim = NULL)$rms), + c(ftime = 2, lon = 1, lat = 4) + ) + + suppressWarnings( + expect_equal( + RMS(exp3, obs3, dat_dim = NULL)$rms[1:6], + c(1.6458118, 0.8860392, 0.8261295, 1.1681939, 2.1693538, 1.3064454), + tolerance = 0.001 + ) +) +suppressWarnings( + expect_equal( + length(which(is.na(RMS(exp3, obs3, dat_dim = NULL)$conf.lower))), + 0 + ) +) +suppressWarnings( + expect_equal( + max(RMS(exp3, obs3, dat_dim = NULL)$conf.lower, na.rm = T), + 1.453144, + tolerance = 0.001 + ) +) +suppressWarnings( + expect_equal( + min(RMS(exp3, obs3, dat_dim = NULL, conf.lev = 0.99)$conf.upper, na.rm = TRUE), + 2.646274, + tolerance = 0.0001 + ) +) +suppressWarnings( + expect_equal( + length(RMS(exp3, obs3, dat_dim = NULL, conf = FALSE)), + 1 + ) +) +}) + +############################################## \ No newline at end of file -- GitLab From 48aebf43d73283e0041c9d41abd12f8d77757a56 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 21 Jul 2022 17:57:11 +0200 Subject: [PATCH 42/47] Changed tests RMS --- tests/testthat/test-RMS.R | 34 +++------------------------------- 1 file changed, 3 insertions(+), 31 deletions(-) diff --git a/tests/testthat/test-RMS.R b/tests/testthat/test-RMS.R index 65edfac..f69cfd8 100644 --- a/tests/testthat/test-RMS.R +++ b/tests/testthat/test-RMS.R @@ -200,39 +200,11 @@ test_that("4. Output checks: dat3", { c(ftime = 2, lon = 1, lat = 4) ) - suppressWarnings( expect_equal( - RMS(exp3, obs3, dat_dim = NULL)$rms[1:6], - c(1.6458118, 0.8860392, 0.8261295, 1.1681939, 2.1693538, 1.3064454), - tolerance = 0.001 - ) -) -suppressWarnings( - expect_equal( - length(which(is.na(RMS(exp3, obs3, dat_dim = NULL)$conf.lower))), - 0 - ) -) -suppressWarnings( - expect_equal( - max(RMS(exp3, obs3, dat_dim = NULL)$conf.lower, na.rm = T), - 1.453144, - tolerance = 0.001 - ) -) -suppressWarnings( - expect_equal( - min(RMS(exp3, obs3, dat_dim = NULL, conf.lev = 0.99)$conf.upper, na.rm = TRUE), - 2.646274, - tolerance = 0.0001 - ) -) -suppressWarnings( - expect_equal( - length(RMS(exp3, obs3, dat_dim = NULL, conf = FALSE)), - 1 + as.vector(RMS(exp3, obs3, dat_dim = NULL)$rms), + c(1.6458118, 0.8860392, 0.8261295, 1.1681939, 2.1693538, 1.3064454, 0.5384229, 1.1215333), + tolerance = 0.00001 ) -) }) ############################################## \ No newline at end of file -- GitLab From bfd7e9dc9a01a8a6fff16871ccc4abc5b8cb4d4e Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Fri, 22 Jul 2022 12:07:14 +0200 Subject: [PATCH 43/47] Revert changes Consist_trend and documentation --- R/Consist_Trend.R | 21 +++++++-------------- man/Consist_Trend.Rd | 6 ++---- 2 files changed, 9 insertions(+), 18 deletions(-) diff --git a/R/Consist_Trend.R b/R/Consist_Trend.R index 84b7699..ee95684 100644 --- a/R/Consist_Trend.R +++ b/R/Consist_Trend.R @@ -14,8 +14,7 @@ #'@param dat_dim A character string indicating the name of the dataset #' dimensions. If data at some point of 'time_dim' are not complete along #' 'dat_dim' in both 'exp' and 'obs', this point in all 'dat_dim' will be -#' discarded. The default value is 'dataset'. If there is no dataset -#' dimension, set NULL. +#' discarded. The default value is 'dataset'. #'@param time_dim A character string indicating the name of dimension along #' which the trend is computed. The default value is 'sdate'. #'@param interval A positive numeric indicating the unit length between two @@ -30,8 +29,7 @@ #' with dimensions c(stats = 2, nexp + nobs, the rest dimensions of 'exp' and #' 'obs' except time_dim), where 'nexp' is the length of 'dat_dim' in 'exp' #' and 'nobs' is the length of 'dat_dim' in 'obs. The 'stats' dimension -#' contains the intercept and the slope. If dat_dim is NULL, nexp and nobs are -#' omitted. +#' contains the intercept and the slope. #'} #'\item{$conf.lower}{ #' A numeric array of the lower limit of 95\% confidence interval with @@ -111,24 +109,19 @@ Consist_Trend <- function(exp, obs, dat_dim = 'dataset', time_dim = 'sdate', int stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.null(dat_dim)) { - if (!is.character(dat_dim) | length(dat_dim) > 1) { - stop("Parameter 'dat_dim' must be a character string or NULL.") - } - if (!dat_dim %in% names(dim(exp)) | !dat_dim %in% names(dim(obs))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", - "Set it as NULL if there is no dataset dimension.") - } + if (!is.character(dat_dim)) { + stop("Parameter 'dat_dim' must be a character string.") + } + if (!all(dat_dim %in% names(dim(exp))) | !all(dat_dim %in% names(dim(obs)))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") } ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - if (!is.null(dat_dim)) { for (i in 1:length(dat_dim)) { name_exp <- name_exp[-which(name_exp == dat_dim[i])] name_obs <- name_obs[-which(name_obs == dat_dim[i])] } - } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have same length of ", "all dimension except 'dat_dim'.")) diff --git a/man/Consist_Trend.Rd b/man/Consist_Trend.Rd index b93dad5..2ac7d42 100644 --- a/man/Consist_Trend.Rd +++ b/man/Consist_Trend.Rd @@ -23,8 +23,7 @@ parameter 'exp' except along 'dat_dim'.} \item{dat_dim}{A character string indicating the name of the dataset dimensions. If data at some point of 'time_dim' are not complete along 'dat_dim' in both 'exp' and 'obs', this point in all 'dat_dim' will be -discarded. The default value is 'dataset'. If there is no dataset -dimension, set NULL.} +discarded. The default value is 'dataset'.} \item{time_dim}{A character string indicating the name of dimension along which the trend is computed. The default value is 'sdate'.} @@ -42,8 +41,7 @@ A list containing: with dimensions c(stats = 2, nexp + nobs, the rest dimensions of 'exp' and 'obs' except time_dim), where 'nexp' is the length of 'dat_dim' in 'exp' and 'nobs' is the length of 'dat_dim' in 'obs. The 'stats' dimension - contains the intercept and the slope. If dat_dim is NULL, nexp and nobs are - omitted. + contains the intercept and the slope. } \item{$conf.lower}{ A numeric array of the lower limit of 95\% confidence interval with -- GitLab From 8cf0bfa6456714b548354d9fac342f49cb211ef1 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Fri, 22 Jul 2022 12:26:29 +0200 Subject: [PATCH 44/47] Update documentation with allow dat_dim NULL for UltimateBrier() --- R/UltimateBrier.R | 9 +++++---- man/UltimateBrier.Rd | 9 +++++---- 2 files changed, 10 insertions(+), 8 deletions(-) diff --git a/R/UltimateBrier.R b/R/UltimateBrier.R index 18c8e5b..d2c4ac9 100644 --- a/R/UltimateBrier.R +++ b/R/UltimateBrier.R @@ -5,14 +5,15 @@ #'to choose. #' #'@param exp A numeric array of forecast anomalies with named dimensions that -#' at least include 'dat_dim', 'memb_dim', and 'time_dim'. It can be provided +#' at least include 'memb_dim', and 'time_dim'. It can be provided #' by \code{Ano()}. #'@param obs A numeric array of observational reference anomalies with named -#' dimensions that at least include 'dat_dim' and 'time_dim'. If it has +#' dimensions that at least include 'time_dim'. If it has #' 'memb_dim', the length must be 1. The dimensions should be consistent with #' 'exp' except 'dat_dim' and 'memb_dim'. It can be provided by \code{Ano()}. #'@param dat_dim A character string indicating the name of the dataset -#' dimension in 'exp' and 'obs'. The default value is 'dataset'. +#' dimension in 'exp' and 'obs'. The default value is 'dataset'. If there is no dataset +#' dimension, set NULL. #'@param memb_dim A character string indicating the name of the member #' dimension in 'exp' (and 'obs') for ensemble mean calculation. The default #' value is 'member'. @@ -55,7 +56,7 @@ #'same dimensions: #'c(nexp, nobs, no. of bins, the rest dimensions of 'exp' except 'time_dim' and #''memb_dim'). 'nexp' and 'nobs' is the length of dataset dimension in 'exp' -#'and 'obs' respectively.\cr +#'and 'obs' respectively. If dat_dim is NULL, nexp and nobs are omitted.\cr #'The list of 4 includes: #' \itemize{ #' \item{$bs: Brier Score} diff --git a/man/UltimateBrier.Rd b/man/UltimateBrier.Rd index 2cad133..0dfa772 100644 --- a/man/UltimateBrier.Rd +++ b/man/UltimateBrier.Rd @@ -19,16 +19,17 @@ UltimateBrier( } \arguments{ \item{exp}{A numeric array of forecast anomalies with named dimensions that -at least include 'dat_dim', 'memb_dim', and 'time_dim'. It can be provided +at least include 'memb_dim', and 'time_dim'. It can be provided by \code{Ano()}.} \item{obs}{A numeric array of observational reference anomalies with named -dimensions that at least include 'dat_dim' and 'time_dim'. If it has +dimensions that at least include 'time_dim'. If it has 'memb_dim', the length must be 1. The dimensions should be consistent with 'exp' except 'dat_dim' and 'memb_dim'. It can be provided by \code{Ano()}.} \item{dat_dim}{A character string indicating the name of the dataset -dimension in 'exp' and 'obs'. The default value is 'dataset'.} +dimension in 'exp' and 'obs'. The default value is 'dataset'. If there is no dataset +dimension, set NULL.} \item{memb_dim}{A character string indicating the name of the member dimension in 'exp' (and 'obs') for ensemble mean calculation. The default @@ -78,7 +79,7 @@ is an array of Brier scores or Brier skill scores. All the arrays have the same dimensions: c(nexp, nobs, no. of bins, the rest dimensions of 'exp' except 'time_dim' and 'memb_dim'). 'nexp' and 'nobs' is the length of dataset dimension in 'exp' -and 'obs' respectively.\cr +and 'obs' respectively. If dat_dim is NULL, nexp and nobs are omitted.\cr The list of 4 includes: \itemize{ \item{$bs: Brier Score} -- GitLab From 1c51947982b1570ff09b4b10cac16f1b7040fc3d Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Mon, 25 Jul 2022 16:32:10 +0200 Subject: [PATCH 45/47] Allow dat_dim NULL in Ano_CrossValid and added tests --- R/Ano_CrossValid.R | 113 +++++++++++++++++++-------- tests/testthat/test-Ano_CrossValid.R | 28 +++++++ 2 files changed, 110 insertions(+), 31 deletions(-) diff --git a/R/Ano_CrossValid.R b/R/Ano_CrossValid.R index 9920502..d1996b9 100644 --- a/R/Ano_CrossValid.R +++ b/R/Ano_CrossValid.R @@ -18,7 +18,8 @@ #'@param dat_dim A character vector indicating the name of the dataset and #' member dimensions. When calculating the climatology, if data at one #' startdate (i.e., 'time_dim') is not complete along 'dat_dim', this startdate -#' along 'dat_dim' will be discarded. The default value is +#' along 'dat_dim' will be discarded. If there is no dataset dimension, it can be NULL. +#' The default value is #' "c('dataset', 'member')". #'@param memb_dim A character string indicating the name of the member #' dimension. Only used when parameter 'memb' is FALSE. It must be one element @@ -83,11 +84,25 @@ Ano_CrossValid <- function(exp, obs, time_dim = 'sdate', dat_dim = c('dataset', stop("Parameter 'time_dim' is not found in 'exp' or 'obs' dimension.") } ## dat_dim - if (!is.character(dat_dim)) { - stop("Parameter 'dat_dim' must be a character vector.") - } - if (!all(dat_dim %in% names(dim(exp))) | !all(dat_dim %in% names(dim(obs)))) { - stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.") + if (!is.null(dat_dim)) { + if (!is.character(dat_dim)) { + stop("Parameter 'dat_dim' must be a character vector.") + } + if (!all(dat_dim %in% names(dim(exp))) | !all(dat_dim %in% names(dim(obs)))) { + stop("Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension.", + " Set it as NULL if there is no dataset dimension.") + } + # If dat_dim is not in obs, add it in + if (any(!dat_dim %in% names(dim(obs)))) { + reset_obs_dim <- TRUE + ori_obs_dim <- dim(obs) + dim(obs) <- c(dim(obs), rep(1, length(dat_dim[which(!dat_dim %in% names(dim(obs)))]))) + names(dim(obs)) <- c(names(ori_obs_dim), dat_dim[which(!dat_dim %in% names(dim(obs)))]) + } else { + reset_obs_dim <- FALSE + } + } else { + reset_obs_dim <- FALSE } ## memb if (!is.logical(memb) | length(memb) > 1) { @@ -115,9 +130,11 @@ Ano_CrossValid <- function(exp, obs, time_dim = 'sdate', dat_dim = c('dataset', ## exp and obs (2) name_exp <- sort(names(dim(exp))) name_obs <- sort(names(dim(obs))) - for (i in 1:length(dat_dim)) { - name_exp <- name_exp[-which(name_exp == dat_dim[i])] - name_obs <- name_obs[-which(name_obs == dat_dim[i])] + if (!is.null(dat_dim)) { + for (i in 1:length(dat_dim)) { + name_exp <- name_exp[-which(name_exp == dat_dim[i])] + name_obs <- name_obs[-which(name_obs == dat_dim[i])] + } } if(!all(dim(exp)[name_exp] == dim(obs)[name_obs])) { stop(paste0("Parameter 'exp' and 'obs' must have the same length of ", @@ -135,36 +152,65 @@ Ano_CrossValid <- function(exp, obs, time_dim = 'sdate', dat_dim = c('dataset', #----------------------------------- # Per-paired method: If any sdate along dat_dim is NA, turn all sdate points along dat_dim into NA. - pos <- rep(0, length(dat_dim)) # dat_dim: [dataset, member] - for (i in 1:length(dat_dim)) { - pos[i] <- which(names(dim(obs)) == dat_dim[i]) - } - outrows_exp <- MeanDims(exp, pos, na.rm = FALSE) + - MeanDims(obs, pos, na.rm = FALSE) - outrows_obs <- outrows_exp - - for (i in 1:length(pos)) { - outrows_exp <- InsertDim(outrows_exp, pos[i], dim(exp)[pos[i]]) - outrows_obs <- InsertDim(outrows_obs, pos[i], dim(obs)[pos[i]]) - } - exp_for_clim <- exp - obs_for_clim <- obs - exp_for_clim[which(is.na(outrows_exp))] <- NA - obs_for_clim[which(is.na(outrows_obs))] <- NA + if (!is.null(dat_dim)) { + pos <- rep(0, length(dat_dim)) # dat_dim: [dataset, member] + for (i in 1:length(dat_dim)) { + pos[i] <- which(names(dim(obs)) == dat_dim[i]) + } + outrows_exp <- MeanDims(exp, pos, na.rm = FALSE) + + MeanDims(obs, pos, na.rm = FALSE) + outrows_obs <- outrows_exp + + for (i in 1:length(pos)) { + outrows_exp <- InsertDim(outrows_exp, pos[i], dim(exp)[pos[i]]) + outrows_obs <- InsertDim(outrows_obs, pos[i], dim(obs)[pos[i]]) + } + exp_for_clim <- exp + obs_for_clim <- obs + exp_for_clim[which(is.na(outrows_exp))] <- NA + obs_for_clim[which(is.na(outrows_obs))] <- NA + } else { + exp_for_clim <- exp + obs_for_clim <- obs + } + #----------------------------------- + res <- Apply(list(exp, obs, exp_for_clim, obs_for_clim), + target_dims = c(time_dim, dat_dim), + fun = .Ano_CrossValid, dat_dim = dat_dim, + memb_dim = memb_dim, memb = memb, + ncores = ncores) - res <- Apply(list(exp, obs, exp_for_clim, obs_for_clim), - target_dims = c(time_dim, dat_dim), - fun = .Ano_CrossValid, - memb_dim = memb_dim, memb = memb, - ncores = ncores) + # Remove dat_dim in obs if obs doesn't have at first place + if (reset_obs_dim) { + res_obs_dim <- ori_obs_dim[-which(names(ori_obs_dim) == time_dim)] + if (!memb & memb_dim %in% names(res_obs_dim)) { + res_obs_dim <- res_obs_dim[-which(names(res_obs_dim) == memb_dim)] + } + if (is.integer(res_obs_dim) & length(res_obs_dim) == 0) { + res$obs <- as.vector(res$obs) + } else { + res$obs <- array(res$obs, dim = res_obs_dim) + } + } return(res) } -.Ano_CrossValid <- function(exp, obs, exp_for_clim, obs_for_clim, +.Ano_CrossValid <- function(exp, obs, exp_for_clim, obs_for_clim, dat_dim = c('dataset', 'member'), memb_dim = 'member', memb = TRUE, ncores = NULL) { + if (is.null(dat_dim)) { + ini_dims_exp <- dim(exp) + ini_dims_obs <- dim(obs) + ini_dims_exp_for_clim <- dim(exp) + ini_dims_obs_for_clim <- dim(exp) + exp <- InsertDim(exp, posdim = 2, lendim = 1, name = 'dataset') + exp_for_clim <- InsertDim(exp_for_clim, posdim = 2, lendim = 1, name = 'dataset') + obs <- InsertDim(obs, posdim = 2, lendim = 1, name = 'dataset') + obs_for_clim <- InsertDim(obs_for_clim, posdim = 2, lendim = 1, name = 'dataset') + } + # exp: [sdate, dat_dim, memb_dim] # obs: [sdate, dat_dim, memb_dim] ano_exp_list <- vector('list', length = dim(exp)[1]) #length: [sdate] @@ -222,5 +268,10 @@ Ano_CrossValid <- function(exp, obs, time_dim = 'sdate', dat_dim = c('dataset', ano_obs <- array(unlist(ano_obs_list), dim = c(dim(obs)[-1], dim(obs)[1])) ano_obs <- Reorder(ano_obs, c(length(dim(obs)), 1:(length(dim(obs)) - 1))) + if (is.null(dat_dim)) { + ano_exp <- array(ano_exp, dim = ini_dims_exp) + ano_obs <- array(ano_obs, dim = ini_dims_obs) + } + return(list(exp = ano_exp, obs = ano_obs)) } diff --git a/tests/testthat/test-Ano_CrossValid.R b/tests/testthat/test-Ano_CrossValid.R index b66fc5f..c50938e 100644 --- a/tests/testthat/test-Ano_CrossValid.R +++ b/tests/testthat/test-Ano_CrossValid.R @@ -13,6 +13,12 @@ exp2 <- array(rnorm(30), dim = c(member = 3, ftime = 2, sdate = 5)) set.seed(2) obs2 <- array(rnorm(20), dim = c(ftime = 2, member = 2, sdate = 5)) +# dat3 +set.seed(1) +exp3 <- array(rnorm(30), dim = c(ftime = 2, sdate = 5)) +set.seed(2) +obs3 <- array(rnorm(20), dim = c(ftime = 2, sdate = 5)) + ############################################## test_that("1. Input checks", { @@ -136,6 +142,28 @@ test_that("3. dat2", { }) +############################################## +test_that("4. dat3", { + expect_equal( + names(Ano_CrossValid(exp3, obs3, dat_dim = NULL)), + c("exp", "obs") + ) + expect_equal( + dim(Ano_CrossValid(exp3, obs3, dat_dim = NULL)$exp), + c(sdate = 5, ftime = 2) + ) + expect_equal( + Ano_CrossValid(exp3, obs3, dat_dim = NULL)$exp[, 2], + c(-0.1182939, 1.6462530, -1.3734335, 0.5750579, -0.7295835), + tolerance = 0.0001 + ) + expect_equal( + Ano_CrossValid(exp3, obs3, dat_dim = NULL)$exp[, 2], + c(-0.1182939, 1.6462530, -1.3734335, 0.5750579, -0.7295835), + tolerance = 0.0001 + ) + +}) -- GitLab From 240c47cd90323ba3b572f69ada7785882e055350 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Mon, 25 Jul 2022 16:41:19 +0200 Subject: [PATCH 46/47] Added text in test-Ano_CrossValid --- tests/testthat/test-Ano_CrossValid.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/testthat/test-Ano_CrossValid.R b/tests/testthat/test-Ano_CrossValid.R index c50938e..2d7c00c 100644 --- a/tests/testthat/test-Ano_CrossValid.R +++ b/tests/testthat/test-Ano_CrossValid.R @@ -57,7 +57,7 @@ test_that("1. Input checks", { ) expect_error( Ano_CrossValid(exp1, obs1, dat_dim = 'dat'), - "Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension." + "Parameter 'dat_dim' is not found in 'exp' or 'obs' dimension. Set it as NULL if there is no dataset dimension." ) # memb expect_error( -- GitLab From 739c51b662396965e955607a8049e314d8a0615b Mon Sep 17 00:00:00 2001 From: aho Date: Thu, 28 Jul 2022 17:08:03 +0200 Subject: [PATCH 47/47] Update documentation --- man/Ano_CrossValid.Rd | 9 +++++---- man/RMS.Rd | 3 +-- man/RMSSS.Rd | 3 +-- 3 files changed, 7 insertions(+), 8 deletions(-) diff --git a/man/Ano_CrossValid.Rd b/man/Ano_CrossValid.Rd index d2234a1..2a82713 100644 --- a/man/Ano_CrossValid.Rd +++ b/man/Ano_CrossValid.Rd @@ -25,10 +25,11 @@ parameter 'exp' except along 'dat_dim'.} The default value is 'sdate'.} \item{dat_dim}{A character vector indicating the name of the dataset and -member dimensions. When calculating the climatology, if data at one -startdate (i.e., 'time_dim') is not complete along 'dat_dim', this startdate -along 'dat_dim' will be discarded. The default value is -"c('dataset', 'member')".} + member dimensions. When calculating the climatology, if data at one + startdate (i.e., 'time_dim') is not complete along 'dat_dim', this startdate + along 'dat_dim' will be discarded. If there is no dataset dimension, it can be NULL. +The default value is + "c('dataset', 'member')".} \item{memb_dim}{A character string indicating the name of the member dimension. Only used when parameter 'memb' is FALSE. It must be one element diff --git a/man/RMS.Rd b/man/RMS.Rd index a84f965..4391df4 100644 --- a/man/RMS.Rd +++ b/man/RMS.Rd @@ -31,8 +31,7 @@ length as 'exp', then the vector will automatically be 'time_dim' and which the correlations are computed. The default value is 'sdate'.} \item{dat_dim}{A character string indicating the name of member (nobs/nexp) -dimension. The default value is 'dataset'. If there is no dataset -dimension, set NULL.} +dimension. The default value is 'dataset'.} \item{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 diff --git a/man/RMSSS.Rd b/man/RMSSS.Rd index 2c23ee4..9ebcf65 100644 --- a/man/RMSSS.Rd +++ b/man/RMSSS.Rd @@ -30,8 +30,7 @@ be 1.} which the RMSSS are computed. The default value is 'sdate'.} \item{dat_dim}{A character string indicating the name of dataset (nobs/nexp) -dimension. The default value is 'dataset'. If there is no dataset -dimension, set NULL.} +dimension. The default value is 'dataset'.} \item{pval}{A logical value indicating whether to compute or not the p-value of the test Ho: RMSSS = 0. If pval = TRUE, the insignificant RMSSS will -- GitLab