model {

          # Likelihood
             for (i in 1 : I) {
          # convolution likelihood: 
          # Y[i] ~ Poisson( sum_ j { mu[i, j] } )
          # where mu[i, j] = pop[i] * lambda[i, j] and
          # pop[i] = (standardised) population (in 100s) in area i<br> 0DX </font><font face="Arial" color="#000100" size="1">&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; 0DX </font><font face="Arial" color="#000100" size="3"># lambda[i, j] = disease rate in area i attributed to jth source, where
          # the sources include both observed and latent covariates.
                Y[i] ~ dpois.conv(mu[i, ])   
  
             for (j in 1 : J) { 
          # rescale kernel matrix
                k1[i, j] <- ker[i, j] / prec
          # jth source (jth latent spatial grid cell)
                lambda[i, j] <- beta2 * gamma[j] * k1[i, j]
             }
          # (J1)th source (overall intercept; represents background rate)
             lambda[i,J1] <- beta0
         # (J2)th source (effect of observed NO2 covariate in area i)
             lambda[i,J2] <- beta1 * NO2[i]
          # Note: if additional covariates have been measured, these 
          # should be included in the same way as NO2, giving terms
          # lambda[i, J3], lambda[i, J4],.....etc. 
             for (j in 1 : J2) { 
                mu[i, j] <- lambda[i, j] * pop[i]
             }
          }
                
          # Priors
          # See Table 22.1 in Best et al. (2000b) for further details.
          
          # Priors for latent spatial grid cell (gamma) parameters:
          # 
          # Assume gamma_j ~ gamma(a, t) 
          # 
          # Prior shape parameter for gamma distn on gamma[j]s will change in proportion to
         # the size of the latent grid cells, to ensure aggregation consistency. The prior inverse
         # scale parameter for gamma distn on gamma[j]s stays constant across different
         # partitions, but depends on the units used to measure the area of the latent grid cells:
         # here we are assuming # area is in sq kilometres, and we take tauG = 1 (per km^2)
         # which corresponds to assuming # that our prior guess at the value of the gamma[j]s<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;# is based on observing about 1 * (area of study region = 30*20km^2) = 600<br> 0DX </font><font face="Arial" color="#000100" size="1">&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; 0DX </font><font face="Arial" color="#808080" size="3"># prior pointsover the study region. 
             
         # Note: area of latent grid cells is read in from data file
             alphaG <- tauG * area
             tauG <- 1 
  
             for (j in 1 : J) {
                gamma[j] ~ dgamma(alphaG, tauG)
             }
 
          # Priors for beta coefficients:
          #
          # Assume priors on each beta are of the form beta ~ gamma(a, t).
          #
          # Here, parameter of the prior are chosen by setting the prior mean to give an equal 
          # number of cases allocated to each of the 3 sources (baseline, NO2 and latent), i.e.
         # prior mean = a / t = a * 3 * Xbar / Ybar,
          # where Ybar is the overall disease rate in number of cases per unit population,
          # and Xbar is the population-weighted mean of covariate X. 
          # 
          # We also take shape parameter a = 0.575, since this gives the ratio of the 90th : 10th
          # percentile of the prior distn = 100, which is suitably diffuse.
          
             YBAR <- sum(Y[])/sum(pop[])   # mean number of cases per unit population
             NO2BAR <- inprod(NO2[], pop[]) / sum(pop[])   # population-weighted mean NO2
             
         # Shape parameter for Intercept (beta0)
             alpha0<- 0.575
         # Inverse scale parameter for Intercept (beta0)
             tau0 <- 3 * alpha0 / YBAR
         # Shape parameter for effect of NO2 (beta1)
             alpha1<- 0.575
         # Inverse scale parameter for effect of NO2 (beta1)
             tau1 <- 3 * alpha1 * NO2BAR / YBAR
         # Shape parameter for Latent coefficient (beta2)
             alpha2 <- 0.575
         # Inverse scale parameter for Latent coefficient (beta2)
             tau2 <- 3 * alpha2   / YBAR

             beta0 ~ dgamma(alpha0, tau0)
             beta1 ~ dgamma(alpha1, tau1)
             beta2 ~ dgamma(alpha2, tau2)

 
          # Summary quantities for posterior inference 
             for(i in 1:I) {
          # Overall disease rate in area i 
                RATE[i] <- sum(mu[i, ])/pop[i]
          # Rate associated with latent spatial factor in area i 
                LATENT[i] <- beta2 * inprod(k1[i,], gamma[])
             }

          # Expected rate (number of cases per unit population) attributed to each source 
          # (see Table 22.2 in Best et al (2000b) )

          # Background (overall baseline) rate per 100 population
          rate.base <- beta0

          # Number of cases per 100 population attributed
          # to NO2 = beta1 * population-weighted mean value
          # of NO2 across the study region. 
             rate.no2 <- beta1 * inprod(pop[], NO2[]) / sum(pop[])
             
          # Number of cases per 100 population attributed to 
          # latent sources = beta2 * population-weighted mean
          # value of the latent spatial factor (sum_j k_ij gamma_j)
             rate.latent <- inprod(pop[], LATENT[]) / sum(pop[])
  
          # Percentage of cases attributed to each source:
             Total <- rate.base + rate.no2 + rate.latent 
         # % cases attributed to baseline risk
            PC.base <- rate.base / Total * 100
         # % cases attributed to NO2 exposure
             PC.no2 <- rate.no2 / Total * 100
         # % cases attributed to latent spatial factors
             PC.latent <- rate.latent / Total * 100

          }