noisy.optimizer {DiceOptim}R Documentation

Optimization of homogenously noisy functions based on Kriging

Description

Sequential optimization of kriging-based criterion conditional on noisy observations, with model update after each evaluation. Eight criteria are proposed to choose the next observation: random search, sequential parameter optimization (SPO), reinterpolation, Expected Improvement (EI) with plugin, Expected Quantile Improvement (EQI), quantile minimization, Augmented Expected Improvement (AEI) and Approximate Knowledge Gradient (AKG). The criterion optimization is based on the package rgenoud.

Usage

noisy.optimizer(
  optim.crit,
  optim.param = NULL,
  model,
  n.ite,
  noise.var = NULL,
  funnoise,
  lower,
  upper,
  parinit = NULL,
  control = NULL,
  CovReEstimate = TRUE,
  NoiseReEstimate = FALSE,
  nugget.LB = 1e-05,
  estim.model = NULL,
  type = "UK"
)

Arguments

optim.crit

String defining the criterion to be optimized at each iteration. Possible values are: "random.search", "SPO", "reinterpolation", "EI.plugin", "EQI", "min.quantile", "AEI", "AKG".

optim.param

List of parameters for the chosen criterion. For "EI.plugin": optim.param$plugin.type is a string defining which plugin is to be used. Possible values are "ytilde", "quantile" and "other". If "quantile" is chosen, optim.param$quantile defines the quantile level. If "other" is chosen, optim.param$plugin directly sets the plugin value.

For "EQI": optim.param$quantile defines the quantile level. If not provided, default value is 0.9.

For "min.quantile": optim.param$quantile defines the quantile level. If not provided, default value is 0.1.

For "AEI": optim.param$quantile defines the quantile level to choose the current best point. If not provided, default value is 0.75.

model

a Kriging model of "km" class

n.ite

Number of iterations

noise.var

Noise variance (scalar). If noiseReEstimate=TRUE, it is an initial guess for the unknown variance (used in optimization).

funnoise

objective (noisy) function

lower

vector containing the lower bounds of the variables to be optimized over

upper

vector containing the upper bounds of the variables to be optimized over

parinit

optional vector of initial values for the variables to be optimized over

control

optional list of control parameters for optimization. One can control "pop.size" (default : [N=3*2^dim for dim<6 and N=32*dim otherwise]]), "max.generations" (N), "wait.generations" (2) and "BFGSburnin" (0) of function "genoud" (see genoud). Numbers into brackets are the default values

CovReEstimate

optional boolean specfiying if the covariance parameters should be re-estimated at every iteration (default value = TRUE)

NoiseReEstimate

optional boolean specfiying if the noise variance should be re-estimated at every iteration (default value = FALSE)

nugget.LB

optional scalar of minimal value for the estimated noise variance. Default value is 1e-5.

estim.model

optional kriging model of "km" class with homogeneous nugget effect (no noise.var). Required if noise variance is reestimated and the initial "model" has heterogenenous noise variances.

type

"SK" or "UK" for Kriging with known or estimated trend

Value

A list with components:

model

the current (last) kriging model of "km" class

best.x

The best design found

best.y

The objective function value at best.x

best.index

The index of best.x in the design of experiments

history.x

The added observations

history.y

The added observation values

history.hyperparam

The history of the kriging parameters

estim.model

If noiseReEstimate=TRUE, the current (last) kriging model of "km" class for estimating the noise variance.

history.noise.var

If noiseReEstimate=TRUE, the history of the noise variance estimate.

Author(s)

Victor Picheny

References

V. Picheny and D. Ginsbourger (2013), Noisy kriging-based optimization methods: A unified implementation within the DiceOptim package, Computational Statistics & Data Analysis

Examples


##########################################################################
### EXAMPLE 1: 3 OPTIMIZATION STEPS USING EQI WITH KNOWN NOISE         ###
### AND KNOWN COVARIANCE PARAMETERS FOR THE BRANIN FUNCTION            ###
##########################################################################

set.seed(10)
library(DiceDesign)
# Set test problem parameters
doe.size <- 9
dim <- 2
test.function <- get("branin2")
lower <- rep(0,1,dim)
upper <- rep(1,1,dim)
noise.var <- 0.1

# Build noisy simulator
funnoise <- function(x)
{     f.new <- test.function(x) + sqrt(noise.var)*rnorm(n=1)
      return(f.new)}

# Generate DOE and response
doe <- as.data.frame(lhsDesign(doe.size, dim)$design)
y.tilde <- funnoise(doe)

# Create kriging model
model <- km(y~1, design=doe, response=data.frame(y=y.tilde),
     covtype="gauss", noise.var=rep(noise.var,1,doe.size), 
     lower=rep(.1,dim), upper=rep(1,dim), control=list(trace=FALSE))

# Optimisation with noisy.optimizer (n.ite can be increased)
n.ite <- 2 
optim.param <- list()
optim.param$quantile <- .9
optim.result <- noisy.optimizer(optim.crit="EQI", optim.param=optim.param, model=model,
		n.ite=n.ite, noise.var=noise.var, funnoise=funnoise, lower=lower, upper=upper,
		NoiseReEstimate=FALSE, CovReEstimate=FALSE)

new.model <- optim.result$model
best.x    <- optim.result$best.x
new.doe   <- optim.result$history.x

## Not run: 
##### DRAW RESULTS #####
# Compute actual function on a grid
n.grid <- 12
x.grid <- y.grid <- seq(0,1,length=n.grid)
design.grid <- expand.grid(x.grid, y.grid)
names(design.grid) <- c("V1","V2")
nt <- nrow(design.grid)
func.grid <- rep(0,1,nt)

for (i in 1:nt)
{ func.grid[i] <- test.function(x=design.grid[i,])}

# Compute initial and final kriging on a grid
pred <- predict(object=model, newdata=design.grid, type="UK", checkNames = FALSE)
mk.grid1 <- pred$m
sk.grid1 <- pred$sd

pred <- predict(object=new.model, newdata=design.grid, type="UK", checkNames = FALSE)
mk.grid2 <- pred$m
sk.grid2 <- pred$sd

# Plot initial kriging mean
z.grid <- matrix(mk.grid1, n.grid, n.grid)
filled.contour(x.grid,y.grid, z.grid, nlevels=50, color = topo.colors,
plot.axes = {title("Initial kriging mean");
points(model@X[,1],model@X[,2],pch=17,col="black"); 
axis(1); axis(2)})

# Plot initial kriging variance
z.grid <- matrix(sk.grid1^2, n.grid, n.grid)
filled.contour(x.grid,y.grid, z.grid, nlevels=50, color = topo.colors,
plot.axes = {title("Initial kriging variance");
points(model@X[,1],model@X[,2],pch=17,col="black"); 
axis(1); axis(2)})

# Plot final kriging mean
z.grid <- matrix(mk.grid2, n.grid, n.grid)
filled.contour(x.grid,y.grid, z.grid, nlevels=50, color = topo.colors,
plot.axes = {title("Final kriging mean");
points(new.model@X[,1],new.model@X[,2],pch=17,col="black"); 
axis(1); axis(2)})

# Plot final kriging variance
z.grid <- matrix(sk.grid2^2, n.grid, n.grid)
filled.contour(x.grid,y.grid, z.grid, nlevels=50, color = topo.colors,
plot.axes = {title("Final kriging variance");
points(new.model@X[,1],new.model@X[,2],pch=17,col="black"); 
axis(1); axis(2)})

# Plot actual function and observations
z.grid <- matrix(func.grid, n.grid, n.grid)
tit <- "Actual function - Black: initial points; red: added points"
filled.contour(x.grid,y.grid, z.grid, nlevels=50, color = topo.colors,
plot.axes = {title(tit);points(model@X[,1],model@X[,2],pch=17,col="black"); 
points(new.doe[1,],new.doe[2,],pch=15,col="red");
axis(1); axis(2)})

## End(Not run)

##########################################################################
### EXAMPLE 2: 3 OPTIMIZATION STEPS USING EQI WITH UNKNOWN NOISE       ###
### AND UNKNOWN COVARIANCE PARAMETERS FOR THE BRANIN FUNCTION          ###
##########################################################################
# Same initial model and parameters as for example 1
n.ite <- 2 # May be changed to a larger value 
res <- noisy.optimizer(optim.crit="min.quantile", 
optim.param=list(type="quantile",quantile=0.01),
model=model, n.ite=n.ite, noise.var=noise.var, funnoise=funnoise, 
lower=lower, upper=upper, 
control=list(print.level=0),CovReEstimate=TRUE, NoiseReEstimate=TRUE)

# Plot actual function and observations
plot(model@X[,1], model@X[,2], pch=17,xlim=c(0,1),ylim=c(0,1))
points(res$history.x[1,], res$history.x[2,], col="blue")

# Restart: requires the output estim.model of the previous run
# to deal with potential repetitions
res2 <- noisy.optimizer(optim.crit="min.quantile", 
optim.param=list(type="quantile",quantile=0.01), 
model=res$model, n.ite=n.ite, noise.var=noise.var, funnoise=funnoise, 
lower=lower, upper=upper, estim.model=res$estim.model,
control=list(print.level=0),CovReEstimate=TRUE, NoiseReEstimate=TRUE)

# Plot new observations
points(res2$history.x[1,], res2$history.x[2,], col="red")



[Package DiceOptim version 2.1.1 Index]