ref.MCMC | R Documentation |
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
)
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
keefe2018ref.ICAR
\insertRefKeefe_2018ref.ICAR
\insertRefFerreira_2021ref.ICAR
#### 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)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.