Convolution priors: Lung cancer in a
London Health Authority
This example has been modified from a London Health Authority annual report. The theme of the report was inequality in health and the aim was to investiagte links between poverty and ill health. To investigate this issue, the Health Authority carried out a series of disease mapping studies at census ward level. There are 44 wards in the Health Authority in this example. The data are simulated observed and exepected counts of lung cancer incidence in males aged 65 and over living in the Health Authority region; a ward level index of socio-economic devprivation is also available.
We fit the following model, allowing a convolution prior for the random effects:
O i ~ Poisson( m i )
log m i = log E i + a + b depriv i + b i + h i
where depriv i is the devprivation covariate, b i are spatial random effects assigned a CAR prior, and we introduce a second set of random effects h i for which we assume an exchangeable Normal prior. The random effect for each area is thus the sum of a spatially structured component b i and an unstructured component h i . This is termed a convolution prior (Besag, York and Mollie 1991; Mollie 1996). Besag, York and Mollie argue that this model is more flexible than assume only CAR random effects, since it allows the data to decide how much of the residual disease risk is due to spatially structured variation, and how much is unstructured over-dispersion.
The code for this model is given below.
Model
model {
for (i in 1 : N) {
# Likelihood
O[i] ~ dpois(mu[i])
log(mu[i]) <- log(E[i]) + alpha + beta * depriv[i] + b[i] + h[i]
# Area-specific relative risk (for maps)
RR[i] <- exp(alpha + beta * depriv[i] + b[i] + h[i])
# Exchangeable prior on unstructured random effects
h[i] ~ dnorm(0, tau.h)
}
# CAR prior distribution for spatial random effects:
b[1 : N] ~ car.normal(adj[], weights[], num[], tau.b)
for(k in 1:sumNumNeigh) {
weights[k] <- 1
}
# Other priors:
alpha ~ dflat()
beta ~ dnorm(0.0, 1.0E-5)
tau.b ~ dgamma(0.5, 0.0005)
sigma.b <- sqrt(1 / tau.b)
tau.h ~ dgamma(0.5, 0.0005)
sigma.h <- sqrt(1 / tau.h)
}
Data click here to open data
Inits click here to open initial values