GLS_chol: GLS Estimate using Cholesky Factor

View source: R/utils.R

GLS_cholR Documentation

GLS Estimate using Cholesky Factor


Computes the GLS estimate using the formula:

μ_{GLS} = (X^\top Σ^{-1} X)^{-1}X^\top Σ^{-1} y.

The computation is done depending on the input class of the Cholesky factor R. It relies on the classical solve or on using forwardsolve and backsolve functions of package spam, see solve. This is much faster than computing the inverse of Σ, especially since we have to compute the Cholesky decomposition of Σ either way.


GLS_chol(R, X, y)

## S3 method for class 'spam.chol.NgPeyton'
GLS_chol(R, X, y)

## S3 method for class 'matrix'
GLS_chol(R, X, y)



(spam.chol.NgPeyton or matrix(n, n))
Cholesky factor of the covariance matrix Σ. If covariance tapering and sparse matrices are used, then the input is of class spam.chol.NgPeyton. Otherwise, R is the output of a standard chol, i.e., a simple matrix


(matrix(n, p))
Data / design matrix.


Response vector


A numeric(p) vector, i.e., the mean effects.


Jakob Dambon


# generate data
n <- 10
X <- cbind(1, 20+1:n)
y <- rnorm(n)
A <- matrix(runif(n^2)*2-1, ncol=n)
Sigma <- t(A) %*% A
# two possibilities
## using standard Cholesky decomposition
R_mat <- chol(Sigma); str(R_mat)
mu_mat <- GLS_chol(R_mat, X, y)
## using spam
R_spam <- chol(spam::as.spam(Sigma)); str(R_spam)
mu_spam <- GLS_chol(R_spam, X, y)
# should be identical to the following
mu <- solve(crossprod(X, solve(Sigma, X))) %*%
      crossprod(X, solve(Sigma, y))
## check
abs(mu - mu_mat)
abs(mu - mu_spam)

varycoef documentation built on Sept. 18, 2022, 1:07 a.m.