## Example taken and adapted from
## http://www.stats.gla.ac.uk/data/NERC/NERC10/Mark/Data%20and%20code%20for%20network/bayes.R

# The rates of lip cancer in 56 counties in Scotland have been analysed by Clayton and
# Kaldor (1987) and Breslow and Clayton (1993). The form of the data includes the ob-
# served and expected cases (expected numbers based on the population and its age and
# sex distribution in the county), a covariate measuring the percentage of the population
# engaged in agriculture, fishing, or forestry, and the “position” of each county expressed as
# a list of adjacent counties.


# Load the INLA package
library(INLA)
# Make the Scottish lip cancer data available for use
data(Scotland)
# Get a description of the data
#?Scotland

# Work out the path for the file of neighbourhood information
scot.neigh <- system.file("demodata/scotland.graph", package="INLA")
# Create a file for the formula, containing only a spatial terms
# and a scaled covariate (divided by 10)
formula <- Counts~f(Region,model="besag", graph = scot.neigh,
    hyper=list(prec=list(prior="loggamma", param=c(0.5,0.0005)))) + I(X/10)
# Fit a Bayesian model using the inla() function
scot.res <- inla(formula,family="poisson",E=E,data=Scotland)

# Look at results
summary(scot.res)
# Look at the estimated marginal posterior distributions
#plot(scot.res)

# Now try plotting some maps to illustrate the INLA results
# Load the SpatialEpi package; if unavailable try 2nd option below:
library(SpatialEpi)
data(scotland)

scotland.map <- scotland$spatial.polygon
# First map - the raw counts
mapvariable(Scotland$Counts/Scotland$E,scotland.map,
    main="Counts, Lip Cancer")
mapvariable(exp(scot.res$summary.random$Region$"0.5quant"),scotland.map,
    main="Smoothed Random Effects, Lip Cancer")