C_GP: C matrix closed form expression for a GP.

In activegp: Gaussian Process Based Design and Analysis for the Active Subspace Method

C matrix closed form expression for a GP.

Description

Computes the integral over the input domain of the outer product of the gradients of a Gaussian process. The corresponding matrix is the C matrix central in active subspace methodology.

Usage

C_GP(
  modelX,
  y,
  measure = "lebesgue",
  xm = NULL,
  xv = NULL,
  S = NULL,
  verbose = TRUE
)

Arguments

modelX	This may be either 1) a homGP or hetGP GP model, see hetGP-package containing, e.g., a vector of thetas, type of covariance ct, an inverse covariance matrix Ki, a design matrix X0, and response vector Z0. 2) A matrix of design locations, in which case a vector of responses must be given as the y argument, and this function will fit a default model for you.
y	A vector of responses corresponding to the design matrix; may be ommited if a GP fit is provided in the modelX argument.
measure	One of c("lebesgue", "gaussian", "trunc_gaussian", "sample", "discrete"), indiciating the probability distribution with respect to which the input points are drawn in the definition of the active subspace. "lebesgue" uses the Lebesgue or Uniform measure over the unit hypercube [0,1]^d. "gaussian" uses a Gaussian or Normal distribution, in which case xm and xv should be specified. "trunc_gaussian" gives a truncated Gaussian or Normal distribution over the unit hypercube [0,1]^d, in which case xm and xv should be specified. "sample" gives the Sample or Empirical measure (dirac deltas located at each design point), which is equivalent to calculating the average expected gradient outer product at the design points. "discrete" gives a measure which puts equal weight at points in the input space specified via the S parameter, which should be a matrix with one row for each atom of the measure.
xm	If measure is "gaussian" or "trunc_gaussian", gives the mean vector.
xv	If measure is "gaussian" or "trunc_gaussian", gives the marginal variance vector. The covariance matrix is assumed to be diagonal.
S	If measure is "discrete", gives the locations of the measure's atoms.
verbose	Should we print progress?

Value

a const_C object with elements

• model: GP model provided or estimated;

• mat: C matrix estimated;

• Wij: list of W matrices, of size number of variables;

• ct: covariance type (1 for "Gaussian", 2 for "Matern3_2", 3 for "Matern5_2").

References

N. Wycoff, M. Binois, S. Wild (2019+), Sequential Learning of Active Subspaces, preprint.

P. Constantine (2015), Active Subspaces, Philadelphia, PA: SIAM.

See Also

print.const_C, plot.const_C

Examples

################################################################################
### Active subspace of a Gaussian process
################################################################################
 
library(hetGP); library(lhs)
set.seed(42)

nvar <- 2
n <- 100

# theta gives the subspace direction
f <- function(x, theta, nugget = 1e-3){
  if(is.null(dim(x))) x <- matrix(x, 1)
  xact <- cos(theta) * x[,1] - sin(theta) * x[,2]
  return(hetGP::f1d(xact) + rnorm(n = nrow(x), sd = rep(nugget, nrow(x))))
}

theta_dir <- pi/6
act_dir <- c(cos(theta_dir), -sin(theta_dir))

# Create design of experiments and initial GP model
design <- X <- matrix(signif(maximinLHS(n, nvar), 2), ncol = nvar)
response <- Y <- apply(design, 1, f, theta = theta_dir)
model <- mleHomGP(design, response, known = list(beta0 = 0))

C_hat <- C_GP(model)

# Subspace distance to true subspace:
print(subspace_dist(C_hat, matrix(act_dir, nrow = nvar), r = 1))
plot(design %*% eigen(C_hat$mat)$vectors[,1], response, 
  main = "Projection along estimated active direction")
plot(design %*% eigen(C_hat$mat)$vectors[,2], response, 
  main = "Projection along estimated inactive direction")
  
# For other plots:
# par(mfrow = c(1, 3)) # uncomment to have all plots together
plot(C_hat)
# par(mfrow = c(1, 1)) # restore graphical window

 

activegp documentation built on June 28, 2022, 1:05 a.m.