covMan: Creator Function for 'covMan' Objects

In kergp: Gaussian Process Laboratory

Creator Function for covMan Objects

Description

Creator function for covMan objects representing a covariance kernel entered manually.

Usage

   covMan(kernel, hasGrad = FALSE, acceptMatrix = FALSE,
          inputs = paste("x", 1:d, sep = ""), 
          d = length(inputs), parNames,
          par = NULL, parLower = NULL, parUpper = NULL,
          label = "covMan", ...) 

Arguments

kernel	A (semi-)positive definite function. This must be an object of class "function" with formal arguments named "x1", "x2" and "par". The first two formal arguments are locations vectors or matrices. The third formal is for the vector \boldsymbol{\theta} of all covariance parameters. An analytical gradient can be computed and returned as an attribute of the result with name "gradient". See Details.
hasGrad	Logical indicating whether the kernel function returns the gradient w.r.t. the vector of parameters as a "gradient" attribute of the result. See Details
acceptMatrix	Logical indicating whether kernel admits matrices as arguments. Default is FALSE. See Examples below.
inputs	Character vector giving the names of the inputs used as arguments of kernel. Optional if d is given.
d	Integer specifying the spatial dimension (equal to the number of inputs). Optional if inputs is given.
parNames	Vector of character strings containing the parameter names.
par, parLower, parUpper	Optional numeric vectors containing the parameter values, lower bounds and upper bounds.
label	Optional character string describing the kernel.
...	Not used at this stage.

Details

The formals and the returned value of the kernel function must be in accordance with the value of acceptMatrix.

• When acceptMatrix is FALSE, the formal arguments x1 and x2 of kernel are numeric vectors with length d. The returned result is a numeric vector of length 1. The attribute named "gradient" of the returned value (if provided in accordance with the value of hasGrad) must then be a numeric vector with length equal to the number of covariance parameters. It must contain the derivative of the kernel value K(\mathbf{x}_1,\,\mathbf{x}_2;\,\boldsymbol{\theta}) with respect to the parameter vector \boldsymbol{\theta}.

• When acceptMatrix is TRUE, the formals x1 and x2 are matrices with d columns and with n_1 and n_2 rows. The result is then a covariance matrix with n_1 rows and n_2 columns. The gradient attribute (if provided in accordance with the value of hasGrad) must be a list with length equal to the number of covariance parameters. The list element \ell must contain a numeric matrix with dimension (n_1, n_2) which is the derivative of the covariance matrix w.r.t. the covariance parameter \theta_\ell.

Note

The kernel function must be symmetric with respect to its first two arguments, and it must be positive definite, which is not checked. If the function returns an object with a "gradient" attribute but hasGrad was set to FALSE, the gradient will not be used in optimization.

The name of the class was motivated by earlier stages in the development.

Examples

myCovMan <- 
      covMan(
         kernel = function(x1, x2, par) { 
         htilde <- (x1 - x2) / par[1]
         SS2 <- sum(htilde^2)
         d2 <- exp(-SS2)
         kern <- par[2] * d2
         d1 <- 2 * kern * SS2 / par[1]            
         attr(kern, "gradient") <- c(theta = d1,  sigma2 = d2)
         return(kern)
      },
      hasGrad = TRUE,    
      d = 1,
      label = "myGauss",
      parLower = c(theta = 0.0, sigma2 = 0.0),
      parUpper = c(theta = Inf, sigma2 = Inf),
      parNames = c("theta", "sigma2"),
      par = c(NA, NA)
      )
      
# Let us now code the same kernel in C
kernCode <- "
       SEXP kern, dkern;
       int nprotect = 0, d;
       double SS2 = 0.0, d2, z, *rkern, *rdkern;

       d = LENGTH(x1);
       PROTECT(kern = allocVector(REALSXP, 1)); nprotect++;
       PROTECT(dkern = allocVector(REALSXP, 2)); nprotect++;
       rkern = REAL(kern);
       rdkern = REAL(dkern);

       for (int i = 0; i < d; i++) {
         z = ( REAL(x1)[i] - REAL(x2)[i] ) / REAL(par)[0];
         SS2 += z * z; 
       }

       d2 = exp(-SS2);
       rkern[0] = REAL(par)[1] * d2;
       rdkern[1] =  d2; 
       rdkern[0] =  2 * rkern[0] * SS2 / REAL(par)[0];

       SET_ATTR(kern, install(\"gradient\"), dkern);
       UNPROTECT(nprotect);
       return kern;
     "
     
myCovMan  

## "inline" the C function into an R function: much more efficient! 

## Not run: 
require(inline)
kernC <- cfunction(sig = signature(x1 = "numeric", x2 = "numeric",
                                   par = "numeric"),
                    body = kernCode)
myCovMan <- covMan(kernel = kernC, hasGrad = TRUE, d = 1,
                   parNames = c("theta", "sigma2"))
myCovMan

## End(Not run)

## A kernel admitting matricial arguments
myCov <- covMan(
    
    kernel = function(x1, x2, par) { 
      # x1 : matrix of size n1xd
      # x2 : matrix of size n2xd
            
      d <- ncol(x1)
            
      SS2 <- 0
      for (j in 1:d){
        Aj <- outer(x1[, j], x2[, j], "-")
        Aj2 <- Aj^2
        SS2 <- SS2 + Aj2 / par[j]^2
      }
      D2 <- exp(-SS2)
      kern <- par[d + 1] * D2
    },
    acceptMatrix = TRUE,
    d = 2,
    label = "myGauss",
    parLower = c(theta1 = 0.0, theta2 = 0.0, sigma2 = 0.0),
    parUpper = c(theta1 = Inf, theta2 = Inf, sigma2 = Inf),
    parNames = c("theta1", "theta2", "sigma2"),
    par = c(NA, NA, NA)

)
      
coef(myCov) <- c(0.5, 1, 4)
show(myCov)

## computing the covariance kernel between two points
X <- matrix(c(0, 0), ncol = 2)
Xnew <- matrix(c(0.5, 1), ncol = 2)
colnames(X) <- colnames(Xnew) <- inputNames(myCov)
covMat(myCov, X)            ## same points
covMat(myCov, X, Xnew)      ## two different points

# computing covariances between sets of given locations
X <- matrix(c(0, 0.5, 0.7, 0, 0.5, 1), ncol = 2)
t <- seq(0, 1, length.out = 3)
Xnew <- as.matrix(expand.grid(t, t))
colnames(X) <- colnames(Xnew) <- inputNames(myCov)
covMat(myCov, X)         ## covariance matrix
covMat(myCov, X, Xnew)   ## covariances between design and new data

kergp documentation built on April 4, 2025, 2:29 a.m.