designL.from.Corr: Computation of “square root” of correlation matrix

Description

This function is somewhat obsolete; see mat_sqrt.

This function is not usually directly called by users, but arguments may be passed to it through higher-level calls (see Examples). For given correlation matrix C, it computes a “design matrix” L such that C = L * t(L). Cholesky factorization (performed by t(base::chol(.)) or some other implementation) is a fast method for this computation, but it is not robust numerically and may even return an error, in which cases more robust methods are used. Matrix roots are not unique (for example, they are lower triangular for t(chol(.)), and symmetric for svd(.). As matrix roots are used to simulate samples under the fitted model (in particular in the parametric bootstrap implemented in fixedLRT), this implies that for given seed of random numbers, these samples will differ with these different methods (although their distribution should be identical).

Usage

designL.from.Corr(m = NULL, symSVD = NULL, try.chol = TRUE, try.eigen = FALSE, 
                  threshold = 1e-06, SVDfix = 1/10)

Arguments

The matrix which 'root' is to be computed. This argument is ignored if symSVD is provided.

symSVD

A list representing the symmetric singular value decomposition of the matrix which 'root' is to be computed. Must have elements $u, a matrix of eigenvectors, and $d, a vector of eigenvalues.

try.chol

If try.chol=TRUE, the Cholesky factorization will be tried.

try.eigen

The default behavior is to try chol, and use svd if chol fails. If try.eigen=TRUE, the sym_eigen factorization will be tried before svd.

threshold

A correction threshold for low eigenvalues is the case and eigensystem or singular-value decomposition are used.

SVDfix

A solution to failures of svd: see Details.

Value

The “square root of the input matrix”. Its rows and columns are labelled according to the columns of the original matrix.

Details

The function may call svd, for singular value decomposition (SVD) of a matrix M. svd may return “error code 1 from Lapack routine 'dgesdd'” (cf. unhelpful discussions on R forums). This can be circumvented by computing the SVD of $(1 - x)$ I $+ x$ M and deducing the singular values of M in a trivial way. The $x$ value to be used in this fix is provided by the SVDfix argument.

svd errors have occurred for correlation matrices that were close to the identity matrix except for a few large non-diagonal elements. Such matrices may occur in particular for the Mat<U+00E9>rn correlation model with low $ν$ , high $ρ$ , and if some samples are spatially close. Then, an alternative fix to the svd problem may be to restrict the $ν$ and/or $ρ$ ranges, using the lower and upper arguments of corrHLfit, although one should make sure that this has no bearing on the inferences.

Examples

Run this code

# NOT RUN {
## try.chol argument passed to designL.from.Corr 
## through the '...' argument of higher-level functions
## such as HLCor, corrHLfit, fixedLRT:
data("scotlip")
opt <- spaMM.options(mat_sqrt_fr="designL.from.Corr")
HLCor(cases~I(prop.ag/10) +adjacency(1|gridcode)+offset(log(expec)),
      ranPars=list(rho=0.174),adjMatrix=Nmatrix,family=poisson(),
      data=scotlip,try.chol=FALSE)
spaMM.options(opt)
# }
# NOT RUN {
# }