Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
#'Multiple Metrics applied in Multiple Model Anomalies
#'
#'@author Mishra Niti, \email{niti.mishra@bsc.es}
#'@description This function calculates correlation (Anomaly Correlation Coefficient; ACC), root mean square error (RMS) and the root mean square error skill score (RMSSS) of individual anomaly models and multi-models mean (if desired) with the observations.
#'
#'@param data an s2dverification object (list) giving as output of \code{Load} function from S2dverification package.
#'@param metric a character string giving the metric for computing the maximum skill. This must be one of the strings 'correlation', 'rms' or 'rmsss.
#'@param multimodel a logical value indicating wether a Multi-Model Mean should be computed.
#'@param names a character vector indicating the names of the models to be compared.
#'
#'@return A list of 3:
#'\itemize{
#'\item\code{$ano} {A list of 2:
#' \itemize{
#' \item\code{$ano_exp} {an array with the same dimensions as the \code{mod} label in the input \code{data}.}
#' \item\code{$ano_exp} {an array with the same dimensions as the \code{obs} label in the input \code{data}.}}}
#'\item\code{$metric} {An array with the same dimensions as the \code{mod} label in the input \code{data} without the 'ftime' dimension containing the correlations coefficients. if \code{multimodel} is TRUE, the greatest first dimension correspons. The third dimension contains the statistics selected: for \code{correlation} and \code{rms} metrics, the third dimension of length four corresponds to the lower limit of the \code{siglev}\% confidence interval, the statistics itselfs, the upper limit of the \code{siglev}\% confidence interval and the \code{siglev}\% significance level, while for the \code{rmsss}, the third dimension of length two corresponds to the statistics itselfs and the p-value of the one-sided Fisher test with Ho: RMSSS = 0.}
#'\item\code{$names} {A vector containing the names of the models.}}
#'@seealso \code{\link[s2dverification]{Corr}} and \code{link[s2dverification]{Load}}
#'@references
#' \url{http://link.springer.com/10.1007/s00382-018-4404-z}
#'
#'@import s2dverification
#'@import stats
#'@examples
#'mod <- 1 : (2 * 3 * 4 * 5 * 6 * 7)
#'dim(mod) <- c(dataset = 2, member = 3, sdate = 4, ftime = 5, lat = 6, lon = 7)
#'obs <- 1 : (1 * 1 * 4 * 5 * 6 * 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)
#'dat <- list(mod = mod, obs = obs, lat = lat, lon = lon)
#'a <- AnoMultiMetric(dat)
#'str(a)
#'@export
AnoMultiMetric <- function(data, metric = "correlation", multimodel = TRUE, names = NULL) {
if (!is.list(data)) {
stop("Parameter 'data' must be a list as output of Load function
from s2dverification package.")
}
if (!(all(c('mod', 'obs', 'lon', 'lat') %in% names(data)) && length(data) >= 4)) {
stop("Parameter 'data' must contain at least 4 elements 'mod', 'obs',
'lon', 'lat'.")
}
if (!is.null(names(dim(data$mod))) & !is.null(names(dim(data$obs)))) {
if (all(names(dim(data$mod)) %in% names(dim(data$obs)))) {
dimnames <- names(dim(data$mod))
} else {
stop("Dimension names of elements 'mod' and 'obs' from parameter 'data'
should have the same name dimmension.")
}
} else {
stop("Elements 'mod' and 'obs' from parameter 'data' should have
dimmension names.")
}
if (is.null(names)) {
if (!is.null(names(data$Datasets$exp))) {
names <- names(data$Datasets$exp)
} else {
names <- paste0('mod', 1 : dim(data$mod)[1])
}
}
if (!is.logical(multimodel)) {
stop("Parameter 'multimodel' must be a logical value.")
}
if (length(multimodel) > 1) {
multimodel <- multimodel[1]
warning("Parameter 'multimodel' has length > 1 and only the first element will be used.")
}
if (length(metric) > 1) {
metric <- metric[1]
warning("Parameter 'multimodel' has length > 1 and only the first element will be used.")
}
ano <- Ano_CrossValid(var_exp = data$mod, var_obs = data$obs)
# seasonal average of anomalies per model
AvgExp <- MeanListDim(ano$ano_exp, narm = T, c(2, 4))
AvgObs <- MeanListDim(ano$ano_obs, narm = T, c(2, 4))
# indv model correlation
if (metric == 'correlation') {
corr <- Corr(AvgExp, AvgObs, posloop = 1, poscor = 2)
} else if (metric == 'rms') {
corr <- RMS(AvgExp, AvgObs, posloop = 1, posRMS = 2)
} else if (metric == 'rmsss') {
corr <- RMSSS(AvgExp, AvgObs, posloop = 1, posRMS = 2)
} else {
stop("Parameter 'metric' must be a character string indicating one of the
options: 'correlation', 'rms' or 'rmse'.")
}
if (multimodel == TRUE) {
# seasonal avg of anomalies for multi-model
AvgExp_MMM <- MeanListDim(AvgExp, narm = T, 1)
AvgObs_MMM <- MeanListDim(AvgObs, narm = T, 1)
# multi model correlation
if (metric == 'correlation') {
corr_MMM <- Corr(var_exp = InsertDim(AvgExp_MMM, 1, 1),
var_obs = InsertDim(AvgObs_MMM, 1, 1),
posloop = 1, poscor = 2)
} else if (metric == 'rms') {
corr_MMM <- RMS(var_exp = InsertDim(AvgExp_MMM, 1, 1),
var_obs = InsertDim(AvgObs_MMM, 1, 1),
posloop = 1, posRMS = 2)
} else if (metric == 'rmsss') {
corr_MMM <- RMSSS(var_exp = InsertDim(AvgExp_MMM, 1, 1),
var_obs = InsertDim(AvgObs_MMM, 1, 1),
posloop = 1, posRMS = 2)
}
corr <- abind::abind(corr, corr_MMM, along = 1)
names <- c(names, 'MMM')
}
names(dim(corr))[3] <- 'statistics'
return(list(ano = ano, metric = corr, names = names))
}