The following R code implements the given recursion for the Cholesky decomposition. The input is assumed to be a real symmetric matrix that is a positive-definite matrix.
mychol <- function(A) {
  n <- nrow(A)
  if (n == 1L) return(matrix(sqrt(A[1L, 1L]), 1L, 1L))
  j <- 2L:n
  L <- matrix(0, n, n)
  L[1L, 1L] <- sqrt(A[1L, 1L])
  b <- A[j, 1L, drop = FALSE] / L[1L, 1L]
  L[j, 1L] <- b
  S <- A[j, j, drop = FALSE] - tcrossprod(b)
  L[j, j] <- mychol(S)
  L
}
Here tcrossprod(b) is the outer product . The drop = FALSE arguments preserve matrix dimensions, including the final one-by-one subproblem, and the base case also returns a matrix.
To verify correctness, partition . The code sets , and recursively factors . This Schur complement is a positive-definite matrix: for ,
(1)
If is the recursive answer, then the returned lower triangular matrix satisfies
(2)
The one-dimensional identity starts this mathematical induction proof. Each recursive call reduces the matrix size by one and therefore terminates.

Ancestors (11)

  1. 3
  2. 2023
  3. Past exam of the simulation and statistical programming course of the University of Oxford
  4. Simulation and statistical programming course of the University of Oxford
  5. Year 2 of the mathematics course of the University of Oxford
  6. Mathematics course of the University of Oxford
  7. Undergraduate course of the University of Oxford
  8. Course of the University of Oxford
  9. University of Oxford
  10. List of universities
  11. Home

View article source

Discussion (0)

New discussion

There are no discussions about this article yet.

Articles by others on the same topic (0)

There are currently no matching articles.
See all articles in the same topic Create my own version