diff --git a/NAMESPACE b/NAMESPACE index 0a7cdfe3151eacd4acffa393cae5d82aa0d76c3f..7072ef735537f9b6ed843d2eaca9ead05693e05a 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -10,7 +10,6 @@ useDynLib(airGR) ##################################### export(Calibration) export(Calibration_Michel) -export(Calibration_optim) export(CreateCalibOptions) export(CreateInputsCrit) export(CreateInputsModel) diff --git a/R/Calibration_optim.R b/R/Calibration_optim.R deleted file mode 100644 index 14fb6eceefb5e2f6940d19fcebc72936cd14ab99..0000000000000000000000000000000000000000 --- a/R/Calibration_optim.R +++ /dev/null @@ -1,152 +0,0 @@ -#************************************************************************************************* -#' Calibration algorithm which minimises the error criterion. \cr -#' \cr -#' The algorithm is based on the "optim" function from the "stats" R-package -#' (using method="L-BFGS-B", i.e. a local optimization quasi-Newton method). -#' -#' To optimise the exploration of the parameter space, transformation functions are used to convert -#' the model parameters. This is done using the TransfoParam functions. -#************************************************************************************************* -#' @title Calibration algorithm which minimises the error criterion using the stats::optim function -#' @author Laurent Coron (August 2013) -#' @example tests/example_Calibration_optim.R -#' @seealso \code{\link{Calibration}}, \code{\link{Calibration_Michel}}, -#' \code{\link{RunModel_GR4J}}, \code{\link{TransfoParam_GR4J}}, \code{\link{ErrorCrit_RMSE}}, -#' \code{\link{CreateInputsModel}}, \code{\link{CreateRunOptions}}, -#' \code{\link{CreateInputsCrit}}, \code{\link{CreateCalibOptions}}. -#' @encoding UTF-8 -#' @export -#_FunctionInputs__________________________________________________________________________________ -#' @param InputsModel [object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details -#' @param RunOptions [object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details -#' @param InputsCrit [object of class \emph{InputsCrit}] see \code{\link{CreateInputsCrit}} for details -#' @param CalibOptions [object of class \emph{CalibOptions}] see \code{\link{CreateCalibOptions}} for details -#' @param FUN_MOD [function] hydrological model function (e.g. RunModel_GR4J, RunModel_CemaNeigeGR4J) -#' @param FUN_CRIT [function] error criterion function (e.g. ErrorCrit_RMSE, ErrorCrit_NSE) -#' @param FUN_TRANSFO (optional) [function] model parameters transformation function, if the FUN_MOD used is native in the package FUN_TRANSFO is automatically defined -#' @param quiet (optional) [boolean] boolean indicating if the function is run in quiet mode or not, default=FALSE -#_FunctionOutputs_________________________________________________________________________________ -#' @return [list] list containing the function outputs organised as follows: -#' \tabular{ll}{ -#' \emph{$ParamFinalR } \tab [numeric] parameter set obtained at the end of the calibration \cr -#' \emph{$CritFinal } \tab [numeric] error criterion obtained at the end of the calibration \cr -#' \emph{$Nruns } \tab [numeric] number of model runs done during the calibration \cr -#' \emph{$CritName } \tab [character] name of the calibration criterion \cr -#' \emph{$CritBestValue} \tab [numeric] theoretical best criterion value \cr -#' } -#************************************************************************************************** -Calibration_optim <- function(InputsModel,RunOptions,InputsCrit,CalibOptions,FUN_MOD,FUN_CRIT,FUN_TRANSFO=NULL,quiet=FALSE){ - - - ##_check_class - if(inherits(InputsModel,"InputsModel")==FALSE){ stop("InputsModel must be of class 'InputsModel' \n"); return(NULL); } - if(inherits(RunOptions,"RunOptions")==FALSE){ stop("RunOptions must be of class 'RunOptions' \n"); return(NULL); } - if(inherits(InputsCrit,"InputsCrit")==FALSE){ stop("InputsCrit must be of class 'InputsCrit' \n"); return(NULL); } - if(inherits(CalibOptions,"CalibOptions")==FALSE){ stop("CalibOptions must be of class 'CalibOptions' \n"); return(NULL); } - if(inherits(CalibOptions,"optim")==FALSE){ stop("CalibOptions must be of class 'optim' if Calibration_optim is used \n"); return(NULL); } - - - ##_check_FUN_TRANSFO - if(is.null(FUN_TRANSFO)){ - if(identical(FUN_MOD,RunModel_GR4H )){ FUN_TRANSFO <- TransfoParam_GR4H ; } - if(identical(FUN_MOD,RunModel_GR4J )){ FUN_TRANSFO <- TransfoParam_GR4J ; } - if(identical(FUN_MOD,RunModel_GR5J )){ FUN_TRANSFO <- TransfoParam_GR5J ; } - if(identical(FUN_MOD,RunModel_GR6J )){ FUN_TRANSFO <- TransfoParam_GR6J ; } - if(identical(FUN_MOD,RunModel_GR2M )){ FUN_TRANSFO <- TransfoParam_GR2M ; } - if(identical(FUN_MOD,RunModel_GR1A )){ FUN_TRANSFO <- TransfoParam_GR1A ; } - if(identical(FUN_MOD,RunModel_CemaNeige )){ FUN_TRANSFO <- TransfoParam_CemaNeige; } - if(identical(FUN_MOD,RunModel_CemaNeigeGR4J) | identical(FUN_MOD,RunModel_CemaNeigeGR5J) | identical(FUN_MOD,RunModel_CemaNeigeGR6J)){ - if(identical(FUN_MOD,RunModel_CemaNeigeGR4J)){ FUN1 <- TransfoParam_GR4J; FUN2 <- TransfoParam_CemaNeige; } - if(identical(FUN_MOD,RunModel_CemaNeigeGR5J)){ FUN1 <- TransfoParam_GR5J; FUN2 <- TransfoParam_CemaNeige; } - if(identical(FUN_MOD,RunModel_CemaNeigeGR6J)){ FUN1 <- TransfoParam_GR6J; FUN2 <- TransfoParam_CemaNeige; } - FUN_TRANSFO <- function(ParamIn,Direction){ - Bool <- is.matrix(ParamIn); - if(Bool==FALSE){ ParamIn <- rbind(ParamIn); } - ParamOut <- NA*ParamIn; - NParam <- ncol(ParamIn); - ParamOut[, 1:(NParam-2)] <- FUN1(ParamIn[, 1:(NParam-2)],Direction); - ParamOut[,(NParam-1):NParam ] <- FUN2(ParamIn[,(NParam-1):NParam ],Direction); - if(Bool==FALSE){ ParamOut <- ParamOut[1,]; } - return(ParamOut); - } - } - if(is.null(FUN_TRANSFO)){ stop("FUN_TRANSFO was not found (in Calibration function) \n"); return(NULL); } - } - - - ##_RunModelAndCrit - RunModelAndCrit <- function(par,InputsModel,RunOptions,InputsCrit,CalibOptions,FUN_MOD,FUN_CRIT,FUN_TRANSFO){ - ParamT <- NA*CalibOptions$FixedParam; - ParamT[CalibOptions$OptimParam] <- par; - Param <- FUN_TRANSFO(ParamIn=ParamT,Direction="TR"); - Param[!CalibOptions$OptimParam] <- CalibOptions$FixedParam[!CalibOptions$OptimParam]; - OutputsModel <- FUN_MOD(InputsModel=InputsModel,RunOptions=RunOptions,Param=Param); - OutputsCrit <- FUN_CRIT(InputsCrit=InputsCrit,OutputsModel=OutputsModel); - return(OutputsCrit$CritValue*OutputsCrit$Multiplier); - } - - - ##_temporary_change_of_Outputs_Sim - RunOptions$Outputs_Sim <- RunOptions$Outputs_Cal; ### this reduces the size of the matrix exchange with fortran and therefore speeds the calibration - ##_screenPrint - if(!quiet){ - cat(paste("\t Calibration in progress (function optim from the stats package) \n",sep="")); - } - - - ##_lower_and_upper_limit_values (transformed) - RangesR <- CalibOptions$SearchRanges; - RangesT <- FUN_TRANSFO(RangesR,"RT"); - lower <- RangesT[1,CalibOptions$OptimParam]; - upper <- RangesT[2,CalibOptions$OptimParam]; - - ##_starting_values (transformed) - ParamStartT <- FUN_TRANSFO(CalibOptions$StartParam,"RT"); - par_start <- ParamStartT[CalibOptions$OptimParam]; - - - ##_calibration - RESULT <- optim(par=par_start,fn=RunModelAndCrit,gr=NULL, - InputsModel,RunOptions,InputsCrit,CalibOptions,FUN_MOD,FUN_CRIT,FUN_TRANSFO, ## arguments for the RunModelAndCrit function (other than par) - method="L-BFGS-B",lower=lower,upper=upper,control=list(),hessian=FALSE) - - - ##_outputs_preparation - ParamFinalT <- NA*ParamStartT; - ParamFinalT[CalibOptions$OptimParam] <- RESULT$par; - ParamFinalR <- FUN_TRANSFO(ParamFinalT,"TR"); - ParamFinalR[!CalibOptions$OptimParam] <- CalibOptions$FixedParam[!CalibOptions$OptimParam]; - CritFinal <- RESULT$value; - - ##_storage_of_crit_info - OutputsModel <- FUN_MOD(InputsModel=InputsModel,RunOptions=RunOptions,Param=ParamFinalR); - OutputsCrit <- FUN_CRIT(InputsCrit=InputsCrit,OutputsModel=OutputsModel); - CritName <- OutputsCrit$CritName; - CritBestValue <- OutputsCrit$CritBestValue; - Multiplier <- OutputsCrit$Multiplier; - - ##_screenPrint - if(!quiet){ - if(RESULT$convergence==0){ - cat(paste("\t Calibration completed: \n",sep="")); - cat(paste("\t Param = ",paste(formatC(ParamFinalR,format="f",width=8,digits=3),collapse=" , "),"\n",sep="")); - cat(paste("\t Crit ",format(CritName,width=12,justify="left")," = ",formatC(CritFinal*Multiplier,format="f",digits=4),"\n",sep="")); - } else { - cat(paste("\t Calibration failed: \n",sep="")); - cat(paste("\t ",RESULT$message,sep="")); - } - } - - - ##_function_output - OutputsCalib <- list(as.double(ParamFinalR),CritFinal*Multiplier,as.integer(RESULT$counts[1]),CritName,CritBestValue); - names(OutputsCalib) <- c("ParamFinalR","CritFinal","NRuns","CritName","CritBestValue"); - class(OutputsCalib) <- c("OutputsCalib","optim"); - return(OutputsCalib); - - -} - - - - diff --git a/man/Calibration_optim.Rd b/man/Calibration_optim.Rd deleted file mode 100644 index d24d91b02a0607324eaef9036a5c8b9007e07867..0000000000000000000000000000000000000000 --- a/man/Calibration_optim.Rd +++ /dev/null @@ -1,109 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/Calibration_optim.R -\encoding{UTF-8} -\name{Calibration_optim} -\alias{Calibration_optim} -\title{Calibration algorithm which optimises the error criterion selected as objective function using the stats::optim function} -\usage{ -Calibration_optim(InputsModel, RunOptions, InputsCrit, CalibOptions, FUN_MOD, - FUN_CRIT, FUN_TRANSFO = NULL, quiet = FALSE) -} -\arguments{ -\item{InputsModel}{[object of class \emph{InputsModel}] see \code{\link{CreateInputsModel}} for details} - -\item{RunOptions}{[object of class \emph{RunOptions}] see \code{\link{CreateRunOptions}} for details} - -\item{InputsCrit}{[object of class \emph{InputsCrit}] see \code{\link{CreateInputsCrit}} for details} - -\item{CalibOptions}{[object of class \emph{CalibOptions}] see \code{\link{CreateCalibOptions}} for details} - -\item{FUN_MOD}{[function] hydrological model function (e.g. RunModel_GR4J, RunModel_CemaNeigeGR4J)} - -\item{FUN_CRIT}{[function] error criterion function (e.g. ErrorCrit_RMSE, ErrorCrit_NSE)} - -\item{FUN_TRANSFO}{(optional) [function] model parameters transformation function, if the FUN_MOD used is native in the package, FUN_TRANSFO is automatically defined} - -\item{quiet}{(optional) [boolean] boolean indicating if the function is run in quiet mode or not, default=FALSE} -} -\value{ -[list] list containing the function outputs organised as follows: - \tabular{ll}{ - \emph{$ParamFinalR } \tab [numeric] parameter set obtained at the end of the calibration \cr - \emph{$CritFinal } \tab [numeric] error criterion selected as objective function obtained at the end of the calibration \cr - \emph{$Nruns } \tab [numeric] number of model runs done during the calibration \cr - \emph{$CritName } \tab [character] name of the calibration criterion selected as objective function \cr - \emph{$CritBestValue} \tab [numeric] theoretical best criterion value \cr - } -} -\description{ -Calibration algorithm which optimises the error criterion selected as objective function. \cr -\cr -The algorithm is based on the \code{\link{optim}} function from the "stats" R-package -(using method="L-BFGS-B", i.e. a local optimization quasi-Newton method). -} -\details{ -To optimise the exploration of the parameter space, transformation functions are used to convert -the model parameters. This is done using the TransfoParam functions. -} -\examples{ -## loading catchment data -library(airGR) -data(L0123001) - -## preparation of InputsModel object -InputsModel <- CreateInputsModel(FUN_MOD = RunModel_GR4J, DatesR = BasinObs$DatesR, - Precip = BasinObs$P, PotEvap = BasinObs$E) - -## calibration period selection -Ind_Run <- seq(which(format(BasinObs$DatesR, format = "\%d/\%m/\%Y \%H:\%M")=="01/01/1990 00:00"), - which(format(BasinObs$DatesR, format="\%d/\%m/\%Y \%H:\%M")=="31/12/1999 00:00")) - -## preparation of RunOptions object -RunOptions <- CreateRunOptions(FUN_MOD = RunModel_GR4J, - InputsModel = InputsModel, IndPeriod_Run = Ind_Run) - -## calibration criterion: preparation of the InputsCrit object -InputsCrit <- CreateInputsCrit(FUN_CRIT = ErrorCrit_NSE, InputsModel = InputsModel, - RunOptions = RunOptions, Qobs = BasinObs$Qmm[Ind_Run]) - -## preparation of CalibOptions object -CalibOptions <- CreateCalibOptions(FUN_MOD = RunModel_GR4J, - FUN_CALIB = Calibration_optim) - -## calibration -OutputsCalib <- Calibration_optim(InputsModel = InputsModel, RunOptions = RunOptions, - InputsCrit = InputsCrit, CalibOptions = CalibOptions, - FUN_MOD = RunModel_GR4J, FUN_CRIT = ErrorCrit_NSE) - -## simulation -Param <- OutputsCalib$ParamFinalR -OutputsModel <- RunModel_GR4J(InputsModel = InputsModel, - RunOptions = RunOptions, Param = Param) - -## results preview -plot_OutputsModel(OutputsModel = OutputsModel, Qobs = BasinObs$Qmm[Ind_Run]) - -## efficiency criterion: Nash-Sutcliffe Efficiency -InputsCrit <- CreateInputsCrit(FUN_CRIT = ErrorCrit_NSE, InputsModel = InputsModel, - RunOptions = RunOptions, Qobs = BasinObs$Qmm[Ind_Run]) -OutputsCrit <- ErrorCrit_NSE(InputsCrit = InputsCrit, OutputsModel = OutputsModel) -cat(paste(" Crit ", OutputsCrit$CritName, " ", - round(OutputsCrit$CritValue, 4), "\\n", sep = "")) - -## efficiency criterion: Kling-Gupta Efficiency -InputsCrit <- CreateInputsCrit(FUN_CRIT = ErrorCrit_KGE, InputsModel = InputsModel, - RunOptions = RunOptions, Qobs = BasinObs$Qmm[Ind_Run]) -OutputsCrit <- ErrorCrit_KGE(InputsCrit = InputsCrit, OutputsModel = OutputsModel) -cat(paste(" Crit ", OutputsCrit$CritName, " ", - round(OutputsCrit$CritValue, 4), "\\n", sep = "")) -} -\author{ -Laurent Coron (August 2013) -} -\seealso{ -\code{\link{Calibration}}, \code{\link{Calibration_Michel}}, - \code{\link{RunModel_GR4J}}, \code{\link{TransfoParam_GR4J}}, \code{\link{ErrorCrit_RMSE}}, - \code{\link{CreateInputsModel}}, \code{\link{CreateRunOptions}}, - \code{\link{CreateInputsCrit}}, \code{\link{CreateCalibOptions}}. -} -