Skip to content
GitLab
Commit 6a61568c authored by unknown's avatar unknown
Browse files

#4410 The package now contains a vignette

parent e0c0f599
dev 115-add-pcict-date-format-managment 118-plot-outputsmodel-remove-case-sensitivity-to-the-which-argument 136-decreased-performance-of-calibration-execution-time 60-migrate-transfoparam-functions-and-createcaliboptions-to-s3-methods 96-parallelize-the-grid-screening-step-during-calibration airGRplus cleanAirGRplus master CRAN_v1.7.0 CRAN_v1.6.12 CRAN_v1.6.11 CRAN_v1.6.10.4 CRAN_v1.6.9.27 CRAN_v1.6.9.24 CRAN_v1.6.9.23 CRAN_v1.6.9.21 CRAN_v1.4.3.65 CRAN_v1.4.3.60 CRAN_v1.4.3.52 CRAN_v1.3.2.42 CRAN_v1.3.2.41 CRAN_v1.3.2.40 CRAN_v1.3.2.39 CRAN_v1.3.2.23 CRAN_v1.3.2.17 CRAN_v1.2.13.16 CRAN_v1.2.13.13 CRAN_v1.0.14.1 CRAN_v1.0.13.19 CRAN_v1.0.10.11 CRAN_v1.0.9.64 CRAN_v1.0.5.12 CRAN_v1.0.5.11 CRAN_v1.0.5.10
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 920 additions and 0 deletions
+920 -0
inst/doc/airGR.R 0 → 100644
+ 54
− 0
View file @ 6a61568c
## ------------------------------------------------------------------------
library(airGR)
## ------------------------------------------------------------------------
data(L0123001)
summary(BasinObs)
## ------------------------------------------------------------------------
InputsModel <- CreateInputsModel(FUN_MOD = RunModel_GR4J, DatesR = BasinObs$DatesR,
Precip = BasinObs$P, PotEvap = BasinObs$E)
str(InputsModel)
## ------------------------------------------------------------------------
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"))
str(Ind_Run)
## ------------------------------------------------------------------------
RunOptions <- CreateRunOptions(FUN_MOD = RunModel_GR4J,
InputsModel = InputsModel, IndPeriod_Run = Ind_Run,
IniStates = NULL, IniResLevels = NULL, IndPeriod_WarmUp = NULL)
str(RunOptions)
## ------------------------------------------------------------------------
InputsCrit <- CreateInputsCrit(FUN_CRIT = ErrorCrit_NSE, InputsModel = InputsModel,
RunOptions = RunOptions, Qobs = BasinObs$Qmm[Ind_Run])
str(InputsCrit)
## ------------------------------------------------------------------------
CalibOptions <- CreateCalibOptions(FUN_MOD = RunModel_GR4J, FUN_CALIB = Calibration_Michel)
str(CalibOptions)
## ------------------------------------------------------------------------
OutputsCalib <- Calibration_Michel(InputsModel = InputsModel, RunOptions = RunOptions,
InputsCrit = InputsCrit, CalibOptions = CalibOptions,
FUN_MOD = RunModel_GR4J, FUN_CRIT = ErrorCrit_NSE)
Param <- OutputsCalib$ParamFinalR
Param
## ------------------------------------------------------------------------
OutputsModel <- RunModel_GR4J(InputsModel = InputsModel, RunOptions = RunOptions, Param = Param)
str(OutputsModel)
## ----eval=F--------------------------------------------------------------
# plot_OutputsModel(OutputsModel = OutputsModel, Qobs = BasinObs$Qmm[Ind_Run])
## ------------------------------------------------------------------------
OutputsCrit <- ErrorCrit_NSE(InputsCrit = InputsCrit, OutputsModel = OutputsModel)
str(OutputsCrit)
## ------------------------------------------------------------------------
OutputsCrit <- ErrorCrit_KGE(InputsCrit = InputsCrit, OutputsModel = OutputsModel)
str(OutputsCrit)
inst/doc/airGR.Rmd 0 → 100644
+ 259
− 0
View file @ 6a61568c
---
title: "airGR -- Overview"
output: rmarkdown::html_vignette
vignette: >
%\VignetteEngine{knitr::rmarkdown}
%\VignetteIndexEntry{airGR}
%\VignetteEncoding{UTF-8}
---
# Introduction
**airGR** is a package which brings into the [**R software**](https://cran.r-project.org/) the hydrological modelling tools used and developed at the [Catchment Hydrology Team](https://webgr.irstea.fr/?lang=en) at [Irstea (France)](http://www.irstea.fr/en/), including the [**GR rainfall-runoff models**](https://webgr.irstea.fr/modeles/?lang=en) and a snowmelt and accumulation model, **CemaNeige**. Each model core is coded in **Fortran** to ensure low computational time. The other package functions (i.e. mainly the calibration algorithm and the efficiency criteria calculation) are coded in **R**.
The **airGR** package has been designed to fulfil two major requirements: to facilitate the use by non-expert users and to allow flexibility regarding the addition of external criteria, models or calibration algorithms. The names of the functions and their arguments were chosen to this end. **airGR** also contains basics plotting facilities.
Six hydrological models and one snowmelt and accumulation model are implemented in **airGR**. The snow model can be used alone or together with the daily hydrological models.
The models can be called within **airGR** using the following functions:
* `RunModel_GR4H()`: four-parameter hourly lumped hydrological model [@mathevet_quels_2005]
* `RunModel_GR4J()`: four-parameter daily lumped hydrological model [@perrin_improvement_2003]
* `RunModel_GR5J()`: five-parameter daily lumped hydrological model [@le_moine_bassin_2008]
* `RunModel_GR6J()`: six-parameter daily lumped hydrological model [@pushpalatha_downward_2011]
* `RunModel_GR2M()`: two-parameter monthly lumped hydrological model [@mouelhi_vers_2003; @mouelhi_stepwise_2006]
* `RunModel_GR1A()`: one-parameter yearly lumped hydrological model [@mouelhi_vers_2003; @mouelhi_linking_2006]
* `RunModel_CemaNeige()`: two-parameter degree-day snowmelt and accumulation model [@valery_as_2014]
* `RunModel_CemaNeigeGR4J()`: combined use of **GR4J** and **CemaNeige**
* `RunModel_CemaNeigeGR5J()`: combined use of **GR5J** and **CemaNeige**
* `RunModel_CemaNeigeGR6J()`: combined use of **GR6J** and **CemaNeige**
The [**GRP**](https://webgr.irstea.fr/modeles/modele-de-prevision-grp/?lang=en) forecasting model and the [**Otamin**](https://webgr.irstea.fr/modeles/otamin/?lang=en) predictive uncertainty tool are not availabe in **airGR**.
# Loading data
In the following example, we use a data sample contained in the package. For real applications, the user has to import its data into **R** and to prepare it with an adequate data.frame format as described below.
First, it is necessary to load the **airGR** package:
```{r}
library(airGR)
```
This is an example of a `data.frame` of hydrometeorological observations time series for a fictional catchment included in the **airGR** package that contains:
* *DatesR*: dates in the POSIXt format
* *P*: average precipitation [mm/day]
* *T*: catchment average air temperature [&#8451;]
* *E*: catchment average potential evapotranspiration [mm/day]
* *Qls*: outlet discharge [l/s]
* *Qmm*: outlet discharge [mm/day]
```{r}
data(L0123001)
summary(BasinObs)
```
The usual functions (e.g. `read.table()`) can be used to load real-case datasets.
# Preparation of functions inputs
To run a model, the functions of the **airGR** package (e.g. the models, calibration and criteria calculation functions) require data and options with specific formats.
To facilitate the use of the package, there are several functions dedicated to the creation of these objects:
* `CreateInputsModel()`: prepares the inputs for the different hydrological models (times series of dates, rainfall, observed streamflow, etc.)
* `CreateRunOptions()`: prepares the options for the hydrological model run (warm-up period, calibration period, etc.)
* `CreateInputsCrit()`: prepares the options in order to compute efficiency criterions (choice of the criterion, choice of the transformation on streamflows: "log", "root", etc.)
* `CreateCalibOptions()`: prepares the options for the hydrological model calibration algorithm (choice of parameters to optimize, predefined values for uncalibrated parameters, etc.)
## InputsModel object
To run a GR hydrologic model, the user has to prepare the input data with the `CreateInputsModel()` function.
As arguments, this function needs the function name corresponding to the model the user wants to run, a vector of dates, a vector of precipitation and a vector of potential evapotranspiration.
In the example below, we already have the potential evapotranspiration. If the user doesn't have these data, it is possible to compute it with the [Oudin's formula](http://dx.doi.org/10.1016/j.jhydrol.2004.08.026) with the `PEdaily_Oudin()` function (this function only needs julian days, daily average air temperature and latitude).
Missing values (`NA`) of precipitation (or potential evapotranspiration) are **not allowed**.
```{r}
InputsModel <- CreateInputsModel(FUN_MOD = RunModel_GR4J, DatesR = BasinObs$DatesR,
Precip = BasinObs$P, PotEvap = BasinObs$E)
str(InputsModel)
```
## RunOptions object
The `CreateRunOptions()` function allows to prepare the options required to the `RunModel*()` functions, which are the actual models functions.
The user must at least define the following arguments:
* `FUN_MOD`: the name of the model function to run
* `InputsModel`: the associated inputs data
* `IndPeriod_Run`: the period on which the model is run
To select a period for which the user wants to run the model, select the corresponding indexes for different time periods (not the POSIXt dates), as follows:
```{r}
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"))
str(Ind_Run)
```
The intialization of hydrological models is of the utmost importance. Indeed, an inaccurate intialisation causes poor quality streamflow simulations during the earliest stages of the running period. For example, in the GR models, the reservoirs levels are by default set to 50 % of their capacity, which may be far from their ideal value. Two solutions are offered to accurately initialize the GR models in **airGR**: manually predefining the initial states or running the models during a warm up period before the actual running period. It is generally advised to set up this warm up period to be equal or superior to one year.
As a consequence, it is possible to define in `CreateRunOptions()` the following arguments:
* `IniStates`: the initial states of the 2 unit hydrographs (20 + 40 = 60 units)
* `IniResLevels`: the initial levels of the production and routing stores
* `IndPeriod_WarmUp`: the warm-up period used to run the model, to be defined in the same format as `IndPeriod_Run`
```{r}
RunOptions <- CreateRunOptions(FUN_MOD = RunModel_GR4J,
InputsModel = InputsModel, IndPeriod_Run = Ind_Run,
IniStates = NULL, IniResLevels = NULL, IndPeriod_WarmUp = NULL)
str(RunOptions)
```
The `CreateRunOptions()` function returns warnings if the default initialization options are used:
* `IniStates` and `IniResLevels` are automatically set to initialize all the model states at 0, except for the production and routing stores which are initialised at 50 % of their capacity
* `IndPeriod_WarmUp` default setting ensures a one-year warm-up using the time steps preceding the `IndPeriod_Run`, if available
## InputsCrit object
The `CreateInputsCrit()` function allows to prepare the input in order to calculate a criterion. It is possible to define the following arguments:
* `FUN_CRIT`: the name of the error criterion function (they are introduced later on)
* `InputsModel`: the inputs of the hydrological model previously prepared by the `CeateInputsModel()` function
* `RunOptions`: the options of the hydrological model previously prepared by the `CreateRunOptions()` function
* `Qobs`: the observed streamflows expressed in *mm/time step*
Missing values (`NA`) are **allowed** for observed streamflows.
```{r}
InputsCrit <- CreateInputsCrit(FUN_CRIT = ErrorCrit_NSE, InputsModel = InputsModel,
RunOptions = RunOptions, Qobs = BasinObs$Qmm[Ind_Run])
str(InputsCrit)
```
## CalibOptions object
The user needs to prepare the calibration options with the `CreateCalibOptions()`function. For that, it is necessary to define the following arguments:
* `FUN_MOD`: the name of the model function
* `FUN_CALIB`: the name of the calibration algorithm
```{r}
CalibOptions <- CreateCalibOptions(FUN_MOD = RunModel_GR4J, FUN_CALIB = Calibration_Michel)
str(CalibOptions)
```
# Criteria
The evaluation of the quality of a simulation is estimated through the calculation of criteria. These criteria can be used both as objective-functions during the calibration of the model, or as a measure for evaluating its control performance.
The package offers the possibility to use different criteria:
* `ErrorCrit_RMSE()`: Root mean square error (RMSE)
* `ErrorCrit_NSE()`: Nash-Sutcliffe model efficiency coefficient (NSE)
* `ErrorCrit_KGE()`: Kling-Gupta efficiency criterion (KGE)
* `ErrorCrit_KGE2()`: modified Kling-Gupta efficiency criterion (KGE')
It is also possible to create user-defined criteria. For doing that, it is only necessary to define the function in **R** following the same syntax that the criteria functions included in **airGR**.
# Calibration
The objective of the calibration algorithm is to identify the model parameters: by comparing the model outputs with observed data, this algorithm determines the combination of parameters that represents the best the behavior of the watershed.
In the **airGR** package, a function called `Calibration_Michel()` is implemented. This functions allows running a calibration with the package models.
The calibration algorithm optimizes the error criterion selected as objective-function. This algorithm works in two steps:
1. a screening of the parameters space is performed using either a rough predefined grid or a user-defined list of parameter sets
2. a simple steepest descent local search algorithm is performed from the best set of parameters found at the first step
```{r}
OutputsCalib <- Calibration_Michel(InputsModel = InputsModel, RunOptions = RunOptions,
InputsCrit = InputsCrit, CalibOptions = CalibOptions,
FUN_MOD = RunModel_GR4J, FUN_CRIT = ErrorCrit_NSE)
Param <- OutputsCalib$ParamFinalR
Param
```
The `Calibration_Michel()` function is the only one implemented in the **airGR** package to calibrate the model, but the user can implement its own calibration function.
This function returns a vector with the parameters of the chosen model, which means that the number of values can differ depending on the model that is used. It is possible to use the `Calibration_Michel()` function with user-implemented hydrological models.
# Validation
This step assess the predictive capacity of the model. Validation is defined as the estimation of the accuracy of the model on datasets that are not used in its construction, and in particular its calibration.
The classical way to perform a validation is to keep data from a period separated from the calibration period. If possible, this control period should correspond to climatic situations rather that differ from those of the calibration period in order to better point out the qualities and weakness of the model. This exercise is necessary for assessing the robustness of the model, that is to say its ability to keep stable performances outside the calibration conditions.
Performing a model validation with **airGR** is similar to running a simulation (see below).
# Simulation
## Simulation run
To run a model, the user has to use the `RunModel*()` functions (`InputsModel`, `RunOptions` and parameters).
All the data needed have already been prepared.
```{r}
OutputsModel <- RunModel_GR4J(InputsModel = InputsModel, RunOptions = RunOptions, Param = Param)
str(OutputsModel)
```
## Results preview
Although it is possible for the user to design its own graphics from the outputs of the `RunModel*()` functions, the **airGR** package offers the possibility to make use of the `plot_OutputsModel()` function. This function returns a dashboard of results including various graphs (depending on the model used):
* time series of total precipitation and simulated streamflows (and observed streamflows if provided)
* interannual median monthly simulated streamflows (and monthly observed streamflows if provided)
* correlation plot between simulated and observed streamflows (if observed streamflows provided)
* cumulative frequency plot for simulated streamflows (and for observed streamflows if provided)
```{r,eval=F}
plot_OutputsModel(OutputsModel = OutputsModel, Qobs = BasinObs$Qmm[Ind_Run])
```
Moreover, if the CemaNeige model is used, the simulated snowpack time series are plotted.
## Efficiency criterion
To evaluate the efficiency of the model, it is possible to use the same criterion as defined at the calibration step or to use an other criterion.
```{r}
OutputsCrit <- ErrorCrit_NSE(InputsCrit = InputsCrit, OutputsModel = OutputsModel)
str(OutputsCrit)
```
```{r}
OutputsCrit <- ErrorCrit_KGE(InputsCrit = InputsCrit, OutputsModel = OutputsModel)
str(OutputsCrit)
```
inst/doc/airGR.html 0 → 100644
+ 348
− 0
View file @ 6a61568c
<!DOCTYPE html>
<html xmlns="http://www.w3.org/1999/xhtml">
<head>
<meta charset="utf-8">
<meta http-equiv="Content-Type" content="text/html; charset=utf-8" />
<meta name="generator" content="pandoc" />
<meta name="viewport" content="width=device-width, initial-scale=1">
<title>airGR – Overview</title>
<style type="text/css">code{white-space: pre;}</style>
<style type="text/css">
div.sourceCode { overflow-x: auto; }
table.sourceCode, tr.sourceCode, td.lineNumbers, td.sourceCode {
margin: 0; padding: 0; vertical-align: baseline; border: none; }
table.sourceCode { width: 100%; line-height: 100%; }
td.lineNumbers { text-align: right; padding-right: 4px; padding-left: 4px; color: #aaaaaa; border-right: 1px solid #aaaaaa; }
td.sourceCode { padding-left: 5px; }
code > span.kw { color: #007020; font-weight: bold; } /* Keyword */
code > span.dt { color: #902000; } /* DataType */
code > span.dv { color: #40a070; } /* DecVal */
code > span.bn { color: #40a070; } /* BaseN */
code > span.fl { color: #40a070; } /* Float */
code > span.ch { color: #4070a0; } /* Char */
code > span.st { color: #4070a0; } /* String */
code > span.co { color: #60a0b0; font-style: italic; } /* Comment */
code > span.ot { color: #007020; } /* Other */
code > span.al { color: #ff0000; font-weight: bold; } /* Alert */
code > span.fu { color: #06287e; } /* Function */
code > span.er { color: #ff0000; font-weight: bold; } /* Error */
code > span.wa { color: #60a0b0; font-weight: bold; font-style: italic; } /* Warning */
code > span.cn { color: #880000; } /* Constant */
code > span.sc { color: #4070a0; } /* SpecialChar */
code > span.vs { color: #4070a0; } /* VerbatimString */
code > span.ss { color: #bb6688; } /* SpecialString */
code > span.im { } /* Import */
code > span.va { color: #19177c; } /* Variable */
code > span.cf { color: #007020; font-weight: bold; } /* ControlFlow */
code > span.op { color: #666666; } /* Operator */
code > span.bu { } /* BuiltIn */
code > span.ex { } /* Extension */
code > span.pp { color: #bc7a00; } /* Preprocessor */
code > span.at { color: #7d9029; } /* Attribute */
code > span.do { color: #ba2121; font-style: italic; } /* Documentation */
code > span.an { color: #60a0b0; font-weight: bold; font-style: italic; } /* Annotation */
code > span.cv { color: #60a0b0; font-weight: bold; font-style: italic; } /* CommentVar */
code > span.in { color: #60a0b0; font-weight: bold; font-style: italic; } /* Information */
</style>
<link href="data:text/css;charset=utf-8,body%20%7B%0Abackground%2Dcolor%3A%20%23fff%3B%0Amargin%3A%201em%20auto%3B%0Amax%2Dwidth%3A%20700px%3B%0Aoverflow%3A%20visible%3B%0Apadding%2Dleft%3A%202em%3B%0Apadding%2Dright%3A%202em%3B%0Afont%2Dfamily%3A%20%22Open%20Sans%22%2C%20%22Helvetica%20Neue%22%2C%20Helvetica%2C%20Arial%2C%20sans%2Dserif%3B%0Afont%2Dsize%3A%2014px%3B%0Aline%2Dheight%3A%201%2E35%3B%0A%7D%0A%23header%20%7B%0Atext%2Dalign%3A%20center%3B%0A%7D%0A%23TOC%20%7B%0Aclear%3A%20both%3B%0Amargin%3A%200%200%2010px%2010px%3B%0Apadding%3A%204px%3B%0Awidth%3A%20400px%3B%0Aborder%3A%201px%20solid%20%23CCCCCC%3B%0Aborder%2Dradius%3A%205px%3B%0Abackground%2Dcolor%3A%20%23f6f6f6%3B%0Afont%2Dsize%3A%2013px%3B%0Aline%2Dheight%3A%201%2E3%3B%0A%7D%0A%23TOC%20%2Etoctitle%20%7B%0Afont%2Dweight%3A%20bold%3B%0Afont%2Dsize%3A%2015px%3B%0Amargin%2Dleft%3A%205px%3B%0A%7D%0A%23TOC%20ul%20%7B%0Apadding%2Dleft%3A%2040px%3B%0Amargin%2Dleft%3A%20%2D1%2E5em%3B%0Amargin%2Dtop%3A%205px%3B%0Amargin%2Dbottom%3A%205px%3B%0A%7D%0A%23TOC%20ul%20ul%20%7B%0Amargin%2Dleft%3A%20%2D2em%3B%0A%7D%0A%23TOC%20li%20%7B%0Aline%2Dheight%3A%2016px%3B%0A%7D%0Atable%20%7B%0Amargin%3A%201em%20auto%3B%0Aborder%2Dwidth%3A%201px%3B%0Aborder%2Dcolor%3A%20%23DDDDDD%3B%0Aborder%2Dstyle%3A%20outset%3B%0Aborder%2Dcollapse%3A%20collapse%3B%0A%7D%0Atable%20th%20%7B%0Aborder%2Dwidth%3A%202px%3B%0Apadding%3A%205px%3B%0Aborder%2Dstyle%3A%20inset%3B%0A%7D%0Atable%20td%20%7B%0Aborder%2Dwidth%3A%201px%3B%0Aborder%2Dstyle%3A%20inset%3B%0Aline%2Dheight%3A%2018px%3B%0Apadding%3A%205px%205px%3B%0A%7D%0Atable%2C%20table%20th%2C%20table%20td%20%7B%0Aborder%2Dleft%2Dstyle%3A%20none%3B%0Aborder%2Dright%2Dstyle%3A%20none%3B%0A%7D%0Atable%20thead%2C%20table%20tr%2Eeven%20%7B%0Abackground%2Dcolor%3A%20%23f7f7f7%3B%0A%7D%0Ap%20%7B%0Amargin%3A%200%2E5em%200%3B%0A%7D%0Ablockquote%20%7B%0Abackground%2Dcolor%3A%20%23f6f6f6%3B%0Apadding%3A%200%2E25em%200%2E75em%3B%0A%7D%0Ahr%20%7B%0Aborder%2Dstyle%3A%20solid%3B%0Aborder%3A%20none%3B%0Aborder%2Dtop%3A%201px%20solid%20%23777%3B%0Amargin%3A%2028px%200%3B%0A%7D%0Adl%20%7B%0Amargin%2Dleft%3A%200%3B%0A%7D%0Adl%20dd%20%7B%0Amargin%2Dbottom%3A%2013px%3B%0Amargin%2Dleft%3A%2013px%3B%0A%7D%0Adl%20dt%20%7B%0Afont%2Dweight%3A%20bold%3B%0A%7D%0Aul%20%7B%0Amargin%2Dtop%3A%200%3B%0A%7D%0Aul%20li%20%7B%0Alist%2Dstyle%3A%20circle%20outside%3B%0A%7D%0Aul%20ul%20%7B%0Amargin%2Dbottom%3A%200%3B%0A%7D%0Apre%2C%20code%20%7B%0Abackground%2Dcolor%3A%20%23f7f7f7%3B%0Aborder%2Dradius%3A%203px%3B%0Acolor%3A%20%23333%3B%0Awhite%2Dspace%3A%20pre%2Dwrap%3B%20%0A%7D%0Apre%20%7B%0Aborder%2Dradius%3A%203px%3B%0Amargin%3A%205px%200px%2010px%200px%3B%0Apadding%3A%2010px%3B%0A%7D%0Apre%3Anot%28%5Bclass%5D%29%20%7B%0Abackground%2Dcolor%3A%20%23f7f7f7%3B%0A%7D%0Acode%20%7B%0Afont%2Dfamily%3A%20Consolas%2C%20Monaco%2C%20%27Courier%20New%27%2C%20monospace%3B%0Afont%2Dsize%3A%2085%25%3B%0A%7D%0Ap%20%3E%20code%2C%20li%20%3E%20code%20%7B%0Apadding%3A%202px%200px%3B%0A%7D%0Adiv%2Efigure%20%7B%0Atext%2Dalign%3A%20center%3B%0A%7D%0Aimg%20%7B%0Abackground%2Dcolor%3A%20%23FFFFFF%3B%0Apadding%3A%202px%3B%0Aborder%3A%201px%20solid%20%23DDDDDD%3B%0Aborder%2Dradius%3A%203px%3B%0Aborder%3A%201px%20solid%20%23CCCCCC%3B%0Amargin%3A%200%205px%3B%0A%7D%0Ah1%20%7B%0Amargin%2Dtop%3A%200%3B%0Afont%2Dsize%3A%2035px%3B%0Aline%2Dheight%3A%2040px%3B%0A%7D%0Ah2%20%7B%0Aborder%2Dbottom%3A%204px%20solid%20%23f7f7f7%3B%0Apadding%2Dtop%3A%2010px%3B%0Apadding%2Dbottom%3A%202px%3B%0Afont%2Dsize%3A%20145%25%3B%0A%7D%0Ah3%20%7B%0Aborder%2Dbottom%3A%202px%20solid%20%23f7f7f7%3B%0Apadding%2Dtop%3A%2010px%3B%0Afont%2Dsize%3A%20120%25%3B%0A%7D%0Ah4%20%7B%0Aborder%2Dbottom%3A%201px%20solid%20%23f7f7f7%3B%0Amargin%2Dleft%3A%208px%3B%0Afont%2Dsize%3A%20105%25%3B%0A%7D%0Ah5%2C%20h6%20%7B%0Aborder%2Dbottom%3A%201px%20solid%20%23ccc%3B%0Afont%2Dsize%3A%20105%25%3B%0A%7D%0Aa%20%7B%0Acolor%3A%20%230033dd%3B%0Atext%2Ddecoration%3A%20none%3B%0A%7D%0Aa%3Ahover%20%7B%0Acolor%3A%20%236666ff%3B%20%7D%0Aa%3Avisited%20%7B%0Acolor%3A%20%23800080%3B%20%7D%0Aa%3Avisited%3Ahover%20%7B%0Acolor%3A%20%23BB00BB%3B%20%7D%0Aa%5Bhref%5E%3D%22http%3A%22%5D%20%7B%0Atext%2Ddecoration%3A%20underline%3B%20%7D%0Aa%5Bhref%5E%3D%22https%3A%22%5D%20%7B%0Atext%2Ddecoration%3A%20underline%3B%20%7D%0A%0Acode%20%3E%20span%2Ekw%20%7B%20color%3A%20%23555%3B%20font%2Dweight%3A%20bold%3B%20%7D%20%0Acode%20%3E%20span%2Edt%20%7B%20color%3A%20%23902000%3B%20%7D%20%0Acode%20%3E%20span%2Edv%20%7B%20color%3A%20%2340a070%3B%20%7D%20%0Acode%20%3E%20span%2Ebn%20%7B%20color%3A%20%23d14%3B%20%7D%20%0Acode%20%3E%20span%2Efl%20%7B%20color%3A%20%23d14%3B%20%7D%20%0Acode%20%3E%20span%2Ech%20%7B%20color%3A%20%23d14%3B%20%7D%20%0Acode%20%3E%20span%2Est%20%7B%20color%3A%20%23d14%3B%20%7D%20%0Acode%20%3E%20span%2Eco%20%7B%20color%3A%20%23888888%3B%20font%2Dstyle%3A%20italic%3B%20%7D%20%0Acode%20%3E%20span%2Eot%20%7B%20color%3A%20%23007020%3B%20%7D%20%0Acode%20%3E%20span%2Eal%20%7B%20color%3A%20%23ff0000%3B%20font%2Dweight%3A%20bold%3B%20%7D%20%0Acode%20%3E%20span%2Efu%20%7B%20color%3A%20%23900%3B%20font%2Dweight%3A%20bold%3B%20%7D%20%20code%20%3E%20span%2Eer%20%7B%20color%3A%20%23a61717%3B%20background%2Dcolor%3A%20%23e3d2d2%3B%20%7D%20%0A" rel="stylesheet" type="text/css" />
</head>
<body>
<h1 class="title toc-ignore">airGR – Overview</h1>
<div id="introduction" class="section level1">
<h1>Introduction</h1>
<p><strong>airGR</strong> is a package which brings into the <a href="https://cran.r-project.org/"><strong>R software</strong></a> the hydrological modelling tools used and developed at the <a href="https://webgr.irstea.fr/?lang=en">Catchment Hydrology Team</a> at <a href="http://www.irstea.fr/en/">Irstea (France)</a>, including the <a href="https://webgr.irstea.fr/modeles/?lang=en"><strong>GR rainfall-runoff models</strong></a> and a snowmelt and accumulation model, <strong>CemaNeige</strong>. Each model core is coded in <strong>Fortran</strong> to ensure low computational time. The other package functions (i.e. mainly the calibration algorithm and the efficiency criteria calculation) are coded in <strong>R</strong>.</p>
<p>The <strong>airGR</strong> package has been designed to fulfil two major requirements: to facilitate the use by non-expert users and to allow flexibility regarding the addition of external criteria, models or calibration algorithms. The names of the functions and their arguments were chosen to this end. <strong>airGR</strong> also contains basics plotting facilities.</p>
<p>Six hydrological models and one snowmelt and accumulation model are implemented in <strong>airGR</strong>. The snow model can be used alone or together with the daily hydrological models.</p>
<p>The models can be called within <strong>airGR</strong> using the following functions:</p>
<ul>
<li><code>RunModel_GR4H()</code>: four-parameter hourly lumped hydrological model <span class="citation">[@mathevet_quels_2005]</span></li>
<li><code>RunModel_GR4J()</code>: four-parameter daily lumped hydrological model <span class="citation">[@perrin_improvement_2003]</span></li>
<li><code>RunModel_GR5J()</code>: five-parameter daily lumped hydrological model <span class="citation">[@le_moine_bassin_2008]</span></li>
<li><code>RunModel_GR6J()</code>: six-parameter daily lumped hydrological model <span class="citation">[@pushpalatha_downward_2011]</span></li>
<li><code>RunModel_GR2M()</code>: two-parameter monthly lumped hydrological model <span class="citation">[@mouelhi_vers_2003; @mouelhi_stepwise_2006]</span></li>
<li><code>RunModel_GR1A()</code>: one-parameter yearly lumped hydrological model <span class="citation">[@mouelhi_vers_2003; @mouelhi_linking_2006]</span></li>
<li><code>RunModel_CemaNeige()</code>: two-parameter degree-day snowmelt and accumulation model <span class="citation">[@valery_as_2014]</span></li>
<li><code>RunModel_CemaNeigeGR4J()</code>: combined use of <strong>GR4J</strong> and <strong>CemaNeige</strong></li>
<li><code>RunModel_CemaNeigeGR5J()</code>: combined use of <strong>GR5J</strong> and <strong>CemaNeige</strong></li>
<li><code>RunModel_CemaNeigeGR6J()</code>: combined use of <strong>GR6J</strong> and <strong>CemaNeige</strong></li>
</ul>
<p>The <a href="https://webgr.irstea.fr/modeles/modele-de-prevision-grp/?lang=en"><strong>GRP</strong></a> forecasting model and the <a href="https://webgr.irstea.fr/modeles/otamin/?lang=en"><strong>Otamin</strong></a> predictive uncertainty tool are not availabe in <strong>airGR</strong>.</p>
</div>
<div id="loading-data" class="section level1">
<h1>Loading data</h1>
<p>In the following example, we use a data sample contained in the package. For real applications, the user has to import its data into <strong>R</strong> and to prepare it with an adequate data.frame format as described below.</p>
<p>First, it is necessary to load the <strong>airGR</strong> package:</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="kw">library</span>(airGR)</code></pre></div>
<p>This is an example of a <code>data.frame</code> of hydrometeorological observations time series for a fictional catchment included in the <strong>airGR</strong> package that contains:</p>
<ul>
<li><em>DatesR</em>: dates in the POSIXt format</li>
<li><em>P</em>: average precipitation [mm/day]</li>
<li><em>T</em>: catchment average air temperature [℃]</li>
<li><em>E</em>: catchment average potential evapotranspiration [mm/day]</li>
<li><em>Qls</em>: outlet discharge [l/s]</li>
<li><em>Qmm</em>: outlet discharge [mm/day]</li>
</ul>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="kw">data</span>(L0123001)
<span class="kw">summary</span>(BasinObs)</code></pre></div>
<pre><code>## DatesR P T E
## Min. :1984-01-01 Min. : 0.000 Min. :-18.700 Min. :0.000
## 1st Qu.:1991-04-02 1st Qu.: 0.000 1st Qu.: 4.100 1st Qu.:0.600
## Median :1998-07-02 Median : 0.300 Median : 9.100 Median :1.400
## Mean :1998-07-02 Mean : 2.915 Mean : 9.147 Mean :1.764
## 3rd Qu.:2005-10-01 3rd Qu.: 3.600 3rd Qu.: 14.500 3rd Qu.:2.900
## Max. :2012-12-31 Max. :66.800 Max. : 28.400 Max. :5.500
##
## Qls Qmm
## Min. : 70 Min. : 0.0168
## 1st Qu.: 1643 1st Qu.: 0.3943
## Median : 4070 Median : 0.9768
## Mean : 6134 Mean : 1.4720
## 3rd Qu.: 7889 3rd Qu.: 1.8933
## Max. :99500 Max. :23.8800
## NA's :755 NA's :755</code></pre>
<p>The usual functions (e.g. <code>read.table()</code>) can be used to load real-case datasets.</p>
</div>
<div id="preparation-of-functions-inputs" class="section level1">
<h1>Preparation of functions inputs</h1>
<p>To run a model, the functions of the <strong>airGR</strong> package (e.g. the models, calibration and criteria calculation functions) require data and options with specific formats.</p>
<p>To facilitate the use of the package, there are several functions dedicated to the creation of these objects:</p>
<ul>
<li><code>CreateInputsModel()</code>: prepares the inputs for the different hydrological models (times series of dates, rainfall, observed streamflow, etc.)</li>
<li><code>CreateRunOptions()</code>: prepares the options for the hydrological model run (warm-up period, calibration period, etc.)</li>
<li><code>CreateInputsCrit()</code>: prepares the options in order to compute efficiency criterions (choice of the criterion, choice of the transformation on streamflows: “log”, “root”, etc.)</li>
<li><code>CreateCalibOptions()</code>: prepares the options for the hydrological model calibration algorithm (choice of parameters to optimize, predefined values for uncalibrated parameters, etc.)</li>
</ul>
<div id="inputsmodel-object" class="section level2">
<h2>InputsModel object</h2>
<p>To run a GR hydrologic model, the user has to prepare the input data with the <code>CreateInputsModel()</code> function. As arguments, this function needs the function name corresponding to the model the user wants to run, a vector of dates, a vector of precipitation and a vector of potential evapotranspiration.</p>
<p>In the example below, we already have the potential evapotranspiration. If the user doesn’t have these data, it is possible to compute it with the <a href="http://dx.doi.org/10.1016/j.jhydrol.2004.08.026">Oudin’s formula</a> with the <code>PEdaily_Oudin()</code> function (this function only needs julian days, daily average air temperature and latitude).</p>
<p>Missing values (<code>NA</code>) of precipitation (or potential evapotranspiration) are <strong>not allowed</strong>.</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">InputsModel &lt;-<span class="st"> </span><span class="kw">CreateInputsModel</span>(<span class="dt">FUN_MOD =</span> RunModel_GR4J, <span class="dt">DatesR =</span> BasinObs$DatesR,
<span class="dt">Precip =</span> BasinObs$P, <span class="dt">PotEvap =</span> BasinObs$E)
<span class="kw">str</span>(InputsModel)</code></pre></div>
<pre><code>## List of 3
## $ DatesR : POSIXlt[1:10593], format: &quot;1984-01-01&quot; &quot;1984-01-02&quot; ...
## $ Precip : num [1:10593] 4.1 15.9 0.8 0 0 0 0 0 2.9 0 ...
## $ PotEvap: num [1:10593] 0.2 0.2 0.3 0.3 0.1 0.3 0.4 0.4 0.5 0.5 ...
## - attr(*, &quot;class&quot;)= chr [1:3] &quot;InputsModel&quot; &quot;daily&quot; &quot;GR&quot;</code></pre>
</div>
<div id="runoptions-object" class="section level2">
<h2>RunOptions object</h2>
<p>The <code>CreateRunOptions()</code> function allows to prepare the options required to the <code>RunModel*()</code> functions, which are the actual models functions.</p>
<p>The user must at least define the following arguments:</p>
<ul>
<li><code>FUN_MOD</code>: the name of the model function to run</li>
<li><code>InputsModel</code>: the associated inputs data</li>
<li><code>IndPeriod_Run</code>: the period on which the model is run</li>
</ul>
<p>To select a period for which the user wants to run the model, select the corresponding indexes for different time periods (not the POSIXt dates), as follows:</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">Ind_Run &lt;-<span class="st"> </span><span class="kw">seq</span>(<span class="kw">which</span>(<span class="kw">format</span>(BasinObs$DatesR, <span class="dt">format =</span> <span class="st">&quot;%d/%m/%Y %H:%M&quot;</span>)==<span class="st">&quot;01/01/1990 00:00&quot;</span>),
<span class="kw">which</span>(<span class="kw">format</span>(BasinObs$DatesR, <span class="dt">format =</span> <span class="st">&quot;%d/%m/%Y %H:%M&quot;</span>)==<span class="st">&quot;31/12/1999 00:00&quot;</span>))
<span class="kw">str</span>(Ind_Run)</code></pre></div>
<pre><code>## int [1:3652] 2193 2194 2195 2196 2197 2198 2199 2200 2201 2202 ...</code></pre>
<p>The intialization of hydrological models is of the utmost importance. Indeed, an inaccurate intialisation causes poor quality streamflow simulations during the earliest stages of the running period. For example, in the GR models, the reservoirs levels are by default set to 50 % of their capacity, which may be far from their ideal value. Two solutions are offered to accurately initialize the GR models in <strong>airGR</strong>: manually predefining the initial states or running the models during a warm up period before the actual running period. It is generally advised to set up this warm up period to be equal or superior to one year.</p>
<p>As a consequence, it is possible to define in <code>CreateRunOptions()</code> the following arguments:</p>
<ul>
<li><code>IniStates</code>: the initial states of the 2 unit hydrographs (20 + 40 = 60 units)</li>
<li><code>IniResLevels</code>: the initial levels of the production and routing stores</li>
<li><code>IndPeriod_WarmUp</code>: the warm-up period used to run the model, to be defined in the same format as <code>IndPeriod_Run</code></li>
</ul>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">RunOptions &lt;-<span class="st"> </span><span class="kw">CreateRunOptions</span>(<span class="dt">FUN_MOD =</span> RunModel_GR4J,
<span class="dt">InputsModel =</span> InputsModel, <span class="dt">IndPeriod_Run =</span> Ind_Run,
<span class="dt">IniStates =</span> <span class="ot">NULL</span>, <span class="dt">IniResLevels =</span> <span class="ot">NULL</span>, <span class="dt">IndPeriod_WarmUp =</span> <span class="ot">NULL</span>)</code></pre></div>
<pre><code>## Warning in CreateRunOptions(FUN_MOD = RunModel_GR4J, InputsModel = InputsModel, : Model warm-up period not defined -&gt; default configuration used
## The year preceding the run period is used</code></pre>
<pre><code>## Warning in CreateRunOptions(FUN_MOD = RunModel_GR4J, InputsModel = InputsModel, : Model states initialisation not defined -&gt; default configuration used</code></pre>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="kw">str</span>(RunOptions)</code></pre></div>
<pre><code>## List of 6
## $ IndPeriod_WarmUp: int [1:365] 1828 1829 1830 1831 1832 1833 1834 1835 1836 1837 ...
## $ IndPeriod_Run : int [1:3652] 2193 2194 2195 2196 2197 2198 2199 2200 2201 2202 ...
## $ IniStates : num [1:67] 0 0 0 0 0 0 0 0 0 0 ...
## $ IniResLevels : num [1:2] 0.3 0.5
## $ Outputs_Cal : chr &quot;Qsim&quot;
## $ Outputs_Sim : chr [1:16] &quot;DatesR&quot; &quot;PotEvap&quot; &quot;Precip&quot; &quot;Prod&quot; ...
## - attr(*, &quot;class&quot;)= chr [1:3] &quot;RunOptions&quot; &quot;GR&quot; &quot;daily&quot;</code></pre>
<p>The <code>CreateRunOptions()</code> function returns warnings if the default initialization options are used:</p>
<ul>
<li><code>IniStates</code> and <code>IniResLevels</code> are automatically set to initialize all the model states at 0, except for the production and routing stores which are initialised at 50 % of their capacity</li>
<li><code>IndPeriod_WarmUp</code> default setting ensures a one-year warm-up using the time steps preceding the <code>IndPeriod_Run</code>, if available</li>
</ul>
</div>
<div id="inputscrit-object" class="section level2">
<h2>InputsCrit object</h2>
<p>The <code>CreateInputsCrit()</code> function allows to prepare the input in order to calculate a criterion. It is possible to define the following arguments:</p>
<ul>
<li><code>FUN_CRIT</code>: the name of the error criterion function (they are introduced later on)</li>
<li><code>InputsModel</code>: the inputs of the hydrological model previously prepared by the <code>CeateInputsModel()</code> function</li>
<li><code>RunOptions</code>: the options of the hydrological model previously prepared by the <code>CreateRunOptions()</code> function</li>
<li><code>Qobs</code>: the observed streamflows expressed in <em>mm/time step</em></li>
</ul>
<p>Missing values (<code>NA</code>) are <strong>allowed</strong> for observed streamflows.</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">InputsCrit &lt;-<span class="st"> </span><span class="kw">CreateInputsCrit</span>(<span class="dt">FUN_CRIT =</span> ErrorCrit_NSE, <span class="dt">InputsModel =</span> InputsModel,
<span class="dt">RunOptions =</span> RunOptions, <span class="dt">Qobs =</span> BasinObs$Qmm[Ind_Run])
<span class="kw">str</span>(InputsCrit)</code></pre></div>
<pre><code>## List of 5
## $ BoolCrit : logi [1:3652] TRUE TRUE TRUE TRUE TRUE TRUE ...
## $ Qobs : num [1:3652] 1.99 1.8 2.86 2.4 3.31 ...
## $ transfo : chr &quot;&quot;
## $ Ind_zeroes: NULL
## $ epsilon : NULL
## - attr(*, &quot;class&quot;)= chr &quot;InputsCrit&quot;</code></pre>
</div>
<div id="caliboptions-object" class="section level2">
<h2>CalibOptions object</h2>
<p>The user needs to prepare the calibration options with the <code>CreateCalibOptions()</code>function. For that, it is necessary to define the following arguments:</p>
<ul>
<li><code>FUN_MOD</code>: the name of the model function</li>
<li><code>FUN_CALIB</code>: the name of the calibration algorithm</li>
</ul>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">CalibOptions &lt;-<span class="st"> </span><span class="kw">CreateCalibOptions</span>(<span class="dt">FUN_MOD =</span> RunModel_GR4J, <span class="dt">FUN_CALIB =</span> Calibration_Michel)
<span class="kw">str</span>(CalibOptions)</code></pre></div>
<pre><code>## List of 4
## $ FixedParam : logi [1:4] NA NA NA NA
## $ SearchRanges : num [1:2, 1:4] 4.59e-05 2.18e+04 -1.09e+04 1.09e+04 4.59e-05 ...
## $ StartParam : num [1:4] 247.151 -0.649 42.098 1.944
## $ StartParamDistrib: num [1:3, 1:4] 169.017 247.151 432.681 -2.376 -0.649 ...
## - attr(*, &quot;class&quot;)= chr [1:3] &quot;CalibOptions&quot; &quot;GR4J&quot; &quot;HBAN&quot;</code></pre>
</div>
</div>
<div id="criteria" class="section level1">
<h1>Criteria</h1>
<p>The evaluation of the quality of a simulation is estimated through the calculation of criteria. These criteria can be used both as objective-functions during the calibration of the model, or as a measure for evaluating its control performance.</p>
<p>The package offers the possibility to use different criteria:</p>
<ul>
<li><code>ErrorCrit_RMSE()</code>: Root mean square error (RMSE)</li>
<li><code>ErrorCrit_NSE()</code>: Nash-Sutcliffe model efficiency coefficient (NSE)</li>
<li><code>ErrorCrit_KGE()</code>: Kling-Gupta efficiency criterion (KGE)</li>
<li><code>ErrorCrit_KGE2()</code>: modified Kling-Gupta efficiency criterion (KGE’)</li>
</ul>
<p>It is also possible to create user-defined criteria. For doing that, it is only necessary to define the function in <strong>R</strong> following the same syntax that the criteria functions included in <strong>airGR</strong>.</p>
</div>
<div id="calibration" class="section level1">
<h1>Calibration</h1>
<p>The objective of the calibration algorithm is to identify the model parameters: by comparing the model outputs with observed data, this algorithm determines the combination of parameters that represents the best the behavior of the watershed.</p>
<p>In the <strong>airGR</strong> package, a function called <code>Calibration_Michel()</code> is implemented. This functions allows running a calibration with the package models. The calibration algorithm optimizes the error criterion selected as objective-function. This algorithm works in two steps:</p>
<ol style="list-style-type: decimal">
<li>a screening of the parameters space is performed using either a rough predefined grid or a user-defined list of parameter sets</li>
<li>a simple steepest descent local search algorithm is performed from the best set of parameters found at the first step</li>
</ol>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">OutputsCalib &lt;-<span class="st"> </span><span class="kw">Calibration_Michel</span>(<span class="dt">InputsModel =</span> InputsModel, <span class="dt">RunOptions =</span> RunOptions,
<span class="dt">InputsCrit =</span> InputsCrit, <span class="dt">CalibOptions =</span> CalibOptions,
<span class="dt">FUN_MOD =</span> RunModel_GR4J, <span class="dt">FUN_CRIT =</span> ErrorCrit_NSE)</code></pre></div>
<pre><code>## Grid-Screening in progress (0% 20% 40% 60% 80% 100%)
## Screening completed (81 runs):
## Param = 247.151 , -0.020 , 83.096 , 2.384
## Crit NSE[Q] = 0.7685
## Steepest-descent local search in progress
## Calibration completed (20 iterations, 226 runs):
## Param = 257.238 , 1.012 , 88.235 , 2.208
## Crit NSE[Q] = 0.7985</code></pre>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">Param &lt;-<span class="st"> </span>OutputsCalib$ParamFinalR
Param</code></pre></div>
<pre><code>## [1] 257.237556 1.012237 88.234673 2.207958</code></pre>
<p>The <code>Calibration_Michel()</code> function is the only one implemented in the <strong>airGR</strong> package to calibrate the model, but the user can implement its own calibration function.</p>
<p>This function returns a vector with the parameters of the chosen model, which means that the number of values can differ depending on the model that is used. It is possible to use the <code>Calibration_Michel()</code> function with user-implemented hydrological models.</p>
</div>
<div id="validation" class="section level1">
<h1>Validation</h1>
<p>This step assess the predictive capacity of the model. Validation is defined as the estimation of the accuracy of the model on datasets that are not used in its construction, and in particular its calibration. The classical way to perform a validation is to keep data from a period separated from the calibration period. If possible, this control period should correspond to climatic situations rather that differ from those of the calibration period in order to better point out the qualities and weakness of the model. This exercise is necessary for assessing the robustness of the model, that is to say its ability to keep stable performances outside the calibration conditions.</p>
<p>Performing a model validation with <strong>airGR</strong> is similar to running a simulation (see below).</p>
</div>
<div id="simulation" class="section level1">
<h1>Simulation</h1>
<div id="simulation-run" class="section level2">
<h2>Simulation run</h2>
<p>To run a model, the user has to use the <code>RunModel*()</code> functions (<code>InputsModel</code>, <code>RunOptions</code> and parameters). All the data needed have already been prepared.</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">OutputsModel &lt;-<span class="st"> </span><span class="kw">RunModel_GR4J</span>(<span class="dt">InputsModel =</span> InputsModel, <span class="dt">RunOptions =</span> RunOptions, <span class="dt">Param =</span> Param)
<span class="kw">str</span>(OutputsModel)</code></pre></div>
<pre><code>## List of 16
## $ DatesR : POSIXlt[1:3652], format: &quot;1990-01-01&quot; &quot;1990-01-02&quot; ...
## $ PotEvap : num [1:3652] 0.3 0.4 0.4 0.3 0.1 0.1 0.1 0.2 0.2 0.3 ...
## $ Precip : num [1:3652] 0 9.3 3.2 7.3 0 0 0 0 0.1 0.2 ...
## $ Prod : num [1:3652] 196 199 199 201 200 ...
## $ AE : num [1:3652] 0.2833 0.4 0.4 0.3 0.0952 ...
## $ Perc : num [1:3652] 0.645 0.696 0.703 0.74 0.725 ...
## $ PR : num [1:3652] 0.645 5.946 2.383 4.992 0.725 ...
## $ Q9 : num [1:3652] 1.78 1.52 3.86 3.17 3.45 ...
## $ Q1 : num [1:3652] 0.2 0.195 0.271 0.387 0.365 ...
## $ Rout : num [1:3652] 53.9 53.6 55.3 56.1 56.9 ...
## $ Exch : num [1:3652] 0.181 0.18 0.177 0.197 0.207 ...
## $ AExch : num [1:3652] 0.362 0.36 0.353 0.393 0.414 ...
## $ QR : num [1:3652] 2.05 1.99 2.36 2.55 2.78 ...
## $ QD : num [1:3652] 0.381 0.375 0.447 0.584 0.572 ...
## $ Qsim : num [1:3652] 2.43 2.37 2.8 3.14 3.35 ...
## $ StateEnd: num [1:67] 188.5 48.9 NA NA NA ...
## - attr(*, &quot;class&quot;)= chr [1:3] &quot;OutputsModel&quot; &quot;daily&quot; &quot;GR&quot;</code></pre>
</div>
<div id="results-preview" class="section level2">
<h2>Results preview</h2>
<p>Although it is possible for the user to design its own graphics from the outputs of the <code>RunModel*()</code> functions, the <strong>airGR</strong> package offers the possibility to make use of the <code>plot_OutputsModel()</code> function. This function returns a dashboard of results including various graphs (depending on the model used):</p>
<ul>
<li>time series of total precipitation and simulated streamflows (and observed streamflows if provided)</li>
<li>interannual median monthly simulated streamflows (and monthly observed streamflows if provided)</li>
<li>correlation plot between simulated and observed streamflows (if observed streamflows provided)</li>
<li>cumulative frequency plot for simulated streamflows (and for observed streamflows if provided)</li>
</ul>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r"><span class="kw">plot_OutputsModel</span>(<span class="dt">OutputsModel =</span> OutputsModel, <span class="dt">Qobs =</span> BasinObs$Qmm[Ind_Run])</code></pre></div>
<p>Moreover, if the CemaNeige model is used, the simulated snowpack time series are plotted.</p>
</div>
<div id="efficiency-criterion" class="section level2">
<h2>Efficiency criterion</h2>
<p>To evaluate the efficiency of the model, it is possible to use the same criterion as defined at the calibration step or to use an other criterion.</p>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">OutputsCrit &lt;-<span class="st"> </span><span class="kw">ErrorCrit_NSE</span>(<span class="dt">InputsCrit =</span> InputsCrit, <span class="dt">OutputsModel =</span> OutputsModel)
<span class="kw">str</span>(OutputsCrit)</code></pre></div>
<pre><code>## List of 5
## $ CritValue : num 0.799
## $ CritName : chr &quot;NSE[Q]&quot;
## $ CritBestValue : num 1
## $ Multiplier : num -1
## $ Ind_notcomputed: int [1:40] 2405 2406 2407 2408 2409 2410 2411 2412 2413 2414 ...</code></pre>
<div class="sourceCode"><pre class="sourceCode r"><code class="sourceCode r">OutputsCrit &lt;-<span class="st"> </span><span class="kw">ErrorCrit_KGE</span>(<span class="dt">InputsCrit =</span> InputsCrit, <span class="dt">OutputsModel =</span> OutputsModel)
<span class="kw">str</span>(OutputsCrit)</code></pre></div>
<pre><code>## List of 7
## $ CritValue : num 0.785
## $ CritName : chr &quot;KGE[Q]&quot;
## $ SubCritValues : num [1:3] 0.898 0.816 1.044
## $ SubCritNames : chr [1:3] &quot;KGE[Q] rPEARSON(sim vs. obs)&quot; &quot;KGE[Q] STDEVsim/STDEVobs&quot; &quot;KGE[Q] MEANsim/MEANobs&quot;
## $ CritBestValue : num 1
## $ Multiplier : num -1
## $ Ind_notcomputed: int [1:40] 2405 2406 2407 2408 2409 2410 2411 2412 2413 2414 ...</code></pre>
</div>
</div>
<!-- dynamically load mathjax for compatibility with self-contained -->
<script>
(function () {
var script = document.createElement("script");
script.type = "text/javascript";
script.src = "https://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML";
document.getElementsByTagName("head")[0].appendChild(script);
})();
</script>
</body>
</html>
vignettes/airGR.Rmd 0 → 100644
+ 259
− 0
View file @ 6a61568c
---
title: "airGR -- Overview"
output: rmarkdown::html_vignette
vignette: >
%\VignetteEngine{knitr::rmarkdown}
%\VignetteIndexEntry{airGR}
%\VignetteEncoding{UTF-8}
---
# Introduction
**airGR** is a package which brings into the [**R software**](https://cran.r-project.org/) the hydrological modelling tools used and developed at the [Catchment Hydrology Team](https://webgr.irstea.fr/?lang=en) at [Irstea (France)](http://www.irstea.fr/en/), including the [**GR rainfall-runoff models**](https://webgr.irstea.fr/modeles/?lang=en) and a snowmelt and accumulation model, **CemaNeige**. Each model core is coded in **Fortran** to ensure low computational time. The other package functions (i.e. mainly the calibration algorithm and the efficiency criteria calculation) are coded in **R**.
The **airGR** package has been designed to fulfil two major requirements: to facilitate the use by non-expert users and to allow flexibility regarding the addition of external criteria, models or calibration algorithms. The names of the functions and their arguments were chosen to this end. **airGR** also contains basics plotting facilities.
Six hydrological models and one snowmelt and accumulation model are implemented in **airGR**. The snow model can be used alone or together with the daily hydrological models.
The models can be called within **airGR** using the following functions:
* `RunModel_GR4H()`: four-parameter hourly lumped hydrological model [@mathevet_quels_2005]
* `RunModel_GR4J()`: four-parameter daily lumped hydrological model [@perrin_improvement_2003]
* `RunModel_GR5J()`: five-parameter daily lumped hydrological model [@le_moine_bassin_2008]
* `RunModel_GR6J()`: six-parameter daily lumped hydrological model [@pushpalatha_downward_2011]
* `RunModel_GR2M()`: two-parameter monthly lumped hydrological model [@mouelhi_vers_2003; @mouelhi_stepwise_2006]
* `RunModel_GR1A()`: one-parameter yearly lumped hydrological model [@mouelhi_vers_2003; @mouelhi_linking_2006]
* `RunModel_CemaNeige()`: two-parameter degree-day snowmelt and accumulation model [@valery_as_2014]
* `RunModel_CemaNeigeGR4J()`: combined use of **GR4J** and **CemaNeige**
* `RunModel_CemaNeigeGR5J()`: combined use of **GR5J** and **CemaNeige**
* `RunModel_CemaNeigeGR6J()`: combined use of **GR6J** and **CemaNeige**
The [**GRP**](https://webgr.irstea.fr/modeles/modele-de-prevision-grp/?lang=en) forecasting model and the [**Otamin**](https://webgr.irstea.fr/modeles/otamin/?lang=en) predictive uncertainty tool are not availabe in **airGR**.
# Loading data
In the following example, we use a data sample contained in the package. For real applications, the user has to import its data into **R** and to prepare it with an adequate data.frame format as described below.
First, it is necessary to load the **airGR** package:
```{r}
library(airGR)
```
This is an example of a `data.frame` of hydrometeorological observations time series for a fictional catchment included in the **airGR** package that contains:
* *DatesR*: dates in the POSIXt format
* *P*: average precipitation [mm/day]
* *T*: catchment average air temperature [&#8451;]
* *E*: catchment average potential evapotranspiration [mm/day]
* *Qls*: outlet discharge [l/s]
* *Qmm*: outlet discharge [mm/day]
```{r}
data(L0123001)
summary(BasinObs)
```
The usual functions (e.g. `read.table()`) can be used to load real-case datasets.
# Preparation of functions inputs
To run a model, the functions of the **airGR** package (e.g. the models, calibration and criteria calculation functions) require data and options with specific formats.
To facilitate the use of the package, there are several functions dedicated to the creation of these objects:
* `CreateInputsModel()`: prepares the inputs for the different hydrological models (times series of dates, rainfall, observed streamflow, etc.)
* `CreateRunOptions()`: prepares the options for the hydrological model run (warm-up period, calibration period, etc.)
* `CreateInputsCrit()`: prepares the options in order to compute efficiency criterions (choice of the criterion, choice of the transformation on streamflows: "log", "root", etc.)
* `CreateCalibOptions()`: prepares the options for the hydrological model calibration algorithm (choice of parameters to optimize, predefined values for uncalibrated parameters, etc.)
## InputsModel object
To run a GR hydrologic model, the user has to prepare the input data with the `CreateInputsModel()` function.
As arguments, this function needs the function name corresponding to the model the user wants to run, a vector of dates, a vector of precipitation and a vector of potential evapotranspiration.
In the example below, we already have the potential evapotranspiration. If the user doesn't have these data, it is possible to compute it with the [Oudin's formula](http://dx.doi.org/10.1016/j.jhydrol.2004.08.026) with the `PEdaily_Oudin()` function (this function only needs julian days, daily average air temperature and latitude).
Missing values (`NA`) of precipitation (or potential evapotranspiration) are **not allowed**.
```{r}
InputsModel <- CreateInputsModel(FUN_MOD = RunModel_GR4J, DatesR = BasinObs$DatesR,
Precip = BasinObs$P, PotEvap = BasinObs$E)
str(InputsModel)
```
## RunOptions object
The `CreateRunOptions()` function allows to prepare the options required to the `RunModel*()` functions, which are the actual models functions.
The user must at least define the following arguments:
* `FUN_MOD`: the name of the model function to run
* `InputsModel`: the associated inputs data
* `IndPeriod_Run`: the period on which the model is run
To select a period for which the user wants to run the model, select the corresponding indexes for different time periods (not the POSIXt dates), as follows:
```{r}
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"))
str(Ind_Run)
```
The intialization of hydrological models is of the utmost importance. Indeed, an inaccurate intialisation causes poor quality streamflow simulations during the earliest stages of the running period. For example, in the GR models, the reservoirs levels are by default set to 50 % of their capacity, which may be far from their ideal value. Two solutions are offered to accurately initialize the GR models in **airGR**: manually predefining the initial states or running the models during a warm up period before the actual running period. It is generally advised to set up this warm up period to be equal or superior to one year.
As a consequence, it is possible to define in `CreateRunOptions()` the following arguments:
* `IniStates`: the initial states of the 2 unit hydrographs (20 + 40 = 60 units)
* `IniResLevels`: the initial levels of the production and routing stores
* `IndPeriod_WarmUp`: the warm-up period used to run the model, to be defined in the same format as `IndPeriod_Run`
```{r}
RunOptions <- CreateRunOptions(FUN_MOD = RunModel_GR4J,
InputsModel = InputsModel, IndPeriod_Run = Ind_Run,
IniStates = NULL, IniResLevels = NULL, IndPeriod_WarmUp = NULL)
str(RunOptions)
```
The `CreateRunOptions()` function returns warnings if the default initialization options are used:
* `IniStates` and `IniResLevels` are automatically set to initialize all the model states at 0, except for the production and routing stores which are initialised at 50 % of their capacity
* `IndPeriod_WarmUp` default setting ensures a one-year warm-up using the time steps preceding the `IndPeriod_Run`, if available
## InputsCrit object
The `CreateInputsCrit()` function allows to prepare the input in order to calculate a criterion. It is possible to define the following arguments:
* `FUN_CRIT`: the name of the error criterion function (they are introduced later on)
* `InputsModel`: the inputs of the hydrological model previously prepared by the `CeateInputsModel()` function
* `RunOptions`: the options of the hydrological model previously prepared by the `CreateRunOptions()` function
* `Qobs`: the observed streamflows expressed in *mm/time step*
Missing values (`NA`) are **allowed** for observed streamflows.
```{r}
InputsCrit <- CreateInputsCrit(FUN_CRIT = ErrorCrit_NSE, InputsModel = InputsModel,
RunOptions = RunOptions, Qobs = BasinObs$Qmm[Ind_Run])
str(InputsCrit)
```
## CalibOptions object
The user needs to prepare the calibration options with the `CreateCalibOptions()`function. For that, it is necessary to define the following arguments:
* `FUN_MOD`: the name of the model function
* `FUN_CALIB`: the name of the calibration algorithm
```{r}
CalibOptions <- CreateCalibOptions(FUN_MOD = RunModel_GR4J, FUN_CALIB = Calibration_Michel)
str(CalibOptions)
```
# Criteria
The evaluation of the quality of a simulation is estimated through the calculation of criteria. These criteria can be used both as objective-functions during the calibration of the model, or as a measure for evaluating its control performance.
The package offers the possibility to use different criteria:
* `ErrorCrit_RMSE()`: Root mean square error (RMSE)
* `ErrorCrit_NSE()`: Nash-Sutcliffe model efficiency coefficient (NSE)
* `ErrorCrit_KGE()`: Kling-Gupta efficiency criterion (KGE)
* `ErrorCrit_KGE2()`: modified Kling-Gupta efficiency criterion (KGE')
It is also possible to create user-defined criteria. For doing that, it is only necessary to define the function in **R** following the same syntax that the criteria functions included in **airGR**.
# Calibration
The objective of the calibration algorithm is to identify the model parameters: by comparing the model outputs with observed data, this algorithm determines the combination of parameters that represents the best the behavior of the watershed.
In the **airGR** package, a function called `Calibration_Michel()` is implemented. This functions allows running a calibration with the package models.
The calibration algorithm optimizes the error criterion selected as objective-function. This algorithm works in two steps:
1. a screening of the parameters space is performed using either a rough predefined grid or a user-defined list of parameter sets
2. a simple steepest descent local search algorithm is performed from the best set of parameters found at the first step
```{r}
OutputsCalib <- Calibration_Michel(InputsModel = InputsModel, RunOptions = RunOptions,
InputsCrit = InputsCrit, CalibOptions = CalibOptions,
FUN_MOD = RunModel_GR4J, FUN_CRIT = ErrorCrit_NSE)
Param <- OutputsCalib$ParamFinalR
Param
```
The `Calibration_Michel()` function is the only one implemented in the **airGR** package to calibrate the model, but the user can implement its own calibration function.
This function returns a vector with the parameters of the chosen model, which means that the number of values can differ depending on the model that is used. It is possible to use the `Calibration_Michel()` function with user-implemented hydrological models.
# Validation
This step assess the predictive capacity of the model. Validation is defined as the estimation of the accuracy of the model on datasets that are not used in its construction, and in particular its calibration.
The classical way to perform a validation is to keep data from a period separated from the calibration period. If possible, this control period should correspond to climatic situations rather that differ from those of the calibration period in order to better point out the qualities and weakness of the model. This exercise is necessary for assessing the robustness of the model, that is to say its ability to keep stable performances outside the calibration conditions.
Performing a model validation with **airGR** is similar to running a simulation (see below).
# Simulation
## Simulation run
To run a model, the user has to use the `RunModel*()` functions (`InputsModel`, `RunOptions` and parameters).
All the data needed have already been prepared.
```{r}
OutputsModel <- RunModel_GR4J(InputsModel = InputsModel, RunOptions = RunOptions, Param = Param)
str(OutputsModel)
```
## Results preview
Although it is possible for the user to design its own graphics from the outputs of the `RunModel*()` functions, the **airGR** package offers the possibility to make use of the `plot_OutputsModel()` function. This function returns a dashboard of results including various graphs (depending on the model used):
* time series of total precipitation and simulated streamflows (and observed streamflows if provided)
* interannual median monthly simulated streamflows (and monthly observed streamflows if provided)
* correlation plot between simulated and observed streamflows (if observed streamflows provided)
* cumulative frequency plot for simulated streamflows (and for observed streamflows if provided)
```{r,eval=F}
plot_OutputsModel(OutputsModel = OutputsModel, Qobs = BasinObs$Qmm[Ind_Run])
```
Moreover, if the CemaNeige model is used, the simulated snowpack time series are plotted.
## Efficiency criterion
To evaluate the efficiency of the model, it is possible to use the same criterion as defined at the calibration step or to use an other criterion.
```{r}
OutputsCrit <- ErrorCrit_NSE(InputsCrit = InputsCrit, OutputsModel = OutputsModel)
str(OutputsCrit)
```
```{r}
OutputsCrit <- ErrorCrit_KGE(InputsCrit = InputsCrit, OutputsModel = OutputsModel)
str(OutputsCrit)
```
Supports Markdown
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