Commit ba3123b3 authored by Delaigue Olivier's avatar Delaigue Olivier
Browse files

docs(vignette): save CaRamel results in a Rda file and make minor code revision

Refs #61
Showing with 29 additions and 17 deletions
+29 -17
......@@ -29,6 +29,7 @@ Imports:
Suggests:
knitr, rmarkdown,
caRamel, coda, DEoptim, dplyr, FME, ggmcmc, hydroPSO, imputeTS, Rmalschains,
ggplot2, GGally,
testthat
Description: Hydrological modelling tools developed at INRAE-Antony (HYCAR Research Unit, France). The package includes several conceptual rainfall-runoff models (GR4H, GR5H, GR4J, GR5J, GR6J, GR2M, GR1A), a snow accumulation and melt model (CemaNeige) and the associated functions for their calibration and evaluation. Use help(airGR) for package description and references.
License: GPL-2
......
......@@ -4,6 +4,10 @@
### 1.6.9.36 Release Notes (2021-01-28)
#### New features
- Added a section 'param_optim' vignette to explain how to manage with multiobjective optimization using the 'CaRamel' package. ([#61](https://gitlab.irstea.fr/HYCAR-Hydro/airgr/-/issues/61))
#### Major user-visible changes
- `Imax()` now returns an error message when `IndPeriod_Run` doesn't select 24 hours by day, instead of `numeric(0)`. ([#92](https://gitlab.irstea.fr/HYCAR-Hydro/airgr/-/issues/92))
......
File added
---
title: "Plugging in new calibration algorithms in airGR"
author: "François Bourgin, Guillaume Thirel"
output: rmarkdown::html_vignette
vignette: >
%\VignetteEngine{knitr::rmarkdown}
......@@ -20,6 +21,7 @@ library(GGally)
# source("airGR.R")
set.seed(321)
load(system.file("vignettesData/vignetteParamOptim.rda", package = "airGR"))
load(system.file("vignettesData/vignetteParamOptimCaramel.rda", package = "airGR"))
```
......@@ -160,7 +162,7 @@ As it can be seen in the table below, the four additional optimization strategie
```{r, warning=FALSE, echo=FALSE, eval=FALSE}
resGLOB <- data.frame(Algo = c("airGR", "PORT", "DE", "PSO", "MA-LS"),
round(rbind(
OutputsCalib$ParamFinalR ,
OutputsCalib$ParamFinalR,
airGR::TransfoParam_GR4J(ParamIn = optPORT$par , Direction = "TR"),
airGR::TransfoParam_GR4J(ParamIn = as.numeric(optDE$optim$bestmem), Direction = "TR"),
airGR::TransfoParam_GR4J(ParamIn = as.numeric(optPSO$par) , Direction = "TR"),
......@@ -189,9 +191,9 @@ InputsCrit_inv <- InputsCrit
InputsCrit_inv$transfo <- "inv"
MOptimGR4J <- function(i) {
if (algo == "caRamel") ParamOptim = x[i, ]
if (algo == "caRamel") {
ParamOptim <- x[i, ]
}
## Transformation of the parameter set to real space
RawParamOptim <- airGR::TransfoParam_GR4J(ParamIn = ParamOptim,
Direction = "TR")
......@@ -224,7 +226,7 @@ optMO <- caRamel::caRamel(nobj = 2,
popsize = 100,
archsize = 100,
maxrun = 15000,
prec = rep(1.e-3,2),
prec = rep(1.e-3, 2),
carallel = FALSE,
graph = FALSE)
```
......@@ -234,34 +236,39 @@ The algorithm returns parameter sets that describe the pareto front, illustratin
```{r, fig.width=6, fig.height=6, warning=FALSE}
ggplot() +
geom_point(aes(optMO$objectives[, 1], optMO$objectives[, 2])) +
coord_equal(xlim = c(0.4,0.9), ylim= c(0.4, 0.9)) +
xlab("KGE") + ylab("KGE_inv") +
coord_equal(xlim = c(0.4, 0.9), ylim = c(0.4, 0.9)) +
xlab("KGE") + ylab("KGE [1/Q]") +
theme_bw()
```
The pameter sets can be viewed in the parameter space, illustrating different populations.
```{r fig.height=6, fig.width=6, message=FALSE, warning=FALSE}
param_optMO = apply(optMO$parameters, 1, FUN = function(x) {
airGR::TransfoParam_GR4J(x, Direction = "TR")})
ggpairs(data.frame(t(param_optMO)), diag = NULL) + theme_bw()
param_optMO <- apply(optMO$parameters, MARGIN = 1, FUN = function(x) {
airGR::TransfoParam_GR4J(x, Direction = "TR")
})
GGally::ggpairs(data.frame(t(param_optMO)), diag = NULL) + theme_bw()
```
```{r fig.height=6, fig.width=12, message=FALSE, warning=FALSE}
RunOptions$Outputs_Sim = "Qsim"
run_optMO = apply(optMO$parameters, 1, FUN = function(x) {
RunOptions$Outputs_Sim <- "Qsim"
run_optMO <- apply(optMO$parameters, MARGIN = 1, FUN = function(x) {
airGR::RunModel_GR4J(InputsModel = InputsModel,
RunOptions = RunOptions,
x)}$Qsim)
run_optMO = data.frame(run_optMO)
Param = x)
}$Qsim)
run_optMO <- data.frame(run_optMO)
ggplot() +
geom_line(aes(x=as.POSIXct(InputsModel$DatesR[Ind_Run]), y=run_optMO$X1)) +
geom_line(aes(x=as.POSIXct(InputsModel$DatesR[Ind_Run]), y=run_optMO$X54), colour="darkred") +
geom_line(aes(x = as.POSIXct(InputsModel$DatesR[Ind_Run]),
y = run_optMO$X1)) +
geom_line(aes(x = as.POSIXct(InputsModel$DatesR[Ind_Run]),
y = run_optMO$X54),
colour = "darkred") +
scale_x_datetime(limits = c(as.POSIXct("1998-01-01"), NA)) +
ylab("Discharge [mm/d]") + xlab("Date") +
scale_y_sqrt() +
theme_classic()
theme_bw()
```
......
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