lineqGPR (version 0.1.1)

simulate.lineqAGP: Simulation Method for the "lineqAGP" S3 Class

Description

Simulation method for the "lineqAGP" S3 class.

Usage

# S3 method for lineqAGP
simulate(object, nsim = 1, seed = NULL, xtest, ...)

Arguments

object

an object with class "lineqAGP".

nsim

the number of simulations.

seed

see simulate.

xtest

a vector (or matrix) with the test input design.

...

further arguments passed to or from other methods.

Value

A "lineqAGP" object with the following elements.

x

a vector (or matrix) with the training input design.

y

the training output vector at x.

xtest

a vector (or matrix) with the test input design.

Phi.test

a matrix corresponding to the hat basis functions evaluated at xtest. The basis functions are indexed by rows.

xi.sim

the posterior sample-path of the finite-dimensional Gaussian vector.

ysim

the posterior sample-path of the observed GP. Note: ysim = Phi.test %*% xi.sim.

Details

The posterior sample-path of the finite-dimensional GP with linear inequality constraints are computed. Here, \(\boldsymbol{\xi}\) is a centred Gaussian vector with covariance \(\boldsymbol{\Gamma}\), s.t. \(\boldsymbol{\Phi} \boldsymbol{\xi} = \boldsymbol{y}\) (interpolation constraints) and \(\boldsymbol{l} \leq \boldsymbol{\Lambda} \boldsymbol{\xi} \leq \boldsymbol{u}\) (inequality constraints).

References

Lopez-Lopera, A. F., Bachoc, F., Durrande, N., and Roustant, O. (2017), "Finite-dimensional Gaussian approximation with linear inequality constraints". ArXiv e-prints [link]

See Also

create.lineqAGP, augment.lineqAGP, predict.lineqAGP

Examples

# NOT RUN {
library(plot3D)
# creating the model
d <- 2
fun1 <- function(x) return(4*(x-0.5)^2)
fun2 <- function(x) return(2*x)
targetFun <- function(x) return(fun1(x[, 1]) + fun2(x[, 2])) 
xgrid <- expand.grid(seq(0, 1, 0.01), seq(0, 1, 0.01))
ygrid <- targetFun(xgrid)
xdesign <- rbind(c(0.5, 0), c(0.5, 0.5), c(0.5, 1), c(0, 0.5), c(1, 0.5))
ydesign <- targetFun(xdesign)
model <- create(class = "lineqAGP", x = xdesign, y = ydesign,
                constrType = c("convexity", "monotonicity"))

# updating and expanding the model
model$localParam$m <- rep(10, d)
model$kernParam[[1]]$type <- "matern52"
model$kernParam[[2]]$type <- "matern52"
model$kernParam[[1]]$par <- c(1, 0.2)
model$kernParam[[2]]$par <- c(1, 0.3)
model$nugget <- 1e-9
model$varnoise <- 1e-5
model <- augment(model)

# sampling from the model
ntest <- 25
xtest  <- cbind(seq(0, 1, length = ntest), seq(0, 1, length = ntest))
ytest <- targetFun(xtest)
sim.model <- simulate(model, nsim = 1e3, seed = 1, xtest = xtest)
persp3D(x = unique(xtest[, 1]), y = unique(xtest[, 2]),
        z = outer(rowMeans(sim.model$ysim[[1]]),
                  rowMeans(sim.model$ysim[[2]]), "+"),
        xlab = "x1", ylab = "x2", zlab = "mode(x1,x2)", zlim = c(0, 3),
        phi = 20, theta = -30, alpha = 1, colkey = FALSE)
points3D(x = xdesign[,1], y = xdesign[,2], z = ydesign, col = "black", pch = 19, add = TRUE)

# }