mgcv (version 1.8-26)

choldrop: Down date a Cholesky factor on dropping a row/col

Description

Given a Cholesky factor, R, of a matrix, A, finds the Cholesky factor of A[-k,-k], where k is an integer.

Usage

choldrop(R,k)

Arguments

R

Cholesky factor of a matrix, A.

k

row and column of A to drop.

Details

If R is upper triangular then t(R[,-k])%*%R[,-k] == A[-k,-k], but R[,-k] has elements on the first sub-diagonal, from its kth column onwards. To get from this to a triangular Cholesky factor of A[-k,-k] we can apply a sequence of Given rotations from the left to eliminate the sub-diagonal elements. The routine does this. If R is a lower triangular factor then Givens rotations from the right are needed to remove the extra elements. If n is the dimension of R then the update has O(n^2) computational cost.

Note that the update is vector oriented, and is hence not susceptible to speed up by use of an optimized BLAS. The update is set up to be relatively Cache friendly, in that in the upper triangular case successive Givens rotations are stored for sequential application columnwise, rather than being applied rowwise as soon as they are computed. Even so, the upper triangular update is slightly slower than the lower triangular update.

Examples

# NOT RUN {
  require(mgcv)
  set.seed(0)
  n <- 6
  A <- crossprod(matrix(runif(n*n),n,n))
  R0 <- chol(A)
  k <- 3
  Rd <- choldrop(R0,k)
  range(Rd-chol(A[-k,-k]))
  Rd;chol(A[-k,-k])
  
  ## same but using lower triangular factor A = LL'
  L <- t(R0)
  Ld <- choldrop(L,k)
  range(Ld-t(chol(A[-k,-k])))
  Ld;t(chol(A[-k,-k]))
# }