10.3 Modeling diseased risk

\[ Y_i| \theta_i \sim Poisson(E_i \times \theta_i) \]

where

  • \(\theta_i\) is the relative risk for county \(i\),

  • \(u_i\) is the structured random effect for county \(i\) modeled with an intrinsic conditional autoregressive (ICAR) model,

\[ u_i|u_{-i} \sim N(\bar{u}_{\delta_i} \frac{1}{\tau_u n_{\delta_i}}) \]

  • \(v_i\) is the random effect for stratum \(i\)

\[ v_i \sim N(0, \frac{1}{\tau_v}) \]

10.3.1 Neighborhood structure

library(spdep)
library(INLA)
nb <- poly2nb(map)
nb2INLA("map.adj", nb)
g <- inla.read.graph(filename = "map.adj")

10.3.2 Model

map$re_u <- 1:nrow(map)
map$re_v <- 1:nrow(map)
formula <- Y ~ smoking +
  f(re_u, 
    model = "besag", 
    graph = g, 
    scale.model = TRUE) +
  f(re_v, model = "iid")
res <- inla(formula, 
            family = "poisson", 
            data = map, 
            E = E,
            control.predictor = list(compute = TRUE),
            control.compute = list(return.marginals.predictor = TRUE))
 res$summary.fixed
                mean     sd 0.025quant 0.5quant
 (Intercept) -0.3235 0.1498   -0.61925  -0.3233
 smoking      1.1546 0.6226   -0.07569   1.1560
             0.975quant    mode       kld
 (Intercept)   -0.02877 -0.3234 3.534e-08
 smoking        2.37845  1.1563 3.545e-08

10.3.3 Relative Risk

 res$summary.fitted.values[1:3, ]
                       mean      sd 0.025quant 0.5quant
 fitted.Predictor.01 0.8781 0.05808     0.7648   0.8778
 fitted.Predictor.02 1.0597 0.02750     1.0072   1.0592
 fitted.Predictor.03 0.9646 0.05089     0.8604   0.9657
                     0.975quant   mode
 fitted.Predictor.01     0.9936 0.8778
 fitted.Predictor.02     1.1150 1.0582
 fitted.Predictor.03     1.0622 0.9681
 
 
 # relative risk
 map$RR <- res$summary.fitted.values[, "mean"]
 
 
 # lower and upper limits 95% CI
 map$LL <- res$summary.fitted.values[, "0.025quant"]
 map$UL <- res$summary.fitted.values[, "0.975quant"]

See the map here:

https://www.paulamoraga.com/book-spatial/disease-risk-modeling.html#mapping-smr