####################################################
#### Run the optmiser on simulated data on a lattice
####################################################
#### Load other libraries required
library(MASS)
#### Set up a square lattice region
x.easting <- 1:10
x.northing <- 1:10
Grid <- expand.grid(x.easting, x.northing)
K <- nrow(Grid)
#### Split the area into two groups between which there will be a boundary.
groups <-rep(1, K)
groups[Grid$Var1>5] <- 2
#### set up distance and neighbourhood (W, based on sharing a common border) matrices
distance <- as.matrix(dist(Grid))
W <-array(0, c(K,K))
W[distance==1] <-1
#### Generate the response data
phi <- mvrnorm(n=1, mu=groups, Sigma=0.05 * exp(-0.1 * distance))
lp <- 4 + phi
Y <- rpois(n=K, lambda =exp(lp))
# Compute the spdata object on the linear predictor scale
spdata <- log(Y) - mean(log(Y))
#### Run W matrix optmiser
W.est <- W.estimate(W, spdata)
Run the code above in your browser using DataLab