From b85a1ce2cbbf7df148fd00a927c531e4d1bb697a Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 13 Oct 2022 17:46:49 +0200 Subject: [PATCH 1/5] Correct memb_dim, documentation and improve code of BiasCorrection. Added tests --- DESCRIPTION | 2 +- R/CST_BiasCorrection.R | 230 ++++++++++++++--------- man/BiasCorrection.Rd | 39 ++-- man/CST_BiasCorrection.Rd | 39 ++-- tests/testthat/test-CST_BiasCorrection.R | 227 ++++++++++++++++------ 5 files changed, 369 insertions(+), 168 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 454397ae..5b1ac599 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -86,4 +86,4 @@ VignetteBuilder: knitr License: Apache License 2.0 Encoding: UTF-8 LazyData: true -RoxygenNote: 7.0.2 +RoxygenNote: 7.2.0 diff --git a/R/CST_BiasCorrection.R b/R/CST_BiasCorrection.R index 263128aa..eb6ca04a 100644 --- a/R/CST_BiasCorrection.R +++ b/R/CST_BiasCorrection.R @@ -1,26 +1,37 @@ #' Bias Correction based on the mean and standard deviation adjustment #' #'@author Verónica Torralba, \email{veronica.torralba@bsc.es} -#'@description This function applies the simple bias adjustment technique described in Torralba et al. (2017). The adjusted forecasts have an equivalent standard deviation and mean to that of the reference dataset. +#'@description This function applies the simple bias adjustment technique +#'described in Torralba et al. (2017). The adjusted forecasts have an equivalent +#'standard deviation and mean to that of the reference dataset. #' -#'@param exp an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the seasonal forecast experiment data in the element named \code{$data} -#'@param obs an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the observed data in the element named \code{$data}. -#'@param exp_cor an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the seasonl forecast experiment to be corrected. If it is NULL, the 'exp' forecast will be corrected. -#'@param na.rm a logical value indicating whether missing values should be stripped before the computation proceeds, by default it is set to FALSE. +#'@param exp An object of class \code{s2dv_cube} as returned by \code{CST_Load} +#' function, containing the seasonal forecast experiment data in the element +#' named \code{$data} +#'@param obs An object of class \code{s2dv_cube} as returned by \code{CST_Load} +#' function, containing the observed data in the element named \code{$data}. +#'@param exp_cor An object of class \code{s2dv_cube} as returned by +#' \code{CST_Load} function, containing the seasonl forecast experiment to be +#' corrected. If it is NULL, the 'exp' forecast will be corrected. +#'@param na.rm A logical value indicating whether missing values should be +#' stripped before the computation proceeds, by default it is set to FALSE. +#'@param memb_dim A character string indicating the name of the member +#' dimension. By default, it is set to 'member'. +#'@param sdate_dim A character string indicating the name of the start date +#' dimension. By default, it is set to 'sdate'. +#'@param ncores An integer that indicates the number of cores for parallel +#' computations using multiApply function. The default value is NULL. +#'@return An object of class \code{s2dv_cube} containing the bias corrected +#'forecasts with the same dimensions of the experimental data. #' -#'@param memb_dim is a character string indicating the name of the member dimension. By default, it is set to 'member'. -#'@param sdate_dim is a character string indicating the name of the start date dimension. By default, it is set to 'sdate'. -#'@param ncores is an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one. -#'@return an object of class \code{s2dv_cube} containing the bias corrected forecasts in the element called \code{$data} with the same dimensions of the experimental data. +#'@references Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. +#'Davis (2017). Seasonal climate prediction: a new source of information for +#'the management of wind energy resources. Journal of Applied Meteorology and +#'Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, +#'EUPORIAS, NEWA, RESILIENCE, SPECS) #' -#'@references Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. Davis (2017). Seasonal climate prediction: a new source of information for the management of wind energy resources. Journal of Applied Meteorology and Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, EUPORIAS, NEWA, RESILIENCE, SPECS) -#' -#'@import multiApply #'@examples -#' #'# Example -#'# Creation of sample s2dverification objects. These are not complete -#'# s2dverification objects though. The Load function returns complete objects. #'mod1 <- 1 : (1 * 3 * 4 * 5 * 6 * 7) #'dim(mod1) <- c(dataset = 1, member = 3, sdate = 4, ftime = 5, lat = 6, lon = 7) #'obs1 <- 1 : (1 * 1 * 4 * 5 * 6 * 7) @@ -32,118 +43,158 @@ #'attr(exp, 'class') <- 's2dv_cube' #'attr(obs, 'class') <- 's2dv_cube' #'a <- CST_BiasCorrection(exp = exp, obs = obs) -#'str(a) +#'@import multiApply #'@export CST_BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE, memb_dim = 'member', sdate_dim = 'sdate', - ncores = 1) { + ncores = NULL) { if (!inherits(exp, 's2dv_cube') || !inherits(obs, 's2dv_cube')) { stop("Parameter 'exp' and 'obs' must be of the class 's2dv_cube', ", "as output by CSTools::CST_Load.") } - if (dim(obs$data)['member'] != 1) { - stop("The length of the dimension 'member' in the component 'data' ", - "of the parameter 'obs' must be equal to 1.") + if (!is.null(exp_cor)) { + if (!inherits(exp_cor, 's2dv_cube')) { + stop("Parameter 'exp_cor' must be of the class 's2dv_cube', ", + "as output by CSTools::CST_Load.") + } + dimnames <- names(dim(exp_cor$data)) + } else { + dimnames <- names(dim(exp$data)) } + + BiasCorrected <- BiasCorrection(exp = exp$data, obs = obs$data, exp_cor = exp_cor$data, + memb_dim = memb_dim, sdate_dim = sdate_dim, + na.rm = na.rm, ncores = ncores) + + pos <- match(dimnames, names(dim(BiasCorrected))) + BiasCorrected <- aperm(BiasCorrected, pos) + if (is.null(exp_cor)) { - dimnames <- names(dim(exp$data)) - BiasCorrected <- BiasCorrection(exp = exp$data, obs = obs$data, na.rm = na.rm, - memb_dim = memb_dim, sdate_dim = sdate_dim, - ncores = ncores) - pos <- match(dimnames, names(dim(BiasCorrected))) - BiasCorrected <- aperm(BiasCorrected, pos) - names(dim(BiasCorrected)) <- dimnames exp$data <- BiasCorrected - exp$Datasets <- c(exp$Datasets, obs$Datasets) - exp$source_files <- c(exp$source_files, obs$source_files) return(exp) - } else { - if (!inherits(exp_cor, 's2dv_cube')) { - stop("Parameter 'var_cor' must be of the class 's2dv_cube'.") - } - dimnames <- names(dim(exp_cor$data)) - BiasCorrected <- BiasCorrection(exp = exp$data, obs = obs$data, - exp_cor = exp_cor$data, na.rm = na.rm) - pos <- match(dimnames, names(dim(BiasCorrected))) - BiasCorrected <- aperm(BiasCorrected, pos) - names(dim(BiasCorrected)) <- dimnames + + } else { exp_cor$data <- BiasCorrected - exp_cor$Datasets <- c(exp_cor$Datasets, exp$Datasets, obs$Datasets) - exp_cor$source_files <- c(exp_cor$source_files, - exp$source_files, obs$source_files) return(exp_cor) } } + #' Bias Correction based on the mean and standard deviation adjustment #' #'@author Verónica Torralba, \email{veronica.torralba@bsc.es} -#'@description This function applies the simple bias adjustment technique described in Torralba et al. (2017). The adjusted forecasts have an equivalent standard deviation and mean to that of the reference dataset. +#'@description This function applies the simple bias adjustment technique +#'described in Torralba et al. (2017). The adjusted forecasts have an equivalent +#'standard deviation and mean to that of the reference dataset. #' -#'@param exp a multidimensional array with named dimensions containing the seasonal forecast experiment data with at least 'member' and 'sdate' dimensions. -#'@param obs a multidimensional array with named dimensions containing the observed data with at least 'sdate' dimension. -#'@param exp_cor a multidimensional array with named dimensions containing the seasonl forecast experiment to be corrected. If it is NULL, the 'exp' forecast will be corrected. -#'@param na.rm a logical value indicating whether missing values should be stripped before the computation proceeds, by default it is set to FALSE. -#'@param memb_dim is a character string indicating the name of the member dimension. By default, it is set to 'member'. -#'@param sdate_dim is a character string indicating the name of the start date dimension. By default, it is set to 'sdate'. -#'@param ncores is an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one. +#'@param exp A multidimensional array with named dimensions containing the +#' seasonal forecast experiment data with at least 'member' and 'sdate' +#' dimensions. +#'@param obs A multidimensional array with named dimensions containing the +#' observed data with at least 'sdate' dimension. +#'@param exp_cor A multidimensional array with named dimensions containing the +#' seasonl forecast experiment to be corrected. If it is NULL, the 'exp' +#' forecast will be corrected. +#'@param na.rm A logical value indicating whether missing values should be +#' stripped before the computation proceeds, by default it is set to FALSE. +#'@param memb_dim A character string indicating the name of the member +#' dimension. By default, it is set to 'member'. +#'@param sdate_dim A character string indicating the name of the start date +#' dimension. By default, it is set to 'sdate'. +#'@param ncores An integer that indicates the number of cores for parallel +#' computations using multiApply function. The default value is NULL. #' -#'@return an object of class \code{s2dv_cube} containing the bias corrected forecasts in the element called \code{$data} with the same dimensions of the experimental data. +#'@return An array containing the bias corrected forecasts with the same +#'dimensions of the experimental data. #' -#'@references Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. Davis (2017). Seasonal climate prediction: a new source of information for the management of wind energy resources. Journal of Applied Meteorology and Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, EUPORIAS, NEWA, RESILIENCE, SPECS) +#'@references Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. +#'Davis (2017). Seasonal climate prediction: a new source of information for the +#'management of wind energy resources. Journal of Applied Meteorology and +#'Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, +#'EUPORIAS, NEWA, RESILIENCE, SPECS) #' -#'@import multiApply #'@examples -#' #'# Example -#'# Creation of sample s2dverification objects. These are not complete -#'# s2dverification objects though. The Load function returns complete objects. #'mod1 <- 1 : (1 * 3 * 4 * 5 * 6 * 7) #'dim(mod1) <- c(dataset = 1, member = 3, sdate = 4, ftime = 5, lat = 6, lon = 7) #'obs1 <- 1 : (1 * 1 * 4 * 5 * 6 * 7) #'dim(obs1) <- c(dataset = 1, member = 1, sdate = 4, ftime = 5, lat = 6, lon = 7) #'a <- BiasCorrection(exp = mod1, obs = obs1) -#'str(a) +#'@import multiApply #'@export -BiasCorrection <- function (exp, obs, exp_cor = NULL, na.rm = FALSE, - memb_dim = 'member', sdate_dim = 'sdate', - ncores = 1) { - - if (!all(c(memb_dim, sdate_dim) %in% names(dim(exp)))) { - stop(paste("Parameter 'exp' must have the dimensions set up in 'memb_dim'", - "and 'sdate_dim' parameters.")) +BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE, + memb_dim = 'member', sdate_dim = 'sdate', + ncores = NULL) { + # Check inputs + ## exp, obs + if (!is.array(exp) || !is.numeric(exp)) { + stop("Parameter 'exp' must be a numeric array.") } - - if (!all(sdate_dim %in% names(dim(obs)))) { - stop("Parameter 'obs' must have the dimension specified in 'sdate_dim'.") + if (!is.array(obs) || !is.numeric(obs)) { + stop("Parameter 'obs' must be a numeric array.") + } + obsdims <- names(dim(obs)) + expdims <- names(dim(exp)) + if (is.null(expdims)) { + stop("Parameter 'exp' must have dimension names.") + } + if (is.null(obsdims)) { + stop("Parameter 'obs' must have dimension names.") } - if (any(is.na(exp))) { warning("Parameter 'exp' contains NA values.") } - if (any(is.na(obs))) { warning("Parameter 'obs' contains NA values.") } - + ## exp_cor + if (!is.null(exp_cor)) { + exp_cordims <- names(dim(exp_cor)) + if (is.null(exp_cordims)) { + stop("Parameter 'exp_cor' must have dimension names.") + } + } + ## sdate_dim, memb_dim + if (!is.character(sdate_dim) || length(sdate_dim) != 1) { + stop("Parameter 'sdate_dim' must be a character string.") + } + if (!sdate_dim %in% expdims || !sdate_dim %in% obsdims) { + stop("Parameter 'sdate_dim' is not found in 'exp' or 'obs' dimension.") + } + if (dim(exp)[sdate_dim] == 1) { + stop("Parameter 'exp' must have dimension length of 'sdate_dim' bigger than 1.") + } + if (!all(c(memb_dim, sdate_dim) %in% expdims)) { + stop("Parameter 'exp' requires 'sdate_dim' and 'member_dim' dimensions.") + } + ## na.rm if (!is.logical(na.rm)) { na.rm <- FALSE warning("Paramater 'na.rm' must be a logical, it has been set to FALSE.") } - - if (length(na.rm)>1) { + if (length(na.rm) > 1) { na.rm <- na.rm[1] warning("Paramter 'na.rm' has length greater than 1, and only the fist element is used.") } + ## ncores + if (!is.null(ncores)) { + if (!is.numeric(ncores) | ncores %% 1 != 0 | ncores <= 0 | + length(ncores) > 1) { + stop("Parameter 'ncores' must be either NULL or a positive integer.") + } + } - target_dims_obs <- sdate_dim - if (memb_dim %in% names(dim(obs))) { - target_dims_obs <- c(memb_dim, target_dims_obs) + if (memb_dim %in% obsdims) { + target_dims_obs <- c(memb_dim, sdate_dim) + } else { + target_dims_obs <- c(sdate_dim) } + if (is.null(exp_cor)) { BiasCorrected <- Apply(data = list(var_obs = obs, var_exp = exp), target_dims = list(target_dims_obs, c(memb_dim, sdate_dim)), - fun = .sbc , na.rm = na.rm, ncores = ncores)$output1 + fun = .sbc, memb_dim = memb_dim, sdate_dim = sdate_dim, + na.rm = na.rm, ncores = ncores)$output1 } else { BiasCorrected <- Apply(data = list(var_obs = obs, var_exp = exp, @@ -151,29 +202,30 @@ BiasCorrection <- function (exp, obs, exp_cor = NULL, na.rm = FALSE, target_dims = list(target_dims_obs, c(memb_dim, sdate_dim), c(memb_dim, sdate_dim)), - fun = .sbc , output_dims = c(memb_dim,sdate_dim), - na.rm = na.rm, ncores = ncores)$output1 + fun = .sbc, + output_dims = c(memb_dim, sdate_dim), + memb_dim = memb_dim, sdate_dim = sdate_dim, + na.rm = na.rm, ncores = ncores)$output1 } return(BiasCorrected) } -.sbc <- function(var_obs, var_exp , var_cor = NULL, na.rm = FALSE) { - nmembers <- dim(var_exp)[1] - ntime <- dim(var_exp)[2][] - +.sbc <- function(var_obs, var_exp , var_cor = NULL, memb_dim = 'member', + sdate_dim = 'sdate', na.rm = FALSE) { + + ntime <- dim(var_exp)[sdate_dim] corrected <- NA * var_exp - if (is.null(var_cor)){ - for (t in 1 : ntime) { + if (is.null(var_cor)) { + for (t in 1:ntime) { # parameters sd_obs <- sd(var_obs[-t], na.rm = na.rm) - sd_exp <- sd(var_exp[ , -t], na.rm = na.rm) - clim_exp <- mean(var_exp[ , -t], na.rm = na.rm) + sd_exp <- sd(var_exp[, -t], na.rm = na.rm) + clim_exp <- mean(var_exp[, -t], na.rm = na.rm) clim_obs <- mean(var_obs[-t], na.rm = na.rm) # bias corrected forecast - corrected[ , t] <- ((var_exp[ , t] - clim_exp) * (sd_obs / sd_exp)) + clim_obs - names(dim(corrected)) <- c('member', 'sdate') + corrected[, t] <- ((var_exp[, t] - clim_exp) * (sd_obs / sd_exp)) + clim_obs } } else { # parameters @@ -184,7 +236,7 @@ BiasCorrection <- function (exp, obs, exp_cor = NULL, na.rm = FALSE, # bias corrected forecast corrected <- ((var_cor - clim_exp) * (sd_obs / sd_exp)) + clim_obs - names(dim(corrected)) <- c('member') } + return(corrected) } diff --git a/man/BiasCorrection.Rd b/man/BiasCorrection.Rd index 3f8d72ce..a00c3461 100644 --- a/man/BiasCorrection.Rd +++ b/man/BiasCorrection.Rd @@ -11,32 +11,43 @@ BiasCorrection( na.rm = FALSE, memb_dim = "member", sdate_dim = "sdate", - ncores = 1 + ncores = NULL ) } \arguments{ -\item{exp}{a multidimensional array with named dimensions containing the seasonal forecast experiment data with at least 'member' and 'sdate' dimensions.} +\item{exp}{A multidimensional array with named dimensions containing the +seasonal forecast experiment data with at least 'member' and 'sdate' +dimensions.} -\item{obs}{a multidimensional array with named dimensions containing the observed data with at least 'sdate' dimension.} +\item{obs}{A multidimensional array with named dimensions containing the +observed data with at least 'sdate' dimension.} -\item{exp_cor}{a multidimensional array with named dimensions containing the seasonl forecast experiment to be corrected. If it is NULL, the 'exp' forecast will be corrected.} +\item{exp_cor}{A multidimensional array with named dimensions containing the +seasonl forecast experiment to be corrected. If it is NULL, the 'exp' +forecast will be corrected.} -\item{na.rm}{a logical value indicating whether missing values should be stripped before the computation proceeds, by default it is set to FALSE.} +\item{na.rm}{A logical value indicating whether missing values should be +stripped before the computation proceeds, by default it is set to FALSE.} -\item{memb_dim}{is a character string indicating the name of the member dimension. By default, it is set to 'member'.} +\item{memb_dim}{A character string indicating the name of the member +dimension. By default, it is set to 'member'.} -\item{sdate_dim}{is a character string indicating the name of the start date dimension. By default, it is set to 'sdate'.} +\item{sdate_dim}{A character string indicating the name of the start date +dimension. By default, it is set to 'sdate'.} -\item{ncores}{is an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one.} +\item{ncores}{An integer that indicates the number of cores for parallel +computations using multiApply function. The default value is NULL.} } \value{ -an object of class \code{s2dv_cube} containing the bias corrected forecasts in the element called \code{$data} with the same dimensions of the experimental data. +An array containing the bias corrected forecasts with the same +dimensions of the experimental data. } \description{ -This function applies the simple bias adjustment technique described in Torralba et al. (2017). The adjusted forecasts have an equivalent standard deviation and mean to that of the reference dataset. +This function applies the simple bias adjustment technique +described in Torralba et al. (2017). The adjusted forecasts have an equivalent +standard deviation and mean to that of the reference dataset. } \examples{ - # Example # Creation of sample s2dverification objects. These are not complete # s2dverification objects though. The Load function returns complete objects. @@ -48,7 +59,11 @@ a <- BiasCorrection(exp = mod1, obs = obs1) str(a) } \references{ -Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. Davis (2017). Seasonal climate prediction: a new source of information for the management of wind energy resources. Journal of Applied Meteorology and Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, EUPORIAS, NEWA, RESILIENCE, SPECS) +Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. +Davis (2017). Seasonal climate prediction: a new source of information for the +management of wind energy resources. Journal of Applied Meteorology and +Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, +EUPORIAS, NEWA, RESILIENCE, SPECS) } \author{ Verónica Torralba, \email{veronica.torralba@bsc.es} diff --git a/man/CST_BiasCorrection.Rd b/man/CST_BiasCorrection.Rd index 202e7651..45014f7d 100644 --- a/man/CST_BiasCorrection.Rd +++ b/man/CST_BiasCorrection.Rd @@ -11,32 +11,43 @@ CST_BiasCorrection( na.rm = FALSE, memb_dim = "member", sdate_dim = "sdate", - ncores = 1 + ncores = NULL ) } \arguments{ -\item{exp}{an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the seasonal forecast experiment data in the element named \code{$data}} +\item{exp}{An object of class \code{s2dv_cube} as returned by \code{CST_Load} +function, containing the seasonal forecast experiment data in the element +named \code{$data}} -\item{obs}{an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the observed data in the element named \code{$data}.} +\item{obs}{An object of class \code{s2dv_cube} as returned by \code{CST_Load} +function, containing the observed data in the element named \code{$data}.} -\item{exp_cor}{an object of class \code{s2dv_cube} as returned by \code{CST_Load} function, containing the seasonl forecast experiment to be corrected. If it is NULL, the 'exp' forecast will be corrected.} +\item{exp_cor}{An object of class \code{s2dv_cube} as returned by +\code{CST_Load} function, containing the seasonl forecast experiment to be +corrected. If it is NULL, the 'exp' forecast will be corrected.} -\item{na.rm}{a logical value indicating whether missing values should be stripped before the computation proceeds, by default it is set to FALSE.} +\item{na.rm}{A logical value indicating whether missing values should be +stripped before the computation proceeds, by default it is set to FALSE.} -\item{memb_dim}{is a character string indicating the name of the member dimension. By default, it is set to 'member'.} +\item{memb_dim}{A character string indicating the name of the member +dimension. By default, it is set to 'member'.} -\item{sdate_dim}{is a character string indicating the name of the start date dimension. By default, it is set to 'sdate'.} +\item{sdate_dim}{A character string indicating the name of the start date +dimension. By default, it is set to 'sdate'.} -\item{ncores}{is an integer that indicates the number of cores for parallel computations using multiApply function. The default value is one.} +\item{ncores}{An integer that indicates the number of cores for parallel +computations using multiApply function. The default value is NULL.} } \value{ -an object of class \code{s2dv_cube} containing the bias corrected forecasts in the element called \code{$data} with the same dimensions of the experimental data. +An object of class \code{s2dv_cube} containing the bias corrected +forecasts with the same dimensions of the experimental data. } \description{ -This function applies the simple bias adjustment technique described in Torralba et al. (2017). The adjusted forecasts have an equivalent standard deviation and mean to that of the reference dataset. +This function applies the simple bias adjustment technique +described in Torralba et al. (2017). The adjusted forecasts have an equivalent +standard deviation and mean to that of the reference dataset. } \examples{ - # Example # Creation of sample s2dverification objects. These are not complete # s2dverification objects though. The Load function returns complete objects. @@ -54,7 +65,11 @@ a <- CST_BiasCorrection(exp = exp, obs = obs) str(a) } \references{ -Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. Davis (2017). Seasonal climate prediction: a new source of information for the management of wind energy resources. Journal of Applied Meteorology and Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, EUPORIAS, NEWA, RESILIENCE, SPECS) +Torralba, V., F.J. Doblas-Reyes, D. MacLeod, I. Christel and M. +Davis (2017). Seasonal climate prediction: a new source of information for +the management of wind energy resources. Journal of Applied Meteorology and +Climatology, 56, 1231-1247, doi:10.1175/JAMC-D-16-0204.1. (CLIM4ENERGY, +EUPORIAS, NEWA, RESILIENCE, SPECS) } \author{ Verónica Torralba, \email{veronica.torralba@bsc.es} diff --git a/tests/testthat/test-CST_BiasCorrection.R b/tests/testthat/test-CST_BiasCorrection.R index 205aef44..1978ad66 100644 --- a/tests/testthat/test-CST_BiasCorrection.R +++ b/tests/testthat/test-CST_BiasCorrection.R @@ -1,60 +1,179 @@ -context("Generic tests") -test_that("Sanity checks", { +context("CSTools::CST_BiasCorrection tests") + +############################################## + +# dat1 +mod <- 1 : (1 * 3 * 4 * 5 * 6 * 7) +obs <- 1 : (1 * 1 * 4 * 5 * 6 * 7) +dim(mod) <- c(dataset = 1, member = 3, sdate = 4, ftime = 5, + lat = 6, lon = 7) +dim(obs) <- c(dataset = 1, member = 1, sdate = 4, ftime = 5, + lat = 6, lon = 7) +lon <- seq(0, 30, 5) +lat <- seq(0, 25, 5) +exp <- list(data = mod, lat = lat, lon = lon) +obs <- list(data = obs, lat = lat, lon = lon) +attr(exp, 'class') <- 's2dv_cube' +attr(obs, 'class') <- 's2dv_cube' + +exp1 <- list(data = array(1:20, dim = c(time = 20))) +class(exp1) <- 's2dv_cube' + +obs1 <- list(data = array(1:20, dim = c(time = 20))) +class(obs1) <- 's2dv_cube' + +exp1_2 <- list(data = array(1:20, dim = c(20))) +class(exp1_2) <- 's2dv_cube' + +obs1_2 <- list(data = array(1:20, dim = c(20))) +class(obs1_2) <- 's2dv_cube' + +exp_cor1 <- list(data = array(1:20, dim = c(20))) +class(exp_cor1) <- 's2dv_cube' + +# dat2 +exp2 <- exp +exp2$data[1, 2, 1, 1, 1, 1] <- NA +obs2 <- obs +obs2$data[1, 1, 2, 1, 1, 1] <- NA + +# dat3 +exp3 <- array(1:6, c(sdate = 3, member = 2)) +obs3 <- array(3:6, c(sdate = 3, member = 1)) +obs3_2 <- array(3:6, c(sdate = 3)) + +# dat4 +exp4 <- array(1:100, dim = c(time = 2, members = 5, lat = 2, lon = 5)) +obs4 <- array(1:200, dim = c(time = 2, members = 5, lat = 2, lon = 5)) +obs4_1 <- obs4 +obs4_1[1,1,1,1] <- NA + +############################################## + +test_that("1. Inpput checks", { + # s2dv_cube expect_error( CST_BiasCorrection(exp = 1), paste0("Parameter 'exp' and 'obs' must be of the class 's2dv_cube', ", - "as output by CSTools::CST_Load.")) - - mod1 <- 1 : (1 * 3 * 4 * 5 * 6 * 7) - obs1 <- 1 : (1 * 1 * 4 * 5 * 6 * 7) - dim(mod1) <- c(dataset = 1, member = 3, sdate = 4, ftime = 5, - lat = 6, lon = 7) - dim(obs1) <- c(dataset = 1, member = 1, sdate = 4, ftime = 5, - lat = 6, lon = 7) - lon <- seq(0, 30, 5) - lat <- seq(0, 25, 5) - exp <- list(data = mod1, lat = lat, lon = lon) - obs <- list(data = obs1, lat = lat, lon = lon) - attr(exp, 'class') <- 's2dv_cube' - attr(obs, 'class') <- 's2dv_cube' - - bc <- CST_BiasCorrection(exp = exp, obs = obs) - expect_equal(length(bc), 3) - expect_equal(dim(bc$data), - c(dataset = 1, member = 3, sdate = 4, ftime = 5, - lat = 6, lon = 7)) - expect_equal(bc$lat, lat) - expect_equal(bc$lon, lon) - - expect_error(CST_BiasCorrection(exp = exp, obs = exp), - paste0("The length of the dimension 'member' in the component 'data' ", - "of the parameter 'obs' must be equal to 1.")) - - exp2 <- exp - exp2$data[1, 2, 1, 1, 1, 1] <- NA - expect_warning(CST_BiasCorrection(exp = exp2, obs = obs), - "Parameter 'exp' contains NA values.") - - obs2 <- obs - obs2$data[1, 1, 2, 1, 1, 1] <- NA - expect_warning(CST_BiasCorrection(exp = exp, obs = obs2), - "Parameter 'obs' contains NA values.") - - expect_warning(CST_BiasCorrection(exp = exp2, obs = obs2), - "Parameter 'obs' contains NA values", "Parameter 'exp' contains NA values.") - - hinc <- array(1:6, c(sdate = 3, member = 2)) - obs <- array(3:6, c(sdate = 3, member = 1)) - res <- round(BiasCorrection(exp = hinc, obs = obs, exp_cor = hinc), 2) - expect_equal(res, array(c(2.66, 4.27, 3.2, 4.8, 3.73, 5.34), c(member = 2, sdate = 3))) - - # if obs doesn't have memb_dim it works the same: - hinc <- array(1:6, c(sdate = 3, member = 2)) - obs <- array(3:6, c(sdate = 3)) - res <- round(BiasCorrection(exp = hinc, obs = obs, exp_cor = hinc), 2) - expect_equal(res, array(c(2.66, 4.27, 3.2, 4.8, 3.73, 5.34), c(member = 2, sdate = 3))) - - + "as output by CSTools::CST_Load.") + ) + expect_error( + CST_BiasCorrection(exp = exp, obs = 1), + paste0("Parameter 'exp' and 'obs' must be of the class 's2dv_cube', ", + "as output by CSTools::CST_Load.") + ) + expect_error( + CST_BiasCorrection(exp = exp1), + 'argument "obs" is missing, with no default' + ) + expect_error( + CST_BiasCorrection(exp = exp1, obs = obs1, exp_cor = 1), + paste0("Parameter 'exp_cor' must be of the class 's2dv_cube', as output ", + "by CSTools::CST_Load.") + ) + # exp and obs + expect_error( + CST_BiasCorrection(exp = exp1_2, obs = obs1), + "Parameter 'exp' must have dimension names." + ) + expect_error( + CST_BiasCorrection(exp = exp1, obs = obs1_2), + "Parameter 'obs' must have dimension names." + ) + expect_warning( + CST_BiasCorrection(exp = exp2, obs = obs2), + "Parameter 'exp' contains NA values." + ) + expect_warning( + CST_BiasCorrection(exp = exp, obs = obs2), + "Parameter 'obs' contains NA values." + ) + expect_warning( + CST_BiasCorrection(exp = exp2, obs = obs2), + "Parameter 'obs' contains NA values", + "Parameter 'exp' contains NA values." + ) + # exp_cor + expect_error( + CST_BiasCorrection(exp = exp1, obs = obs1, exp_cor = exp_cor1, sdate_dim = 'time'), + "Parameter 'exp_cor' must have dimension names." + ) + # sdate_dim, member_dim + expect_error( + CST_BiasCorrection(exp = exp1, obs = obs1, sdate_dim = 1), + paste0("Parameter 'sdate_dim' must be a character string.") + ) + expect_error( + CST_BiasCorrection(exp = exp, obs = obs, sdate_dim = 'time'), + paste0("Parameter 'sdate_dim' is not found in 'exp' or 'obs' dimension.") + ) + expect_error( + BiasCorrection(exp = array(1:20, dim = c(time = 1, member = 1)), + obs = array(1:20, dim = c(time = 2, member = 1)), sdate_dim = 'time'), + paste0("Parameter 'exp' must have dimension length of 'sdate_dim' bigger than 1.") + ) + expect_error( + CST_BiasCorrection(exp = exp1, obs = obs1, sdate_dim = 'time'), + paste0("Parameter 'exp' requires 'sdate_dim' and 'member_dim' dimensions.") + ) + ## na.rm + expect_warning( + CST_BiasCorrection(exp = exp, obs = obs, na.rm = 1), + "Paramater 'na.rm' must be a logical, it has been set to FALSE." + ) + expect_warning( + CST_BiasCorrection(exp = exp, obs = obs, na.rm = c(T,F)), + "Paramter 'na.rm' has length greater than 1, and only the fist element is used." + ) + # ncores + expect_error( + CST_BiasCorrection(exp = exp, obs = obs, ncores = TRUE), + "Parameter 'ncores' must be either NULL or a positive integer." + ) +}) +############################################## +test_that("2. Output checks: dat1", { + bc <- CST_BiasCorrection(exp = exp, obs = obs) + expect_equal( + length(bc), + 3 + ) + expect_equal( + dim(bc$data), + c(dataset = 1, member = 3, sdate = 4, ftime = 5, lat = 6, lon = 7) + ) + expect_equal( + bc$lat, + lat + ) + expect_equal( + bc$lon, + lon + ) + expect_equal( + round(BiasCorrection(exp = exp3, obs = obs3, exp_cor = exp3), 2), + array(c(2.66, 4.27, 3.2, 4.8, 3.73, 5.34), c(member = 2, sdate = 3)) + ) + expect_equal( + round(BiasCorrection(exp = exp3, obs = obs3_2, exp_cor = exp3), 2), + array(c(2.66, 4.27, 3.2, 4.8, 3.73, 5.34), c(member = 2, sdate = 3)) + ) + expect_equal( + dim(BiasCorrection(exp = exp4, obs = obs4, sdate_dim = 'time', memb_dim = 'members')), + c(members = 5, time = 2, lat = 2, lon = 5) + ) + suppressWarnings( + expect_equal( + sum(is.na(BiasCorrection(exp = exp4, obs = obs4_1, sdate_dim = 'time', memb_dim = 'members', na.rm = TRUE))), + 0 + ) + ) + suppressWarnings( + expect_equal( + sum(is.na(BiasCorrection(exp = exp4, obs = obs4_1, sdate_dim = 'time', memb_dim = 'members', na.rm = FALSE))), + 5 + ) + ) }) -- GitLab From 844dcc13b07c03dfaf86c1d38f35203e708a8bab Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Thu, 13 Oct 2022 18:30:02 +0200 Subject: [PATCH 2/5] Correct output CST_BiasCorrection() --- R/CST_BiasCorrection.R | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/R/CST_BiasCorrection.R b/R/CST_BiasCorrection.R index eb6ca04a..e176640d 100644 --- a/R/CST_BiasCorrection.R +++ b/R/CST_BiasCorrection.R @@ -71,10 +71,16 @@ CST_BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE, if (is.null(exp_cor)) { exp$data <- BiasCorrected + exp$Datasets <- c(exp$Datasets, obs$Datasets) + exp$source_files <- c(exp$source_files, obs$source_files) + return(exp) } else { exp_cor$data <- BiasCorrected + exp_cor$Datasets <- c(exp_cor$Datasets, exp$Datasets, obs$Datasets) + exp_cor$source_files <- c(exp_cor$source_files, exp$source_files, obs$source_files) + return(exp_cor) } } -- GitLab From 7fed642e5e58afbfe2bb9c2fb7f718e34cf3541a Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Fri, 14 Oct 2022 18:34:17 +0200 Subject: [PATCH 3/5] Remove sdate_dim and memb_dim from .sbc and add restriction for members --- R/CST_BiasCorrection.R | 25 +++++++++++++----------- tests/testthat/test-CST_BiasCorrection.R | 21 ++++++++++++-------- 2 files changed, 27 insertions(+), 19 deletions(-) diff --git a/R/CST_BiasCorrection.R b/R/CST_BiasCorrection.R index e176640d..3d2c5d2e 100644 --- a/R/CST_BiasCorrection.R +++ b/R/CST_BiasCorrection.R @@ -170,7 +170,13 @@ BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE, stop("Parameter 'exp' must have dimension length of 'sdate_dim' bigger than 1.") } if (!all(c(memb_dim, sdate_dim) %in% expdims)) { - stop("Parameter 'exp' requires 'sdate_dim' and 'member_dim' dimensions.") + stop("Parameter 'exp' requires 'sdate_dim' and 'memb_dim' dimensions.") + } + if (memb_dim %in% obsdims) { + if (dim(obs)[memb_dim] != 1) { + stop("The length of the dimension 'member' in the component 'data' ", + "of the parameter 'obs' must be equal to 1.") + } } ## na.rm if (!is.logical(na.rm)) { @@ -189,17 +195,16 @@ BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE, } } - if (memb_dim %in% obsdims) { - target_dims_obs <- c(memb_dim, sdate_dim) - } else { - target_dims_obs <- c(sdate_dim) + target_dims_obs <- sdate_dim + if (memb_dim %in% names(dim(obs))) { + target_dims_obs <- c(memb_dim, target_dims_obs) } if (is.null(exp_cor)) { BiasCorrected <- Apply(data = list(var_obs = obs, var_exp = exp), target_dims = list(target_dims_obs, c(memb_dim, sdate_dim)), - fun = .sbc, memb_dim = memb_dim, sdate_dim = sdate_dim, + fun = .sbc, na.rm = na.rm, ncores = ncores)$output1 } else { BiasCorrected <- Apply(data = list(var_obs = obs, @@ -209,17 +214,15 @@ BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE, c(memb_dim, sdate_dim), c(memb_dim, sdate_dim)), fun = .sbc, - output_dims = c(memb_dim, sdate_dim), - memb_dim = memb_dim, sdate_dim = sdate_dim, + output_dims = c(memb_dim, sdate_dim), na.rm = na.rm, ncores = ncores)$output1 } return(BiasCorrected) } -.sbc <- function(var_obs, var_exp , var_cor = NULL, memb_dim = 'member', - sdate_dim = 'sdate', na.rm = FALSE) { +.sbc <- function(var_obs, var_exp , var_cor = NULL, na.rm = FALSE) { - ntime <- dim(var_exp)[sdate_dim] + ntime <- dim(var_exp)[2] corrected <- NA * var_exp if (is.null(var_cor)) { diff --git a/tests/testthat/test-CST_BiasCorrection.R b/tests/testthat/test-CST_BiasCorrection.R index 1978ad66..5b40fb3d 100644 --- a/tests/testthat/test-CST_BiasCorrection.R +++ b/tests/testthat/test-CST_BiasCorrection.R @@ -39,12 +39,13 @@ obs2$data[1, 1, 2, 1, 1, 1] <- NA # dat3 exp3 <- array(1:6, c(sdate = 3, member = 2)) -obs3 <- array(3:6, c(sdate = 3, member = 1)) -obs3_2 <- array(3:6, c(sdate = 3)) +obs3 <- array(3:6, c(sdate = 3, member = 1)) +obs3_2 <- array(3:6, c(sdate = 3)) +obs3_3 <- array(3:6, c(sdate = 3, member = 2)) # dat4 -exp4 <- array(1:100, dim = c(time = 2, members = 5, lat = 2, lon = 5)) -obs4 <- array(1:200, dim = c(time = 2, members = 5, lat = 2, lon = 5)) +exp4 <- array(1:100, dim = c(time = 5, members = 5, lat = 2, lon = 5)) +obs4 <- array(1:200, dim = c(time = 5, members = 1, lat = 2, lon = 5)) obs4_1 <- obs4 obs4_1[1,1,1,1] <- NA @@ -98,7 +99,7 @@ test_that("1. Inpput checks", { CST_BiasCorrection(exp = exp1, obs = obs1, exp_cor = exp_cor1, sdate_dim = 'time'), "Parameter 'exp_cor' must have dimension names." ) - # sdate_dim, member_dim + # sdate_dim, memb_dim expect_error( CST_BiasCorrection(exp = exp1, obs = obs1, sdate_dim = 1), paste0("Parameter 'sdate_dim' must be a character string.") @@ -114,7 +115,11 @@ test_that("1. Inpput checks", { ) expect_error( CST_BiasCorrection(exp = exp1, obs = obs1, sdate_dim = 'time'), - paste0("Parameter 'exp' requires 'sdate_dim' and 'member_dim' dimensions.") + paste0("Parameter 'exp' requires 'sdate_dim' and 'memb_dim' dimensions.") + ) + expect_error( + BiasCorrection(exp = exp3, obs = obs3_3), + paste0("The length of the dimension 'member' in the component 'data' of the parameter 'obs' must be equal to 1.") ) ## na.rm expect_warning( @@ -162,7 +167,7 @@ test_that("2. Output checks: dat1", { ) expect_equal( dim(BiasCorrection(exp = exp4, obs = obs4, sdate_dim = 'time', memb_dim = 'members')), - c(members = 5, time = 2, lat = 2, lon = 5) + c(members = 5, time = 5, lat = 2, lon = 5) ) suppressWarnings( expect_equal( @@ -173,7 +178,7 @@ test_that("2. Output checks: dat1", { suppressWarnings( expect_equal( sum(is.na(BiasCorrection(exp = exp4, obs = obs4_1, sdate_dim = 'time', memb_dim = 'members', na.rm = FALSE))), - 5 + 20 ) ) }) -- GitLab From 0bbc6c51f1b1d5f43df9dbca9662a92aeaf55450 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Fri, 14 Oct 2022 18:36:45 +0200 Subject: [PATCH 4/5] Correct word --- R/CST_BiasCorrection.R | 4 ++-- tests/testthat/test-CST_BiasCorrection.R | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/R/CST_BiasCorrection.R b/R/CST_BiasCorrection.R index 3d2c5d2e..e829f8ba 100644 --- a/R/CST_BiasCorrection.R +++ b/R/CST_BiasCorrection.R @@ -174,8 +174,8 @@ BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE, } if (memb_dim %in% obsdims) { if (dim(obs)[memb_dim] != 1) { - stop("The length of the dimension 'member' in the component 'data' ", - "of the parameter 'obs' must be equal to 1.") + stop("The length of the dimension 'memb_dim' in the component 'data' ", + "of the parameter 'obs' must be equal to 1.") } } ## na.rm diff --git a/tests/testthat/test-CST_BiasCorrection.R b/tests/testthat/test-CST_BiasCorrection.R index 5b40fb3d..2ded6252 100644 --- a/tests/testthat/test-CST_BiasCorrection.R +++ b/tests/testthat/test-CST_BiasCorrection.R @@ -119,7 +119,7 @@ test_that("1. Inpput checks", { ) expect_error( BiasCorrection(exp = exp3, obs = obs3_3), - paste0("The length of the dimension 'member' in the component 'data' of the parameter 'obs' must be equal to 1.") + paste0("The length of the dimension 'memb_dim' in the component 'data' of the parameter 'obs' must be equal to 1.") ) ## na.rm expect_warning( -- GitLab From eed2b7ac1dab837e1db555151a703b97e8a5a642 Mon Sep 17 00:00:00 2001 From: Eva Rifa Date: Mon, 17 Oct 2022 09:27:15 +0200 Subject: [PATCH 5/5] Correct output text for memb_dim in obs --- R/CST_BiasCorrection.R | 3 +-- tests/testthat/test-CST_BiasCorrection.R | 2 +- 2 files changed, 2 insertions(+), 3 deletions(-) diff --git a/R/CST_BiasCorrection.R b/R/CST_BiasCorrection.R index e829f8ba..c3097070 100644 --- a/R/CST_BiasCorrection.R +++ b/R/CST_BiasCorrection.R @@ -174,8 +174,7 @@ BiasCorrection <- function(exp, obs, exp_cor = NULL, na.rm = FALSE, } if (memb_dim %in% obsdims) { if (dim(obs)[memb_dim] != 1) { - stop("The length of the dimension 'memb_dim' in the component 'data' ", - "of the parameter 'obs' must be equal to 1.") + stop("If parameter 'obs' has dimension 'memb_dim' its length must be equal to 1.") } } ## na.rm diff --git a/tests/testthat/test-CST_BiasCorrection.R b/tests/testthat/test-CST_BiasCorrection.R index 2ded6252..826fcf12 100644 --- a/tests/testthat/test-CST_BiasCorrection.R +++ b/tests/testthat/test-CST_BiasCorrection.R @@ -119,7 +119,7 @@ test_that("1. Inpput checks", { ) expect_error( BiasCorrection(exp = exp3, obs = obs3_3), - paste0("The length of the dimension 'memb_dim' in the component 'data' of the parameter 'obs' must be equal to 1.") + paste0("If parameter 'obs' has dimension 'memb_dim' its length must be equal to 1.") ) ## na.rm expect_warning( -- GitLab