spaMM (version 1.4.1)

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). 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 other methods (eigen or svd) are used. Note that these methods differ by much more than simply numerical accuracy, and in particular that this will affect the simulation of samples for the parametric bootstrap implemented in fixedLRT.

Usage

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

Arguments

m

The matrix which 'root' is to be computed.

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 factori

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

Singular value decomposition (SVD) of a matrix M using svd may encounter problems resulting in 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 is here 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 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

## 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)