From a0b60f071bb9a1fd3a78659798e351f37bebaf6e Mon Sep 17 00:00:00 2001
From: Mauricio Zambrano-Bigiarini <hzambran@users.noreply.github.com>
Date: Tue, 27 Nov 2012 11:10:09 +0000
Subject: [PATCH] hydroPSO.R: small bug fixed in IPSO. updated vignette
---
R/PSO_v2012.R | 15 ++++++----
inst/vignette/hydroPSO_vignette.Rnw | 45 +++++++++++++++--------------
2 files changed, 33 insertions(+), 27 deletions(-)
diff --git a/R/PSO_v2012.R b/R/PSO_v2012.R
index dbce3d6..1849492 100755
--- a/R/PSO_v2012.R
+++ b/R/PSO_v2012.R
@@ -259,13 +259,16 @@ compute.veloc <- function(x, v, w, c1, c2, CF, Pbest, part.index, gbest,
if(MinMax == "min") {
c2i <- c2 * ( (1/ngbest.fit)/sum(1/ngbest.fit) )
} else c2i <- c2 * ( ngbest.fit/sum(ngbest.fit) )
+
+ nan.index <- which(is.nan(c2i))
+ if (length(nan.index) > 0) c2i[nan.index] <- c2
# transforming 'x' into a matrix, with the same values in each row, in
# order to be able to substract 'x' from 'ngest'
X <- matrix(rep(x, nngbest), ncol=n, byrow=TRUE)
# computing the velocity
- vn <- CF * ( w*v + r1*c1*(pbest-x) + colSums(R2*c2i*(ngbest-X) ) )
+ vn <- CF * ( w*v + r1*c1*(pbest-x) + colSums(R2*c2i*(ngbest-X) ) )
} else if ( method=="fips" ) {
@@ -1085,7 +1088,7 @@ ComputeSwarmRadiusAndDiameter <- function(x, gbest, Lmax) {
radius <- numeric(npart)
gbest <- matrix(rep(gbest, npart), nrow=npart, byrow=TRUE)
- radius <- sqrt( rowSums( (x-gbest)^2 ) )
+ radius <- sqrt( rowSums( (x-gbest)^2, na.rm=TRUE) )
#swarm.radius <- max(radius, na.rm=TRUE)
swarm.radius <- median(radius, na.rm=TRUE)
@@ -1224,7 +1227,7 @@ RegroupingSwarm <- function(x,
# Relative change achieved in each dimension
#rel.change <- (RangeNew-RangeO)/RangeO
#names(rel.change) <- param.IDs
-
+
out <- list(3)
out[[1]] <- x
out[[2]] <- v
@@ -1326,7 +1329,7 @@ hydromod.eval <- function(part, Particles, iter, npart, maxit,
# 14-Jun-2012 ; 15-Jun-2012 ; 03-Jul-2012 ; 06-Jul-2012 #
# 11-Jul-2012 ; 17-Jul-2012 ; 18-Jul-2012 ; 13-Sep-2012 ; 14-Sep-2012 #
# 17-Sep-2012 ; 23-Sep-2012 ; 15-Oct-2012 ; 25-Oct-2012 ; 28-Oct-2012 #
-# 08-Nov-2012 ; 26-Nov-2012 #
+# 08-Nov-2012 ; 26-Nov-2012 ; 27-Nov-2012 #
################################################################################
# 'lower' : minimum possible value for each parameter
# 'upper' : maximum possible value for each parameter
@@ -2576,7 +2579,7 @@ hydroPSO <- function(
if ( (topology=="lbest") & (iter <= iter.ini) ) {
ltopology <- "gbest"
} else ltopology <- topology
-
+
V[j,] <- compute.veloc(
x= X[j, ],
v= V[j, ],
@@ -2597,7 +2600,7 @@ hydroPSO <- function(
ngbest=X.best.part[ngbest.pos, ], # topology="ipso"
lpbest.fit= pbest.fit[X.neighbours[j, ]] # method="wfips"
)
-
+
V[j,] <- velocity.boundary.treatment(v= V[j,], vmax=Vmax)
########################################################################
diff --git a/inst/vignette/hydroPSO_vignette.Rnw b/inst/vignette/hydroPSO_vignette.Rnw
index df3f769..5cbda5b 100644
--- a/inst/vignette/hydroPSO_vignette.Rnw
+++ b/inst/vignette/hydroPSO_vignette.Rnw
@@ -66,11 +66,11 @@ options(width=80)
\addcontentsline{toc}{section}{Preface}
The next tutorial shows how to use the \emph{hydroPSO} R package to calibrate different model codes. Two of the main properties of \emph{hydroPSO} are the independence from the model to be calibrated and the simple interfacing between the model code and the calibration-engine, i.e. \emph{hydroPSO}. These two properties are illustrated through the calibration of real-world case studies, which should be the starting point for a user planning to implement \emph{hydroPSO} for her/his own model calibration exercise.
-We show here how to use \emph{hydroPSO} to calibrate a semi-distributed hydrological model and a steady-state groundwater flow model. The first is fully implemented in SWAT-2005 \citep{neitsch+al2005a} whereas the second is implemented in MODFLOW-2005 (MF2005) \citep{mf2005}. Both model codes are widely-used programs to solve surface water and groundwater flow problems. At the same time, for illustrative purposes we use the program ZONEBUDGET (ZB) \citep{zonebudget} to obtain groundwater flows at particular aquifer cross sections. The purpose of including ZB is to illustrate how to interface \emph{hydroPSO} with a sequential modelling application. In this case, however, only MF2005 is calibrated against observed data and including other observations can be easily implemented in a similar way as described here.
+First, we illustrate the main features of \emph{hydroPSO} when used to calibrate any objective function fully defined within the R environment. Second, we show how to use \emph{hydroPSO} to calibrate a semi-distributed hydrological model and a steady-state groundwater flow model, SWAT-2005 \citep{neitsch+al2005a} and MODFLOW-2005 (MF2005) \citep{mf2005}, respectively. Both model codes are widely-used programs to simulate surface water and groundwater flow problems. At the same time, for illustrative purposes we use the program ZONEBUDGET (ZB) \citep{zonebudget} to obtain groundwater flows at particular aquifer cross sections. The purpose of including ZB is to illustrate how to interface \emph{hydroPSO} with a modelling application that involves more than one executable code. In the latter case, however, only MF2005 is calibrated against observed data, and including other observations can be easily implemented in a similar way as described here.
Both programs SWAT-2005 and MF2005/ZB are public domain and accessible from \url{http://swatmodel.tamu.edu/software/swat-model/} and \url{http://water.usgs.gov/nrp/gwsoftware/modflow.html}, respectively. For simplicity, we must assume some basic knowledge from part of the reader about the file naming convention used in SWAT-2005, MF2005/ZB, as well as the setting up of the internal options for these programs.
-As an integral part of this tutorial, the reader can download from \url{http://www.rforge.net/hydroPSO/} two files \href{http://www.rforge.net/hydroPSO/files/MF2005.zip}{MF2005.zip} and \href{http://www.rforge.net/hydroPSO/files/SWAT2005.zip}{SWAT2005.zip} including all the required files to run SWAT-2005 and MF2005, sample R scripts to interface \emph{hydroPSO} with SWAT-2005 and MF2005, as well as several auxiliary files.
+As an integral part of this tutorial, the reader can download from \url{http://www.rforge.net/hydroPSO/} two files \href{http://www.rforge.net/hydroPSO/files/MF2005.zip}{MF2005.zip} and \href{http://www.rforge.net/hydroPSO/files/SWAT2005.zip}{SWAT2005.zip} including all the required files to run the SWAT-2005 and MF2005 applications, sample R scripts to interface \emph{hydroPSO} with SWAT-2005 and MF2005, as well as several auxiliary files.
For assistance, bugs report, comments, and suggestions please contact the authors of the \emph{hydroPSO} package at:
\href{mailto:mzb.devel@gmail.com}{mzb.devel@gmail.com} and/or \href{mailto:Rodrigo.RojasMujica@gmail.com}{Rodrigo.RojasMujica@gmail.com}.
@@ -78,20 +78,20 @@ For assistance, bugs report, comments, and suggestions please contact the author
\clearpage
\section{Introduction\label{sec:intro}}
-\emph{hydroPSO} is an R package implementing the Particle Swarm Optimisation (PSO) algorithm \citep{kennedyeberhart1995,eberhartkennedy1995}. PSO is a population-based stochastic optimisation technique inspired by social behaviour of bird flocking, which shares few similarities with other evolutionary optimisation techniques such as Genetic Algorithms (GA) \citep{poli+al2007}. In PSO, however, the solution-space is explored on the basis of individual and neighbourhood best-known ``particle positions'' with no presence of evolution operators (e.g. mutation or crossover) as in GA. PSO has recently received a surge of attention in literature given its flexibility, simplicity of implementation (programming), low memory and computational requirements, low number of adjustable parameters, and efficiency \citep[see, e.g.,][]{eberhartshi1998,shieberhart1999,eberhartshi2001,poli+al2007}.
+\emph{hydroPSO} is an R package implementing the Particle Swarm Optimisation (PSO) algorithm \citep{kennedyeberhart1995,eberhartkennedy1995}. PSO is a population-based global optimisation technique inspired by social behaviour of bird flocking, which shares few similarities with other evolutionary optimisation techniques such as Genetic Algorithms (GA) \citep{poli+al2007}. In PSO, however, the solution-space is explored on the basis of individual and neighbourhood best-known ``particle positions'' with no presence of evolution operators (e.g. mutation or crossover) as in GA. PSO has recently received a surge of attention in literature given its flexibility, simplicity of implementation (programming), low memory and computational requirements, low number of adjustable parameters, and efficiency \citep[see, e.g.,][]{eberhartshi1998,shieberhart1999,eberhartshi2001,poli+al2007}.
The main motivations for developing the \emph{hydroPSO} package are:
\begin{enumerate}
\item bring this powerful optimisation/calibration technique to the attention of practitioners and scientists working on environmental modelling,
\item provide a model-independent software package allowing the user to calibrate (different) environmental models without having to invest considerable programming efforts in customizing the calibration engine to the environmental model, and
-\item allow the PSO research community to explore alternative advancements and configurations of the standard PSO using a versatile single-package software.
+\item allow the PSO research community to explore alternative configurations and recent developments of the standard PSO using a versatile single-package software.
\end{enumerate}
-Unlike other R packages recently developed for similar purposes (see e.g. \emph{hydromad}, \citealt{andrews+al2011}, \emph{R-SWAT-FME}, \citealt{wuliu2012}, and \emph{pso}, \citealt{psoptim2012}), \emph{hydroPSO} is not restricted to a limited number of hard-coded models, can be interfaced with any model with a relatively low programming effort, is fully compatible with calibration tools employing PEST-like template files \citep[see][]{pest2010}, and is easily parallelizable.
+Unlike other R packages recently developed for similar purposes (see e.g. \emph{hydromad}, \citealt{andrews+al2011}, \emph{R-SWAT-FME} and \citealt{wuliu2012}, \emph{hydroPSO} is not restricted to a limited number of hard-coded models, can be interfaced with any model with a relatively low programming effort, is fully compatible with calibration tools employing PEST-like template files \citep[see][]{pest2010}, and is easily parallelizable.
-\emph{hydroPSO} is also capable of performing sensitivity analysis using the Latin Hypercube One-factor-At-a-Time (LH-OAT) method \citep[see][]{vangriensven+al2006} and provides detailed information about the evolution of PSO's performance. In addition, advanced plotting capabilities and a complete family of built-in functions contained in the R language \citep{R2011} are available to visualize and summarize the calibration results, respectively. At the same time, \emph{hydroPSO} features a suite of controlling options and PSO variants for fine-tuning and improving the performance of PSO, thus, allowing the user to adapt the calibration/optimisation engine to different modelling problems. For a full description of the PSO enhancements included in the \emph{hydroPSO} package the reader is referred to Section~\ref{sec:pso}.
+\emph{hydroPSO} is also capable of performing sensitivity analysis using the Latin Hypercube One-factor-At-a-Time (LH-OAT) method \citep[see][]{vangriensven+al2006} and provides detailed information about the evolution of PSO's performance. In addition, advanced plotting capabilities and a complete family of built-in functions contained in the R language \citep{R2011} are available to visualize, analize and summarize the calibration results. At the same time, \emph{hydroPSO} features a suite of controlling options and PSO variants for fine-tuning and improving the performance of PSO, thus, allowing the user to adapt the calibration/optimisation engine to different modelling problems. For a full description of the PSO enhancements included in the \emph{hydroPSO} package the reader is referred to Section~\ref{sec:pso}.
-In principle, \emph{hydroPSO} only needs to know \emph{which} model parameters need to be calibrated, \emph{where} they need to be written, and \emph{from where} and \emph{how} to read the main model output. Then, it will take control of the model(s) to be calibrated until either a maximum number of iterations or an error tolerance in the objective function are reached: these two being problem-specific and user-defined. The basic interaction between \emph{hydroPSO} and the model to be calibrated is shown in Figure~\ref{fig:f01vign}. For models with numerous input and output files or cascading models (i.e. serial modelling applications), several I/O functions can be combined to interface \emph{hydroPSO} and the model.
+In principle, \emph{hydroPSO} only needs to know \emph{which} model parameters need to be calibrated, \emph{where} they need to be written, and \emph{from where} and \emph{how} to read the main model output. Then, it will take control of the model(s) to be calibrated until either a maximum number of iterations or an error tolerance in the objective function are reached: these two being problem-specific and user-defined. The basic interaction between \emph{hydroPSO} and the model to be calibrated is shown in Figure~\ref{fig:f01vign}. For models with numerous input and output files or cascading models (i.e. applications involving more than one model code), several I/O functions can be combined to interface \emph{hydroPSO} and the model.
\begin{figure}
\centering
@@ -100,7 +100,7 @@ In principle, \emph{hydroPSO} only needs to know \emph{which} model parameters n
\label{fig:f01vign}
\end{figure}
-In this tutorial we show first how to interface \emph{hydroPSO} with SWAT-2005 \citep{neitsch+al2005a} to calibrate a semi-distributed hydrological model for the Ega River basin in Spain. That section illustrates three main points: 1) basic interfacing between \emph{hydroPSO} and the model through the definition of ASCII files \emph{ParamRanges.txt} and \emph{ParamFiles.txt}, 2) advanced sensitivity analysis using LH-OAT and the use of \emph{hydroPSO} pre-defined goodness-of-fit measures as objective functions, and 3) the set up for calibrating a model output variable in time and static in space (transient application). Second, we show how to interface \emph{hydroPSO} with MF2005 \citep{mf2005} to calibrate a groundwater flow model for the regional aquifer of the Pampa del Tamarugal basin in Chile. That section, in turn, illustrates: 1) advanced interfacing between \emph{hydroPSO} and the model through the implementation of user-defined I/O wrapper functions written in R, 2) how to define/use a customized goodness-of-fit measure as objective function, and 3) how to set up \emph{hydroPSO} for a model output fixed in time but variable in space (steady-state application). Details about both applications are summarized in Table~\ref{tab:table1}.
+In this tutorial we show first how to interface \emph{hydroPSO} with SWAT-2005 \citep{neitsch+al2005a} to calibrate a semi-distributed hydrological model for the Ega River basin in Spain. That section illustrates three main points: 1) basic interfacing between \emph{hydroPSO} and the model through the definition of ASCII files \emph{ParamRanges.txt} and \emph{ParamFiles.txt}, used to define which model parameters need to be calibrated, and where they need to be written; 2) advanced sensitivity analysis using LH-OAT and the use of \emph{hydroPSO} with a pre-defined goodness-of-fit measure as objective function, and 3) the set up for calibrating a model output variable in time and static in space (transient application). Second, we show how to interface \emph{hydroPSO} with MF2005 \citep{mf2005} to calibrate a groundwater flow model for the regional aquifer of the Pampa del Tamarugal basin in Chile. That section, in turn, illustrates: 1) advanced interfacing between \emph{hydroPSO} and the model through the implementation of user-defined I/O wrapper functions written in R, 2) how to define/use a customized goodness-of-fit measure as objective function, and 3) how to set up \emph{hydroPSO} for a model output fixed in time but variable in space (steady-state application). Details about both applications are summarized in Table~\ref{tab:table1}.
%%%%%%% TABLE1
\begin{table}
@@ -124,17 +124,17 @@ Goodness-of-fit measure & Pre-defined Nash-Sutcliffe efficiency & U
\end{table}
%%% END TABLE1
-We therefore believe that this tutorial cover two standard-problems commonly found by the modelling community. We must note, however, that these tutorials deal with surface water and groundwater modelling applications as these two are the authors' areas of expertise. However, based on the flexibility of the \emph{hydroPSO} package and the benefits added by programming it in R, we believe this package can be implemented to a wider range of environmental models requiring some form of parameter estimation. To the date of writing this tutorial, we have provided support to users calibrating models such as: SWAT-2005, MF2005, LISFLOOD, AGNPS, HBV-96, and SWIM.
+We therefore believe that this tutorial cover two standard-problems commonly found by the modelling community. We must note, however, that these tutorials deal with surface water and groundwater modelling applications as these two are the authors' areas of expertise. However, based on the flexibility of the \emph{hydroPSO} package and the benefits added by programming it in R, we believe this package can be implemented to a wider range of environmental models requiring some form of parameter estimation. To the date of writing this tutorial, we have provided support to users calibrating models such as: SWAT-2005, MF2005, LISFLOOD, AGNPS and HBV-96. %(and SWIM)
-In Section~\ref{sec:pso} we present a brief description of PSO, including the main algorithm, discussion on topologies and PSO variants and fine-tuning options available in \emph{hydroPSO}. Section~\ref{sec:package} briefly describes the main package functions as well as the basic requirements to interface \emph{hydroPSO} with a given model code. An introductory application of \emph{hydroPSO} to optimise test functions commonly used to assess the performance of optimisation algorithms is described in Section~\ref{sec:application}. The procedure to interface \emph{hydroPSO} with SWAT-2005 and MF2005 as well as the corresponding calibration results, are reported in Sections~\ref{sec:swat2005} and~\ref{sec:mf2005}, respectively.
+The outline of this document is as follow. In Section~\ref{sec:pso} we present a brief description of PSO, including the main algorithm, discussion on topologies, PSO variants and fine-tuning options available in \emph{hydroPSO}. Section~\ref{sec:package} briefly describes the main package functions as well as the basic requirements to interface \emph{hydroPSO} with a given model code. An introductory application of \emph{hydroPSO} to optimise test functions commonly used to assess the performance of optimisation algorithms is described in Section~\ref{sec:application}. The procedure to interface \emph{hydroPSO} with SWAT-2005 and MF2005 as well as the corresponding calibration results, are reported in Sections~\ref{sec:swat2005} and~\ref{sec:mf2005}, respectively.
\clearpage
\section{Particle Swarm Optimisation (PSO)\label{sec:pso}}
\subsection{Canonical PSO Algorithm\label{sec:psoalg}}
-Particle Swarm Optimisation (PSO) is a population-based search algorithm developed by \citet{kennedyeberhart1995} and \citet{eberhartkennedy1995}. A feature distinguishing PSO from evolutionary algorithms is the lack of genetic operators, instead each individual of the population, termed \textit{particles} in PSO terminology, adjusts its \textit{flying trajectory} around the multi-dimensional search-space according to its own flying experience (\textit{local search}) and the one of all particles in the swarm (\textit{global search}) \citep{eberhartshi1998}.
+Particle Swarm Optimisation (PSO) is a population-based global search algorithm developed by \citet{kennedyeberhart1995} and \citet{eberhartkennedy1995}. A feature distinguishing PSO from evolutionary algorithms is the lack of genetic operators, instead each individual of the population, termed \textit{particles} in PSO terminology, adjusts its \textit{flying trajectory} around the multi-dimensional search-space according to its own flying experience and the one of all neighbouring particles in the swarm \citep{eberhartshi1998}.
-Considering a $D$-dimensional search-space, \textit{position} and \textit{velocity} for the $i$-th particle are represented by $\vec{X}_i=x_{i1},x_{i2},\ldots,x_{iD}$ and $\vec{V}_i=v_{i1},v_{i2},\ldots,v_{iD}$, respectively. The performance of each particle is assessed through a problem-specific ``goodness-of-fit'' measure, which is the basis for calculating $\vec{X}_i$, thus, reflecting the quality of the particle's position. The best-known position of the $i$-th particle (known as \textit{personal best}) is represented by $\vec{P}_i=p_{i1},p_{i2},\ldots, p_{iD}$, whereas the best-known position for the whole swarm (known as \textit{local best}) is represented by $\vec{G}=g_{1},g_{2},\ldots, g_{D}$. Velocity and position of the $i$-th particle are updated according to the following equations,
+Considering a $D$-dimensional search-space, \textit{position} and \textit{velocity} for the $i$-th particle are represented by $\vec{X}_i=x_{i1},x_{i2},\ldots,x_{iD}$ and $\vec{V}_i=v_{i1},v_{i2},\ldots,v_{iD}$, respectively. The performance of each particle is assessed through a problem-specific ``goodness-of-fit'' measure, which is the basis for calculating $\vec{X}_i$, thus, reflecting the quality of the particle's position. The best-known position of the $i$-th particle (known as \textit{personal best}) is represented by $\vec{P}_i=p_{i1},p_{i2},\ldots, p_{iD}$, whereas the best-known position for the neighbouring particles (known as \textit{local best}) is represented by $\vec{G}=g_{1},g_{2},\ldots, g_{D}$. Velocity and position of the $i$-th particle are updated according to the following equations,
\begin{equation}\label{eq:update1}
\vec{V}_{i}^{t+1} = \omega \vec{V}_{i}^t + c_1 \vec{U}_1^t \otimes (\vec{P}_{i}^t-\vec{X}_{i}^t) + c_2 \vec{U}_2^t \otimes (\vec{G}^t-\vec{X}_{i}^t)
@@ -168,7 +168,7 @@ v_{id} = \left\{
\end{array} \right.
\end{equation}
-Equations~\ref{eq:update1} to~\ref{eq:PSOvmax} can be considered as the canonical PSO algorithm \citep[see][]{kennedy2006}, which is summarised as follows
+Equations~\ref{eq:update1} to~\ref{eq:PSOvmax} can be considered as the \emph{canonical} PSO algorithm \citep[see][]{kennedy2006}, considering the whole swarm as a single neighbourhood, which is summarised as follows
\begin{algorithm}[H]
\caption{Canonical PSO algorithm.}
@@ -235,14 +235,14 @@ In this section several options used to tailor the calibration engine to specifi
\begin{enumerate}
\item PSO Parameters \\
-The PSO algorithm includes a few parameters that have to be carefully selected in order to avoid swarm explosion and to ensure proper convergence of the algorithm. These parameters include: a) inertia weight $(\omega)$, b) Clerc's constriction factor $(\chi)$, c) the cognitive coefficient $(c_1)$, and d) the social coefficient $(c_2)$.
+The PSO algorithm includes a few parameters that have to be carefully selected in order to avoid swarm explosion and to ensure proper convergence of the algorithm. These parameters include: a) inertia weight $(\omega)$ or the Clerc's constriction factor $(\chi)$, b) the cognitive acceleration coefficient $(c_1)$, and c) the social coefficient $(c_2)$.
Several definitions for the \emph{inertia weight} ($\omega$) have been extensively investigated and, hence, will not be further developed here. \citet{shieberhart1998b} initially propose a linearly decreasing variation for $\omega$. This strategy forces an extensive exploration of the search-space at initial iterations, which gradually shifts to a more dissipative (``exploitative'') search at later iterations \citep{poli+al2007}. Alternative $\omega$ definitions implemented in \emph{hydroPSO} include: a) linear \citep{shieberhart1998b}, b) non-linear \citep{chatterjeesiarry2006}, c) adaptive factor \citep{liu+al2005}, d) global-local best ratio \citep{arumugamrao2008}, and e) fully random \citep{eberhartshi2001}. These alternative definitions for $\omega$ allow the modeller to fine-tune the ``exploitation/exploration'' properties of the algorithm for different calibration problems.
Depending on the nature of the optimisation/calibration problem, the definition of the acceleration coefficients $c_1$ and $c_2$ may have a drastic impact on the algorithm's performance. \citet{ratnaweera+al2004} suggest that time-varying acceleration coefficients (TVAC) may reduce premature convergence at the first stages of the search and improve convergence at later iterations. For example, a large value of $c_1$ and a small value of $c_2$ at initial iterations allow the particles to ``fly'' around the search-space instead of moving towards their personal best. For later iterations, in turn, a small value of $c_1$ and a large value of $c_2$ allow the particles to converge to the global attraction zone. Based on this, linear, non-linear, and global-local best ratio (only for $c_1$) TVAC are implemented in \emph{hydroPSO} to fine-tune the coefficients $c_1$ and $c_2$.
\item Boundary Conditions \\
-In many optimisation/calibration problems it is a good practice to confine particle position to physically meaningful parameter ranges. In \emph{hydroPSO} this is done by using a $D$-dimensional vector defining the range of the search-space $(\Psi=\{(\psi_{1}^{min},\psi_{1}^{max}),(\psi_{2}^{min},\psi_{2}^{max}),\ldots,(\psi_{D}^{min},\psi_{D}^{max})\})$. The behaviour of the particles at these boundaries may influence the exploration of the search-space \citep{robinsonrahmatsamii2004}. Generally, four types of boundary condition may be imposed on a particle reaching the boundary of the search-space: a) \emph{absorbing}, where the particle's position is set to the boundary value and the velocity is set to zero; b) \emph{invisible}, where the particle is allowed to ``fly'' beyond the boundaries without limiting either its position or velocity, but its goodness-of-fit is not evaluated, c) \emph{reflecting}, where the sign of the velocity in the overshot dimension is changed and the particle is reflected back to the search-space; and d) \emph{damping}, where part (1-$f$) of the velocity is absorbed by the boundary and part ($f$) is reflected back, with $f$ as a random damping factor in [0,1] \citep[see][]{robinsonrahmatsamii2004,huangmohan2005}.
+In many optimisation/calibration problems it is a good practice to confine particle position to physically meaningful parameter ranges. In \emph{hydroPSO} this is done by using a $D$-dimensional vector defining the range of the search-space $(\Psi=\{(\psi_{1}^{min},\psi_{1}^{max}),(\psi_{2}^{min},\psi_{2}^{max}),\ldots,(\psi_{D}^{min},\psi_{D}^{max})\})$. The behaviour of the particles at these boundaries may influence the exploration of the search-space \citep{robinsonrahmatsamii2004}. Generally, four types of boundary condition may be imposed on a particle reaching the boundary of the search-space: a) \emph{absorbing}, where the particle's position is set to the boundary value and the velocity is set to zero; b) \emph{invisible}, where the particle is allowed to ``fly'' beyond the boundaries without limiting either its position or velocity, but its goodness-of-fit is not evaluated, c) \emph{reflecting}, where the sign of the velocity in the overshot dimension is changed and the particle is reflected back to the search-space; and d) \emph{damping}, where part (1-$f$) of the velocity is absorbed by the boundary and part ($f$) is reflected back, with $f$ as a random damping factor in [0,1] \citep[see][]{robinsonrahmatsamii2004,huangmohan2005,clerc2012}.
\item Regrouping Strategy \\
In the case of premature convergence (i.e. when the proposed solution approximates a local rather than a global optimum) or to validate the global attraction zone at final iterations, a suitable mechanism allowing particles to escape from sub-optimal solutions is required. \citet{eversghalia2009} propose a strategy that avoids stagnation of particles by automatically triggering swarm regrouping (see Table~\ref{tab:table2} for details). At each iteration the \emph{swarm radius} ($\delta^t$) is computed as the maximum Euclidean distance between the position of the swarm optimum ($\vec{G}$) and any given particle ($\vec{X}_i$). Then, regrouping is triggered when the \emph{normalized swarm radius} ($\delta^t_{norm}$) is less than a user-defined \emph{stagnation threshold} ($\epsilon$). The swarm is regrouped around the swarm optimum by using an updated range ($range_{d}(\Psi^{r})$ where $r$ is the \emph{swarm regrouping index}) calculated as the minimum between (i) the original range of the search-space on dimension $d$, and (ii) the product of the \emph{regrouping factor} ($\rho$) with the maximum distance along dimension $d$ between any particle and the swarm optimum. This regrouping method defines a reduced search-space ($\Psi^{r}$), which is small enough for an efficient search and large enough to allow particles to escape from sub-optimal solutions.
@@ -254,7 +254,7 @@ To enhance a faster propagation of the best solutions found throughout the swarm
For later stage iterations, oscillation of the particles around the local best may be restricted by using a time-varying $\vec{V}_{max}$. This reduces the search-space for latter stage iterations and, thus, complements the effect of ($\omega$). This is suggested by the authors of \emph{hydroPSO} following \citet{chatterjeesiarry2006}.
-Finally, \emph{hydroPSO} offers parallel capabilities to speed up the optimization process, which are easily configured.
+Finally, \emph{hydroPSO} (>=0.3-0) offers \textbf{parallel capabilities} to speed up the optimization process in multi-core machines or cluster networks, which are easily configured.
\end{enumerate}
\subsection{Standard PSO 2011 (SPSO 2011)\label{ch:}}
@@ -278,7 +278,7 @@ The most recent benchmark proposed, SPSO 2011, starts by defining a centre of gr
\noindent where $r$ is a random number in $[0,1]$. The particle's position is updated following equation~\eqref{eq:PositionUpdate}.
-Despite the fact that SPSO 2011 is an improvement to the previous SPSO 2007 in terms of accounting for rotational invariance, SPSO 2007 could perform better than SPSO 2011 for separable and (least likely) for non-separable functions (Maurice Clerc, developer SPSO 2007 \& 2011, personal communication, 2012). This is the main reason for including both SPSO 2007 \& SPSO 2011 in the \emph{hydroPSO} R package.
+Despite the fact that SPSO 2011 is an improvement to the previous SPSO 2007 in terms of accounting for rotational invariance, SPSO 2007 could perform better than SPSO 2011 for separable and (less likely) for non-separable functions (Maurice Clerc, developer SPSO 2007 \& 2011, personal communication, 2012). This is the main reason for including both SPSO 2007 \& SPSO 2011 in the \emph{hydroPSO} R package.
\singlespacing
%%%%%% TABLE2
@@ -319,9 +319,10 @@ Feature & \emph{hydroPSO} argument & Description & References \\
& \Verb+TVc2.type="non-linear"+ & $c2_{iter}=\left[ \frac{iter}{iter_{max}} (c2_{fin}-c2_{ini}) \right]+c2_{ini}$~\tnote[e] & \citet{ratnaweera+al2004} \\
& & & \\
3. Boundary conditions & \Verb+boundary.wall="absorbing"+ & 3.1. Absorbing & \citet{robinsonrahmatsamii2004} \\
- definition & \Verb+boundary.wall="invisible"+ & 3.2. Invisible & \citet{robinsonrahmatsamii2004} \\
- & \Verb+boundary.wall="reflecting"+ & 3.3. Reflecting & \citet{robinsonrahmatsamii2004} \\
- & \Verb+boundary.wall="damping"+ & 3.4. Damping & \citet{huangmohan2005} \\
+ & \Verb+boundary.wall="absorbing"+ & 3.1. Absorbing & \citet{robinsonrahmatsamii2004} \\
+ & \Verb+boundary.wall="invisible"+ & 3.2. Invisible & \citet{robinsonrahmatsamii2004} \\
+ & \Verb+boundary.wall="reflecting"+ & 3.3. Reflecting & \citet{robinsonrahmatsamii2004} \\
+ & \Verb+boundary.wall="damping"+ & 3.4. Damping & \citet{huangmohan2005} \\
& & & \\
4. Regrouping & \Verb+use.RG=TRUE+ & 4.1. $swarm~radius=\delta^t=max \Arrowvert \overrightarrow{X}_{i}^t-\overrightarrow{G}^t\Arrowvert$~\tnote[f] & \citet{eversghalia2009} \\
& \Verb+RG.thr=+$\epsilon$ & $normalized~swarm~radius=\delta^t_{norm}=\frac{\delta^t}{diam(\Psi)}<\epsilon$ & \\
@@ -906,7 +907,9 @@ The interaction between SWAT-2005 and \emph{hydroPSO} is shown in Figure~\ref{fi
\end{figure}
\subsection{Definition of \emph{ParamFiles.txt} and~\emph{ParamRanges.txt} files}
-The basic interaction between \emph{hydroPSO} and SWAT-2005 is defined through the specification of, first, the names and location of the parameters to be calibrated and, second, the definition of meaningful parameter ranges. This information is provided by two (problem-specific) ASCII files, namely, \emph{ParamFiles.txt} and \emph{ParamRanges.txt}, which must be contained in the \Verb+./PSO.in+ folder. In general, if a sensitivity analysis is implemented resulting in the identification of (in)sensitive parameters, the \emph{hydroPSO} user will have to define two sets of interfacing files, one for the sensitivity analysis (most likely including all parameters considered) and other for the calibration stage (most likely including only sensitive parameters). In this study case, \emph{ParamFiles-Sens.txt} and \emph{ParamRanges-Sens.txt} are defined as the interfacing files used for the sensitivity analysis, which correspond to the full set of hydrology-related parameters listed in Table~\ref{tab:swatparams}.
+\emph{hydroPSO} needs to know \emph{which} model parameters need to be calibrated and \emph{where} they need to be written. This information is provided by two (problem-specific) ASCII interfacing files, namely, \emph{ParamFiles.txt} and \emph{ParamRanges.txt}, which must be contained in the \Verb+./PSO.in+ folder. \emph{ParamFiles.txt} contains the names and location of the parameters to be calibrated, while \emph{ParamRanges.txt} stores physically meaningful parameter ranges to be used for the calibration of each parameter.
+
+In general, if a sensitivity analysis is going to be implemented for identifying (in)sensitive parameters, the \emph{hydroPSO} user will have to define two sets of interfacing files, one for the sensitivity analysis (most likely including all parameters considered) and other for the calibration stage (most likely including only sensitive parameters). In the SWAT-2005 study case, \emph{ParamFiles-Sens.txt} and \emph{ParamRanges-Sens.txt} are defined as the interfacing files used for the sensitivity analysis, which correspond to the full set of hydrology-related parameters listed in Table~\ref{tab:swatparams}: while \emph{ParamFiles.txt} and \emph{ParamRanges.txt} are the interfacing files used during calibration.
\begin{Verbatim}[frame=single,label=ParamFiles-Sens.txt,fontsize=\scriptsize]
ParameterNmbr ParameterName Filename Row.Number Col.Start Col.End DecimalPlaces
--
GitLab