Cholesky decomposition
chol0(x)
A numeric matrix.
chol0
is implemented as (t o chol o t) becauase the Cholesky decomposition
in R (chol
)
returns an upper-triangle matrix;
uses only the upper half of the matrix when the matrix is real.
This does not match with the usual math notation of A = LL^T, where L is a lower triangular matrix. (Note that the Cholesky-Banachiewicz algorithm for example only uses the lower triangular part of A to compute L, so one should expect d L / d A to look like something you get by differentiating a lower triangular matrix w.r.t. a lower triangular matrix, i.e. the resulting Jacobian matrix is also lower triangular.)
To convert the R version of chol to the usual math version of chol, we define
chol0 <- function(x) { t(chol(t(x))) }
The first transpose ensures the output is lower-triangular, and the second
ensures the lower-triangular part of the input is used. Now, chol0
returns a lower-trianglar matrix
uses only the lower half of the matrix when the matrix is real
and this is what we want. (Additional note: finite-differencing with Cholesky is not too accurate due to the many floating point operations involved.)