Title: | Objective Bayes Intrinsic Conditional Autoregressive Model for Areal Data |
---|---|
Description: | Implements an objective Bayes intrinsic conditional autoregressive prior. This model provides an objective Bayesian approach for modeling spatially correlated areal data using an intrinsic conditional autoregressive prior on a vector of spatial random effects. |
Authors: | Erica M. Porter [aut, cre], Matthew J. Keefe [aut], Christopher T. Franck [aut], Marco A.R. Ferreira [aut] |
Maintainer: | Erica M. Porter <[email protected]> |
License: | MIT + file LICENSE |
Version: | 2.0.2 |
Built: | 2025-02-21 07:06:17 UTC |
Source: | CRAN |
Performs simultaneous selection of covariates and spatial model structure for areal data.
probs.icar( Y, X, H, H.spectral = NULL, Sig_phi = NULL, b = 0.05, verbose = FALSE )
probs.icar( Y, X, H, H.spectral = NULL, Sig_phi = NULL, b = 0.05, verbose = FALSE )
Y |
A vector of responses. |
X |
A matrix of covariates, which should include a column of 1's for models with a non-zero intercept |
H |
Neighborhood matrix for spatial subregions. |
H.spectral |
Spectral decomposition of neighborhood matrix, if user wants to pre-compute it to save time. |
Sig_phi |
Pseudo inverse of the neighborhood matrix, if user wants to pre-compute it to save time. |
b |
Training fraction for the fractional Bayes factor (FBF) approach. |
verbose |
If FALSE, marginal likelihood progress is not printed. |
A list containing a data frame with all posterior model probabilities and other selection information.
probs.mat |
Data frame containing posterior model probabilities for all candidate OLMs and ICAR models from the data. |
mod.prior |
Vector of model priors used to obtain the posterior model probabilities. |
logmargin.all |
Vector of all (log) fractional integrated likelihoods. |
base.model |
Maximum (log) fractional integrated likelihood among all candidate models. All fractional Bayes factors are obtained with respect to this model. |
BF.vec |
Vector of fractional Bayes factors for all candidate models. |
Erica M. Porter, Christopher T. Franck, and Marco A.R. Ferreira
Porter EM, Franck CT, Ferreira MAR (2024). “Objective Bayesian model selection for spatial hierarchical models with intrinsic conditional autoregressive priors.” Bayesian Analysis, 19(4), 985–1011. doi:10.1214/23-BA1375.
Performs analysis on a geographical areal data set using the objective prior for intrinsic conditional autoregressive (ICAR) random effects (Keefe et al. 2019). It takes a shapefile, data, and region names to construct a neighborhood matrix and perform Markov chain Monte Carlo sampling on the unstructured and spatial random effects. Finally, the function obtains regional estimates and performs posterior inference on the model parameters.
ref.analysis( shape.file, X, y, x.reg.names, y.reg.names, shp.reg.names = NULL, iters = 10000, burnin = 5000, verbose = TRUE, tauc.start = 1, beta.start = 1, sigma2.start = 1, step.tauc = 0.5, step.sigma2 = 0.5 )
ref.analysis( shape.file, X, y, x.reg.names, y.reg.names, shp.reg.names = NULL, iters = 10000, burnin = 5000, verbose = TRUE, tauc.start = 1, beta.start = 1, sigma2.start = 1, step.tauc = 0.5, step.sigma2 = 0.5 )
shape.file |
A shapefile corresponding to the regions for analysis. |
X |
A matrix of covariates, which should include a column of 1's for models with a non-zero intercept |
y |
A vector of responses. |
x.reg.names |
A vector specifying the order of region names contained in |
y.reg.names |
A vector specifying the order of region names contained in |
shp.reg.names |
A vector specifying the order of region names contained in the shapefile, if there is not a NAME column in the file. |
iters |
Number of MCMC iterations to perform. Defaults to 10,000. |
burnin |
Number of MCMC iterations to discard as burn-in. Defaults to 5,000. |
verbose |
If FALSE, MCMC progress is not printed. |
tauc.start |
Starting MCMC value for the spatial dependence parameter. |
beta.start |
Starting MCMC value for the fixed effect regression coefficients. |
sigma2.start |
Starting MCMC value for the variance of the unstructured random effects. |
step.tauc |
Step size for the spatial dependence parameter. |
step.sigma2 |
Step size for the variance of the unstructured random effects. |
A list containing H
, MCMC chains, parameter summaries, fitted regional values,
and regional summaries.
H |
The neighborhood matrix. |
MCMC |
Matrix of MCMC chains for all model parameters. |
beta.median |
Posterior medians of the fixed effect regression coefficients. |
beta.hpd |
Highest Posterior Density intervals for the fixed effect regression coefficients. |
tauc.median |
Posterior median of the spatial dependence parameter. |
tauc.hpd |
Highest Posterior Density interval for the spatial dependence parameter. |
sigma2.median |
Posterior median of the unstructured random effects variance. |
sigma2.hpd |
Highest Posterior Density interval for the unstructured random effects variance. |
tauc.accept |
Final acceptance rate for the spatial dependence parameter. |
sigma2.accept |
Final acceptance rate for the unstructured random effects variance. |
fit.dist |
Matrix of fitted posterior values for each region in the data. |
reg.medians |
Vector of posterior medians for fitted response by region. |
reg.hpd |
Data frame of Highest Posterior Density intervals by region. |
Erica M. Porter, Matthew J. Keefe, Christopher T. Franck, and Marco A.R. Ferreira
## Refer to the vignette attached to the package.
## Refer to the vignette attached to the package.
Implements the Metropolis-within-Gibbs sampling algorithm proposed by Ferreira et al. (2021), to perform posterior inference for the intrinsic conditional autoregressive model with spatial random effects. This algorithm uses the spectral domain for the hierarchical model to create the Spectral Gibbs Sampler (SGS), which provides notable speedups to the MCMC algorithm proposed by Keefe et al (2019).
ref.MCMC( y, X, H, iters = 10000, burnin = 5000, verbose = TRUE, tauc.start = 1, beta.start = 1, sigma2.start = 1, step.tauc = 0.5, step.sigma2 = 0.5 )
ref.MCMC( y, X, H, iters = 10000, burnin = 5000, verbose = TRUE, tauc.start = 1, beta.start = 1, sigma2.start = 1, step.tauc = 0.5, step.sigma2 = 0.5 )
y |
Vector of responses. |
X |
Matrix of covariates. This should include a column of 1's for models with a non-zero intercept. |
H |
The neighborhood matrix. |
iters |
Number of MCMC iterations to perform. Defaults to 10,000. |
burnin |
Number of MCMC iterations to discard as burn-in. Defaults to 5,000. |
verbose |
If FALSE, MCMC progress is not printed. |
tauc.start |
Starting value for the spatial dependence parameter. |
beta.start |
Starting value for the vector of fixed effect regression coefficients. |
sigma2.start |
Starting value for the variance of the unstructured random effects. |
step.tauc |
Step size for the spatial dependence parameter |
step.sigma2 |
Step size for the variance of the unstructured random effects. |
A list containing MCMC chains and parameter summaries.
MCMCchain |
Matrix of MCMC chains. |
tauc.MCMC |
MCMC chains for the spatial dependence parameter. |
sigma2.MCMC |
MCMC chains for the variance of the unstructured random effects. |
phi.MCMC |
MCMC chains for the spatial random effects. |
beta.MCMC |
MCMC chains for the fixed effect regression coefficients. |
accept.sigma2 |
Final acceptance number for variance of the unstructured random effects. |
accept.tauc |
Final acceptance number for spatial dependence parameter. |
accept.phi |
Final acceptance number for spatial random effects. |
Erica M. Porter, Matthew J. Keefe, Christopher T. Franck, and Marco A.R. Ferreira
Keefe MJ, Ferreira MAR, Franck CT (2019). “Objective Bayesian analysis for Gaussian hierarchical models with intrinsic conditional autoregressive priors.” Bayesian Analysis, 14(1), 181 – 209. doi:10.1214/18-BA1107.
Keefe MJ, Ferreira MAR, Franck CT (2018). “On the formal specification of sum-zero constrained intrinsic conditional autoregressive models.” Spatial Statistics, 24, 54–65. doi:10.1016/j.spasta.2018.03.007.
Ferreira MAR, Porter EM, Franck CT (2021). “Fast and scalable computations for Gaussian hierarchical models with intrinsic conditional autoregressive spatial random effects.” Computational Statistics and Data Analysis, 162, 107264. ISSN 0167-9473, doi:10.1016/j.csda.2021.107264.
#### Fit the model for simulated areal data on a grid #### ### Load extra libraries library(sp) library(methods) library(spdep) library(mvtnorm) ### Generate areal data on a grid rows=5; cols=5 tauc=1 sigma2=2; beta=c(1,5) ### Create grid grid <- GridTopology(c(1,1), c(1,1), c(cols,rows)) polys <- as(grid, "SpatialPolygons") spgrid <- SpatialPolygonsDataFrame(polys,data=data.frame(row.names=row.names(polys))) ### Create neighborhood matrix grid.nb <- poly2nb(spgrid,queen=FALSE) W <- nb2mat(grid.nb, style="B") ### Put spatially correlated data in grid p <- length(beta) num.reg <- (rows*cols) if(p>1){x1<-rmvnorm(n=num.reg,mean=rep(0,p-1),sigma=diag(p-1))} else{x1<-NULL} X <- cbind(rep(1,num.reg),x1) Dmat <- diag(apply(W,1,sum)) H <- Dmat - W row.names(H) <- NULL ### Obtain true response vector theta_true <- rnorm(num.reg,mean=0,sd=sqrt(sigma2)) Q <- eigen(H,symmetric=TRUE)$vectors eigH <- eigen(H,symmetric=TRUE)$values D <- diag(eigH) Qmat <- Q[,1:(num.reg-1)] phimat <- diag(1/sqrt(eigH[1:(num.reg-1)])) z <- t(rmvnorm(1,mean=rep(0,num.reg-1),sigma=diag(num.reg-1))) phi_true <- sqrt((1/tauc)*sigma2)*(Qmat%*%phimat%*%z) Y <- X%*%beta + theta_true + phi_true ### Fit the model set.seed(5432) model <- ref.MCMC(y=Y,X=X,H=H,iters=15000,burnin=5000,verbose=TRUE,tauc.start=.1,beta.start=-1, sigma2.start=.1,step.tauc=0.5, step.sigma2=0.5) #### Small example for checking model <- ref.MCMC(y=Y,X=X,H=H,iters=1000,burnin=50,verbose=TRUE,tauc.start=.1,beta.start=-1, sigma2.start=.1,step.tauc=0.5, step.sigma2=0.5)
#### Fit the model for simulated areal data on a grid #### ### Load extra libraries library(sp) library(methods) library(spdep) library(mvtnorm) ### Generate areal data on a grid rows=5; cols=5 tauc=1 sigma2=2; beta=c(1,5) ### Create grid grid <- GridTopology(c(1,1), c(1,1), c(cols,rows)) polys <- as(grid, "SpatialPolygons") spgrid <- SpatialPolygonsDataFrame(polys,data=data.frame(row.names=row.names(polys))) ### Create neighborhood matrix grid.nb <- poly2nb(spgrid,queen=FALSE) W <- nb2mat(grid.nb, style="B") ### Put spatially correlated data in grid p <- length(beta) num.reg <- (rows*cols) if(p>1){x1<-rmvnorm(n=num.reg,mean=rep(0,p-1),sigma=diag(p-1))} else{x1<-NULL} X <- cbind(rep(1,num.reg),x1) Dmat <- diag(apply(W,1,sum)) H <- Dmat - W row.names(H) <- NULL ### Obtain true response vector theta_true <- rnorm(num.reg,mean=0,sd=sqrt(sigma2)) Q <- eigen(H,symmetric=TRUE)$vectors eigH <- eigen(H,symmetric=TRUE)$values D <- diag(eigH) Qmat <- Q[,1:(num.reg-1)] phimat <- diag(1/sqrt(eigH[1:(num.reg-1)])) z <- t(rmvnorm(1,mean=rep(0,num.reg-1),sigma=diag(num.reg-1))) phi_true <- sqrt((1/tauc)*sigma2)*(Qmat%*%phimat%*%z) Y <- X%*%beta + theta_true + phi_true ### Fit the model set.seed(5432) model <- ref.MCMC(y=Y,X=X,H=H,iters=15000,burnin=5000,verbose=TRUE,tauc.start=.1,beta.start=-1, sigma2.start=.1,step.tauc=0.5, step.sigma2=0.5) #### Small example for checking model <- ref.MCMC(y=Y,X=X,H=H,iters=1000,burnin=50,verbose=TRUE,tauc.start=.1,beta.start=-1, sigma2.start=.1,step.tauc=0.5, step.sigma2=0.5)
This function creates trace plots for the parameters in the ICAR reference prior model (Keefe et al. 2019).
ref.plot(MCMCchain, X, burnin, num.reg)
ref.plot(MCMCchain, X, burnin, num.reg)
MCMCchain |
Matrix of MCMC chains for the model parameters. |
X |
Matrix of covariates. |
burnin |
Number of MCMC iterations from |
num.reg |
Number of regions in the areal data set. |
Trace plots for the fixed effect regression coefficients, the precision parameter, and the unstructured random effects variance.
Erica M. Porter, Matthew J. Keefe, Christopher T. Franck, and Marco A.R. Ferreira
## Refer to the vignette attached to the package.
## Refer to the vignette attached to the package.
Takes a matrix of MCMC chains, iterations, and acceptance values to return posterior summaries of the parameters, including posterior medians, intervals, and acceptance rates.
ref.summary( MCMCchain, tauc.MCMC, sigma2.MCMC, beta.MCMC, phi.MCMC, accept.phi, accept.sigma2, accept.tauc, iters = 10000, burnin = 5000 )
ref.summary( MCMCchain, tauc.MCMC, sigma2.MCMC, beta.MCMC, phi.MCMC, accept.phi, accept.sigma2, accept.tauc, iters = 10000, burnin = 5000 )
MCMCchain |
Matrix of MCMC chains for the ICAR model parameters. |
tauc.MCMC |
MCMC chains for the spatial dependence parameter. |
sigma2.MCMC |
MCMC chains for the variance of the unstructured random effects. |
beta.MCMC |
MCMC chains for the fixed effect regression coefficients. |
phi.MCMC |
MCMC chains for the spatial random effects. |
accept.phi |
Final acceptance number for spatial random effects. |
accept.sigma2 |
Final acceptance number for variance of the unstructured random effects. |
accept.tauc |
Final acceptance number for the spatial dependence parameter. |
iters |
Number of MCMC iterations in |
burnin |
Number of MCMC iterations discarded as burn-in for |
Parameter summaries
beta.median |
Posterior medians of the fixed effect regression coefficients. |
beta.hpd |
Highest Posterior Density intervals for the fixed effect regression coefficients. |
tauc.median |
Posterior median of the spatial dependence parameter. |
tauc.hpd |
Highest Posterior Density interval for the spatial dependence parameter. |
sigma2.median |
Posterior median of the unstructured random effects variance. |
sigma2.hpd |
Highest Posterior Density interval for the unstructured random effects variance. |
tauc.accept |
Final acceptance rate for the spatial dependence parameter. |
sigma2.accept |
Final acceptance rate for the unstructured random effects variance. |
Erica M. Porter, Matthew J. Keefe, Christopher T. Franck, and Marco A.R. Ferreira
## Refer to the vignette attached to the package.
## Refer to the vignette attached to the package.
This function takes data and sampled MCMC chains for an areal data set and gives fitted posterior values and summaries by region using the model by (Keefe et al. 2019).
reg.summary(MCMCchain, X, Y, burnin)
reg.summary(MCMCchain, X, Y, burnin)
MCMCchain |
Matrix of MCMC chains, using the sampling from (Keefe et al. 2019). |
X |
Matrix of covariates. |
Y |
Vector of responses. |
burnin |
Number of MCMC iterations discarded as burn-in in |
A list of the fitted distributions by region, and medians and credible intervals by region.
fit.dist |
Matrix of fitted posterior values for each region in the data. |
reg.medians |
Vector of posterior medians for fitted response by region. |
reg.cred |
Data frame of credbile intervals by region. |
Erica M. Porter, Matthew J. Keefe, Christopher T. Franck, and Marco A.R. Ferreira
## Refer to the vignette attached to the package.
## Refer to the vignette attached to the package.
Takes a path to a shape file and creates a neighborhood matrix. This neighborhood matrix can be used with the objective ICAR model (Keefe et al. 2018).
shape.H(shape.file)
shape.H(shape.file)
shape.file |
File path to a shapefile. |
A list containing a neighborhood matrix and the SpatialPolygonsDataFrame object corresponding to the shape file.
H |
A neighborhood matrix. |
map |
SpatialPolygonsDataFrame object from the provided shapefile. |
Erica M. Porter, Matthew J. Keefe, Christopher T. Franck, and Marco A.R. Ferreira
#### Load extra libraries library(sp) library(sf) ### Read in a shapefile of the contiguous U.S. from package data system.path <- system.file("extdata", "us.shape48.shp", package = "ref.ICAR", mustWork = TRUE) shp.layer <- gsub('.shp','',basename(system.path)) shp.path <- dirname(system.path) us.shape48 <- st_read(dsn = path.expand(shp.path), layer = shp.layer) shp.data <- shape.H(system.path) names(shp.data)
#### Load extra libraries library(sp) library(sf) ### Read in a shapefile of the contiguous U.S. from package data system.path <- system.file("extdata", "us.shape48.shp", package = "ref.ICAR", mustWork = TRUE) shp.layer <- gsub('.shp','',basename(system.path)) shp.path <- dirname(system.path) us.shape48 <- st_read(dsn = path.expand(shp.path), layer = shp.layer) shp.data <- shape.H(system.path) names(shp.data)