laGP (version 1.0)

predGP: GP Prediction/Kriging

Description

Perform Gaussian processes prediction at new XX locations using a GP object stored on the C side

Usage

predGP(gpi, XX, lite = FALSE)

Arguments

gpi

a C-side GP object identifier (positive integer); e.g., as returned by newGP

XX

a matrix or data.frame containing a design of predictive locations

lite

a scalar logical indicating whether (lite = FALSE, default) or not (lite = TRUE) a full predictive covariance matrix should be returned

Value

• The output is a list with the following components.
• mean
  a vector of predictive means of length nrow(Xref)
• Sigma
  a Sigma containing the covariance matrix of for a multivariate Student-t distribution; or if lite = TRUE, then the field sigma2 contains the diagonal of this matrix
• df
  a Student-t degrees of freedom scalar (applies to all XX)

Details

Returns the parameters of Student-t predictive equations. By default, these include a full predictive covariance matrix between all XX locations. However, this can be slow when nrow(XX) is large, so a lite options is provided, which only returns the diagonal of that matrix.

GP prediction is sometimes called kriging, especially in the spatial statistics literature. So this function could also be described as returning the kriging equations

References

for standard GP prediction, refer to any graduate text, e.g., Rasmussen & Williams Gaussian Processes for Machine Learning

See Also

newGP, mleGP, jmleGP, predGP

Examples

## a "computer experiment" -- a much smaller version than the one shown
## in the aGP docs

## Simple 2-d test function used in Gramacy & Apley (2013);
## 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)
  }

## design with N=441
x <- seq(-2, 2, length=11)
X <- as.matrix(expand.grid(x, x))
Z <- f2d(X)

## fit a GP
gpi <- newGP(X, Z, d=0.35, g=1/1000)

## predictive grid with NN=400
xx <- seq(-1.9, 1.9, length=20)
XX <- as.matrix(expand.grid(xx, xx))
ZZ <- f2d(XX)

## predict
p <- predGP(gpi, XX, lite=TRUE)
## RMSE: compare to similar experiment in aGP docs
sqrt(mean((p$mean - ZZ)^2))

## visualize the result
par(mfrow=c(1,2))
image(xx, xx, matrix(p$mean, nrow=length(xx)), col=heat.colors(128),
      xlab="x1", ylab="x2", main="predictive mean")
image(xx, xx, matrix(p$mean-ZZ, nrow=length(xx)), col=heat.colors(128),
      xlab="x1", ylab="x2", main="bas")

## clean up
deleteGP(gpi)

## see the newGP and mleGP docs for examples using lite = FALSE for
## sampling from the joint predictive distribution