## examining a particular laGP call from the larger analysis provided
## in the aGP documentation
## Simple 2-d test function used in Gramacy & Apley (2014);
## thanks to Lee, Gramacy, Taddy, and others who have used it before
f2d <- function(x, y=NULL)
{
if(is.null(y)) {
if(!is.matrix(x)) x <- matrix(x, ncol=2)
y <- x[,2]; x <- x[,1]
}
g <- function(z)
return(exp(-(z-1)^2) + exp(-0.8*(z+1)^2) - 0.05*sin(8*(z+0.1)))
z <- -g(x)*g(y)
}
## build up a design with N=~40K locations
x <- seq(-2, 2, by=0.02)
X <- as.matrix(expand.grid(x, x))
Z <- f2d(X)
## local analysis, first pass
Xref <- matrix(c(-1.725, 1.725), nrow=TRUE)
out <- laGP(Xref, 6, 50, X, Z, method="nn")
## second and pass via ALC, MSPE, and ALC-ray respectively
out2 <- laGP(Xref, 6, 50, X, Z, d=out$mle$d)
out2.mspe <- laGP(Xref, 6, 50, X, Z, d=out2$mle$d, method="mspe")
out2.alcray <- laGP(Xref, 6, 50, X, Z, d=out2$mle$d, method="alcray")
## look at the different designs
plot(rbind(X[out2$Xi,], X[out2.mspe$Xi,]), type="n",
xlab="x1", ylab="x2", main="comparing local designs")
points(Xref[1], Xref[2], col=2, cex=0.5)
text(X[out2$Xi,], labels=1:50, cex=0.6)
text(X[out2.mspe$Xi,], labels=1:50, cex=0.6, col="green")
text(X[out2.alcray$Xi,], labels=1:50, cex=0.6, col="red")
legend("topright", c("pass 2 ALC", "pass 3 MSPE", "pass 3 ALCRAY"),
text.col=c("black", "green", "red"), bty="n")
## compare times
data.frame(nn=out$time, alc=out2$time,
mspe=out2.mspe$time, alcray=out2.alcray$time)
## Here is the example from the Gramacy & Haaland (2014) paper;
## the lower lengthscale (d) setting generates more spread
out <- laGP(Xref, 6, 50, X, Z, d=0.1, method="alc")
out2 <- laGP(Xref, 6, 50, X, Z, d=0.1, method="alcray")
plot(X[out$Xi,], xlab="x1", ylab="x2", type="n")
text(X[out$Xi,], labels=1:length(out$Xi), cex=0.7)
text(X[out2$Xi,], labels=1:length(out2$Xi), cex=0.7, col=3)
points(Xref[1], Xref[2], pch=19, col=2)
legend("topright", c("exhaustive", "via rays"), text.col=c(1,3), bty="n")
Run the code above in your browser using DataLab