getSparseInv: Entries of a sparse inverse

In timflutre/rutilstimflutre: Timothee Flutre's personal R code

Entries of a sparse inverse

Description

Efficiently computes the entries of the sparse inverse (C) of a square, symmetric matrix (W) from its square-root-free Cholesky decomposition (W = L D L') as implemented in LDLT. This is notably used in quantitative genetics as described initially by Misztal and Perez-Enciso (1993). See the section 15.7 "Computing entries of the inverse of a sparse matrix" from Duff et al (2017) for the correct equation.

Usage

getSparseInv(vecD, L, reprSparse = NULL, verbose = FALSE)

Arguments

vecD	vector
L	lower triangular matrix
reprSparse	if not NULL, the output matrix will be a sparseMatrix of the given representation ("C", "R" or "T")
verbose	verbosity level

Value

matrix

Author(s)

Timothee Flutre

See Also

LDLT

Examples

## Not run: ## example matrix from Misztal and Perez-Enciso (1993):
W <- matrix(NA, nrow=5, ncol=5)
W[1,] <- c(5, 0, 0, 3, 0)
W[2,-1] <- c(3, 0, 0, 1)
W[3,-c(1:2)] <- c(4, 1, 1)
W[4,-c(1:3)] <- c(4, 0)
W[5,-c(1:4)] <- 2
W
W[lower.tri(W)] <- t(W)[lower.tri(W)]
W

## LDLT decomposition (a.k.a. square-root-free Cholesky):
decomp <- LDLT(W, returnFormat="list")

## C, inverse of W:
getSparseInv(decomp$vecD, decomp$L)
getSparseInv(decomp$vecD, decomp$L, "T")

## End(Not run)

timflutre/rutilstimflutre documentation built on Feb. 12, 2025, 11:35 p.m.