# The following script was inspired by a discussion at http://gis.stackexchange.com/q/19522/4490
# (later merged with http://gis.stackexchange.com/q/4484/4490)
#
# The script performs permutations to compare two marked point patterns. A simple ANN statistic
# is used in this example.
#
library(spatstat)
# Generate random points
Z = as.im(function(x,y){ 5 * exp(2 * x - 1) }, owin()) # Generates a sloped density surface
P1 = rpoint(100, Z)
P2 = rpoint(50, Z)
# Combine points into a single markes PPP
P = superimpose(A=P1,B=P2) # Assign label A to P1 and B to P2
# Plot points
plot.owin(P$window,col="#DDDDDD")
plot(split(P)$B,pch=16,cols="white",add=T)
plot(split(P)$A,pch=16,cols="blue",add=T)
# Permutate labels (do not change point locations)
nn.sim = vector()
P.r = P
for(i in 1:1999){
marks(P.r) = sample(P$marks) # Randomly reassign labels
nn.sim[i] = mean(nncross(split(P.r)$A,split(P.r)$B)$dist)
}
# Plot distribution and original statistic (using red line)
hist(nn.sim, breaks=30)
abline(v=mean(nncross(split(P)$A,split(P)$B)$dist),col="red")
# Compute empirical cumulative distribution
nn.sim.ecdf = ecdf(nn.sim)
# See how the original stat compares to the simulated distribution</pre>
nn.sim.ecdf(mean(nncross(split(P)$A,split(P)$B)$dist))