From d22415e0e4e45b40608cbb43327fbdf6840b8356 Mon Sep 17 00:00:00 2001 From: Mauricio Zambrano-Bigiarini <hzambran@users.noreply.github.com> Date: Thu, 22 Nov 2012 17:06:51 +0000 Subject: [PATCH] minor changes --- DESCRIPTION | 7 +++-- NEWS | 3 +- R/hydromod.R | 13 +++++--- R/test_functions.R | 38 +++++++++++++++++----- man/hydroPSO.Rd | 73 +++++++++++++++++-------------------------- man/hydromod.Rd | 5 ++- man/lhoat.Rd | 5 +-- man/pest2hydroPSO.Rd | 13 ++++---- man/test_functions.Rd | 28 ++++++++--------- 9 files changed, 100 insertions(+), 85 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index bc33681..0b55da7 100755 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,9 +1,9 @@ Package: hydroPSO Type: Package Title: Model-Independent Particle Swarm Optimisation for Environmental Models -Version: 0.1-58-28 -Date: 2012-11-21 -Authors@R: c(person("Mauricio", "Zambrano-Bigiarini", email="mzb.devel@gmail.com", role=c("aut","cre"), comment="feel free to write me in English, Spanish or Italian"), person("Rodrigo", "Rojas", email="Rodrigo.RojasMujica@gmail.com", role=c("ctb")) ) +Version: 0.1-58-29 +Date: 2012-11-22 +Authors@R: c(person("Mauricio", "Zambrano-Bigiarini", email="mzb.devel@gmail.com", role=c("aut","cre")), person("Rodrigo", "Rojas", email="Rodrigo.RojasMujica@gmail.com", role=c("ctb")) ) Maintainer: Mauricio Zambrano-Bigiarini <mzb.devel@gmail.com> Description: This package implements a state-of-the-art version of the Particle Swarm Optimisation (PSO) algorithm (SPSO-2011 and SPSO-2007 capable), with a special focus on the calibration of environmental models. hydroPSO is model-independent, allowing the user to easily interface any model code with the calibration engine (PSO). It includes a series of controlling options and PSO variants to fine-tune the performance of the calibration engine to different calibration problems. An advanced sensitivity analysis function together with user-friendly plotting summaries facilitate the interpretation and assessment of the calibration results. Bugs reports/comments/questions are very welcomed. License: GPL (>=2) @@ -14,3 +14,4 @@ URL: http://www.rforge.net/hydroPSO, http://cran.r-project.org/web/packages/hydr Classification: optimisation, optimization, calibration, environment, environmental sciences, hydrology, PSO LazyLoad: yes ByteCompile: TRUE +X-CRAN-Comment: feel free to write the authors in English, Spanish or Italian diff --git a/NEWS b/NEWS index 8672155..06ca073 100755 --- a/NEWS +++ b/NEWS @@ -82,7 +82,8 @@ NEWS/ChangeLog for hydroPSO o 'hydromod' : -) now it is able to handle sub-daily models -) if the 'hydroGOF' package is not available, a simple correlation plot is produced between observations and the best simulation - -) improved documentation + -) improved documentation + -) removed default value for 'exe.name' (previously './swat.out') o 'test_functions' : -) new benchmark functions: 'schwefel' : Schwefel function diff --git a/R/hydromod.R b/R/hydromod.R index a68e877..1a722ce 100755 --- a/R/hydromod.R +++ b/R/hydromod.R @@ -22,13 +22,13 @@ # Updates: 20-Dec-2010 # # 19-Jan-2011 ; 22-Jan-2011 ; 02-Feb-2011 ; 11-May-2011 # # 13-Jan-2012 ; 16-Jan-2012 ; 23-Jan-2012 ; 02-May-2012 ; 12-Oct-2012 # -# 15-Oct-2012 # +# 15-Oct-2012 ; 22-Nov-2012 # ################################################################################ hydromod <- function( param.values, # Numeric vector with the paramter values that will be used in the input files of the hydrological model param.files="ParamFiles.txt", # Character, with the name of the file (with full path) that stores the name of the files that have to be modified for each parameter model.drty=getwd(), # Character, with the path of the directory that stores the exe file of the hydrological model and ALL the input files required for the simulation - exe.fname="./swat2005.out", + exe.fname, stdout=FALSE, # a logical (not NA) indicating whether messages written to ‘stdout’ should be sent or not. See '?system2' stderr="", # a logical (not NA) indicating whether messages written to ‘stderr’ should be sent or not. See '?system2' verbose= FALSE, # logical, indicating if progress messages have to be printed during the simulations. If \code{verbose=TRUE}, the following messages will appear: i)parameter values for each particle; (ii) model execution; iii) extraction of simulated values; and iv) computation of the goodness-of-fit measure. @@ -79,12 +79,15 @@ hydromod <- function( # 0) Checkings # ############################################################################## - if (!file.exists(param.files)) - stop( "Invalid argument: the file '", param.files, "' doesn't exist!" ) - # Verifying 'param.values' #if (missing(param.values)) # stop( "Missing argument: 'param.values' has to be given !") + + if (!file.exists(param.files)) + stop( "Invalid argument: the file '", param.files, "' doesn't exist!" ) + + if ( missing(exe.fname) ) + stop( "Missing argument: 'exe.fname'" ) if ( missing(out.FUN) ) { stop( "Missing argument: 'out.FUN'" ) diff --git a/R/test_functions.R b/R/test_functions.R index 9deed4a..1f3d83a 100755 --- a/R/test_functions.R +++ b/R/test_functions.R @@ -5,7 +5,7 @@ # All these function were started on 2008, with updates on: # # 13-Dec-2010 ; 20-Dec-2010; 21-Dec-2010 # # 24-Jan-2011 ; 02-Feb-2011 # -# 14-Nov-2011 ; 21-Sep-2012 ; 25-Sep-2012 +# 14-Nov-2011 ; 21-Sep-2012 ; 25-Sep-2012 ; 21-Nov-2012 ; 22-Nov-2012 # # MZB, 21-Jun-2011 # 3D sinc function: f(1,..,1)=1. Maximization @@ -34,7 +34,12 @@ rosenbrock <- function(x) { # MZB, RR, 21-Jun-2011; 21-Nov-2012 # Sphere function: f(1,..,1)=0. Minimization. In [-100, 100]^n. AcceptableError < 0.01 # Properties : Unimodal, additively separable -# Description: The Sphere test function is one of the most simple test functions available in the specialized literature. This unimodal and separable test function can be scaled up to any number of variables. It belongs to a family of functions called quadratic functions and only has one optimum in the point o = (0,...,0). The search range commonly used for the Sphere function is [−100, 100] for each decision variable. +# Description: The Sphere test function is one of the most simple test functions +# available in the specialized literature. This unimodal and separable +# test function can be scaled up to any number of variables. +# It belongs to a family of functions called quadratic functions and +# only has one optimum in the point o = (0,...,0). The search range +# commonly used for the Sphere function is [−100, 100] for each decision variable. # Reference : http://www.it.lut.fi/ip/evo/functions/node2.html sphere <- function(x) { return(sum(x^2)) @@ -51,7 +56,12 @@ rastrigrin <- function(x) { # MZB, RR, 17-Jul-2012. 21-Nov-2012. The correct name of the function is 'Rastrigin' and NOT 'Rastrigrin' !!! # Rastrigin function: f(0,..,0)=0. Minimization. In [-5.12, 5.12]^n. AcceptableError < 100 # Properties : Multimodal, additively separable -# Description: The generalized Rastrigin test function is non-convex and multimodal. It has several local optima arranged in a regular lattice, but it only has one global optimum located at the point \preformatted{o=(0,...,0)}. The search range for the Rastrigin function is [-5.12, 5.12] in each variable. This function is a fairly difficult problem due to its large search space and its large number of local minima +# Description: The generalized Rastrigin test function is non-convex and multimodal. +# It has several local optima arranged in a regular lattice, but it +# only has one global optimum located at the point \preformatted{o=(0,...,0)}. +# The search range for the Rastrigin function is [-5.12, 5.12] in each variable. +# This function is a fairly difficult problem due to its large search +# space and its large number of local minima # Reference : http://www.it.lut.fi/ip/evo/functions/node6.html, http://en.wikipedia.org/wiki/Rastrigin_function rastrigin <- function(x) { n <- length(x) @@ -62,17 +72,28 @@ rastrigin <- function(x) { # MZB, RR, 21-Jun-2011 ; 21-Nov-2012 # Griewank function: f(0,..,0)=0. Minimization. In [-600, 600]^n. AcceptableError < 0.05 # Properties : Multimodal, Non-separable -# Description: The Griewank test function is multimodal and non-separable, with has several local optima within the search region defined by [-600, 600]. It is similar to the Rastrigin function, but the number of local optima is larger in this case. It only has one global optimum located at the point \kbd{o=(0,...,0)}. While this function has an exponentially increasing number of local minima as its dimension increases, it turns out that a simple multistart algorithm is able to detect its global minimum more and more easily as the dimension increases (Locatelli, 2003) +# Description: The Griewank test function is multimodal and non-separable, with +# several local optima within the search region defined by [-600, 600]. +# It is similar to the Rastrigin function, but the number of local +# optima is larger in this case. It only has one global optimum +# located at the point \kbd{o=(0,...,0)}. While this function has +# an exponentially increasing number of local minima as its dimension +# increases, it turns out that a simple multistart algorithm is able +# to detect its global minimum more and more easily as the dimension +# increases (Locatelli, 2003) # Reference : http://www.geatbx.com/docu/fcnindex-01.html -# Locatelli, M. 2003. A note on the griewank test function, Journal of Global Optimization, 25 (2), 169-174, doi:10.1023/A:1021956306041 +# Locatelli, M. 2003. A note on the griewank test function, +# Journal of Global Optimization, 25 (2), 169-174, doi:10.1023/A:1021956306041 griewank <- function(x) { n <- length(x) return( 1 + (1/4000)*sum( x^2 ) - prod( cos( x/sqrt(seq(1:n)) ) ) ) } # 'griewank' END -# MZB, RR, 21-Jun-2011, 14-Nov-2011, 13-Sep-2012 +# MZB, RR, 21-Jun-2011, 14-Nov-2011, 13-Sep-2012 ; 22-Nov-2012 # Schaffer's f6 function: f(0,..,0)=0. Minimization. In [-100, 100]^n. AcceptableError < 1E-4 +# Reference: Xiaohong Qiu, Jun Liu. 2009. A Novel Adaptive PSO Algorithm on Schaffer's F6 Function. +# vol. 2, pp.94-98. Ninth International Conference on Hybrid Intelligent Systems schafferF6 <- function(x) { return( 0.5 + ( ( sin( sqrt( sum( x^2 ) ) ) )^2 - 0.5) / ( ( 1 + 0.001*sum(x^2) )^2 ) ) } # 'schafferF6' END @@ -81,7 +102,10 @@ schafferF6 <- function(x) { # MZB, RR, 14-Nov-2011, 21-Nov-2012 # Ackley function: f(0,..,0)=0. Minimization. In [-32.768, 32.768]^n. AcceptableError < 0.01, a=20 ; b=0.2 ; c=2*pi # Properties : Multimodal, Separable -# Description: The Ackley test function is multimodal and separable, with several local optima that, for the search range [-32, 32], look more like noise, although they are located at regular intervals. The Ackley function only has one global optimum located at the point o=(0,...,0). +# Description: The Ackley test function is multimodal and separable, with several +# local optima that, for the search range [-32, 32], look more like noise, +# although they are located at regular intervals. The Ackley function +# only has one global optimum located at the point o=(0,...,0). # Reference : http://www.it.lut.fi/ip/evo/functions/node14.html ackley <- function(x) { n <- length(x) diff --git a/man/hydroPSO.Rd b/man/hydroPSO.Rd index 3472cd0..fb1c9d5 100755 --- a/man/hydroPSO.Rd +++ b/man/hydroPSO.Rd @@ -26,8 +26,8 @@ OPTIONAL. numeric with a first guess for the parameters to be optimised (\env{n} All the particles are randomly initialised according to the value of \code{Xini.type}. If the user provides \env{m} parameter sets for \code{par}, they are used to overwrite the first \env{m} parameter sets randomly defined according to the value of \code{Xini.type}. If some elements in \code{par} are non finite (lower than \code{lower} or larger than \code{upper}) they are ignored } \item{fn}{ -character, name of a valid R function to be optimised (minimized or maximized) or the character value \sQuote{'hydromod'} \cr --) When \code{fn!='hydromod'}, the first argument of \code{fn} has to be a vector of parameters over which optimisation is to take place. It should return a scalar result. When \code{fn!='hydromod'} the algorithm uses the value(s) returned by \code{fn} as both model output and its corresponding goodness-of-fit measure \cr +character, name of a valid R function to be optimised (minimized or maximized) or the character value \sQuote{hydromod} \cr +-) When \code{fn!='hydromod'}, the first argument of \code{fn} has to be a vector of parameters over which optimisation is going to take place. It should return a scalar result. When \code{fn!='hydromod'} the algorithm uses the value(s) returned by \code{fn} as both model output and its corresponding goodness-of-fit measure \cr -) When \code{fn=='hydromod'} the algorithm will optimise the model defined by \code{model.FUN} and \code{model.args}, which are used to extract the values simulated by the model and to compute its corresponding goodness-of-fit measure } \item{\dots}{ @@ -35,7 +35,7 @@ OPTIONAL. Only used when \code{fn!='hydromod'}. \cr further arguments to be passed to \code{fn}. } \item{method}{ -character, variant of the PSO algorithm to be used. By default \code{method='spso2011'}, while valid values are \sQuote{spso2011}, \sQuote{spso2007}, \sQuote{ipso}, \sQuote{fips}, \sQuote{wfips}, \sQuote{canonical}: \cr +character, variant of the PSO algorithm to be used. By default \code{method='spso2011'}, while valid values are \sQuote{spso2011}, \sQuote{spso2007}, \sQuote{ipso}, \sQuote{fips}, \sQuote{wfips}, \sQuote{canonical}:\cr \kbd{spso2011}: At each iteration particles are attracted to its own best-known \sQuote{personal} and to the best-known position in its \sQuote{local} neighbourhood, which depens on the value of \code{topology}. In addition, values of the PSO engine are set to the values defined in the Standard PSO 2011 (SPSO 2011, see Clerc 2012) \cr @@ -148,9 +148,9 @@ By default \code{Xini.type=lhs} \item{Vini.type}{ character, indicates how to initialise the particles' velocities in the swarm. Valid values are: \cr -) \kbd{random2011}: random initialisation of velocities within \code{lower-Xini} and \code{upper-Xini}, as defined in SPSO 2011 (\samp{Vini=U(lower-Xini, upper-Xini)}) (see Clerc, 2012, 2010) \cr --) \kbd{lhs2011}: same as in \kbd{random2011}, but using a Latin Hypercube initialisation with \code{npart} number of strata instead of a random uniform distribution for each parameter (\samp{Vini=LHS(lower-Xini, upper-Xini)}). \bold{It requires the \pkg{lhs} package} \cr +-) \kbd{lhs2011}: same as in \kbd{random2011}, but using a Latin Hypercube initialisation with \code{npart} number of strata instead of a random uniform distribution for each parameter. \bold{It requires the \pkg{lhs} package} \cr -) \kbd{random2007}: random initialisation of velocities within \code{lower} and \code{upper} using the \sQuote{half-diff} method defined in SPSO 2007 (\samp{Vini=[U(lower, upper)-Xini]/2}) (see Clerc, 2012, 2010) \cr --) \kbd{lhs2007}: same as in \kbd{random2007}, but using a Latin Hypercube initialisation with \code{npart} number of strata instead of a random uniform distribution for each parameter (\samp{Vini=[LHS(lower, upper)-Xini]/2}) (see Clerc, 2010). \bold{It requires the \pkg{lhs} package} \cr +-) \kbd{lhs2007}: same as in \kbd{random2007}, but using a Latin Hypercube initialisation with \code{npart} number of strata instead of a random uniform distribution for each parameter. \bold{It requires the \pkg{lhs} package} \cr -) \kbd{zero}: all the particles are initialised with zero velocity \cr By default \code{Vini.type=NA}, which means that \code{Vini.type} depends on the value of \code{method}: when \kbd{method='spso2007'} \code{Vini.type='random2007'}, or \code{Vini.type='random2011'} otherwise } @@ -165,7 +165,7 @@ logical, if \code{TRUE} the particles are processed in random order to update th By default \code{random.update=TRUE} } \item{boundary.wall}{ -character, indicates the type of boundary condition to be applied during optimisation. Valid values are in \code{c(NA, 'absorbing2011', 'absorbing2007', 'reflecting', 'damping', 'invisible')} \cr +character, indicates the type of boundary condition to be applied during optimisation. Valid values are: \code{NA}, \sQuote{absorbing2011}, \sQuote{absorbing2007}, \sQuote{reflecting}, \sQuote{damping}, \sQuote{invisible} \cr By default \code{boundary.wall=NA}, which means that \code{boundary.wall} depends on the value of \code{method}: when \kbd{method='spso2007'} \code{boundary.wall='aborbing2007'}, or \code{boundary.wall='aborbing2011'} otherwise\cr Experience has shown that Clerc's constriction factor and the inertia weights do not always confine the particles within the solution space. To address this problem, Robinson and Rahmat-Samii (2004) and Huang and Mohan (2005) propose different boundary conditions, namely, \kbd{reflecting}, \kbd{damping}, \kbd{absorbing} and \kbd{invisible} to define how particles are treated when reaching the boundary of the searching space (see Robinson and Rahmat-Samii (2004) and Huang and Mohan (2005) for further details) @@ -182,8 +182,7 @@ By default \code{topology='random'} OPTIONAL. Only used when \code{topology} is in \code{c(lbest, vonNeumann, random)} \cr numeric, neighbourhood size, i.e. the number of informants for each particle (including the particle itself) to be considered in the computation of their personal best \cr When \code{topology=lbest} \code{K} MUST BE an even number in order to consider the same amount of neighbours to the left and the right of each particle \cr -As special case, \code{K} could be equal to \code{npart} \cr -By default \code{K=3} +As special case, \code{K} could be equal to \code{npart}. By default \code{K=3} } \item{iter.ini}{ OPTIONAL. Only used when \code{topology=='lbest'} \cr @@ -196,10 +195,9 @@ OPTIONAL. Only used when \code{method=='ipso'} \cr numeric, number of particles considered in the computation of the global best \cr By default \code{ngbest=4} (see Zhao, 2006) } - \item{randomise}{ + \item{normalise}{ logical, indicates whether the parameter values have to be normalised to the [0,1] interval during the optimisation or not\cr -This option appears in the C and Matlab version of SPSO-2011 (See \url{http://www.particleswarm.info/standard_pso_2011_c.zip}) and there it is recommended to use this option when the search space is not an hypercube. If the search space is an hypercube, it is better not normalise (there is a small difference between the position without any normalisation and the de-normalised one)\cr -By default \code{randomise=FALSE} +This option appears in the C and Matlab version of SPSO-2011 (See \url{http://www.particleswarm.info/standard_pso_2011_c.zip}) and there it is recommended to use this option when the search space is not an hypercube. If the search space is an hypercube, it is better not normalise (there is a small difference between the position without any normalisation and the de-normalised one). By default \code{normalise=FALSE} } \item{IW.type}{ OPTIONAL. Used only when \code{use.IW= TRUE} \& \code{length(IW.w)>1}\cr @@ -214,11 +212,10 @@ By default \code{IW.type='linear'} \item{IW.exp}{ OPTIONAL. Used only when \code{use.IW= TRUE} AND \code{IW.type= 'non-linear'} \cr numeric, non-linear modulation index (see Chatterjee and Siarry, 2006) \cr -When \code{IW.type='linear'}, \code{IW.exp} is set to 1 +When \code{IW.type='linear'}, \code{IW.exp} is set to 1. By default \code{IW.exp=1} } \item{use.TVc1}{ -logical, indicates if the cognitive acceleration coefficient \code{c1} will have a time-varying value instead of a constant one provided by the user (see Ratnaweera et al. 2004) \cr -By default \code{use.TVc1=FALSE} +logical, indicates if the cognitive acceleration coefficient \code{c1} will have a time-varying value instead of a constant one provided by the user (see Ratnaweera et al. 2004). By default \code{use.TVc1=FALSE} } \item{TVc1.type}{ character, required only when \code{use.TVc1 = TRUE}. Valid values are: \cr @@ -236,11 +233,10 @@ By default \code{TVc1.rng=c(1.28, 1.05)} OPTIONAL. Used only when \code{use.TVc1= TRUE} AND \code{TVc1.type= 'non-linear'} \cr numeric, non-linear modulation index \cr When \code{TVc1.exp} is equal to 1, \code{TVc1} corresponds to the improvement proposed by Ratnaweera et al., (2004), whereas when \code{TVc1.exp} is different from one, no reference has been found in literature by the authors, but it was included as an option based on the work of Chatterjee and Siarry (2006) for the inertia weight \cr -When \code{TVc1.type='linear'} \code{TVc1.exp} is automatically set to 1 +When \code{TVc1.type='linear'} \code{TVc1.exp} is automatically set to 1. By default \code{TVc1.exp=1} } \item{use.TVc2}{ -logical, indicates whether the social acceleration coefficient \code{c2} will have a time-varying value or a constant one provided by the user (see Ratnaweera et al. 2004) \cr -By default \code{use.TVc2=FALSE} +logical, indicates whether the social acceleration coefficient \code{c2} will have a time-varying value or a constant one provided by the user (see Ratnaweera et al. 2004). By default \code{use.TVc2=FALSE} } \item{TVc2.type}{ character, required only when \code{use.TVc2=TRUE}. Valid values are: \cr @@ -257,7 +253,7 @@ By default \code{TVc2.rng=c(1.05, 1.28)} OPTIONAL. Used only when \code{use.TVc2= TRUE} AND \code{TVc2.type='non-linear'} \cr numeric, non-linear modulation index \cr When \code{TVc2.exp} is equal to 1, \code{TVc2} corresponds to the improvement proposed by Ratnaweera et al., 2004, whereas when \code{TVc2.exp} is different from one, no reference has been found in literature by the authors, but it was included as an option based on the work of Chatterjee and Siarry (2006) for the inertia weight \cr -When \code{TVc2.type= linear} \code{TVc2.exp} is automatically set to 1 +When \code{TVc2.type= linear} \code{TVc2.exp} is automatically set to 1. By default \code{TVc2.exp=1} } \item{use.TVlambda}{ logical, indicates whether the percentage to limit the maximum velocity \code{lambda} will have a time-varying value or a constant value provided by the user. Proposed by the authors of hydroPSO based on the work of Chatterjee and Siarry (2006) for the inertia weight\cr @@ -277,7 +273,7 @@ By default \code{TVlambda.rng=c(1, 0.25)} \item{TVlambda.exp}{ OPTIONAL. only required when \code{use.TVlambda= TRUE} AND \code{TVlambda.type='non-linear'} \cr numeric, non-linear modulation index \cr -When \code{TVlambda.type='linear'} \code{TVlambda.exp} is automatically set to 1 +When \code{TVlambda.type='linear'} \code{TVlambda.exp} is automatically set to 1. By default \code{TVlambda.exp=1} } \item{use.RG}{ logical, indicates if the swarm should be regrouped when premature convergence is detected. By default \code{use.RG=FALSE} \cr @@ -291,13 +287,11 @@ There are 4 differences wrt Evers and Ghalia 2009: \cr \item{RG.thr}{ ONLY required when \code{use.RG=TRUE} \cr numeric, positive number representing the \kbd{stagnation threshold} used to decide whether the swarm has to be regrouped or not. See Evers and Galia (2009) for further details \cr -Regrouping occurs when the \kbd{normalised swarm radius} is less than \code{RG.thr}\cr -By default \code{RG.thr=1E-5} +Regrouping occurs when the \emph{normalised swarm radius} is less than \code{RG.thr}. By default \code{RG.thr=1E-5} } \item{RG.r}{ ONLY required when \code{use.RG=TRUE}. \cr -numeric, positive number representing the \kbd{regrouping factor}, which is used to regroup the swarm in a search space centred around the current global best (see Evers and Galia, 2009 for further details)\cr -By default \code{RG.thr=2} +numeric, positive number representing the \kbd{regrouping factor}, which is used to regroup the swarm in a search space centred around the current global best (see Evers and Galia, 2009 for further details). By default \code{RG.thr=2} } \item{RG.miniter}{ ONLY required when \code{use.RG=TRUE} \cr @@ -405,6 +399,8 @@ integer code where \code{0} indicates that the algorithm terminated by reaching \cite{Yong-Ling Zheng; Long-Hua Ma; Li-Yan Zhang; Ji-Xin Qian. On the convergence analysis and parameter selection in particle swarm optimization. Machine Learning and Cybernetics, 2003 International Conference on , vol.3, no., pp. 1802-1807 Vol.3, 2-5 Nov. 2003. doi: 10.1109/ICMLC.2003.1259789} +\cite{Zambrano-Bigiarini, M., and R. Rojas. 2012. hydroPSO: A model-independent particle swarm optimization software for calibration of environmental models, Environmental Modelling \& Software, (under-review)} + \cite{Zhao, B. An Improved Particle Swarm Optimization Algorithm for Global Numerical Optimization. In Proceedings of International Conference on Computational Science (1). 2006, 657-664} \cite{Neighborhood Topologies, \url{http://tracer.uc3m.es/tws/pso/neighborhood.html}. Last visited [15-Feb-2012] } @@ -425,7 +421,7 @@ Note for \code{\link[stats]{optim}} users: \cr } \examples{ # Number of dimensions of the optimisation problem (for all the examples) -nparam <- 3 +nparam <- 5 # Boundaries of the search space (Rastrigin function) lower <- rep(-5.12, nparam) @@ -444,11 +440,11 @@ set.seed(100) # Basic use 1. Rastrigin function (non-linear and multimodal with many local minima) # Without saving results to the hard disk -hydroPSO(fn=rastrigrin, lower=lower, upper=upper, control=list(write2disk=FALSE) ) +hydroPSO(fn=rastrigin, lower=lower, upper=upper, control=list(write2disk=FALSE) ) # Basic use 2. Rastrigin function (non-linear and multimodal with many local minima) # Saving results to the hard disk. Slower than before -hydroPSO(fn=rastrigrin, lower=lower, upper=upper ) +hydroPSO(fn=rastrigin, lower=lower, upper=upper ) # Plotting the results plot_results(MinMax="min") @@ -466,7 +462,7 @@ set.seed(100) # Defining the swarm size ('npart'), the relative tolerance ('reltol'), # the frequency of report messages printed to the screen ('REPORT'), and no # output files ('write2disk') -hydroPSO( fn=rastrigrin, lower=lower, upper=upper, +hydroPSO( fn=rastrigin, lower=lower, upper=upper, control=list(npart=15, reltol=1e-15, REPORT=10, write2disk=FALSE) ) @@ -479,9 +475,8 @@ hydroPSO( fn=rastrigrin, lower=lower, upper=upper, set.seed(100) # Same as Example 2, but setting the topology to global gest ('topology="gbest"') -hydroPSO( fn=rastrigrin, lower=lower, upper=upper, - control=list(npart=15, reltol=1e-15, REPORT=10, write2disk=FALSE, - topology="gbest") ) +hydroPSO( fn=rastrigin, lower=lower, upper=upper, + control=list(reltol=1e-15,REPORT=10,write2disk=FALSE,topology="gbest") ) @@ -493,19 +488,9 @@ hydroPSO( fn=rastrigrin, lower=lower, upper=upper, set.seed(100) # Same as Example 3 ('topology="gbest"') but using regrouping -hydroPSO( fn=rastrigrin, lower=lower, upper=upper, - control=list(npart=15, reltol=1e-15, REPORT=10, write2disk=FALSE, - topology="gbest", - use.RG=TRUE, RG.thr=1e-4, RG.r=3, RG.miniter=50) ) - - -# Setting the seed (for reproducible results) -set.seed(100) - -# Same as Example 2 ('topology="random"') but using regrouping -hydroPSO( fn=rastrigrin, lower=lower, upper=upper, - control=list(npart=15, reltol=1e-15, REPORT=10, write2disk=FALSE, - use.RG=TRUE, RG.thr=1e-4, RG.r=3, RG.miniter=50) ) +hydroPSO( fn=rastrigin, lower=lower, upper=upper, + control=list(reltol=1e-15,REPORT=10,write2disk=FALSE,topology="gbest", + use.RG=TRUE, RG.thr=1e-5, RG.r=3, RG.miniter=50) ) ################################ @@ -517,7 +502,7 @@ set.seed(100) # Same as Example 3, but using asynchronus update of previus and local best # ('best.update="async"') -hydroPSO( fn=rastrigrin, lower=lower, upper=upper, +hydroPSO( fn=rastrigin, lower=lower, upper=upper, control=list(npart=15, reltol=1e-15, abstol=1e-14, REPORT=10, topology="gbest", best.update="async") ) diff --git a/man/hydromod.Rd b/man/hydromod.Rd index 475798a..f08a3d9 100755 --- a/man/hydromod.Rd +++ b/man/hydromod.Rd @@ -15,8 +15,8 @@ It runs the user-defined model to be calibrated/optimised and returns a goodness } \usage{ hydromod(param.values, param.files = "ParamFiles.txt", model.drty = getwd(), - exe.fname = "./swat2005.out", stdout= FALSE, stderr="", - verbose = FALSE, out.FUN, out.FUN.args, gof.FUN, gof.FUN.args=list(), + exe.fname, stdout= FALSE, stderr="", verbose = FALSE, + out.FUN, out.FUN.args, gof.FUN, gof.FUN.args=list(), gof.Ini, gof.Fin, date.fmt = "\%Y-\%m-\%d", obs, do.png=FALSE, png.fname, width = 1200, height = 600, res=90, main, leg.cex=1.2, tick.tstep= "auto", lab.tstep= "auto", @@ -35,7 +35,6 @@ character, path to the executable file of the model specified in \code{exe.fname } \item{exe.fname}{ character, model command line arguments to be entered through a prompted string to execute the user-defined model - } \item{stdout, stderr}{ where output to \sQuote{stdout} or \sQuote{stderr} should be sent. Possible values are \code{FALSE} (discard output, the default), \code{""}, to the R console. See \code{\link[base]{system2}} diff --git a/man/lhoat.Rd b/man/lhoat.Rd index 4c5a579..f29f10e 100755 --- a/man/lhoat.Rd +++ b/man/lhoat.Rd @@ -20,8 +20,9 @@ lhoat(fn="hydromod", lower=-Inf, upper=Inf, control=list(), %- maybe also 'usage' for other objects documented here. \arguments{ \item{fn}{ -character, name of a valid R function to be optimised or character value \kbd{'hydromod'}. When \code{fn='hydromod'} the algorithm uses \code{model.FUN} and \code{model.FUN.args} to extract the values simulated by the model and to compute its corresponding goodness-of-fit function. When \code{fn!='hydromod'} the algorithm uses the value(s) returned by \code{fn} as both model output and its corresponding goodness-of-fit \cr -When \code{fn='hydromod'} the algorithm will optimise the model defined by \code{model.FUN} and \code{model.args} +character, name of a valid R function to be analysed or the character value \sQuote{hydromod}. \cr +-) When \code{fn!='hydromod'}, the first argument of \code{fn} has to be a vector of parameters over which the analysis is going to take place. It should return a scalar result. When \code{fn!='hydromod'} the algorithm uses the value(s) returned by \code{fn} as both model output and its corresponding goodness-of-fit measure \cr +-) When \code{fn=='hydromod'} the algorithm will analyse the model defined by \code{model.FUN} and \code{model.args}, which are used to extract the values simulated by the model and to compute its corresponding goodness-of-fit measure } %%\item{\dots}{ %%further arguments to be passed to \code{fn} diff --git a/man/pest2hydroPSO.Rd b/man/pest2hydroPSO.Rd index ce8506c..c62561b 100755 --- a/man/pest2hydroPSO.Rd +++ b/man/pest2hydroPSO.Rd @@ -8,10 +8,10 @@ \alias{pest2hydroPSO} %- Also NEED an '\alias' for EACH other topic documented here. \title{ -Export hydroPSO input files to PEST +Import PEST input files to hydroPSO } \description{ -This function exports the content of the hydroPSO input files \sQuote{ParamRanges.txt} and \sQuote{ParamFiles.txt} to PEST, into a single \sQuote{.pst} files with its corresponding \sQuote{.tpl} and \sQuote{.ins} files +This function imports the PEST input files (a master \sQuote{.pst} and its corresponding \sQuote{.tpl} and \sQuote{.ins}) into hydroPSO (\sQuote{ParamRanges.txt} and \sQuote{ParamFiles.txt}) } \usage{ pest2hydroPSO(pst.fname, drty.pest=NULL, drty.model=NULL, drty.out="PSO.in", @@ -40,7 +40,7 @@ character, name of the output file that will store the location and names of the By default this file is called \sQuote{ParamFiles.txt} and -if the full path it is not specified- it is searched for within the \sQuote{PSO.in} subdirectory of \code{drty.model} } \item{param.ranges}{ -character, name of the output file defining the minimum and maximum boundary values for each one of the parameters to be calibrated by hydroPSO. \cr +character, name of the output file defining the minimum and maximum boundary values for each one of the parameters to be calibrated with hydroPSO. \cr By default this file is called \sQuote{ParamRanges.txt} and -if the full path it is not specified- it is searched for within the \sQuote{PSO.in} subdirectory of \code{drty.model} } \item{decimals}{ @@ -59,9 +59,10 @@ logical, indicates if progress messages are to be printed. By default \code{verb %% ~~ If necessary, more details than the description above ~~ %%} \value{ -A list of two elements: -\item{sim}{numeric, with the simulated values obtained by running the model} -\item{GoF}{numeric, goodness-of-fit value representing how close each one of the simulated values in \code{sim} are to their observed counterparts, by using the USER-DEFINED \code{gof.FUN} function} +Two input files for \code{\link{hydroPSO}}: +\item{param.files}{plain text file with the location and names of the model files that have to be modified for each parameter subject to calibration with hydroPSO} +\item{param.ranges}{plain text file defining the minimum and maximum boundary values for each one of the parameters to be calibrated with hydroPSO} + } %%\references{ %% ~put references to the literature/web site here ~ diff --git a/man/test_functions.Rd b/man/test_functions.Rd index 68eb02d..9662bf6 100755 --- a/man/test_functions.Rd +++ b/man/test_functions.Rd @@ -56,7 +56,7 @@ numeric with the bias to be imposed } \details{ The \bold{Ackley} test function is multimodal and separable, with several local optima that, for the search range [-32, 32], look more like noise, although they are located at regular intervals. The Ackley function only has one global optimum located at the point \kbd{o=(0,...,0)}. It is defined by: -\deqn{ ackley = 20+\exp(1)-20\exp\left( -0.2\sqrt{\frac{1}{n}\sum_{i=1}^{n}x_{i}^2} \right)-\exp\left(\frac{1}{\textcolor{red}{n}}\sum_{i=1}^{n}\cos(2\pi x_{i})\right) ; -32 \leq x_i \leq 32 \ ; \ i=1,2,\ldots,n } +\deqn{ ackley = 20+\exp(1)-20\exp\left( -0.2\sqrt{\frac{1}{n}\sum_{i=1}^{n}x_{i}^2} \right)-\exp\left(\frac{1}{n}\sum_{i=1}^{n}\cos(2\pi x_{i})\right) ; -32 \leq x_i \leq 32 \ ; \ i=1,2,\ldots,n } The generalized \bold{Rastrigin} test function is non-convex, multimodal and additively separable. It has several local optima arranged in a regular lattice, but it only has one global optimum located at the point \kbd{o=(0,...,0)}. The search range for the Rastrigin function is [-5.12, 5.12] in each variable. This function is a fairly difficult problem due to its large search space and its large number of local minima. It is defined by: @@ -64,20 +64,20 @@ The generalized \bold{Rastrigin} test function is non-convex, multimodal and add \deqn{ rastrigin = 10n+\sum_{i=1}^{n}\left[x_{i}^{2}-10\cos(2\pi x_{i})\right] \ ; \ -5.12 \leq x_i \leq 5.12 \ ; \ i=1,2,\ldots,n } -The \bold{Griewank} test function is multimodal and non-separable, with has several local optima within the search region defined by [-600, 600]. It is similar to the Rastrigin function, but the number of local optima is larger in this case. It only has one global optimum located at the point \kbd{o=(0,...,0)}. The function interpretation changes with the scale; the general overview suggests convex function, medium-scale view suggests existence of local extremum, and finally zoom on the details indicates complex structure of numerous local minima. While this function has an exponentially increasing number of local minima as its dimension increases, it turns out that a simple multistart algorithm is able to detect its global minimum more and more easily as the dimension increases (Locatelli, 2003). It is defined by: +The \bold{Griewank} test function is multimodal and non-separable, with several local optima within the search region defined by [-600, 600]. It is similar to the Rastrigin function, but the number of local optima is larger in this case. It only has one global optimum located at the point \kbd{o=(0,...,0)}. The function interpretation changes with the scale; the general overview suggests convex function, medium-scale view suggests existence of local minima, and finally zoom on the details indicates complex structure of numerous local minima. While this function has an exponentially increasing number of local minima as its dimension increases, it turns out that a simple multistart algorithm is able to detect its global minimum more and more easily as the dimension increases (Locatelli, 2003). It is defined by: \deqn{ griewank = \frac{1}{4000}\sum_{i=1}^{n}x_{i}^{2}-\prod_{i=1}^{n}\cos\left(\frac{x_i}{\sqrt{i}}\right)+1 \ ; \ -600 \leq x_i \leq 600 \ ; \ i=1,2,\ldots,n } The \bold{Rosenbrock} function is non-convex, unimodal and non-separable. It is also known as \emph{Rosenbrock's valley} or \emph{Rosenbrock's banana} function. The global minimum is inside a long, narrow, parabolic shaped flat valley. To find the valley is trivial. To converge to the global minimum, however, is difficult. It only has one optimum located at the point \kbd{o=(1,...,1)}. It is a quadratic function, and its search range is [-30, 30] for each variable. It is defined by: -\deqn{ rosenbrock = \sum_{i=1}^{n-1}\left[100(x_{i+1}\textcolor{red}{-}x_{i}^{2})^{2}+(1-x_{i})^{2}\right] \ ; \ -30 \leq x_i \leq 30 \ ; \ i=1,2,\ldots,n } +\deqn{ rosenbrock = \sum_{i=1}^{n-1}\left[100(x_{i+1}-x_{i}^{2})^{2}+(1-x_{i})^{2}\right] \ ; \ -30 \leq x_i \leq 30 \ ; \ i=1,2,\ldots,n } -The main difficulty of the \bold{Schaffer's F6} test function is that the size of the potential maxima that need to be overcomed to get to a minimum increases the closer +The main difficulty of the \bold{Schaffer's F6} test function is that the size of the potential maxima that need to be overcome to get to a minimum increases the closer one gets to the global minimum. It is defined by: -\deqn{ schafferF6 = 0.5+\frac{\sin^{2}\sqrt{\textcolor{red}{\sum_{i=1}^{n}}x_{i}^{2}}-0.5}{(1+0.001\textcolor{red}{\sum_{i=1}^{n}}x_{i}^{2})^{2}} \ ; \ -100 \leq x_i \leq 100 \ ; \ i=1,2,\ldots,n } \cr +\deqn{ schafferF6 = 0.5+\frac{\sin^{2}\sqrt{\sum_{i=1}^{n}x_{i}^{2}}-0.5}{(1+0.001\sum_{i=1}^{n}x_{i}^{2})^{2}} \ ; \ -100 \leq x_i \leq 100 \ ; \ i=1,2,\ldots,n } \cr The \emph{first function of De Jong's} or \bold{Sphere} function is one of the most simple test functions available in the specialized literature. This continuous, convex, unimodal and additively separable test function can be scaled up to any number of variables. It belongs to a family of functions called quadratic functions and only has one optimum in the point \kbd{o=(0,...,0)}. The search range commonly used for the Sphere function is [-100, 100] for each decision variable. It is defined by: @@ -127,23 +127,23 @@ Test functions for optimization needs: \cite{\url{http://www.zsd.ict.pwr.wroc.pl \bold{Web pages}: -GEATbx: Example Functions (single and multi-objective functions). \cite{\url{http://www.geatbx.com/docu/fcnindex-01.html}}\cr +GEATbx: Example Functions (single and multi-objective functions). \cite{\url{http://www.geatbx.com/docu/fcnindex-01.html}} -Benchmark Problems \cite{\url{http://www.cs.cmu.edu/afs/cs/project/jair/pub/volume24/ortizboyer05a-html/node6.html}}\cr +Benchmark Problems \cite{\url{http://www.cs.cmu.edu/afs/cs/project/jair/pub/volume24/ortizboyer05a-html/node6.html}} -Test Functions for Unconstrained Global Optimization \cite{\url{http://www-optima.amp.i.kyoto-u.ac.jp/member/student/hedar/Hedar_files/TestGO_files/Page364.htm}}\cr +Test Functions for Unconstrained Global Optimization \cite{\url{http://www-optima.amp.i.kyoto-u.ac.jp/member/student/hedar/Hedar_files/TestGO_files/Page364.htm}} -Rosenbrock: \cite{\url{http://www.it.lut.fi/ip/evo/functions/node5.html}}, \cite{\url{http://en.wikipedia.org/wiki/Rosenbrock_function}}\cr +Rosenbrock: \cite{\url{http://www.it.lut.fi/ip/evo/functions/node5.html}}, \cite{\url{http://en.wikipedia.org/wiki/Rosenbrock_function}} -Sphere: \cite{\url{http://www.it.lut.fi/ip/evo/functions/node2.html}}\cr +Sphere: \cite{\url{http://www.it.lut.fi/ip/evo/functions/node2.html}} -Rastrigin: \cite{\url{http://www.it.lut.fi/ip/evo/functions/node6.html}}, \cite{\url{http://en.wikipedia.org/wiki/Rastrigin_function}}\cr +Rastrigin: \cite{\url{http://www.it.lut.fi/ip/evo/functions/node6.html}}, \cite{\url{http://en.wikipedia.org/wiki/Rastrigin_function}} -Ackley: \cite{\url{http://www.it.lut.fi/ip/evo/functions/node14.html}}\cr +Ackley: \cite{\url{http://www.it.lut.fi/ip/evo/functions/node14.html}} -Griewank: \cite{Locatelli, M. 2003. A note on the griewank test function, Journal of Global Optimization, 25 (2), 169-174, doi:10.1023/A:1021956306041}\cr +Griewank: \cite{Locatelli, M. 2003. A note on the griewank test function, Journal of Global Optimization, 25 (2), 169-174, doi:10.1023/A:1021956306041} -Schaffer's F6 \cite{Xiaohong Qiu, Jun Liu. 2009. A Novel Adaptive PSO Algorithm on Schaffer's F6 Function. Hybrid Intelligent Systems, International Conference on, pp. 94-98, 2009 Ninth International Conference on Hybrid Intelligent Systems} \cr +Schaffer's F6 \cite{Xiaohong Qiu, Jun Liu. 2009. A Novel Adaptive PSO Algorithm on Schaffer's F6 Function. Hybrid Intelligent Systems, International Conference on, pp. 94-98, 2009 Ninth International Conference on Hybrid Intelligent Systems} Schwefel: \cite{\url{http://www.geatbx.com/docu/fcnindex-01.html\#P150_6749}} } -- GitLab