1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247
|
\name{gosolnp}
\alias{gosolnp}
\title{
Random Initialization and Multiple Restarts of the solnp solver.
}
\description{
When the objective function is non-smooth or has many local minima, it is hard
to judge the optimality of the solution, and this usually depends critically on
the starting parameters. This function enables the generation of a set of
randomly chosen parameters from which to initialize multiple restarts of the
solver (see note for details).
}
\usage{
gosolnp(pars = NULL, fixed = NULL, fun, eqfun = NULL, eqB = NULL, ineqfun = NULL,
ineqLB = NULL, ineqUB = NULL, LB = NULL, UB = NULL, control = list(),
distr = rep(1, length(LB)), distr.opt = list(), n.restarts = 1, n.sim = 20000,
cluster = NULL, rseed = NULL, ...)
}
\arguments{
\item{pars}{
The starting parameter vector. This is not required unless the fixed option is
also used.
}
\item{fixed}{
The numeric index which indicates those parameters which should stay fixed
instead of being randomly generated.
}
\item{fun}{
The main function which takes as first argument the parameter vector and returns
a single value.
}
\item{eqfun}{
(Optional) The equality constraint function returning the vector of evaluated
equality constraints.
}
\item{eqB}{
(Optional) The equality constraints.
}
\item{ineqfun}{
(Optional) The inequality constraint function returning the vector of evaluated
inequality constraints.
}
\item{ineqLB}{
(Optional) The lower bound of the inequality constraints.
}
\item{ineqUB}{
(Optional) The upper bound of the inequality constraints.
}
\item{LB}{
The lower bound on the parameters. This is not optional in this function.
}
\item{UB}{
The upper bound on the parameters. This is not optional in this function.
}
\item{control}{
(Optional) The control list of optimization parameters. The \code{eval.type}
option in this control list denotes whether to evaluate the function as is and
exclude inequality violations in the final ranking (default, value = 1),
else whether to evaluate a penalty barrier function comprised of the objective
and all constraints (value = 2). See \code{solnp} function documentation for
details of the remaining control options.
}
\item{distr}{
A numeric vector of length equal to the number of parameters, indicating the
choice of distribution to use for the random parameter generation. Choices are
uniform (1), truncated normal (2), and normal (3).
}
\item{distr.opt}{
If any choice in \code{distr} was anything other than uniform (1), this is a
list equal to the length of the parameters with sub-components for the mean and
sd, which are required in the truncated normal and normal distributions.
}
\item{n.restarts}{
The number of solver restarts required.
}
\item{n.sim}{
The number of random parameters to generate for every restart of the solver.
Note that there will always be significant rejections if inequality bounds are
present. Also, this choice should also be motivated by the width of the upper
and lower bounds.
}
\item{cluster}{
If you want to make use of parallel functionality, initialize and pass a cluster
object from the parallel package (see details), and remember to terminate it!
}
\item{rseed}{
(Optional) A seed to initiate the random number generator, else system time will
be used.
}
\item{\dots}{
(Optional) Additional parameters passed to the main, equality or inequality
functions
}
}
\details{
Given a set of lower and upper bounds, the function generates, for those
parameters not set as fixed, random values from one of the 3 chosen
distributions. Depending on the \code{eval.type} option of the \code{control}
argument, the function is either directly evaluated for those points not
violating any inequality constraints, or indirectly via a penalty barrier
function jointly comprising the objective and constraints. The resulting values
are then sorted, and the best N (N = random.restart) parameter vectors
(corresponding to the best N objective function values) chosen in order to
initialize the solver. Since version 1.14, it is up to the user to prepare and
pass a cluster object from the parallel package for use with gosolnp, after
which the parLapply function is used. If your function makes use of additional
packages, or functions, then make sure to export them via the \code{clusterExport}
function of the parallel package. Additional arguments passed to the solver via the
\dots option are evaluated and exported by gosolnp to the cluster.
}
\value{
A list containing the following values:
\item{pars}{Optimal Parameters.}
\item{convergence }{Indicates whether the solver has converged (0) or not (1).}
\item{values}{Vector of function values during optimization with last one the
value at the optimal.}
\item{lagrange}{The vector of Lagrange multipliers.}
\item{hessian}{The Hessian at the optimal solution.}
\item{ineqx0}{The estimated optimal inequality vector of slack variables used
for transforming the inequality into an equality constraint.}
\item{nfuneval}{The number of function evaluations.}
\item{elapsed}{Time taken to compute solution.}
\item{start.pars}{The parameter vector used to start the solver}
}
\references{
Y.Ye, \emph{Interior algorithms for linear, quadratic, and linearly constrained
non linear programming}, PhD Thesis, Department of EES Stanford University,
Stanford CA.\cr
Hu, X. and Shonkwiler, R. and Spruill, M.C. \emph{Random Restarts in Global
Optimization}, 1994, Georgia Institute of technology, Atlanta.
}
\author{
Alexios Ghalanos and Stefan Theussl\cr
Y.Ye (original matlab version of solnp)
}
\note{
The choice of which distribution to use for randomly sampling the parameter
space should be driven by the user's knowledge of the problem and confidence or
lack thereof of the parameter distribution. The uniform distribution indicates
a lack of confidence in the location or dispersion of the parameter, while the
truncated normal indicates a more confident choice in both the location and
dispersion. On the other hand, the normal indicates perhaps a lack of knowledge
in the upper or lower bounds, but some confidence in the location and dispersion
of the parameter. In using choices (2) and (3) for \code{distr},
the \code{distr.opt} list must be supplied with \code{mean} and \code{sd} as
subcomponents for those parameters not using the uniform (the examples section
hopefully clarifies the usage).
}
\examples{
\dontrun{
# [Example 1]
# Distributions of Electrons on a Sphere Problem:
# Given n electrons, find the equilibrium state distribution (of minimal Coulomb
# potential) of the electrons positioned on a conducting sphere. This model is
# from the COPS benchmarking suite. See http://www-unix.mcs.anl.gov/~more/cops/.
gofn = function(dat, n)
{
x = dat[1:n]
y = dat[(n+1):(2*n)]
z = dat[(2*n+1):(3*n)]
ii = matrix(1:n, ncol = n, nrow = n, byrow = TRUE)
jj = matrix(1:n, ncol = n, nrow = n)
ij = which(ii<jj, arr.ind = TRUE)
i = ij[,1]
j = ij[,2]
# Coulomb potential
potential = sum(1.0/sqrt((x[i]-x[j])^2 + (y[i]-y[j])^2 + (z[i]-z[j])^2))
potential
}
goeqfn = function(dat, n)
{
x = dat[1:n]
y = dat[(n+1):(2*n)]
z = dat[(2*n+1):(3*n)]
apply(cbind(x^2, y^2, z^2), 1, "sum")
}
n = 25
LB = rep(-1, 3*n)
UB = rep(1, 3*n)
eqB = rep(1, n)
ans = gosolnp(pars = NULL, fixed = NULL, fun = gofn, eqfun = goeqfn, eqB = eqB,
LB = LB, UB = UB, control = list(outer.iter = 100, trace = 1),
distr = rep(1, length(LB)), distr.opt = list(), n.restarts = 2, n.sim = 20000,
rseed = 443, n = 25)
# should get a function value around 243.813
# [Example 2]
# Parallel functionality for solving the Upper to Lower CVaR problem (not properly
# formulated...for illustration purposes only).
mu =c(1.607464e-04, 1.686867e-04, 3.057877e-04, 1.149289e-04, 7.956294e-05)
sigma = c(0.02307198,0.02307127,0.01953382,0.02414608,0.02736053)
R = matrix(c(1, 0.408, 0.356, 0.347, 0.378, 0.408, 1, 0.385, 0.565, 0.578, 0.356,
0.385, 1, 0.315, 0.332, 0.347, 0.565, 0.315, 1, 0.662, 0.378, 0.578,
0.332, 0.662, 1), 5,5, byrow=TRUE)
# Generate Random deviates from the multivariate Student distribution
set.seed(1101)
v = sqrt(rchisq(10000, 5)/5)
S = chol(R)
S = matrix(rnorm(10000 * 5), 10000) \%*\% S
ret = S/v
RT = as.matrix(t(apply(ret, 1, FUN = function(x) x*sigma+mu)))
# setup the functions
.VaR = function(x, alpha = 0.05)
{
VaR = quantile(x, probs = alpha, type = 1)
VaR
}
.CVaR = function(x, alpha = 0.05)
{
VaR = .VaR(x, alpha)
X = as.vector(x[, 1])
CVaR = VaR - 0.5 * mean(((VaR-X) + abs(VaR-X))) / alpha
CVaR
}
.fn1 = function(x,ret)
{
port=ret\%*\%x
obj=-.CVaR(-port)/.CVaR(port)
return(obj)
}
# abs(sum) of weights ==1
.eqn1 = function(x,ret)
{
sum(abs(x))
}
LB=rep(0,5)
UB=rep(1,5)
pars=rep(1/5,5)
ctrl = list(delta = 1e-10, tol = 1e-8, trace = 0)
cl = makePSOCKcluster(2)
# export the auxilliary functions which are used and cannot be seen by gosolnp
clusterExport(cl, c(".CVaR", ".VaR"))
ans = gosolnp(pars, fun = .fn1, eqfun = .eqn1, eqB = 1, LB = LB, UB = UB,
n.restarts = 2, n.sim=500, cluster = cl, ret = RT)
ans
# don't forget to stop the cluster!
stopCluster(cl)
}
}
\keyword{optimize}
|