Optimize (minimize) a noisy blackbox function (i.e., a function which
may not be differentiable, and may not execute deterministically).
A b*
tgp model is used as a surrogate for
adaptive sampling via improvement (and other) statistics. Note that
this function is intended as a skeleton to be tailored as required
for a particular application
1 2 3 4 5  optim.step.tgp(f, rect, model = btgp, prev = NULL, X = NULL,
Z = NULL, NN = 20 * length(rect), improv = c(1, 5), cands = c("lhs", "tdopt"),
method = c("LBFGSB", "NelderMead", "BFGS", "CG", "SANN", "optimize"), ...)
optim.ptgpf(start, rect, tgp.obj,
method=c("LBFGSB", "NelderMead", "BFGS", "CG", "SANN", "optimize"))

f 
A function to be optimized, having only one free argument 
rect 

model 
The 
prev 
The output from a previous call to 
X 

Z 
Vector of current output responses 
NN 
Number of candidate locations ( 
improv 
Indicates the 
cands 
The type of candidates ( 
method 
A method from 
... 
Further arguments to the 
start 
A starting value for optimization of the MAP predictive
(kriging) surface of a 
tgp.obj 
A 
optim.step.tgp
executes one step in a search for the global
optimum (minimum) of a noisy function (Z~f(X)
) in a bounded
rectangle (rect
). The procedure essentially fits a tgp
model
and samples from
the posterior distribution of improvement
statistics at NN+1
candidates locations.
NN
of the candidates come from
cands
, i.e., "lhs"
or "tdopt"
, plus one which
is the location of the minima found in a previous run via
prev
by using optim
(with a particular
method
or optimize
instead) on the MAP
model
predictive surface using the "tgp"
class object
contained therein.
The improv[2]
with the the highest expected improvement are
recommended for adding into the design on output.
optim.ptgpf
is the subroutine used by
optim.step.tgp
to find optimize on the MAP (surrogate)
predictive surface for the "tgp"
class object contained in
prev
.
Please see vignette("tgp2")
for a detailed illustration
The list
return has the following components.
X 
A 
progress 
A onerow 
obj 
the 
The ellipses (...
) argument is used differently here, as
compared to optim
, and optimize
. It
allows further arguments to be passed to the b*
model
function, whereas for optim
it would describe
further (static) arguments to the function f
to be optimized.
If static arguments need to be set for f
, then we recommend
setting defaults via the formals
of f
Robert B. Gramacy, rbgramacy@chicagobooth.edu, and Matt Taddy, taddy@chicagobooth.edu
Matthew Taddy, Herbert K.H. Lee, Genetha A. Gray, and Joshua D. Griffin. (2009) Bayesian guided pattern search for robust local optimization. Technometrics, to appear.
http://bobby.gramacy.com/r_packages/tgp
btgp
, etc., optim
,
optimize
, tgp.design
,
predict.tgp
, dopt.gp
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  ## optimize the simple exponential function
f < function(x) { exp2d.Z(x)$Z }
## create the initial design with Doptimal candidates
rect < rbind(c(2,6), c(2,6))
Xcand < lhs(500, rect)
X < dopt.gp(50, X=NULL, Xcand)$XX
Z < f(X)
## do 10 rounds of adaptive sampling
out < progress < NULL
for(i in 1:10) {
## get recommendations for the next point to sample
out < optim.step.tgp(f, X=X, Z=Z, rect=rect, prev=out)
## add in the inputs, and newly sampled outputs
X < rbind(X, out$X)
Z < c(Z, f(out$X))
## keep track of progress and best optimum
progress < rbind(progress, out$progress)
print(progress[i,])
}
## plot the progress so far
par(mfrow=c(2,2))
plot(out$obj, layout="surf")
plot(out$obj, layout="as", as="improv")
matplot(progress[,1:nrow(rect)], main="optim results",
xlab="rounds", ylab="x[,1:2]", type="l", lwd=2)
plot(log(progress$improv), type="l", main="max log improv",
xlab="rounds", ylab="max log(improv)")

Questions? Problems? Suggestions? Tweet to @rdrrHQ or email at ian@mutexlabs.com.
All documentation is copyright its authors; we didn't write any of that.