crit_EHI | R Documentation |
Multi-objective Expected Hypervolume Improvement with respect to the current Pareto front. With two objectives the analytical formula is used, while Sample Average Approximation (SAA) is used with more objectives. To avoid numerical instabilities, the new point is penalized if it is too close to an existing observation.
crit_EHI( x, model, paretoFront = NULL, critcontrol = list(nb.samp = 50, seed = 42), type = "UK" )
x |
a vector representing the input for which one wishes to calculate |
model |
list of objects of class |
paretoFront |
(optional) matrix corresponding to the Pareto front of size |
critcontrol |
optional list with arguments:
Options for the |
type |
" |
The computation of the analytical formula with two objectives is adapted from the Matlab source code by Michael Emmerich and Andre Deutz, LIACS, Leiden University, 2010 available here : http://liacs.leidenuniv.nl/~csmoda/code/HV_based_expected_improvement.zip.
The Expected Hypervolume Improvement at x
.
J. D. Svenson (2011), Computer Experiments: Multiobjective Optimization and Sensitivity Analysis, Ohio State University, PhD thesis.
M. T. Emmerich, A. H. Deutz, J. W. Klinkenberg (2011), Hypervolume-based expected improvement: Monotonicity properties and exact computation,
Evolutionary Computation (CEC), 2147-2154.
EI
from package DiceOptim, crit_EMI
, crit_SUR
, crit_SMS
.
#--------------------------------------------------------------------------- # Expected Hypervolume Improvement surface associated with the "P1" problem at a 15 points design #--------------------------------------------------------------------------- set.seed(25468) library(DiceDesign) n_var <- 2 f_name <- "P1" n.grid <- 26 test.grid <- expand.grid(seq(0, 1, length.out = n.grid), seq(0, 1, length.out = n.grid)) n_appr <- 15 design.grid <- round(maximinESE_LHS(lhsDesign(n_appr, n_var, seed = 42)$design)$design, 1) response.grid <- t(apply(design.grid, 1, f_name)) Front_Pareto <- t(nondominated_points(t(response.grid))) mf1 <- km(~., design = design.grid, response = response.grid[,1]) mf2 <- km(~., design = design.grid, response = response.grid[,2]) EHI_grid <- crit_EHI(x = as.matrix(test.grid), model = list(mf1, mf2), critcontrol = list(refPoint = c(300,0))) filled.contour(seq(0, 1, length.out = n.grid), seq(0, 1, length.out = n.grid), nlevels = 50, matrix(EHI_grid, n.grid), main = "Expected Hypervolume Improvement", xlab = expression(x[1]), ylab = expression(x[2]), color = terrain.colors, plot.axes = {axis(1); axis(2); points(design.grid[,1], design.grid[,2], pch = 21, bg = "white") } )
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.