Return an adjacency matrix for a square lattice.

Description

Return an adjacency matrix for a square lattice.

Usage

adjacency.matrix(m, n = NULL)

Arguments

Details

This function builds the adjacency matrix for the m by n square lattice.

Value

A matrix A of 0s and 1s, where A_{ij} is equal to 1 iff vertices i and j are adjacent.

Fit copCAR model to discrete areal data.

Description

Fit the copCAR model to Bernoulli, negative binomial, or Poisson observations.

Usage

copCAR(formula, family, data, offset, A, method = c("CML", "DT", "CE"),
  confint = c("none", "bootstrap", "asymptotic"), model = TRUE, x = FALSE,
  y = TRUE, verbose = FALSE, control = list())

Arguments

Details

This function performs frequentist inference for the copCAR model proposed by Hughes (2015). copCAR is a copula-based areal regression model that employs the proper conditional autoregression (CAR) introduced by Besag, York, and Mollié (1991). Specifically, copCAR uses the CAR copula, a Caussian copula based on the proper CAR.

The spatial dependence parameter \rho \in [0, 1), regression coefficients \beta = (\beta_1, \dots, \beta_p)' \in R^p, and, for negative binomial margins, dispersion parameter \theta>0 can be estimated using the continous extension (CE) (Madsen, 2009), distributional transform (DT) (Kazianka and Pilz, 2010), or composite marginal likelihood (CML) (Varin, 2008) approaches.

The CE approach transforms the discrete observations to continous outcomes by convolving them with independent standard uniforms (Denuit and Lambert, 2005). The true likelihood for the discrete outcomes is the expected likelihood for the transformed outcomes. An estimate (sample mean) of the expected likelihood is optimized to estimate the copCAR parameters. The number of standard uniform vectors, m, can be chosen by the user. The default value is 1,000. The CE approach is exact up to Monte Carlo standard error but is computationally intensive (the computational burden grows rapidly with increasing m). The CE approach tends to perform poorly when applied to Bernoulli outcomes, and so that option is not permitted.

The distributional transform stochastically "smoothes" the jumps of a discrete distribution function (Ferguson, 1967). The DT-based approximation (Kazianka and Pilz, 2010) for copCAR performs well for Poisson and negative binomial marginals but, like the CE approach, tends to perform poorly for Bernoulli outcomes.

The CML approach optimizes a composite marginal likelihood formed as the product of pairwise likelihoods of adjacent observations. This approach performs well for Bernoulli, negative binomial, and Poisson outcomes.

In the CE and DT approaches, the CAR variances are approximated. The quality of the approximation is determined by the values of control parameters \epsilon > 0 and \rho^{\max} \in [0, 1). The default values are 0.01 and 0.999, respectively.

When confint = "bootstrap", a parametric bootstrap is carried out, and confidence intervals are computed using the quantile method. Monte Carlo standard errors (Flegal et al., 2008) of the quantile estimators are also provided.

When confint = "asymptotic", confidence intervals are computed using an estimate of the asymptotic covariance matrix of the estimator. For the CE method, the inverse of the observed Fisher information matrix is used. For the CML and DT methods, the objective function is misspecified, and so the asymptotic covariance matrix is the inverse of the Godambe information matrix (Godambe, 1960), which has a sandwich form. The "bread" is the inverse of the Fisher information matrix, and the "meat" is the covariance matrix of the score function. The former is estimated using the inverse of the observed Fisher information matrix. The latter is estimated using a parametric bootstrap.

Value

copCAR returns an object of class "copCAR", which is a list containing the following components:

References

Besag, J., York, J., and Mollié, A. (1991) Bayesian image restoration, with two applications in spatial statistics. Annals of the Institute of Statistical Mathematics, 43(1), 1–20.

Denuit, M. and Lambert, P. (2005) Constraints on concordance measures in bivariate discrete data. Journal of Multivariate Analysis, 93, 40–57.

Ferguson, T. (1967) Mathematical statistics: a decision theoretic approach, New York: Academic Press.

Flegal, J., Haran, M., and Jones, G. (2008) Markov Chain Monte Carlo: can we trust the third significant figure? Statistical Science, 23(2), 250–260.

Godambe, V. (1960) An optimum property of regular maximum likelihood estimation. The Annals of Mathmatical Statistics, 31(4), 1208–1211.

Hughes, J. (2015) copCAR: A flexible regression model for areal data. Journal of Computational and Graphical Statistics, 24(3), 733–755.

Kazianka, H. and Pilz, J. (2010) Copula-based geostatistical modeling of continuous and discrete data including covariates. Stochastic Environmental Research and Risk Assessment, 24(5), 661–673.

Madsen, L. (2009) Maximum likelihood estimation of regression parameters with spatially dependent discrete data. Journal of Agricultural, Biological, and Environmental Statistics, 14(4), 375–391.

Varin, C. (2008) On composite marginal likelihoods. Advances in Statistical Analysis, 92(1), 1–28.

Examples

## Not run: 
# Simulate data and fit copCAR to them.

# Use the 20 x 20 square lattice as the underlying graph.

m = 20
A = adjacency.matrix(m)

# Create a design matrix by assigning coordinates to each vertex
# such that the coordinates are restricted to the unit square.

x = rep(0:(m - 1) / (m - 1), times = m)
y = rep(0:(m - 1) / (m - 1), each = m)
X = cbind(x, y)

# Set the dependence parameter, regression coefficients, and dispersion parameter.

rho = 0.995      # strong dependence
beta = c(1, 1)   # the mean surface increases in the direction of (1, 1)
theta = 2        # dispersion parameter

# Simulate negative binomial data from the model.

z = rcopCAR(rho, beta, X, A, family = negbinomial(theta))

# Fit the copCAR model using the continous extension, and compute 95% (default)
# asymptotic confidence intervals. Give theta the initial value of 1. Use m equal to 100.

fit.ce = copCAR(z ~ X - 1, A = A, family = negbinomial(1), method = "CE", confint = "asymptotic",
                control = list(m = 100))
summary(fit.ce)

# Fit the copCAR model using the DT approximation, and compute 90% confidence
# intervals. Bootstrap the intervals, based on a bootstrap sample of size 100.
# Do the bootstrap in parallel, using ten nodes.

fit.dt = copCAR(z ~ X - 1, A = A, family = negbinomial(1), method = "DT", confint = "bootstrap",
                control = list(bootit = 100, nodes = 10))
summary(fit.dt, alpha = 0.9)

# Fit the copCAR model using the composite marginal likelihood approach.
# Do not compute confidence intervals.

fit.cml = copCAR(z ~ X - 1, A = A, family = negbinomial(1), method = "CML", confint = "none")
summary(fit.cml)

## End(Not run)

Family function for negative binomial GLMs.

Description

Provides the information required to apply copCAR with negative binomial marginal distributions.

Usage

negbinomial(theta = stop("'theta' must be specified."), link = "log")

Arguments

Value

An object of class “family”, a list of functions and expressions needed to fit a negative binomial GLM.

Simulate areal data.

Description

rcopCAR simulates areal data from the copCAR model.

Usage

rcopCAR(rho, beta, X, A, family)

Arguments

Details

This function simulates data from the copCAR model with the given spatial dependence parameter \rho, regression coefficients \beta, design matrix X, and adjacency structure A. For negative binomial marginal distributions, a value for the dispersion parameter \theta>0 is also required; this value must be passed to the negbinomial family function. For more details on the copCAR model, see copCAR.

Value

A vector of length n distributed according to the specified copCAR model.

Examples

# Use the 20 x 20 square lattice as the underlying graph.

m = 20
A = adjacency.matrix(m)

# Create a design matrix by assigning coordinates to each vertex
# such that the coordinates are restricted to the unit square.

x = rep(0:(m - 1) / (m - 1), times = m)
y = rep(0:(m - 1) / (m - 1), each = m)
X = cbind(x, y)

# Set the dependence parameter and regression coefficients.

rho = 0.995      # strong dependence
beta = c(1, 1)   # the mean surface increases in the direction of (1, 1)

# Simulate Poisson data from the corresponding copCAR model.

z = rcopCAR(rho, beta, X, A, family = poisson(link = "log"))

# Simulate Bernoulli outcomes.

z = rcopCAR(rho, beta, X, A, family = binomial(link = "logit"))

# Set the dispersion parameter.

theta = 10

# Simulate negative binomial outcomes.

z = rcopCAR(rho, beta, X, A, family = negbinomial(theta))

Extract model residuals.

Description

Extract model residuals.

Usage

## S3 method for class 'copCAR'
residuals(object, type = c("deviance", "pearson",
  "response"), ...)

Arguments

Value

A vector of residuals.

See Also

copCAR, residuals.glm

Print a summary of a copCAR model fit.

Description

Print a summary of a copCAR model fit.

Usage

## S3 method for class 'copCAR'
summary(object, alpha = 0.05, digits = 4, ...)

Arguments

Details

This function displays (1) the call to copCAR, (2) the values of the control parameters, (3) a table of estimates, and (when applicable) (4) confidence intervals and (5) Monte Carlo standard errors.

Each row of the table of estimates shows a parameter estimate and (when applicable) the confidence interval for the parameter. If copCAR was called with confint = "bootstrap", Monte Carlo standard errors are provided.

References

Flegal, J., Haran, M., and Jones, G. (2008) Markov Chain Monte Carlo: can we trust the third significant figure? Statistical Science, 23(2), 250–260.

See Also

copCAR

Return the estimated covariance matrix for a copCAR model object.

Description

Return the estimated covariance matrix for a copCAR model object.

Usage

## S3 method for class 'copCAR'
vcov(object, ...)

Arguments

Details

Unless copCAR was called with confint = "none", this function returns an estimate of the covariance matrix of the CE/CML/DT estimator of the parameters. If confint = "bootstrap", cov is applied to the bootstrap sample to compute the estimate. If confint = "asymptotic", an estimate of the asymptotic covariance matrix is returned; this is an estimate of the inverse Fisher information matrix if method = "CE", or an estimate of the inverse of the Godambe information matrix if method = "CML" or method = "DT". Note that the entries involving the spatial dependence parameter are for \gamma=\Phi^{-1}(\rho) rather than for \rho (Hughes, 2015).

Value

An estimate of the covariance matrix of the CE/CML/DT estimator of the parameters.

References

Hughes, J. (2015) copCAR: A flexible regression model for areal data. Journal of Computational and Graphical Statistics, 24(3), 733–755.