Skip to content

GitLab

Commit 595a2e0d authored by Delaigue Olivier's avatar Delaigue Olivier
Browse files

v1.0.13.1 ErrorCrit_KGE cleaned

parent 0a005c7c
Show whitespace changes
Inline Side-by-side
DESCRIPTION
View file @ 595a2e0d
Package: airGR Package: airGR
Type: Package Type: Package
Title: Suite of GR Hydrological Models for Precipitation-Runoff Modelling Title: Suite of GR Hydrological Models for Precipitation-Runoff Modelling
Version: 1.0.13.0 Version: 1.0.13.1
Date: 2018-08-29 Date: 2018-08-29
Authors@R: c( Authors@R: c(
person("Laurent", "Coron", role = c("aut", "trl"), comment = c(ORCID = "0000-0002-1503-6204")), person("Laurent", "Coron", role = c("aut", "trl"), comment = c(ORCID = "0000-0002-1503-6204")),
......
NEWS.rmd
View file @ 595a2e0d
...@@ -14,7 +14,7 @@ output: ...@@ -14,7 +14,7 @@ output:
### 1.0.13.0 Release Notes (2018-08-29) ### 1.0.13.1 Release Notes (2018-08-29)
#### Deprectated and defunct #### Deprectated and defunct
......
R/ErrorCrit_KGE.R
View file @ 595a2e0d
ErrorCrit_KGE <- function(InputsCrit,OutputsModel, warnings = TRUE, verbose = TRUE){ ErrorCrit_KGE <- function(InputsCrit, OutputsModel, warnings = TRUE, verbose = TRUE) {
##Arguments_check________________________________ ##Arguments_check________________________________
if(inherits(InputsCrit,"InputsCrit")==FALSE){ stop("InputsCrit must be of class 'InputsCrit' \n"); return(NULL); } if (inherits(InputsCrit, "InputsCrit") == FALSE) {
if(inherits(OutputsModel,"OutputsModel")==FALSE){ stop("OutputsModel must be of class 'OutputsModel' \n"); return(NULL); } stop("InputsCrit must be of class 'InputsCrit' \n")
return(NULL)
}
if (inherits(OutputsModel, "OutputsModel") == FALSE) {
stop("OutputsModel must be of class 'OutputsModel' \n")
return(NULL)
}
##Initialisation_________________________________ ##Initialisation_________________________________
CritName <- NA; CritName <- NA
if(InputsCrit$transfo=="" ){ CritName <- "KGE[Q]" ; } if (InputsCrit$transfo == "") {
if(InputsCrit$transfo=="sqrt"){ CritName <- "KGE[sqrt(Q)]"; } CritName <- "KGE[Q]"
if(InputsCrit$transfo=="log" ){ CritName <- "KGE[log(Q)]" ; } }
if(InputsCrit$transfo=="inv" ){ CritName <- "KGE[1/Q]" ; } if (InputsCrit$transfo == "sqrt") {
if(InputsCrit$transfo=="sort"){ CritName <- "KGE[sort(Q)]"; } CritName <- "KGE[sqrt(Q)]"
CritValue <- NA; }
CritBestValue <- +1; if (InputsCrit$transfo == "log") {
Multiplier <- -1; ### must be equal to -1 or +1 only CritName <- "KGE[log(Q)]"
}
if (InputsCrit$transfo == "inv") {
CritName <- "KGE[1/Q]"
}
if (InputsCrit$transfo == "sort") {
CritName <- "KGE[sort(Q)]"
}
CritValue <- NA
CritBestValue <- +1
Multiplier <- -1
### must be equal to -1 or +1 only
##Data_preparation_______________________________
VarObs <- InputsCrit$Qobs
VarObs[!InputsCrit$BoolCrit] <- NA
VarSim <- OutputsModel$Qsim
VarSim[!InputsCrit$BoolCrit] <- NA
##Data_preparation_______________________________
VarObs <- InputsCrit$Qobs ; VarObs[!InputsCrit$BoolCrit] <- NA;
VarSim <- OutputsModel$Qsim; VarSim[!InputsCrit$BoolCrit] <- NA;
##Data_transformation ##Data_transformation
if("Ind_zeroes" %in% names(InputsCrit) & "epsilon" %in% names(InputsCrit)){ if(length(InputsCrit$Ind_zeroes)>0){ if ("Ind_zeroes" %in% names(InputsCrit) & "epsilon" %in% names(InputsCrit)) {
VarObs <- VarObs + InputsCrit$epsilon; if (length(InputsCrit$Ind_zeroes) > 0) {
VarSim <- VarSim + InputsCrit$epsilon; VarObs <- VarObs + InputsCrit$epsilon
} } VarSim <- VarSim + InputsCrit$epsilon
if(InputsCrit$transfo=="sqrt"){ VarObs <- sqrt(VarObs); VarSim <- sqrt(VarSim); } }
if(InputsCrit$transfo=="log" ){ VarObs <- log(VarObs) ; VarSim <- log(VarSim) ; VarSim[VarSim < -1E100] <- NA; } }
if(InputsCrit$transfo=="inv" ){ VarObs <- 1/VarObs ; VarSim <- 1/VarSim ; VarSim[abs(VarSim) > 1E+100] <- NA; } if (InputsCrit$transfo == "sqrt") {
if(InputsCrit$transfo=="sort"){ VarObs <- sqrt(VarObs)
VarSim <- sqrt(VarSim)
}
if (InputsCrit$transfo == "log") {
VarObs <- log(VarObs)
VarSim <- log(VarSim)
VarSim[VarSim < -1e100] <- NA
}
if (InputsCrit$transfo == "inv") {
VarObs <- 1 / VarObs
VarSim <- 1 / VarSim
VarSim[abs(VarSim) > 1e+100] <- NA
}
if (InputsCrit$transfo == "sort") {
VarSim[is.na(VarObs)] <- NA VarSim[is.na(VarObs)] <- NA
VarSim <- sort(VarSim, na.last = TRUE) VarSim <- sort(VarSim, na.last = TRUE)
VarObs <- sort(VarObs, na.last = TRUE) VarObs <- sort(VarObs, na.last = TRUE)
InputsCrit$BoolCrit <- sort(InputsCrit$BoolCrit, decreasing = TRUE) InputsCrit$BoolCrit <- sort(InputsCrit$BoolCrit, decreasing = TRUE)
} }
##TS_ignore ##TS_ignore
TS_ignore <- !is.finite(VarObs) | !is.finite(VarSim) | !InputsCrit$BoolCrit ; TS_ignore <- !is.finite(VarObs) | !is.finite(VarSim) | !InputsCrit$BoolCrit
Ind_TS_ignore <- which(TS_ignore); if(length(Ind_TS_ignore)==0){ Ind_TS_ignore <- NULL; } Ind_TS_ignore <- which(TS_ignore)
if(sum(!TS_ignore)==0){ OutputsCrit <- list(NA); names(OutputsCrit) <- c("CritValue"); return(OutputsCrit); } if (length(Ind_TS_ignore) == 0) {
if(sum(!TS_ignore)==1){ OutputsCrit <- list(NA); names(OutputsCrit) <- c("CritValue"); return(OutputsCrit); } ### to avoid a problem in standard deviation computation Ind_TS_ignore <- NULL
if(inherits(OutputsModel,"hourly" )){ WarningTS <- 365; } }
if(inherits(OutputsModel,"daily" )){ WarningTS <- 365; } if (sum(!TS_ignore) == 0) {
if(inherits(OutputsModel,"monthly")){ WarningTS <- 12; } OutputsCrit <- list(NA)
if(inherits(OutputsModel,"yearly" )){ WarningTS <- 3; } names(OutputsCrit) <- c("CritValue")
if(sum(!TS_ignore)<WarningTS & warnings){ warning("\t criterion computed on less than ", WarningTS, " time-steps ") } return(OutputsCrit)
}
if (sum(!TS_ignore) == 1) {
OutputsCrit <- list(NA)
names(OutputsCrit) <- c("CritValue")
return(OutputsCrit)
} ### to avoid a problem in standard deviation computation
if (inherits(OutputsModel, "hourly")) {
WarningTS <- 365
}
if (inherits(OutputsModel, "daily")) {
WarningTS <- 365
}
if (inherits(OutputsModel, "monthly")) {
WarningTS <- 12
}
if (inherits(OutputsModel, "yearly")) {
WarningTS <- 3
}
if (sum(!TS_ignore) < WarningTS & warnings) {
warning("\t criterion computed on less than ", WarningTS, " time-steps ")
}
##Other_variables_preparation ##Other_variables_preparation
meanVarObs <- mean(VarObs[!TS_ignore]); meanVarObs <- mean(VarObs[!TS_ignore])
meanVarSim <- mean(VarSim[!TS_ignore]); meanVarSim <- mean(VarSim[!TS_ignore])
iCrit <- 0 iCrit <- 0
SubCritPrint <- NULL SubCritPrint <- NULL
SubCritNames <- NULL SubCritNames <- NULL
SubCritValues <- NULL SubCritValues <- NULL
##SubErrorCrit_____KGE_rPearson__________________
##SubErrorCrit_____KGE_rPearson__________________ iCrit <- iCrit + 1
iCrit <- iCrit+1; SubCritPrint[iCrit] <- paste(CritName, " cor(sim, obs, \"pearson\") =", sep = "")
SubCritPrint[iCrit] <- paste(CritName," cor(sim, obs, \"pearson\") =", sep = "") SubCritValues[iCrit] <- NA
SubCritValues[iCrit] <- NA;
SubCritNames[iCrit] <- "r" SubCritNames[iCrit] <- "r"
Numer <- sum( (VarObs[!TS_ignore]-meanVarObs)*(VarSim[!TS_ignore]-meanVarSim) );
Deno1 <- sqrt( sum((VarObs[!TS_ignore]-meanVarObs)^2) ); Numer <- sum((VarObs[!TS_ignore] - meanVarObs) * (VarSim[!TS_ignore] - meanVarSim))
Deno2 <- sqrt( sum((VarSim[!TS_ignore]-meanVarSim)^2) ); Deno1 <- sqrt(sum((VarObs[!TS_ignore] - meanVarObs) ^ 2))
if(Numer==0){ if(Deno1==0 & Deno2==0){ Crit <- 1; } else { Crit <- 0; } Deno2 <- sqrt(sum((VarSim[!TS_ignore] - meanVarSim) ^ 2))
} else { Crit <- Numer/(Deno1*Deno2); }
if(is.numeric(Crit) & is.finite(Crit)){ SubCritValues[iCrit] <- Crit; } if (Numer == 0) {
if (Deno1 == 0 & Deno2 == 0) {
Crit <- 1
##SubErrorCrit_____KGE_alpha_____________________ } else {
iCrit <- iCrit+1; Crit <- 0
SubCritPrint[iCrit] <- paste(CritName," sd(sim)/sd(obs) =", sep = "") }
SubCritValues[iCrit] <- NA; } else {
Crit <- Numer / (Deno1 * Deno2)
}
if (is.numeric(Crit) & is.finite(Crit)) {
SubCritValues[iCrit] <- Crit
}
##SubErrorCrit_____KGE_alpha_____________________
iCrit <- iCrit + 1
SubCritPrint[iCrit] <- paste(CritName, " sd(sim)/sd(obs) =", sep = "")
SubCritValues[iCrit] <- NA
SubCritNames[iCrit] <- "alpha" SubCritNames[iCrit] <- "alpha"
Numer <- sd(VarSim[!TS_ignore]);
Denom <- sd(VarObs[!TS_ignore]);
if(Numer==0 & Denom==0){ Crit <- 1; } else { Crit <- Numer/Denom ; }
if(is.numeric(Crit) & is.finite(Crit)){ SubCritValues[iCrit] <- Crit; }
Numer <- sd(VarSim[!TS_ignore])
Denom <- sd(VarObs[!TS_ignore])
##SubErrorCrit_____KGE_beta______________________ if (Numer == 0 & Denom == 0) {
iCrit <- iCrit+1; Crit <- 1
SubCritPrint[iCrit] <- paste(CritName," mean(sim)/mean(obs) =", sep = "") } else {
SubCritValues[iCrit] <- NA; Crit <- Numer / Denom
}
if (is.numeric(Crit) & is.finite(Crit)) {
SubCritValues[iCrit] <- Crit
}
##SubErrorCrit_____KGE_beta______________________
iCrit <- iCrit + 1
SubCritPrint[iCrit] <- paste(CritName, " mean(sim)/mean(obs) =", sep = "")
SubCritValues[iCrit] <- NA
SubCritNames[iCrit] <- "beta" SubCritNames[iCrit] <- "beta"
if(meanVarSim==0 & meanVarObs==0){ Crit <- 1; } else { Crit <- meanVarSim/meanVarObs ; }
if(is.numeric(Crit) & is.finite(Crit)){ SubCritValues[iCrit] <- Crit; } if (meanVarSim == 0 & meanVarObs == 0) {
Crit <- 1
} else {
Crit <- meanVarSim / meanVarObs
}
if (is.numeric(Crit) & is.finite(Crit)) {
SubCritValues[iCrit] <- Crit
}
##ErrorCrit______________________________________ ##ErrorCrit______________________________________
if(sum(is.na(SubCritValues))==0){ if (sum(is.na(SubCritValues)) == 0) {
CritValue <- ( 1 - sqrt( (SubCritValues[1]-1)^2 + (SubCritValues[2]-1)^2 + (SubCritValues[3]-1)^2 ) ); CritValue <- (1 - sqrt((SubCritValues[1] - 1)^2 + (SubCritValues[2] - 1)^2 + (SubCritValues[3] - 1)^2))
} }
##Verbose______________________________________ ##Verbose______________________________________
if(verbose) { if (verbose) {
message("Crit. ", CritName, " = ", sprintf("%.4f", CritValue)) message("Crit. ", CritName, " = ", sprintf("%.4f", CritValue))
message(paste("\tSubCrit.", SubCritPrint, sprintf("%.4f", SubCritValues), "\n", sep = " ")) message(paste("\tSubCrit.", SubCritPrint, sprintf("%.4f", SubCritValues), "\n", sep = " "))
} }
##Output_________________________________________ ##Output_________________________________________
OutputsCrit <- list(CritValue = CritValue, CritName = CritName, OutputsCrit <- list(CritValue = CritValue,
SubCritValues = SubCritValues, SubCritNames = SubCritNames, CritBestValue = CritBestValue, CritName = CritName,
Multiplier = Multiplier, Ind_notcomputed = Ind_TS_ignore) SubCritValues = SubCritValues,
SubCritNames = SubCritNames,
CritBestValue = CritBestValue,
Multiplier = Multiplier,
Ind_notcomputed = Ind_TS_ignore
)
return(OutputsCrit) return(OutputsCrit)
} }
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment