mbo_model_fail_rate.R

library(mlrMBO)
library(data.table)
library(ggplot2)
configureMlr(on.par.without.desc = "warn")
set.seed(1)
# Define objective function
fn = makeRosenbrockFunction(5)

# define mbo control object
ctrl = makeMBOControl(propose.points = 4, store.model.at = 0:20)
ctrl = setMBOControlInfill(ctrl, crit = "ei")
ctrl = setMBOControlMultiPoint(ctrl, method = "cl", cl.lie = max)
ctrl = setMBOControlTermination(ctrl, iters = 6)

# surrogate learner
lrn = makeLearner("regr.km", nugget.stability = 10^-8, covtype = "matern5_2", predict.type = "se", multistart = 1)

# initial design n = sqrt(d) * 10
des = generateDesign(par.set = smoof::getParamSet(fn), n = floor(10 * sqrt(getParamLengths(smoof::getParamSet(fn)))))

# start mbo
res = mbo(fn, design = des, learner = lrn, control = ctrl)

# evaluation
op.dt = as.data.table(res$opt.path)

# generates surrogate model for a given mbo-iteration (i.dob) and for a specific lie (i.lie)
# i.lie = 0 generates the model without any liar
reconstructModel = function(op.dt, i.dob = 3, i.lie = 1) {
  base = op.dt[cumsum(dob == i.dob) <= i.lie + lie.length & dob <= i.dob | dob == 0]
  base[cumsum(dob == i.dob) <= i.lie & dob == i.dob, y := replicate(i.lie, ctrl$multipoint.cl.lie(base[dob <= i.dob, y]))]
  tsk.data = as.data.frame(base[, c(ctrl$y.name, getParamIds(smoof::getParamSet(fn), with.nr = TRUE, repeated = TRUE)) , with = FALSE])
  tsk = makeRegrTask(id = "base", data = tsk.data, target = ctrl$y.name)
  mod = train(lrn, tsk)
  mod
}

# function to simply check if the model only proposas constant values
checkModelConstant = function(mod) {
  grid = generateDesign(n = 200, par.set = smoof::getParamSet(fn))
  #grid = generateGridDesign(par.set = smoof::getParamSet(fn), resolution = 100)
  pred = predict(mod, newdata = grid)
  sd(getPredictionResponse(pred)) < .Machine$double.eps ^ 0.5
}

# generate all dob and lie combinations
dob.lie = expand.grid(lie = seq_len(ctrl$propose.points)-1, dob = unique(op.dt$dob)[-1])
# add first model from init design
dob.lie = rbind(c(0,0), dob.lie)
was.constant = Map(function(lie, dob) {
    mean(replicate(20, checkModelConstant(reconstructModel(op.dt, dob, lie))))
  }, 
  lie = dob.lie$lie, dob = dob.lie$dob)
dob.lie$fail = unlist(was.constant)
g = ggplot(dob.lie, aes (x = sprintf("%i;%i", dob, lie), y = fail, fill = lie==0))
g + geom_bar(stat="identity") + xlab("dob; lie") + ggtitle("Kriging Failrate")

Sample Output

The same but with multistart = 10:

cl.lie = min and multistart = 1

cl.lie = function(y) runif(1, min = min(y), max = max(y)) and multistart = 1