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

docs(vignette): speed up the parameters optimization

- vignettes :
  - param_optim
  - param_mcmc
- RunModel outputs: return only Qsim
Reviewed-by: @francois.bourgin
Refs #101
Showing with 27 additions and 6 deletions
+27 -6
...@@ -45,6 +45,13 @@ Please note that the calibration period is defined in the `CreateRunOptions()` f ...@@ -45,6 +45,13 @@ Please note that the calibration period is defined in the `CreateRunOptions()` f
example("Calibration_Michel") example("Calibration_Michel")
``` ```
In order for the `RunModel_*()` functions to run faster during the parameter estimation process, it is recommended that the outputs contain only the simulated flows (see the `Outputs_Sim` argument in the `CreateRunOptions()` help page).
```{r, results='hide', eval=FALSE}
RunOptions <- airGR::CreateRunOptions(FUN_MOD = RunModel_GR4J, InputsModel = InputsModel,
IndPeriod_Run = Ind_Run,
Outputs_Sim = "Qsim")
```
Regarding the different optimization strategies presented here, we refer to each package for in-depth information about the description of the methods used. Regarding the different optimization strategies presented here, we refer to each package for in-depth information about the description of the methods used.
...@@ -58,7 +65,6 @@ Here we choose to minimize the root mean square error. ...@@ -58,7 +65,6 @@ Here we choose to minimize the root mean square error.
The change of the repository from the "real" parameter space to a "transformed" space ensures homogeneity of displacement in the different dimensions of the parameter space during the step-by-step procedure of the calibration algorithm of the model. The change of the repository from the "real" parameter space to a "transformed" space ensures homogeneity of displacement in the different dimensions of the parameter space during the step-by-step procedure of the calibration algorithm of the model.
```{r, warning=FALSE, results='hide', eval=FALSE} ```{r, warning=FALSE, results='hide', eval=FALSE}
OptimGR4J <- function(ParamOptim) { OptimGR4J <- function(ParamOptim) {
## Transformation of the parameter set to real space ## Transformation of the parameter set to real space
...@@ -142,6 +148,7 @@ Here we use the following R implementation of some popular strategies: ...@@ -142,6 +148,7 @@ Here we use the following R implementation of some popular strategies:
* [Rmalschains: memetic algorithms](https://cran.r-project.org/package=Rmalschains) * [Rmalschains: memetic algorithms](https://cran.r-project.org/package=Rmalschains)
## Differential Evolution ## Differential Evolution
```{r, warning=FALSE, results='hide', eval=FALSE} ```{r, warning=FALSE, results='hide', eval=FALSE}
optDE <- DEoptim::DEoptim(fn = OptimGR4J, optDE <- DEoptim::DEoptim(fn = OptimGR4J,
lower = lowerGR4J, upper = upperGR4J, lower = lowerGR4J, upper = upperGR4J,
...@@ -150,6 +157,7 @@ optDE <- DEoptim::DEoptim(fn = OptimGR4J, ...@@ -150,6 +157,7 @@ optDE <- DEoptim::DEoptim(fn = OptimGR4J,
## Particle Swarm ## Particle Swarm
```{r, warning=FALSE, results='hide', message=FALSE, eval=FALSE} ```{r, warning=FALSE, results='hide', message=FALSE, eval=FALSE}
optPSO <- hydroPSO::hydroPSO(fn = OptimGR4J, optPSO <- hydroPSO::hydroPSO(fn = OptimGR4J,
lower = lowerGR4J, upper = upperGR4J, lower = lowerGR4J, upper = upperGR4J,
...@@ -157,6 +165,7 @@ optPSO <- hydroPSO::hydroPSO(fn = OptimGR4J, ...@@ -157,6 +165,7 @@ optPSO <- hydroPSO::hydroPSO(fn = OptimGR4J,
``` ```
## MA-LS-Chains ## MA-LS-Chains
```{r, warning=FALSE, results='hide', eval=FALSE} ```{r, warning=FALSE, results='hide', eval=FALSE}
optMALS <- Rmalschains::malschains(fn = OptimGR4J, optMALS <- Rmalschains::malschains(fn = OptimGR4J,
lower = lowerGR4J, upper = upperGR4J, lower = lowerGR4J, upper = upperGR4J,
...@@ -222,6 +231,7 @@ MOptimGR4J <- function(i) { ...@@ -222,6 +231,7 @@ MOptimGR4J <- function(i) {
``` ```
## caRamel ## caRamel
caRamel is a multiobjective evolutionary algorithm combining the MEAS algorithm and the NGSA-II algorithm. caRamel is a multiobjective evolutionary algorithm combining the MEAS algorithm and the NGSA-II algorithm.
```{r, warning=FALSE, results='hide', eval=FALSE} ```{r, warning=FALSE, results='hide', eval=FALSE}
......
...@@ -38,6 +38,14 @@ We use the **GR4J** model and we assume that the R global environment contains d ...@@ -38,6 +38,14 @@ We use the **GR4J** model and we assume that the R global environment contains d
example("Calibration_Michel") example("Calibration_Michel")
``` ```
In order for the `RunModel_*()` functions to run faster during the parameter estimation process, it is recommended that the outputs contain only the simulated flows (see the `Outputs_Sim` argument in the `CreateRunOptions()` help page).
```{r, results='hide', eval=FALSE}
RunOptions <- airGR::CreateRunOptions(FUN_MOD = RunModel_GR4J, InputsModel = InputsModel,
IndPeriod_Run = Ind_Run,
Outputs_Sim = "Qsim")
```
Please refer to the [2.1 Plugging in new calibration](V02.1_param_optim.html) vignette for explanations on how to plug in a parameter estimation algorithm to **airGR**. Please refer to the [2.1 Plugging in new calibration](V02.1_param_optim.html) vignette for explanations on how to plug in a parameter estimation algorithm to **airGR**.
...@@ -49,6 +57,7 @@ We show how to use the DRAM algorithm for SLS Bayesian inference, with the `modM ...@@ -49,6 +57,7 @@ We show how to use the DRAM algorithm for SLS Bayesian inference, with the `modM
First, we need to define a function that returns twice the opposite of the log-likelihood for a given parameter set. First, we need to define a function that returns twice the opposite of the log-likelihood for a given parameter set.
Nota: in the `LogLikeGR4J()` function, the computation of the log-likelihood is simplified in order to ensure a good computing performance. It corresponds to a translation of the two following lines. Nota: in the `LogLikeGR4J()` function, the computation of the log-likelihood is simplified in order to ensure a good computing performance. It corresponds to a translation of the two following lines.
```{r, echo=TRUE, eval=FALSE, purl=FALSE} ```{r, echo=TRUE, eval=FALSE, purl=FALSE}
Likelihood <- sum((ObsY - ModY)^2, na.rm = TRUE)^(-sum(!is.na(ObsY)) / 2) Likelihood <- sum((ObsY - ModY)^2, na.rm = TRUE)^(-sum(!is.na(ObsY)) / 2)
LogLike <- -2 * log(Likelihood) LogLike <- -2 * log(Likelihood)
...@@ -96,8 +105,10 @@ Nota: in this example, there are relatively few iterations (2000), in order to l ...@@ -96,8 +105,10 @@ Nota: in this example, there are relatively few iterations (2000), in order to l
With the DRAM algorithm, the covariance of the proposal is updated every 100 runs and delayed rejection is applied. With the DRAM algorithm, the covariance of the proposal is updated every 100 runs and delayed rejection is applied.
```{r, results='hide', eval=FALSE} ```{r, results='hide', eval=FALSE}
iniParPORT <- data.frame(Chain1 = iniParPORT, Chain2 = iniParPORT, Chain3 = iniParPORT, iniParPORT <- data.frame(Chain1 = iniParPORT,
row.names = paste0("X", 1:4)) Chain2 = iniParPORT,
Chain3 = iniParPORT,
row.names = paste0("X", 1:4))
iniParPORT <- sweep(iniParPORT, MARGIN = 2, STATS = c(1, 0.9, 1.1), FUN = "*") iniParPORT <- sweep(iniParPORT, MARGIN = 2, STATS = c(1, 0.9, 1.1), FUN = "*")
iniParPORT[iniParPORT < -9.9] <- -9.9 iniParPORT[iniParPORT < -9.9] <- -9.9
iniParPORT[iniParPORT > +9.9] <- +9.9 iniParPORT[iniParPORT > +9.9] <- +9.9
...@@ -111,8 +122,8 @@ mcmcDRAM <- apply(iniParPORT, MARGIN = 2, FUN = function(iIniParPORT) { ...@@ -111,8 +122,8 @@ mcmcDRAM <- apply(iniParPORT, MARGIN = 2, FUN = function(iIniParPORT) {
jump = 0.01, jump = 0.01,
outputlength = 2000, outputlength = 2000,
burninlength = 0, burninlength = 0,
updatecov = 100, ## Adaptative Metropolis updatecov = 100, ## adaptative Metropolis (AM)
ntrydr = 2) ## Delayed Rejection ntrydr = 2) ## delayed rejection (RD)
}) })
``` ```
...@@ -143,7 +154,7 @@ In addition, graphical tools can be used, with for example the [ggmcmc](https:// ...@@ -143,7 +154,7 @@ In addition, graphical tools can be used, with for example the [ggmcmc](https://
First, the evolution of the Markov chains can be seen with a traceplot: First, the evolution of the Markov chains can be seen with a traceplot:
```{r, fig.width=6, fig.height=9, warning=FALSE} ```{r, fig.width=6, fig.height=9, warning=FALSE}
parDRAM <- ggmcmc::ggs(multDRAM) ## to convert objet for using by all ggs_* graphical functions parDRAM <- ggmcmc::ggs(multDRAM) ## to convert object for using by all ggs_* graphical functions
ggmcmc::ggs_traceplot(parDRAM) ggmcmc::ggs_traceplot(parDRAM)
``` ```
......
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