spaMM (version 1.9.16)

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

Description

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). t(chol(.)) (Cholesky factorization) 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 (eigen or svd) 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, symSVD=NULL, try.chol=TRUE, try.eigen=FALSE, threshold=1e-06, debug=FALSE, SVDfix=1/10)

Arguments

m
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 eigen factorization will be tried before svd. eigen is a compromise between speed and accuracy, but in our experience it may *hang* so by default it is not tried.
threshold
A correction threshold for low eigenvalues is the case and eigensystem or singular-value decomposition are used.
debug
Not documented, only for development purposes.
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érn correlation model with low $\nu$, high $\rho$, and if some samples are spatially close. Then, an alternative fix to the svd problem may be to restrict the $\nu$ and/or $\rho$ 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)
# HLCor(cases~I(prop.ag/10) +adjacency(1|gridcode)+offset(log(scotlip$expec)),
#       ranPars=list(rho=0.174),adjMatrix=Nmatrix,family=poisson(),
#       data=scotlip,try.chol=FALSE)
# ## End(Not run)

Run the code above in your browser using DataCamp Workspace