# NOT RUN {
set.seed(123)
# a 9-points factorial design, and the corresponding response
d <- 2; n <- 9
design.fact <- expand.grid(seq(0,1,length=3), seq(0,1,length=3))
names(design.fact)<-c("x1", "x2")
design.fact <- data.frame(design.fact)
names(design.fact)<-c("x1", "x2")
response.branin <- apply(design.fact, 1, branin)
response.branin <- data.frame(response.branin)
names(response.branin) <- "y"
# model identification
fitted.model1 <- km(~1, design=design.fact, response=response.branin,
covtype="gauss", control=list(pop.size=50,trace=FALSE), parinit=c(0.5, 0.5))
# graphics
n.grid <- 9 # Increase to 50 for a nicer picture
x.grid <- y.grid <- seq(0,1,length=n.grid)
design.grid <- expand.grid(x.grid, y.grid)
#response.grid <- apply(design.grid, 1, branin)
EI.grid <- apply(design.grid, 1, EI,fitted.model1)
#EI.grid <- apply(design.grid, 1, EI.plot,fitted.model1, gr=TRUE)
z.grid <- matrix(EI.grid, n.grid, n.grid)
contour(x.grid,y.grid,z.grid,20)
title("Expected Improvement for the Branin function known at 9 points")
points(design.fact[,1], design.fact[,2], pch=17, col="blue")
# graphics
n.gridx <- 5 # increase to 15 for nicer picture
n.gridy <- 5 # increase to 15 for nicer picture
x.grid2 <- seq(0,1,length=n.gridx)
y.grid2 <- seq(0,1,length=n.gridy)
design.grid2 <- expand.grid(x.grid2, y.grid2)
EI.envir <- new.env()
environment(EI) <- environment(EI.grad) <- EI.envir
for(i in seq(1, nrow(design.grid2)) )
{
x <- design.grid2[i,]
ei <- EI(x, model=fitted.model1, envir=EI.envir)
eigrad <- EI.grad(x , model=fitted.model1, envir=EI.envir)
if(!(is.null(ei)))
{
suppressWarnings(arrows(x$Var1,x$Var2,
x$Var1 + eigrad[1]*2.2*10e-5, x$Var2 + eigrad[2]*2.2*10e-5,
length = 0.04, code=2, col="orange", lwd=2))
}
}
# }
Run the code above in your browser using DataLab