Ickstadt and Wolpert (1998) and Wolpert and Ickstadt (1998) analyse the spatial distribution of hickory trees (a subdominant species) in a 140 metre x 140 metre research plot in Duke Forest, North Carolina, USA, with the aim of investigating the regularity of the distribution of trees: spatial clustering would reveal that the subdominant species is receding from, or encroaching into, a stand recovering from a disturbance, while regularity would suggest a more stable and undisturbed recent history.

The data comprise counts of trees, Y

Y

μ

γ

Ickstadt and Wolpert (1998) consider two alternative kernel functions: an adjacency-based kernel and a distance-based kernel. The adjacency based kernel is defined as:

k

k

neighbours considered) and n

k

The distance based kernel is defined as:

k

For both kernels, the parameter exp(θ

Ickstadt and Wolpert (1998) elicit expert prior opinion concerning the unknown hyperparameters θ

This model may be implemented in BUGS using the pois.conv distribution. The code is given below.

Modelmodel {

# likelihood for (i in 1:N) {

Y[i] ~ dpois.conv(mu[i,])

for (j in 1:J) {

mu[i, j] <- A[i] * lambda[i, j]

lambda[i, j] <- k[i, j] * gamma[j]

}

}

# priors (see Ickstadt and Wolpert (1998) for details of prior elicitation)

for (j in 1:J) {

gamma[j] ~ dgamma(alpha, beta)

}

alpha <- exp(theta0)

beta <- exp(-theta1)

theta0 ~ dnorm(-0.383, 1)

theta1 ~ dnorm(-5.190, 1)

theta2 ~ dnorm(-1.775, 1) # prior on theta2 for adjacency-based kernel

# theta2 ~ dnorm(1.280, 1) # prior on theta2 for distance-based kernel

# compute adjacency-based kernel # Note N = J in this example (necessary for adjacency-based kernel) for (i in 1:N) {

k[i, i] <- 1

for (j in 1:J) {

# distance between areas i and j

d[i, j] <- sqrt((x[i] - x[j])*(x[i] - x[j]) + (y[i] - y[j])*(y[i] - y[j]))

# (only needed to help compute nearest neighbour weights;

# alternatively, read matrix w from file)

w[i, j] <- step(30.1 - d[i, j]) # nearest neighbour weights:

# areas are 30 sq m, so any pair of areas with centroids > 30m apart are not

# nearest neighbours (0.1 added to account for numeric imprecision in d)

}

for (j in (i+1):J) {

k[i, j] <- w[i, j] * exp(theta2) / (sum(w[i,]) - 1)

k[j, i] <- w[j, i] * exp(theta2) / (sum(w[j,]) - 1)

}

}

# alternatively, compute distance-based kernel # for (i in 1:N) {

# k[i, i] <- 1

# for (j in 1:J) {

# distance between areas i and j

# d[i, j] <- sqrt((x[i] - x[j])*(x[i] - x[j]) + (y[i] - y[j])*(y[i] - y[j]))

# }

# for (j in (i+1):J) {

# k[i, j] <- exp(-pow(d[i, j]/exp(theta2), 2))

# k[j, i] <- exp(-pow(d[j, i]/exp(theta2), 2))

# }

# }

# summary quantities for posterior inference

for (i in 1:N) {

# smoothed density of hickory trees (per sq metre) in area i

density[i] <- sum(lambda[i, ])

# observed density of hickory trees (per sq metre) in area i

obs.density[i] <- Y[i]/A[i]

}

# large values indicate strong spatial dependence;

# spatial.effect -> 0 indicates spatial independence

spatial.effect <- exp(theta2)

}

```
list(N = 16, J = 16,
```

Y = c(4, 3, 0, 2, 3, 2, 0, 0, 2, 1, 2, 3, 17, 3, 2, 13),

A = c(900, 900, 900, 900, 900, 900, 900, 900, 900, 900, 900, 900, 900, 900, 900, 900),

x = c(15, 45, 75, 105, 15, 45, 75, 105, 15, 45, 75, 105, 15, 45, 75, 105),

y = c(15, 15, 15, 15, 45, 45, 45, 45, 75, 75, 75, 75, 105, 105, 105, 105)

)

```
list(theta0 = -1, theta1 = -3, theta2 = 0,
```

gamma = c(0.00406, 0.00518, 0.00916, 0.0117, 0.0101, 0.00249, 0.0178, 0.00112, 0.000703,

0.00427, 0.0112, 0.00599, 0.00131, 0.00599, 0.000738, 0.00259))

```
list(theta0 = -0.1, theta1 = -0.3, theta2 = 1.0,
```

gamma = c(0.0406, 0.0518, 0.0916, 0.117, 0.101, 0.0249, 0.178, 0.0112, 0.00703,

0.0427, 0.112, 0.0599, 0.0131, 0.0599, 0.00738, 0.0259))

Results

Ickstadt and Wolpert report a posterior mean (sd) of 0.281 (0.156) for the spatial effect, exp(θ

Assuming distance-based kernel (equivalent to Ickstadt and Wolpert's model M

Ickstadt and Wolpert report a posterior mean (sd) of 7.449 (6.764) for the spatial effect, exp(θ

Note however that convergence of this distance-based model is not very good, and the results reported above are not particularly stable (the adjacency-based model appears to converge much better). This may be partly due to the rather strong prior information about theta2 for the distance-based model which, as noted by Ickstadt and Wolpert, appears to conflict with the data.