,
Silvia Liverani [aut] ,
Peter Congdon [aut] | Title: | CAR-MM Modelling in Stan |
|---|---|
| Description: | 'Stan' based functions to estimate CAR-MM models. These models allow to estimate Generalised Linear Models with CAR (conditional autoregressive) spatial random effects for spatially and temporally misaligned data, provided a suitable Multiple Membership matrix. The main references are Gramatica, Liverani and Congdon (2023) <doi:10.1214/23-BA1370>, Petrof, Neyens, Nuyts, Nackaerts, Nemery and Faes (2020) <doi:10.1002/sim.8697> and Gramatica, Congdon and Liverani <doi:10.1111/rssc.12480>. |
| Authors: | Marco Gramatica [aut, cre] ,
Silvia Liverani [aut] ,
Peter Congdon [aut] |
| Maintainer: | Marco Gramatica <[email protected]> |
| License: | GPL (>=3) |
| Version: | 0.1.1 |
| Built: | 2025-03-11 06:03:46 UTC |
| Source: | https://github.com/marcogramatica/carme |
CAR-MM modelling in Stan
Stan Development Team (2023). RStan: the R interface to Stan. R package version 2.26.11. https://mc-stan.org
Marco Gramatica. Silvia Liverani. Peter Congdon. Structure Induced by a Multiple Membership Transformation on the Conditional Autoregressive Model. Bayesian Analysis Advance Publication 1 - 25, 2023. https://doi.org/10.1214/23-BA1370
Petrof, O, Neyens, T, Nuyts, V, Nackaerts, K, Nemery, B, Faes, C. On the impact of residential history in the spatial analysis of diseases with a long latency period: A study of mesothelioma in Belgium. Statistics in Medicine. 2020; 39: 3840– 3866. https://doi.org/10.1002/sim.8697
Marco Gramatica, Peter Congdon, Silvia Liverani, Bayesian Modelling for Spatially Misaligned Health Areal Data: A Multiple Membership Approach, Journal of the Royal Statistical Society Series C: Applied Statistics, Volume 70, Issue 3, June 2021, Pages 645–666, https://doi.org/10.1111/rssc.12480
CAR-MM prior model
car_mm(d_list, ...)car_mm(d_list, ...)
d_list |
List of data inputs for the stan model. |
... |
Arguments passed to |
An object of class stanfit returned by rstan::sampling
Marco Gramatica. Silvia Liverani. Peter Congdon. "Structure Induced by a Multiple Membership Transformation on the Conditional Autoregressive Model." Bayesian Analysis Advance Publication 1 - 25, 2023. https://doi.org/10.1214/23-BA1370
Petrof, O, Neyens, T, Nuyts, V, Nackaerts, K, Nemery, B, Faes, C. On the impact of residential history in the spatial analysis of diseases with a long latency period: A study of mesothelioma in Belgium. Statistics in Medicine. 2020; 39: 3840– 3866. https://doi.org/10.1002/sim.8697
set.seed(455) #---- Load data data(W_sel) ## Number of areas n <- nrow(W_sel) ## Number of memberships m <- 153 #---- Simulate covariates X <- cbind(rnorm(nrow(W_sel)), rnorm(nrow(W_sel))) ## Min-max normalisation X_cent <- apply(X, 2, function(x) (x - min(x))/diff(range(x))) #---- Simulate MM matrix w_ord <- c(.5, .35, .15) # Weight of each neighbours orders ord <- length(w_ord) - 1 # Order of neighbours to include H_sel_sim <- sim_MM_matrix( W = W_sel, m = m, ord = ord, w_ord = w_ord, id_vec = rep(1, nrow(W_sel)) ) #---- Simulate outcomes ## Linear term parameters gamma <- -.5 # Intercept beta <- c(1, .5) # Covariates coefficients ## CAR random effects phi_car <- sim_car(W = W_sel, alpha = .9, tau = 5) # Areal log relative risks l_RR <- X_cent %*% beta + phi_car ## Membership log relative risks l_RR_mm <- as.numeric(apply(H_sel_sim, 1, function(x) x %*% l_RR)) ## Expected rates exp_rates <- rpois(m, lambda = 20) ## Outcomes y <- rpois(m, lambda = exp_rates*exp(l_RR_mm)) #---- Create dataset for stan function d_sel <- list( # Number of areas n = nrow(W_sel), # Covariates k = ncol(X_cent), X_cov = X_cent, # Adjacency W_n = sum(W_sel) / 2, # Number of neighbour pairs W = W_sel, # Memberships m = nrow(H_sel_sim), H = H_sel_sim, # Outcomes y = y, log_offset = log(exp_rates), # Prior parameters ## Intercept (mean and sd of normal prior) mu_gamma = 0, sigma_gamma = 1, ## Covariates (mean and sd of normal prior) mu_beta = 0, sigma_beta = 1, ## Marginal precision gamma prior tau_shape = 2, tau_rate = 0.2 ) #---- HMC parameters niter <- 1E4 nchains <- 4 #---- Stan sampling fit <- car_mm( d_list = d_sel, # arguments passed to sampling iter = niter, chains = nchains, refresh = 500, control = list(adapt_delta = .99, max_treedepth = 15) )set.seed(455) #---- Load data data(W_sel) ## Number of areas n <- nrow(W_sel) ## Number of memberships m <- 153 #---- Simulate covariates X <- cbind(rnorm(nrow(W_sel)), rnorm(nrow(W_sel))) ## Min-max normalisation X_cent <- apply(X, 2, function(x) (x - min(x))/diff(range(x))) #---- Simulate MM matrix w_ord <- c(.5, .35, .15) # Weight of each neighbours orders ord <- length(w_ord) - 1 # Order of neighbours to include H_sel_sim <- sim_MM_matrix( W = W_sel, m = m, ord = ord, w_ord = w_ord, id_vec = rep(1, nrow(W_sel)) ) #---- Simulate outcomes ## Linear term parameters gamma <- -.5 # Intercept beta <- c(1, .5) # Covariates coefficients ## CAR random effects phi_car <- sim_car(W = W_sel, alpha = .9, tau = 5) # Areal log relative risks l_RR <- X_cent %*% beta + phi_car ## Membership log relative risks l_RR_mm <- as.numeric(apply(H_sel_sim, 1, function(x) x %*% l_RR)) ## Expected rates exp_rates <- rpois(m, lambda = 20) ## Outcomes y <- rpois(m, lambda = exp_rates*exp(l_RR_mm)) #---- Create dataset for stan function d_sel <- list( # Number of areas n = nrow(W_sel), # Covariates k = ncol(X_cent), X_cov = X_cent, # Adjacency W_n = sum(W_sel) / 2, # Number of neighbour pairs W = W_sel, # Memberships m = nrow(H_sel_sim), H = H_sel_sim, # Outcomes y = y, log_offset = log(exp_rates), # Prior parameters ## Intercept (mean and sd of normal prior) mu_gamma = 0, sigma_gamma = 1, ## Covariates (mean and sd of normal prior) mu_beta = 0, sigma_beta = 1, ## Marginal precision gamma prior tau_shape = 2, tau_rate = 0.2 ) #---- HMC parameters niter <- 1E4 nchains <- 4 #---- Stan sampling fit <- car_mm( d_list = d_sel, # arguments passed to sampling iter = niter, chains = nchains, refresh = 500, control = list(adapt_delta = .99, max_treedepth = 15) )
sim_car returns a vector of CAR distributed random effects
sim_car(W, alpha = 0.5, tau = 5)sim_car(W, alpha = 0.5, tau = 5)
W |
Symmetric adjacency matrix of size |
alpha |
properness parameter between 0 and 1. Defaults to 0.5 |
tau |
marginal precision. Defaults to 5 |
a vector of length n
Jin, X., Carlin, B.P. and Banerjee, S. (2005), Generalized Hierarchical Multivariate CAR Models for Areal Data. Biometrics, 61: 950-961. https://doi.org/10.1111/j.1541-0420.2005.00359.x
data(W_sel) sim_car(W = W_sel, alpha = .9, tau = 5)data(W_sel) sim_car(W = W_sel, alpha = .9, tau = 5)
sim_MM_matrix returns a multiple membership matrix simulated based on an
adjacency matrix according to the method described in
sim_MM_matrix(W, m, ord = 3, w_ord, id_vec, excess_areas = FALSE, red_areas)sim_MM_matrix(W, m, ord = 3, w_ord, id_vec, excess_areas = FALSE, red_areas)
W |
Symmetric adjacency matrix of size |
m |
Integer. Number of membership to simulate |
ord |
Integer. Maximum order of neighbours to be used to simulate the
memberships based on the adjacency matrix |
w_ord |
A vector of length |
id_vec |
Vector of zeros and ones of length |
excess_areas |
if different from FALSE it indicates the indices of the
areas to reuse in simulating memberships, whenever |
red_areas |
vector of indices of areas to use if |
an m x n matrix of weights
Marco Gramatica. Silvia Liverani. Peter Congdon. "Structure Induced by a Multiple Membership Transformation on the Conditional Autoregressive Model." Bayesian Analysis Advance Publication 1 - 25, 2023. https://doi.org/10.1214/23-BA1370
set.seed(455) #---- Load data data(W_sel) ## Number of areas n <- nrow(W_sel) ## Number of memberships m <- 153 #---- Simulate MM matrix w_ord <- c(.5, .35, .15) # Weight of each neighbours orders ord <- length(w_ord) - 1 # Order of neighbours to include H_sel_sim <- sim_MM_matrix( W = W_sel, m = m, ord = ord, w_ord = w_ord, id_vec = rep(1, nrow(W_sel)) )set.seed(455) #---- Load data data(W_sel) ## Number of areas n <- nrow(W_sel) ## Number of memberships m <- 153 #---- Simulate MM matrix w_ord <- c(.5, .35, .15) # Weight of each neighbours orders ord <- length(w_ord) - 1 # Order of neighbours to include H_sel_sim <- sim_MM_matrix( W = W_sel, m = m, ord = ord, w_ord = w_ord, id_vec = rep(1, nrow(W_sel)) )
Adjacency matrix of 152 MSOAs in South East London, used for the data analysis in the paper "Structure induced by a multiple membership transformation on the Conditional Autoregressive model". Column and rows names indicate the MSOA code.
data(W_sel)data(W_sel)
A 152x152 symmetric matrix
Marco Gramatica. Silvia Liverani. Peter Congdon. "Structure Induced by a Multiple Membership Transformation on the Conditional Autoregressive Model." Bayesian Analysis Advance Publication 1 - 25, 2023. https://doi.org/10.1214/23-BA1370