sqrt_methods: Matrix factorization functions

In watanabe-j/eigvaldisp: Statistics of Eigenvalue Dispersion Indices

Matrix factorization functions

Description

Assuming A is a symmetric, non-negative definite matrix:

matsqrt(A) returns a matrix square root of A with spectral or singular value decomposition

chol_piv(A) conducts Cholesky factorization of A with “pivoting”. crossprod(chol_piv(A)) will be equal to A, unless multiple zero eigenvalues exist (in which case the results can be spurious; a warning results).

chol_qr(A) conducts Cholesky factorization of A with QR factorization. crossprod(chol_qr(A)) is equal to A, even with multiple zero eigenvalues. This is more reliable, but slower than chol_piv(). matsqrt(), which is called internally, usually suffices.

chol_qrsvd(A) is essentially chol_qr(A, method = "svd")

Usage

matsqrt(A, method = c("svd", "eigen"))

chol_piv(A)

chol_qr(A, method = c("svd", "eigen"))

chol_qrsvd(A)

Arguments

A	Numeric matrix, assumed to be symmetric
method	Either "svd" (default) or "eigen" to specify the function to be used for eigendecomposition. Results should be identical, but svd would be faster in most environments.

Details

Cholesky factors are useful in simulating multivariate normal variates, but chol(A) by default returns an error unless A is strictly positive definite (i.e., singular matrices are not allowed). One way to avoid this is to use chol(A, pivot = TRUE) and arrange the columns appropriately (see the documentation of chol). chol_piv(A) is a utility to conduct this.

Unexpectedly, however, this can still fail in certain situations, e.g., when more than one eigenvalues are 0. This is despite that a Cholesky factorization can be defined for any nonnegative definite matrix, as can be confirmed via QR factorization, (e.g., Schott, 2017; although this factor may not be unique). chol_piv() returns a warning in this case. chol_qr() handles this factorization via a naive implementation; first taking a matrix square root via eigendecomposition (matsqrt(A)) and then conducting QR factorization of this (with pivoting). This allows for a correct Cholesky factorization of any positive semidefinite matrix, although for most simulation purposes the matrix square root by matsqrt() usually suffices.

Value

matsqrt(A): Symmetric matrix which is a square root of A

chol_piv(A)/chol_qr(A): Upper triangular matrix which is a Cholesky factor of A

References

Schott, J. R. (2017) Matrix Analysis for Statistics, 3rd edition. John Wiley & Sons, Hoboken, New Jersey.

See Also

chol, qr

Examples

(A <- diag(c(2, 1, 1, 0)))
## Not run: chol(A) # This returns an error because of singularity
cA <- eigvaldisp:::chol_piv(A)
all.equal(A, crossprod(cA))
# TRUE, as expected

B <- matrix(1, 4, 4)
B[1:2, 3:4] <- B[3:4, 1:2] <- 0
print(B)
cB1 <- eigvaldisp:::chol_piv(B)
all.equal(B, crossprod(cB1))
# not TRUE! (though perhaps environment-dependent)
crossprod(cB1)

cB2 <- eigvaldisp:::chol_qr(B)
all.equal(B, crossprod(cB2))
# TRUE, as it should be
crossprod(cB2)

watanabe-j/eigvaldisp documentation built on Dec. 8, 2023, 4:38 a.m.