vcov.ppm: Variance-Covariance Matrix for a Fitted Point Process Model

In spatstat.model: Parametric Statistical Modelling and Inference for the 'spatstat' Family

Variance-Covariance Matrix for a Fitted Point Process Model

Description

Returns the variance-covariance matrix of the estimates of the parameters of a fitted point process model.

Usage

  ## S3 method for class 'ppm'
vcov(object, ...,
                    what = c("vcov", "corr", "fisher"),
                    verbose = TRUE,
                    fine=FALSE,
                    gam.action=c("warn", "fatal", "silent"),
                    matrix.action=c("warn", "fatal", "silent"),
                    logi.action=c("warn", "fatal", "silent"),
                    nacoef.action=c("warn", "fatal", "silent"),
                    hessian=FALSE)

Arguments

object	A fitted point process model (an object of class "ppm".)
...	Ignored.
what	Character string (partially-matched) that specifies what matrix is returned. Options are "vcov" for the variance-covariance matrix, "corr" for the correlation matrix, and "fisher" or "Fisher" for the Fisher information matrix.
fine	Logical value indicating whether to use a quick estimate (fine=FALSE, the default) or a slower, more accurate estimate (fine=TRUE).
verbose	Logical. If TRUE, a message will be printed if various minor problems are encountered.
gam.action	String indicating what to do if object was fitted by gam.
matrix.action	String indicating what to do if the matrix is ill-conditioned (so that its inverse cannot be calculated).
logi.action	String indicating what to do if object was fitted via the logistic regression approximation using a non-standard dummy point process.
nacoef.action	String indicating what to do if some of the fitted coefficients are NA (so that variance cannot be calculated).
hessian	Logical. Use the negative Hessian matrix of the log pseudolikelihood instead of the Fisher information.

Details

This function computes the asymptotic variance-covariance matrix of the estimates of the canonical parameters in the point process model object. It is a method for the generic function vcov.

object should be an object of class "ppm", typically produced by ppm.

The canonical parameters of the fitted model object are the quantities returned by coef.ppm(object). The function vcov calculates the variance-covariance matrix for these parameters.

The argument what provides three options:

what="vcov"

return the variance-covariance matrix of the parameter estimates

what="corr"

return the correlation matrix of the parameter estimates

what="fisher"

return the observed Fisher information matrix.

In all three cases, the result is a square matrix. The rows and columns of the matrix correspond to the canonical parameters given by coef.ppm(object). The row and column names of the matrix are also identical to the names in coef.ppm(object).

For models fitted by the Berman-Turner approximation (Berman and Turner, 1992; Baddeley and Turner, 2000) to the maximum pseudolikelihood (using the default method="mpl" in the call to ppm), the implementation works as follows.

• If the fitted model object is a Poisson process, the calculations are based on standard asymptotic theory for the maximum likelihood estimator (Kutoyants, 1998). The observed Fisher information matrix of the fitted model object is first computed, by summing over the Berman-Turner quadrature points in the fitted model. The asymptotic variance-covariance matrix is calculated as the inverse of the observed Fisher information. The correlation matrix is then obtained by normalising.

• If the fitted model is not a Poisson process (i.e. it is some other Gibbs point process) then the calculations are based on Coeurjolly and Rubak (2012). A consistent estimator of the variance-covariance matrix is computed by summing terms over all pairs of data points. If required, the Fisher information is calculated as the inverse of the variance-covariance matrix.

For models fitted by the Huang-Ogata method (method="ho" in the call to ppm), the implementation uses the Monte Carlo estimate of the Fisher information matrix that was computed when the original model was fitted.

For models fitted by the logistic regression approximation to the maximum pseudolikelihood (method="logi" in the call to ppm),

• Calculations are based on Baddeley et al. (2013). A consistent estimator of the variance-covariance matrix is computed by summing terms over all pairs of data points. If required, the Fisher information is calculated as the inverse of the variance-covariance matrix.

• The calculations depend on the type of dummy pattern used when the model was fitted:

  • currently only the dummy types "stratrand" (the default), "binomial" and "poisson" as generated by quadscheme.logi are supported.

  • For other dummy types the behavior depends on the argument logi.action. If logi.action="fatal" an error is produced. Otherwise, for dummy types "grid" and "transgrid" the formulas for "stratrand" are used which in many cases should be conservative. For an arbitrary, user-specified dummy pattern (type "given"), the formulas for "poisson" are used which in many cases should be conservative. If logi.action="warn" a warning is issued, otherwise the calculation proceeds without a warning.

• The result of the calculation is random (i.e. not deterministic) when dummy type is "stratrand" (the default) because some of the variance terms are estimated by random sampling. This can be avoided by specifying dummytype='poisson' or dummytype='binomial' in the call to ppm when the model is fitted.

The argument verbose makes it possible to suppress some diagnostic messages.

The asymptotic theory is not correct if the model was fitted using gam (by calling ppm with use.gam=TRUE). The argument gam.action determines what to do in this case. If gam.action="fatal", an error is generated. If gam.action="warn", a warning is issued and the calculation proceeds using the incorrect theory for the parametric case, which is probably a reasonable approximation in many applications. If gam.action="silent", the calculation proceeds without a warning.

If hessian=TRUE then the negative Hessian (second derivative) matrix of the log pseudolikelihood, and its inverse, will be computed. For non-Poisson models, this is not a valid estimate of variance, but is useful for other calculations.

Note that standard errors and 95% confidence intervals for the coefficients can also be obtained using confint(object) or coef(summary(object)).

Value

A square matrix.

Error messages

An error message that reports system is computationally singular indicates that the determinant of the Fisher information matrix was either too large or too small for reliable numerical calculation.

If this message occurs, try repeating the calculation using fine=TRUE.

Singularity can occur because of numerical overflow or collinearity in the covariates. To check this, rescale the coordinates of the data points and refit the model. See the Examples.

In a Gibbs model, a singular matrix may also occur if the fitted model is a hard core process: this is a feature of the variance estimator.

Author(s)

Original code for Poisson point process was written by \adrian and \rolf. New code for stationary Gibbs point processes was generously contributed by \ege and Jean-\Francois Coeurjolly. New code for generic Gibbs process written by \adrian. New code for logistic method written by \ege.

References

Baddeley, A., Coeurjolly, J.-F., Rubak, E. and Waagepetersen, R. (2014) Logistic regression for spatial Gibbs point processes. Biometrika 101 (2) 377–392.

Coeurjolly, J.-F. and Rubak, E. (2013) Fast covariance estimation for innovations computed from a spatial Gibbs point process. Scandinavian Journal of Statistics 40 669–684.

Kutoyants, Y.A. (1998) Statistical Inference for Spatial Poisson Processes, Lecture Notes in Statistics 134. New York: Springer 1998.

See Also

vcov for the generic,

ppm for information about fitted models,

confint for confidence intervals.

Examples

  X <- rpoispp(42)
  fit <- ppm(X ~ x + y)
  vcov(fit)
  vcov(fit, what="Fish")

  # example of singular system
  m <- ppm(demopat ~polynom(x,y,2))
  
    try(v <- vcov(m))
  
  # rescale x, y coordinates to range [0,1] x [0,1] approximately
  demopatScale <- rescale(demopat, 10000)
  m <- ppm(demopatScale ~ polynom(x,y,2))
  v <- vcov(m)

  # Gibbs example
  fitS <- ppm(swedishpines ~1, Strauss(9))
  coef(fitS)
  sqrt(diag(vcov(fitS)))

spatstat.model documentation built on Sept. 30, 2024, 9:26 a.m.