VGAM (version 1.0-5)

calibrate.qrrvglm: Calibration for CQO and CAO models

Description

Performs maximum likelihood calibration for constrained and unconstrained quadratic and additive ordination models (CQO and CAO models are better known as QRR-VGLMs and RR-VGAMs respectively).

Usage

calibrate.qrrvglm(object, newdata = NULL,
        type = c("latvar", "predictors", "response", "vcov", "all3or4"),
        se.type = c("asbefore", "wald"),
        initial.vals = NULL, ...)

Arguments

object

The fitted CQO/CAO model.

newdata

A data frame with new response data (usually new species data). The default is to use the original data used to fit the model; however, the calibration may take a long time to compute because the computations are expensive.

type

What type of result is to be returned. The first are the calibrated latent variables or site scores. This must be computed always. The "predictors" are the linear/quadratic or additive predictors evaluated at the calibrated latent variables or site scores. The "response" are the fitted means evaluated at the calibrated latent variables or site scores. The "vcov" are the estimated variance-covariance matrices of the calibrated latent variables or site scores. The "all3or4" is for all of them, i.e., all types. For CAO models, "vcov" is unavailable, so all 3 are returned. For CQO models, "vcov" is available, so all 4 are returned.

se.type

What type of standard errors are to be returned. The choice "asbefore" comes from a long time ago. The choice "wald" is based on work by David Zucker.

initial.vals

Initial values for the search. For rank-1 models, this should be a vector of length nrow(newdata), and for rank 2 models this should be a two column matrix with the number of rows equalling the number of rows in newdata. The default is a grid defined by arguments in calibrate.qrrvglm.control.

…

Arguments that are fed into calibrate.qrrvglm.control.

Value

The argument type determines what is returned. If type = "all3or4" then all the type values are returned in a list, with the following components. Each component has length nrow(newdata).

latvar

Calibrated latent variables or site scores.

predictors

linear/quadratic or additive predictors. For example, for Poisson families, this will be on a log scale, and for binomial families, this will be on a logit scale.

response

Fitted values of the response, evaluated at the calibrated latent variables or site scores.

vcov

Estimated variance-covariance matrix of the calibrated latent variables or site scores. Actually, these are stored in an array whose last dimension is nrow(newdata).

Warning

This function is computationally expensive. Setting trace = TRUE to get a running log is a good idea.

Details

Given a fitted regression CQO/CAO model, maximum likelihood calibration is theoretically easy and elegant. However, the method assumes that all species are independent, which is not really true in practice. More details and references are given in Yee (2018) and ch.6 of Yee (2015).

The function optim is used to search for the maximum likelihood solution. Good initial values are needed, and calibrate.qrrvglm.control allows the user some control over the choice of these.

References

Yee, T. W. (2018) On constrained and unconstrained quadratic ordination. Manuscript in preparation.

ter Braak, C. J. F. 1995. Calibration. In: Data Analysis in Community and Landscape Ecology by Jongman, R. H. G., ter Braak, C. J. F. and van Tongeren, O. F. R. (Eds.) Cambridge University Press, Cambridge.

See Also

calibrate.qrrvglm.control, calibrate.rrvglm, calibrate, cqo, cao.

Examples

# NOT RUN {
hspider[, 1:6] <- scale(hspider[, 1:6])  # Stdze environmental variables
set.seed(123)
siteNos <- c(1, 5)  # Calibrate these sites
pet1 <- cqo(cbind(Pardlugu, Pardmont, Pardnigr, Pardpull, Zoraspin) ~
          WaterCon + BareSand + FallTwig + CoveMoss + CoveHerb + ReflLux,
          data = hspider[-siteNos, ],  # Sites not in fitted model
          family = poissonff, I.toler = TRUE, Crow1positive = TRUE)
y0 <- hspider[siteNos, colnames(depvar(pet1))]  # Species counts
cpet1 <- calibrate(pet1, trace = TRUE, newdata = data.frame(y0))
# }
# NOT RUN {
# }
# NOT RUN {
# Graphically compare the actual site scores with their calibrated values
persp(pet1, main = "Site scores: solid=actual, dashed=calibrated",
      label = TRUE, col = "blue", las = 1)
# Actual site scores:
xvars <- rownames(concoef(pet1))  # Variables comprising the latvar
est.latvar <- as.matrix(hspider[siteNos, xvars]) <!-- %*% concoef(pet1) -->
abline(v = est.latvar, lty = 1, col = 1:length(siteNos))
abline(v = cpet1, lty = 2, col = 1:length(siteNos))  # Calibrated values
# }