msqrt {npreg}

Matrix (Inverse?) Square Root

Description

Stable computation of the square root (or inverse square root) of a positive semi-definite matrix.

Usage

msqrt(x, inverse = FALSE, symmetric = FALSE, 
      tol = .Machine$double.eps, checkx = TRUE)

Arguments

x

positive semi-definite matrix

inverse

compute inverse square root?

symmetric

does the square root need to be symmetric? See Details.

tol

tolerance for detecting linear dependencies in x

checkx

should x be checked for symmetry using isSymmetric?

Details

If symmetric = FALSE, this function computes the matrix z such that x = tcrossprod(z)

If symmetric = TRUE, this function computes the matrix z such that x = crossprod(z) = tcrossprod(z)

If inverse = TRUE, the matrix x is replaced by the pseudo-inverse of x in these equations (see psolve)

Value

The matrix z that gives the (inverse?) square root of x. See Details.

Note

The matrix (inverse?) square root is calculated by (inverting and) square rooting the eigenvalues that are greater than the first value multiplied by tol * nrow(x)

Author(s)

Nathaniel E. Helwig <[e-mail]>

See Also

psolve

Examples

# generate x
set.seed(0)
x <- crossprod(matrix(rnorm(100), 20, 5))

# asymmetric square root (default)
xsqrt <- msqrt(x)
mean(( x - crossprod(xsqrt) )^2)
mean(( x - tcrossprod(xsqrt) )^2)

# symmetric square root
xsqrt <- msqrt(x, symmetric = TRUE)
mean(( x - crossprod(xsqrt) )^2)
mean(( x - tcrossprod(xsqrt) )^2)

# asymmetric inverse square root (default)
xsqrt <- msqrt(x, inverse = TRUE)
mean(( solve(x) - crossprod(xsqrt) )^2)
mean(( solve(x) - tcrossprod(xsqrt) )^2)

# symmetric inverse square root
xsqrt <- msqrt(x, inverse = TRUE, symmetric = TRUE)
mean(( solve(x) - crossprod(xsqrt) )^2)
mean(( solve(x) - tcrossprod(xsqrt) )^2)