-
Notifications
You must be signed in to change notification settings - Fork 2
/
verifyRMSE.R
72 lines (62 loc) · 2.48 KB
/
verifyRMSE.R
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
# "`-''-/").___..--''"`-._
# (`6_ 6 ) `-. ( ).`-.__.`) WE ARE ...
# (_Y_.)' ._ ) `._ `. ``-..-' PENN STATE!
# _ ..`--'_..-_/ /--'_.' ,'
# (il),-'' (li),' ((!.-'
#
#
# Author: Guido Cervone (cervone@psu.edu), Martina Calovi (mxc895@psu.edu), Laura Clemente-Harding (laura@psu.edu)
# Geoinformatics and Earth Observation Laboratory (http://geolab.psu.edu)
# Department of Geography and Institute for CyberScience
# The Pennsylvania State University
#
#' RAnEnExtra::verifyRMSE
#'
#' RAnEnExtra::verifyRMSE calculates root mean square errors.
#'
#' RMSE ^ 2 = CRMSE ^ 2 + Bias ^ 2
#'
#' To set the number of cores to use when parallel is used,
#' `options(mc.cores = 8)`.
#'
#' @details Bootstrap confidence interval is defaulted to 0.95.
#' To change this, use `options(RAnEnExtra_boot_conf = 0.9)`.
#'
#' @author Guido Cervone \email{cervone@@psu.edu}
#' @author Martina Calovi \email{mxc895@@psu.edu}
#' @author Laura Clemente-Harding \email{laura@@psu.edu}
#'
#' @param anen.ver A 4-dimensional array. This array is usually created from the `value` column of
#' the `analogs` member in the results of `RAnEn::generateAnalogs`. The dimensions should be
#' `[stations, times, lead times, members]`.
#' @param obs.ver A 3-dimensional array. The dimensions should be `[stations, times, lead times]`.
#' You can generate the array using `RAnEn::alignObservations`.
#' @param boot Whether to use bootstrap.
#' @param R The number of bootstrap replicates. Used by the function `boot::boot`.
#' @param na.rm Whether to remove NA values.
#' @param parallel Whether to turn on parallel processing.
#'
#' @md
#' @export
verifyRMSE <- function(anen.ver, obs.ver, boot=F, R=1000, na.rm=T, parallel = F) {
stopifnot(length(dim(anen.ver)) == 4)
stopifnot(length(dim(obs.ver)) == 3)
if ( !identical(dim(anen.ver)[1:3], dim(obs.ver)[1:3]) ) {
cat("Error: Observations and Forecasts have incompatible dimensions.\n")
return(NULL)
}
# dim(obs.ver) [1] 105 362 288
# dim(anen.ver) [1] 105 362 288 21
obs <- matrix(obs.ver,ncol=dim(obs.ver)[3])
anen <- anen.mean(anen.ver, na.rm, parallel = parallel)
err <- (anen - obs)^2
if ( boot == F) {
err.tot <- sqrt( mean(err, na.rm=na.rm) )
err.mean <- sqrt( apply( err, 2, mean, na.rm=na.rm) )
return(list(mean=err.tot, flt=err.mean, mat=err))
} else {
boot <- sqrt( apply( err, 2, boot.fun.ver, R))
boot.tot <- mean(boot[1,],na.rm=na.rm)
return(list(mean=boot.tot, flt=boot[1,], mat=boot))
}
}