Commit 138dc72d authored by Delaigue Olivier's avatar Delaigue Olivier
Browse files

v1.0.14.0 NEW: add Vignette_Param datasets to reduce runtime of vignettes to pass CRAN check

Showing with 46 additions and 16 deletions
+46 -16
Package: airGR
Type: Package
Title: Suite of GR Hydrological Models for Precipitation-Runoff Modelling
Version: 1.0.13.19
Date: 2018-09-27
Version: 1.0.14.0
Date: 2018-09-28
Authors@R: c(
person("Laurent", "Coron", role = c("aut", "trl"), comment = c(ORCID = "0000-0002-1503-6204")),
person("Charles", "Perrin", role = c("aut", "ths"), comment = c(ORCID = "0000-0001-8552-1881")),
......
......@@ -14,7 +14,7 @@ output:
### 1.0.13.19 Release Notes (2018-09-27)
### 1.0.14.0 Release Notes (2018-09-28)
#### Deprecated and defunct
......@@ -57,6 +57,8 @@ output:
- As recomanded by CRAN managers, the NEWS file is now at the text format and is no more just a link to the airGR Website
- Added the <code>Vignette_Param.</code> datasets in order to reduce runtime during the re-building of vignettes. It contains different objects needed for param_optim and param_mcmc vignettes.
____________________________________________________________________________________
......
File added
\docType{data}
\encoding{UTF-8}
\name{Vignette_Param}
\alias{IniParam}
\alias{ListIniParam}
\alias{list_opt}
\alias{mcmcDRAM}
\alias{optDE}
\alias{optMALS}
\alias{optPORT}
\alias{optPORT_}
\alias{optPSO}
\alias{startGR4J}
\title{Datasets needed to run param_optim or param_mcmc vignettes}
\description{Datasets needed to run param_optim or param_mcmc vignettes.}
......@@ -17,6 +17,7 @@ library(hydroPSO)
library(Rmalschains)
# source("airGR.R")
set.seed(321)
data("Vignette_Param")
```
......@@ -80,9 +81,8 @@ upperGR4J <- rep(+9.99, times = 4)
# Local optimization
We start with a local optimization strategy by using the PORT routines (using the `nlminb()` of the `stats` package) and by setting a starting point in the transformed parameter space:
```{r, warning=FALSE, results='hide'}
startGR4J <- c(4.1, 3.9, -0.9, -8.7)
optPORT <- stats::nlminb(start = startGR4J,
```{r, warning=FALSE, results='hide', eval=FALSE}
optPORT <- stats::nlminb(start = c(4.1, 3.9, -0.9, -8.7),
objective = OptimGR4J,
lower = lowerGR4J, upper = upperGR4J,
control = list(trace = 1))
......@@ -92,7 +92,7 @@ The RMSE value reaches a local minimum value after 35 iterations.
We can also try a multi-start approach to test the consistency of the local optimization.
Here we use the same grid used for the filtering step of the Michel's calibration strategy (`Calibration_Michel()` function).
For each starting point, a local optimization is performed.
```{r, warning=FALSE, results='hide'}
```{r, warning=FALSE, results='hide', eval=FALSE}
startGR4J <- expand.grid(data.frame(CalibOptions$StartParamDistrib))
optPORT_ <- function(x) {
opt <- stats::nlminb(start = x,
......@@ -128,7 +128,7 @@ Here we use the following R implementation of some popular strategies:
* [Rmalschains: memetic algorithms](https://cran.r-project.org/package=Rmalschains)
## Differential Evolution
```{r, warning=FALSE, results='hide'}
```{r, warning=FALSE, results='hide', eval=FALSE}
optDE <- DEoptim::DEoptim(fn = OptimGR4J,
lower = lowerGR4J, upper = upperGR4J,
control = DEoptim::DEoptim.control(NP = 40, trace = 10))
......@@ -136,14 +136,14 @@ optDE <- DEoptim::DEoptim(fn = OptimGR4J,
## Particle Swarm
```{r, warning=FALSE, results='hide', message=FALSE}
```{r, warning=FALSE, results='hide', message=FALSE, eval=FALSE}
optPSO <- hydroPSO::hydroPSO(fn = OptimGR4J,
lower = lowerGR4J, upper = upperGR4J,
control = list(write2disk = FALSE, verbose = FALSE))
```
## MA-LS-Chains
```{r, warning=FALSE, results='hide'}
```{r, warning=FALSE, results='hide', eval=FALSE}
optMALS <- Rmalschains::malschains(fn = OptimGR4J,
lower = lowerGR4J, upper = upperGR4J,
maxEvals = 2000)
......
......@@ -18,6 +18,7 @@ library(ggmcmc)
library(dplyr)
# source("airGR.R")
set.seed(123)
data("Vignette_Param")
```
......@@ -49,9 +50,8 @@ First, we need to define a function that returns twice the opposite of the log-l
Nota: in the `RunAirGR4J()` 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}
Likelihood <- sum((ObsY - ModY)^2, na.rm = TRUE)^(-sum(!is.na(ObsY)) / 2)
LL <- -2 * log(Likelihood)
Likelihood <- sum((ObsY - ModY)^2, na.rm = TRUE)^(-sum(!is.na(ObsY)) / 2)
LL <- -2 * log(Likelihood)
```
In our simplified setting of Gaussian likelihood with measurement error integrated out, the log of the sum of squared error is related to the log-likelihood.
......@@ -78,7 +78,7 @@ RunAirGR4J <- function(Param_Optim) {
## Estimation of the best-fit parameters as a starting point
We start by using the PORT optimization routine to estimate the best-fit parameters.
```{r, results='hide'}
```{r, results='hide', eval=FALSE}
optPORT <- stats::nlminb(start = c(4.1, 3.9, -0.9, -8.7),
objective = RunAirGR4J,
lower = rep(-9.9, times = 4), upper = rep(9.9, times = 4),
......@@ -95,7 +95,7 @@ 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.
```{r, results='hide'}
```{r, results='hide', eval=FALSE}
ListIniParam <- data.frame(Chain1 = IniParam, Chain2 = IniParam, Chain3 = IniParam,
row.names = paste0("X", 1:4))
ListIniParam <- sweep(ListIniParam, MARGIN = 2, STATS = c(1, 0.9, 1.1), FUN = "*")
......@@ -143,7 +143,6 @@ First, the evolution of the Markov chains can be seen with a traceplot:
```{r, fig.width=6, fig.height=9, warning=FALSE}
ParamDRAM <- ggmcmc::ggs(MultDRAM) ## to convert objet for using by all ggs_* graphical functions
ggmcmc::ggs_traceplot(ParamDRAM)
```
The posterior density for each parameter can then be visualised:
......
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