Simulation envelopes for Foxall J not containing the line y=1 (fully reproducible example)

Question

During some exploratory analysis of my data I noticed that when Jfox was calculated for a point pattern relative to some polygons, I frequently obtained that the pooled envelopes area did not include the Theoretical line (at least in some regions). This was not the case when Jfox was computed between point patterns only.

More specifically, my dataset consisted of different replicate transects in the same geographical area, and for each transect I had (at least):

• one point pattern of interest
• a set of polygons representing landscape features

I was thus computing envelops for the Jfox function in each transect, and then pooling the envelopes together to see the overall effect (at least that was my intention).

In a previous post I was asking why this was happening and Prof. Adrian Baddeley suggested that it could have been an issue due to the fact that:

Pooling of objects (such as simulation envelopes) is only justified when the individual objects (envelopes) were created under identical conditions. In your example you appear to be pooling several envelopes of the Foxall J function Jfox(X,Y) that were computed for different datasets. Each envelope was computed by taking a particular dataset (X,Y), holding Y fixed, and generating random point patterns Xsim with the same average intensity as X. These are not identical conditions between the different envelopes, so it is not valid to pool them, and paradoxes can occur, such as you have found.

My interpretation of this answer is that, the problem is the pooling of the envelopes of the different transect.
I thus plotted the JFox envelopes of each replicate transect of that dataset (available here) and obtained the following:

library(spatstat)
#> Warning: package 'spatstat' was built under R version 4.1.3
#> Loading required package: spatstat.data
#> Warning: package 'spatstat.data' was built under R version 4.1.3
#> Loading required package: spatstat.geom
#> Warning: package 'spatstat.geom' was built under R version 4.1.3
#> spatstat.geom 3.0-3
#> Loading required package: spatstat.random
#> Warning: package 'spatstat.random' was built under R version 4.1.3
#> spatstat.random 3.0-1
#> Loading required package: spatstat.core
#> Warning: package 'spatstat.core' was built under R version 4.1.3
#> Loading required package: nlme
#> Loading required package: rpart
#> spatstat.core 2.4-4
#> Loading required package: spatstat.linnet
#> Warning: package 'spatstat.linnet' was built under R version 4.1.3
#> spatstat.linnet 2.3-2
#> 
#> spatstat 2.3-4       (nickname: 'Watch this space') 
#> For an introduction to spatstat, type 'beginner'

# Load the Hyperframe create for the prevoius question (envelopes nsims=6)
load('Q575569-envelopes-hyperframe.Rdata')
hypCV_06<-hypCV

This is how the dataset looks like:

plot(hypCV_06, quote({plot(rkaIM);plot(arub, pch=16, cols='white', cex=0.4,add=T)}))

plot(hypCV_06, quote({plot(arub_LSsubs_jfoxPTW)}))

The problem seems to persist without the pooling.

To verify if the “problem” persists in simulated data I do the following:

# Create some artificial data

# create random points
set.seed(120)
pts<-  runifpoint(20, win=owin(c(0,10),c(0,10)))
# create 2 windows
sq1<-square(c(2,3))
ds1<-disc(radius=1, centre=c(4,7))
ds2<-disc(radius=0.3, centre=c(7,3))
poly<-union.owin(sq1,ds1,ds2)


plot(pts)
plot(poly, add=T)

To use the Jfox function betwween a ppp and an owin object, I follow the approach suggested by A. Baddeley in his answere here.

# Convert owin object (representing the polygon to be used in Jfox) into a binary mask 
EPS<-0.01
winIM<- as.im(poly,  value=TRUE, na.replace=FALSE,eps=EPS)
winIM<- as.im(winIM,W=Window(pts),na.replace=FALSE,eps=EPS)

plot(winIM)
plot(pts,add=T)

I then compute the envelopes for Jfox on the simulated data using different amounts of simulations.
For consistency with my previuous question I used 6 simulation but checked what happens with less(i.e., 3), and more (i.e., 39 and 199).
Also for consistency, I calculate the envelopes of Jfox computed only on points.

# JFox only on ppp

envJpp_3<-envelope(Y=pts, fun=Jfox, funargs=list(Y=pts), nsim=3, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS) 
#> Generating 3 simulations of CSR with fixed number of points  ...
#> 1, 2,  3.
#> 
#> Done.
envJpp_6<-envelope(Y=pts, fun=Jfox, funargs=list(Y=pts), nsim=6, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS)
#> Generating 6 simulations of CSR with fixed number of points  ...
#> 1, 2, 3, 4, 5,  6.
#> 
#> Done.

# JFox of ppp relative to polygons

envJpw_3<-envelope(Y=pts, fun=Jfox, funargs=list(Y=winIM),  nsim=3, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS) 
#> Generating 3 simulations of CSR with fixed number of points  ...
#> 1, 2,  3.
#> 
#> Done.
envJpw_6<-envelope(Y=pts, fun=Jfox, funargs=list(Y=winIM),  nsim=6, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS) 
#> Generating 6 simulations of CSR with fixed number of points  ...
#> 1, 2, 3, 4, 5,  6.
#> 
#> Done.
envJpw_39<-envelope(Y=pts, fun=Jfox, funargs=list(Y=winIM),  nsim=39, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS) 
#> Generating 39 simulations of CSR with fixed number of points  ...
#> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38,  39.
#> 
#> Done.
envJpw_199<-envelope(Y=pts, fun=Jfox, funargs=list(Y=winIM),  nsim=199, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS) 
#> Generating 199 simulations of CSR with fixed number of points  ...
#> 1, 2,  [etd 3:43] 3, 4 [etd 3:49] .6 [etd 3:56] .8
#>  [etd 3:59] .10 [etd 3:57] .12 [etd 4:02] .14 [etd 4:00] .16
#>  [etd 3:55] .18 [etd 3:52] .20 [etd 3:49] .22 [etd 3:47] .24
#>  [etd 3:43] .26 [etd 3:40] .28 [etd 3:36] .30 [etd 3:33] .32
#>  ..... similar output cut to save space ...
#>  ]  199.
#> 
#> Done.

# Visualize results
par(mfrow=c(2,3))
# only points
plot(envJpp_3, xlim=c(0,1), main='envJpp_3')
plot(envJpp_6, xlim=c(0,1), main='envJpp_6')
# points vs polygons
plot(envJpw_3, main='envJpw_3')
plot(envJpw_6, main='envJpw_6')
plot(envJpw_39, main='envJpw_39')
plot(envJpw_199, main='envJpw_199')

Take home message here:

1. pointwise envelopes do not always include the Theoretical line for Jfox between points and polygons
2. This is particularly evident when few (<=6) simulations are used
3. the problem seems to desappear when the number of simulation increase
4. even for Jfox only between points the problem seems to appear if very few simulation are used. I want to check if the behaviour is the same when working with irregular polygons.
   Once again, data are available here)

# Load irregularly-shaped polygons 
load('TEST-stack_exchange-poly1.Rdata')
load('TEST-stack_exchange-poly2.Rdata')
load('TEST-stack_exchange-poly2IM.Rdata')
load('TEST-stack_exchange-pts2.Rdata')

This is how this second set of artificial data looks like:

plot(poly2IM)
plot(pts2,add=T, cols='white', pch=16)
plot(poly1,add=T, border=2)
plot(poly2,add=T)

Calulate the envelope in a similar way as did before

poly2_3<-envelope(Y=pts2, fun=Jfox, funargs=list(Y=poly2IM),  nsim=3, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS) 
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> Generating 3 simulations of CSR with fixed number of points  ...
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> 1,
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> 2,
#>  ..... similar output cut to save space ...
#> Done.
poly2_6<-envelope(Y=pts2, fun=Jfox, funargs=list(Y=poly2IM),  nsim=6, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS) 
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> Generating 6 simulations of CSR with fixed number of points  ...
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> 1,
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> 2,
#>  ..... similar output cut to save space ...
#> 
#> Done.
poly2_39<-envelope(Y=pts2, fun=Jfox, funargs=list(Y=poly2IM),  nsim=39, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS) 
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> Generating 39 simulations of CSR with fixed number of points  ...
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> 1,
#>  ..... similar output cut to save space ...
#> Done.
poly2_199<-envelope(Y=pts2, fun=Jfox, funargs=list(Y=poly2IM),  nsim=199, fix.n=T, savefuns=T, global=F, nrank=1, correction='best', eps=EPS) 
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> Generating 199 simulations of CSR with fixed number of points  ...
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> 1,
#> Warning in resolve.foxall.window(X, Y, W, isTRUE(warn.trim)): Trimming the
#> window of X to be a subset of the observation window of Y
#> 2,
#>  ..... similar output cut to save space ...
#> 
#> Done.

par(mfrow=c(2,2))
plot(poly2_3, main='poly2_3')
plot(poly2_6, main='poly2_6')
plot(poly2_39, main='poly2_39')
plot(poly2_199, main='poly2_199')

Once again, it seems that the envelope area can fail to include the theoretical line when few simulation are used.

If the problem is due to the number of simulation used, then I should solve it by increasing the simulations used for Jfox in the real dataset.
I thus recalclulated the envelopes for each transect in the real dataset suing 39 simulations:

# Load the real hyperframe created for the [prevoius question](https://stats.stackexchange.com/questions/575511/pointwise-envelopes-not-including-theoretical-line-foxall-j/575569#575569) and updated with 39 simulations 
load('Q575569-envelopes-hyperframe-sim39.Rdata')
plot(hypCV_39, quote({plot(arub_LSsubs_jfoxPTW)}))

The problem seems to persist. And note that in this last plot, no pooling has been done.

Questions

1. Is this a problem then?
2. Is this a problem of Jfox or of envelopes?
3. If this is not a problem, am I doing something wrong?

^{Created on 2022-11-24 by the reprex package (v2.0.1)}

Answer 1

1. No.
2. No.
3. No.

A simulation envelope is constructed by finding the largest and smallest values of the simulated summary functions. The simulated function values are random, so there is no guarantee that the "theoretical" version of the summary function will lie in between these two extremes.

For any choice of summary function, it can occur and does sometimes occur that the significance band does not include the theoretical curve.

This is more likely to occur when the number of simulations nsim is small, the number of points in the data point pattern is small, and the summary function has high variance. These conditions all hold in your example.

Try plot(envelope(cells[1:5], Kest, nsim=3, fix.n=TRUE)) a few times.

For a simple analogy, suppose that, at a particular distance, the theoretical value is equal to 0.5, and the simulated values are random numbers uniformly distributed on the interval from 0 to 1. The theoretical mean of the simulated numbers is equal to the theoretical value 0.5. But if you just pick nsim=3 numbers at random between 0 and 1, there is a chance that they will all lie above 0.5 (probability 1/8) or will all lie below 0.5 (probability 1/8) so that there is a 1/4 chance that the "envelope" will not contain the theoretical value. However if we pick nsim=19 or 39 numbers at random, the "envelope" is almost certain to contain the theoretical value.

This phenomenon is not specific to Jfox. It is more likely to occur with J functions than with K functions, because the J functions have higher variance.

In cases where this phenomenon occurs, the hypothesis test associated with the envelope is still valid (i.e. still has the quoted significance level), but typically has very low power, so is not trustworthy.