spaMM: corrHLfit – R documentation

Become an expert in R — Interactive courses, Cheat Sheets, certificates and more!

corrHLfit

Fits a mixed model, typically a spatial GLMM.

Description

This was the first function for fitting all spatial models in spaMM, and is still fully functional, but it is recommended to use fitme which has different defaults and generally selects more efficient fitting methods, and will handle all classes of models that spaMM can fit, including non-spatial ones. corrHLfit performs the joint estimation of correlation parameters, fixed effect and dispersion parameters.

Usage

corrHLfit(formula, data, init.corrHLfit = list(), 
          init.HLfit = list(), ranFix = list(), lower = list(), 
          upper = list(), 
          objective = NULL, resid.model = ~1, 
          control.dist = list(), control.corrHLfit = list(),
          processed = NULL, family = gaussian(), method="REML",
          nb_cores = NULL, ...)

Arguments

`formula`	Either a linear model `formula` (as handled by various fitting functions) or a `predictor`, i.e. a formula with attributes (see `Predictor` and examples below). See Details in `spaMM` for allowed terms in the formula.
`data`	A data frame containing the variables in the response and the model formula.
`init.corrHLfit`	An optional list of initial values for correlation and/or dispersion parameters, e.g. `list(rho=1,nu=1,lambda=1,phi=1)` where `rho` and `nu` are parameters of the Matérn family (see `Matern`), and `lambda` and `phi` are dispersion parameters (see Details in `spaMM` for the meaning of these parameters). All are optional, but giving values for a dispersion parameter changes the ways it is estimated (see Details). `rho` may be a vector (see `make_scaled_dist`) and, in that case, it is possible that some or all of its elements are `NA`, for which `corrHLfit` substitutes automatically determined values.
`init.HLfit`	See identically named `HLfit` argument.
`ranFix`	A list similar to `init.corrHLfit`, but specifying fixed values of the parameters not estimated. See `ranFix` for further information.
`lower`	An optional (sub)list of values of the parameters specified through `init.corrHLfit`, in the same format as `init.corrHLfit`, used as lower values in calls to `optim`. See Details for default values.
`upper`	Same as `lower`, but for upper values.
`objective`	For development purpose, not documented (this had a distinct use in the first version of spaMM, but has been deprecated as such).
`resid.model`	See identically named `HLfit` argument.
`control.dist`	See `control.dist` in `HLCor`
`control.corrHLfit`	This may be used control the optimizer. See `spaMM.options`

for default values.

`processed`	For programming purposes, not documented.
`family`	Either a `family` or a `multi` value.
`method`	Character: the fitting method to be used, such as `"ML"`, `"REML"` or `"PQL/L"`. `"REML"` is the default. Other possible values of `HLfit`'s `method` argument are handled.
`nb_cores`	Not yet operative, only for development purposes. Number of cores to use for parallel computations.
`...`	Optional arguments passed to `HLCor`, `HLfit` or `mat_sqrt`, for example the `distMatrix` argument of `HLCor`, or the `verbose` argument of `HLfit`. Arguments that do not fit within these functions are detected and a warning is issued. In a `corrHLfit` call, the `verbose` vector of booleans may include a `TRACE=TRUE` element, in which case information is displayed for each set of correlation and dispersion parameter values considered by the optimiser. Non-boolean values of `TRACE` are meaningful, but the source code of `spaMM:::.do_TRACE` should be consulted for their meaning.

Details

For approximations of likelihood, see method. For the possible structures of random effects, see random-effects,

By default corrHLfit will estimate correlation parameters by maximizing the objective value returned by HLCor calls wherein the dispersion parameters are estimated jointly with fixed effects for given correlation parameters. If dispersion parameters are specified in init.corrHLfit, they will also be estimated by maximizing the objective value, and HLCor calls will not estimate them jointly with fixed effects. This means that in general the fixed effect estimates may vary depending on init.corrHLfit when any form of REML correction is applied.

Correctly using corrHLfit for likelihood ratio tests of fixed effects may then be tricky. It is safe to perform full ML fits of all parameters (using method="ML") for such tests (see Examples). The higher level function fixedLRT is a safe interface for likelihood ratio tests using some form of REML estimation in corrHLfit.

attr(<fitted object>,"optimInfo")$lower and ...$upper gives the lower and upper bounds for optimization of correlation parameters. These are the default values if the user did not provide explicit values. For the adjacency model, the default values are the inverse of the maximum and minimum eigenvalues of the adjMatrix. For the Matérn model, the default values are not so easily summarized: they are intended to cover the range of values for which there is statistical information to distinguish among them.

Value

The return value of an HLCor call, with additional attributes. The HLCor call is evaluated at the estimated correlation parameter values. These values are included in the return object as its $corrPars member. The attributes added by corrHLfit include the original call of the function (which can be retrived by getCall(<fitted object>), and information about the optimization call within corrHLfit.

Examples

# Example with an adjacency matrix (autoregressive model):
if (spaMM.getOption("example_maxtime")>0.7) {          
  corrHLfit(cases~I(prop.ag/10) +adjacency(1|gridcode)+offset(log(expec)),
          adjMatrix=Nmatrix,family=poisson(),data=scotlip,method="ML") 
}

#### Examples with Matern correlations
## A likelihood ratio test based on the ML fits of a full and of a null model.
if (spaMM.getOption("example_maxtime")>1.4) {
 data("blackcap")
 (fullfit <- corrHLfit(migStatus ~ means+ Matern(1|longitude+latitude),data=blackcap,
                    method="ML") )
 (nullfit <- corrHLfit(migStatus ~ 1 + Matern(1|longitude+latitude),data=blackcap,
                    method="ML",init.corrHLfit=list(phi=1e-6))) 
 ## p-value:
 1-pchisq(2*(logLik(fullfit)-logLik(nullfit)),df=1)
}

spaMM

Mixed-Effect Models, with or without Spatial Random Effects

v3.10.0

CeCILL-2

Authors

François Rousset [aut, cre, cph] (<https://orcid.org/0000-0003-4670-0371>), Jean-Baptiste Ferdy [aut, cph], Alexandre Courtiol [aut] (<https://orcid.org/0000-0003-0637-2959>), GSL authors [ctb] (src/gsl_bessel.*)

Initial release

2022-02-06