GPareto-package: Package GPareto

Description

Multi-objective optimization and quantification of uncertainty on Pareto fronts, using Gaussian process models.

Arguments

Author

Mickael Binois, Victor Picheny

Details

Important functions:
GParetoptim
easyGParetoptim
crit_optimizer
plotGPareto
CPF

References

M. Binois, D. Ginsbourger and O. Roustant (2015), Quantifying Uncertainty on Pareto Fronts with Gaussian process conditional simulations, European Journal of Operational Research, 243(2), 386-394.

O. Roustant, D. Ginsbourger and Yves Deville (2012), DiceKriging, DiceOptim: Two R Packages for the Analysis of Computer Experiments by Kriging-Based Metamodeling and Optimization, Journal of Statistical Software, 51(1), 1-55, tools:::Rd_expr_doi("10.18637/jss.v051.i01").

M. T. Emmerich, A. H. Deutz, J. W. Klinkenberg (2011), Hypervolume-based expected improvement: Monotonicity properties and exact computation, Evolutionary Computation (CEC), 2147-2154.

V. Picheny (2015), Multiobjective optimization using Gaussian process emulators via stepwise uncertainty reduction, Statistics and Computing, 25(6), 1265-1280.

T. Wagner, M. Emmerich, A. Deutz, W. Ponweiser (2010), On expected-improvement criteria for model-based multi-objective optimization, Parallel Problem Solving from Nature, 718-727, Springer, Berlin.

J. D. Svenson (2011), Computer Experiments: Multiobjective Optimization and Sensitivity Analysis, Ohio State University, PhD thesis.

C. Chevalier (2013), Fast uncertainty reduction strategies relying on Gaussian process models, University of Bern, PhD thesis.

M. Binois, V. Picheny (2019), GPareto: An R Package for Gaussian-Process-Based Multi-Objective Optimization and Analysis, Journal of Statistical Software, 89(8), 1-30, tools:::Rd_expr_doi("10.18637/jss.v089.i08").

Examples

Run this code

if (FALSE) {
#------------------------------------------------------------
# Example 1 : Surrogate-based multi-objective Optimization with postprocessing
#------------------------------------------------------------
set.seed(25468)

d <- 2 
fname <- P2

plotParetoGrid(P2) # For comparison

# Optimization
budget <- 25 
lower <- rep(0, d) 
upper <- rep(1, d)
     
omEGO <- easyGParetoptim(fn = fname, budget = budget, lower = lower, upper = upper)

# Postprocessing
plotGPareto(omEGO, add= FALSE, UQ_PF = TRUE, UQ_PS = TRUE, UQ_dens = TRUE)

}
#------------------------------------------------------------
# Example 2 : Surrogate-based multi-objective Optimization including a cheap function
#------------------------------------------------------------
set.seed(42)
library(DiceDesign)

d <- 2 

fname <- P1
n.grid <- 19
test.grid <- expand.grid(seq(0, 1, length.out = n.grid), seq(0, 1, length.out = n.grid))
nappr <- 15 
design.grid <- maximinESE_LHS(lhsDesign(nappr, d, seed = 42)$design)$design
response.grid <- t(apply(design.grid, 1, fname))

mf1 <- km(~., design = design.grid, response = response.grid[,1])
mf2 <- km(~., design = design.grid, response = response.grid[,2])
model <- list(mf1, mf2)

nsteps <- 1 
lower <- rep(0, d)
upper <- rep(1, d)

# Optimization with fastfun: hypervolume with discrete search

optimcontrol <- list(method = "discrete", candidate.points = test.grid)
omEGO2 <- GParetoptim(model = model, fn = fname, cheapfn = branin, crit = "SMS",
                      nsteps = nsteps, lower = lower, upper = upper,
                      optimcontrol = optimcontrol)
print(omEGO2$par)
print(omEGO2$values)

if (FALSE) { 
plotGPareto(omEGO2)

#------------------------------------------------------------
# Example 3 : Surrogate-based multi-objective Optimization (4 objectives)
#------------------------------------------------------------
set.seed(42)
library(DiceDesign)

d <- 5 

fname <- DTLZ3
nappr <- 25
design.grid <- maximinESE_LHS(lhsDesign(nappr, d, seed = 42)$design)$design
response.grid <- t(apply(design.grid, 1, fname, nobj = 4))
mf1 <- km(~., design = design.grid, response = response.grid[,1])
mf2 <- km(~., design = design.grid, response = response.grid[,2])
mf3 <- km(~., design = design.grid, response = response.grid[,3])
mf4 <- km(~., design = design.grid, response = response.grid[,4])

# Optimization
nsteps <- 5 
lower <- rep(0, d) 
upper <- rep(1, d)     
omEGO3 <- GParetoptim(model = list(mf1, mf2, mf3, mf4), fn = fname, crit = "EMI",
                      nsteps = nsteps, lower = lower, upper = upper, nobj = 4)
print(omEGO3$par)
print(omEGO3$values)
plotGPareto(omEGO3)

#------------------------------------------------------------
# Example 4 : quantification of uncertainty on Pareto front
#------------------------------------------------------------
library(DiceDesign)
set.seed(42)

nvar <- 2

# Test function P1
fname <- "P1"

# Initial design
nappr <- 10
design.grid <- maximinESE_LHS(lhsDesign(nappr, nvar, seed = 42)$design)$design
response.grid <- t(apply(design.grid, 1, fname))

PF <- t(nondominated_points(t(response.grid)))

# kriging models : matern5_2 covariance structure, linear trend, no nugget effect
mf1 <- km(~., design = design.grid, response = response.grid[,1])
mf2 <- km(~., design = design.grid, response = response.grid[,2])

# Conditional simulations generation with random sampling points 
nsim <- 100 # increase for better results
npointssim <- 1000 # increase for better results
Simu_f1 <- matrix(0, nrow = nsim, ncol = npointssim)
Simu_f2 <- matrix(0, nrow = nsim, ncol = npointssim)
design.sim <- array(0, dim = c(npointssim, nvar, nsim))

for(i in 1:nsim){
  design.sim[,,i] <- matrix(runif(nvar*npointssim), nrow = npointssim, ncol = nvar)
  Simu_f1[i,] <- simulate(mf1, nsim = 1, newdata = design.sim[,,i], cond = TRUE,
                         checkNames = FALSE, nugget.sim = 10^-8)
  Simu_f2[i,] <- simulate(mf2, nsim = 1, newdata = design.sim[,,i], cond = TRUE,
                         checkNames = FALSE, nugget.sim = 10^-8)
}

# Computation of the attainment function and Vorob'ev Expectation
CPF1 <- CPF(Simu_f1, Simu_f2, response.grid, paretoFront = PF)

summary(CPF1)

plot(CPF1)

# Display of the symmetric deviation function
plotSymDevFun(CPF1)
}