C_GP {activegp} R Documentation

## C matrix closed form expression for a GP.

### Description

Computes the integral (uniform measure) over the [0,1] box domain of the outer product of the gradients of a Gaussian process. The corresponding matrix is the C matrix central in the active subspace methodology.

### Usage

```C_GP(model)
```

### Arguments

 `model` `homGP` or `hetGP` GP model, see `hetGP-package` containing, e.g., a vector of `theta`s, type of covariance `ct`, an inverse covariance matrix `Ki`, a design matrix `X0`, and response vector `Z0`.

### Value

a `const_C` object with elements

• `model`: GP model provided;

• `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.

`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

```

[Package activegp version 1.0.6 Index]