# Estimating the 𝐽 function without edge correctionAdrian Baddeley and Katja Schladitz were supported by a grant from the Australian Research Council. Martin Kerscher was supported by the “Sonderforschungsbereich SFB 375 für Astroteilchenphysik der Deutschen Forschungsgemeinschaft”. Bryan Scott was supported by an Australian Postgraduate Award and CSIRO Postgraduate Studentship.

## 1 Introduction

A spatial point pattern is often studied by estimating point process characteristics such as the empty space function $`F`$, the nearest-neighbour distance distribution function $`G`$, and the $`K`$-function. Here we consider

$$J(r)=\frac{1G(r)}{1F(r)}$$
(1)   
as advocated by van Lieshout and Baddeley (1996). The $`J`$ function is identically equal to 1 for a Poisson process, and values of $`J(r)`$ less than or greater than 1 are suggestive of clustering or regularity, respectively. (Note that it is of course possible to find non-Poisson point processes for which $`J(r)=1`$, as Bedford and van den Berg (1997) have shown).

In practice, observation of the point process is usually restricted to some bounded window $`W`$. As a consequence, estimation of the summary functions, which is based on the measurement of various distances of the point process, is hampered by the “edge effects” (bias and censoring) introduced by restricting observation of these distances to $`W`$. In order to counter edge effects it is necessary to apply some form of edge correction to the empirical estimates of the summary functions. For further details see Baddeley (1999); Stoyan et al. (1995); Cressie (1991); Ripley (1988).

Unbiasedness is highly desirable when a summary function estimate is to be compared directly to the corresponding theoretical value for a point process model. However, as Diggle argues in the discussion of Ripley (1977) and in Diggle (1983), unbiasedness is not essential when using a summary function estimator as the test statistic in a hypothesis test, since the bias will be accounted for in the null distribution of the test statistic.

This paper studies the uncorrected estimator of $`J`$ obtained by ignoring edge effects and computing $`\widehat{J}=(1\widehat{G})/(1\widehat{F})`$ from the uncorrected, empirical distributions $`\widehat{G}`$ and $`\widehat{F}`$ of distances observed in a compact window. It was prompted by the accidental discovery that this uncorrected estimator remarkably still yields values $`\widehat{J}(r)`$ approximately equal to 1 for the Poisson process. An intuitive explanation is that the relative bias due to edge effects is roughly equal for the estimates of $`1G`$ and $`1F`$, so that these biases approximately cancel in the ratio estimator of $`J`$. It follows that the uncorrected estimate of $`J`$ could be used for the direct visual assessment of deviations from the Poisson process, something not possible with the uncorrected estimates of $`F`$, $`G`$ or $`K`$. Our aim is to formalise this uncorrected estimator of $`J`$, and to investigate its use as a summary statistic and as a test statistic in point pattern analysis.

The paper is organised as follows. Section 2 outlines and generalises some fundamental point process ideas. In Section 3 we define the $`J_W`$ function estimated by the uncorrected procedure described above and derives some of its properties. In Section 4 we verify that the natural estimator of $`J_W`$ is the uncorrected estimate $`\widehat{J}_W`$. Finally, Section 5 presents the results of a computational experiment to compare the power of Monte Carlo Tests constructed from estimates of $`J`$ and $`J_W`$ as well as estimation results from the simulations of various point process models in square and rectangular windows in $`^2`$ and in a cubic window in $`^3`$.

## 2 Background

Let $`X`$ be a stationary point process in $`^d`$ with intensity $`\lambda `$. (For details of the theory of point processes see Daley and Vere-Jones, 1988; Cressie, 1991; Stoyan et al., 1995). The empty space function $`F(r)`$ is the probability of finding a point of the process within a radius $`r`$ of an arbitrary fixed point:

$$F(r)=(XB(0,r)\mathrm{}),$$
(2)   
where $`B(x,r)`$ denotes the ball of radius $`r`$ centred at $`x`$. The nearest neighbour distance distribution function $`G(r)`$ is the probability of finding another point of the process in the ball of radius $`r`$ centred at a “typical” point of the process:

$$G(r)=^{0!}(XB(0,r)\mathrm{}),$$
(3)   
where $`^{0!}`$ denotes the reduced Palm distribution at the origin $`0`$. Roughly speaking $`^{0!}`$ is the distribution of rest of the process $`X\{0\}`$ given there is a point of the process at the origin (see Daley and Vere-Jones, 1988).

Let $`W`$ be a compact observation window in $`^d`$ with nonempty interior. The construction of estimators for $`F`$ is based on the stationarity of $`X`$ yielding

$$F(r)=\frac{1}{|W|}\underset{W}{}(XB(x,r)\mathrm{})𝑑x,$$
(4)   
where $`|W|`$ denotes the volume of $`W`$. For $`G`$ the Campbell-Mecke formula (Stoyan et al., 1995, (4.4.3)) gives

$$G(r)=\frac{1}{\lambda |W|}𝔼\underset{xXW}{}\mathrm{𝟙}\{XB(x,r)\{x\}\mathrm{}\}.$$
(5)   
In both cases we need information about $`XB(x,r)`$ for $`xW`$, whereas we only observe $`XB(x,r)W`$. Usually this edge effect problem is countered by restricting the integration in (4) and summation in (5) to those points $`x`$ for which $`B(x,r)W`$ (the “border method”) or by weighting the contributions to the integral and sum so as to correct for the bias (see for example Baddeley, 1999).

The uncorrected, empirical distributions of distances observed in the window $`W`$ correspond to simply replacing $`X`$ by $`XW`$ in (4) and (5). In order to investigate the effect of this, we extend $`F`$ and $`G`$ to functionals as follows.

###### Definition 1.

For every compact set $`K^d`$ containing the origin define

$`𝔽(K)`$ $`:=`$ $`(XK\mathrm{})`$   
$`𝔾(K)`$ $`:=`$ $`^{0!}(XK\mathrm{})`$   
$`𝕁(K)`$ $`:=`$ $`{\displaystyle \frac{1𝔾(K)}{1𝔽(K)}}.`$   
(Note that the empty space functional bears some relation to the contact distribution function (Stoyan et al., 1995, p. 105) in that $`H_B(r)=𝔽(rB)`$). ¿From these we are able to define the window based $`J`$ function.

## 3 The $`J_W`$ function

###### Definition 2.

For every compact set $`W^d`$ with nonempty interior let

$$J_W(r):=\frac{\underset{W}{}[1𝔾(B(0,r)W_x)]𝑑x}{\underset{W}{}[1𝔽(B(0,r)W_x)]𝑑x}$$
(6)   
be the *window based $`J`$ function*, where $`W_x=\{yx:yW\}`$ is the translate of $`W`$ by $`x^d`$.

If $`X`$ is a stationary Poisson process, then by Slivnyak’s Theorem (Stoyan et al., 1995, (4.4.7)) $`^{0!}`$. Thus $`𝔽𝔾`$ and we arrive at the following proposition.

###### Proposition 1.

Let $`X`$ be a stationary Poisson process. Then

$$J_W(r)1\text{for all}W\text{and}r0.$$

Explicit evaluation of $`J_W`$ for other point process models seems difficult. However, we can show that $`J_W`$ behaves similarly to the $`J`$ function for ordered and clustered processes, suggesting it can also be interpreted in the same way as the $`J`$ function. For some processes it is also possible to demonstrate that the $`J_W`$ function exhibits less deviation from the Poisson hypothesis than the equivalent $`J`$ function.

###### Proposition 2.

Suppose $`X`$ is a process which is “ordered” in the sense that its $`J`$ functional is non-decreasing, that is $`K_1K_2`$ implies $`𝕁(K_1)𝕁(K_2)`$. Then

$$1J_W(r)J(r)\text{for all }r.$$

Similarly if a process is “clustered”, $`K_1K_2`$ implies $`𝕁(K_1)𝕁(K_2)`$, then $`1J_W(r)J(r)`$ for all $`r`$.

Proof: Observe that (6) can be rewritten

$$J_W(r)=_W𝕁(B(0,r)W_x)h_{W,r}(x)𝑑x$$
(7)   
where

$$h_{W,r}(x)=\frac{\left(1𝔽(B(0,r)W_x)\right)}{_W\left(1𝔽(B(0,r)W_y)\right)dy}$$

satisfies $`h_{W,r}(x)0`$ for all $`x^d`$ and $`_Wh_{W,r}(x)𝑑x=1`$. Hence

$$\underset{W}{\mathrm{min}}𝕁(B(0,r)W_x)J_W(r)\underset{W}{\mathrm{max}}𝕁(B(0,r)W_x)$$

and since $`𝕁`$ is nondecreasing

$$1=𝕁(\mathrm{})J_W(r)𝕁(B(0,r))=J(r).\text{ }\mathrm{}$$

The latter result can be strengthened to strict inequality for specific examples.

###### Proposition 3.

Let $`X`$ be a Neyman-Scott cluster process with mean number of points per cluster greater than $`1`$. Assuming the support of the distribution of the cluster points contains a neighbourhood of the origin, then

$$J(r)<J_W(r)<1\text{for all}W\text{and}r>0.$$

Examples of processes satisfying the conditions of Proposition 3 are Matérn’s cluster process and the modified Thomas process described, for example, in Stoyan et al. (1995).

Proof: The Palm distribution of a Neyman-Scott process is the convolution of the original distribution $`P`$ of the process and the Palm distribution $`c_0`$ of the representative cluster $`N`$ (Stoyan et al., 1995, (5.3.2)). Thus, for every compact set $`K`$, we have

$`1𝔾(K)`$ $`=`$ $`{\displaystyle \mathrm{𝟙}\{(\phi \psi )K\{0\}=\mathrm{}\}c_0(d\psi )P(d\phi )}`$   
$`=`$ $`{\displaystyle \mathrm{𝟙}\{\phi K\{0\}=\mathrm{}\}\mathrm{𝟙}\{\psi K\{0\}=\mathrm{}\}c_0(d\psi )P(d\phi )}`$   
$`=`$ $`(XK\{0\}=\mathrm{})c_0(NK\{0\}=\mathrm{})`$   
$`=`$ $`(1𝔽(K))c_0(NK\{0\}=\mathrm{}).`$   
Hence

$$𝕁(K)=c_0(NK\{0\}=\mathrm{}).$$

Now the assumption on the cluster distribution ensures $`c_0(NK\{0\}=\mathrm{})<1`$ for all $`K`$ containing a neighbourhood of the origin. The conclusion then follows since $`B(0,r)W_x`$ contains a neighbourhood of $`0`$ whenever $`x`$ is an interior point of $`W`$.  $`\mathrm{}`$

###### Proposition 4.

Let $`X`$ be a hard-core process with hard-core radius $`R`$. Then

$$1<J_W(r)<J(r)\text{for all}W\text{and}0<r<R,$$

and $`J_W(r)`$ is non-decreasing in $`r`$ for all $`r<R`$.

Proof: Trivially we have $`𝔾(K)=0`$ for all $`KB(0,R)`$, while for any point process the empty space functional $`𝔽`$ is nondecreasing. Therefore the $`J`$ functional becomes

$$𝕁(K)=\frac{1}{1𝔽(K)}\text{for all}KB(0,R),$$

which is also non-decreasing and the result follows by proposition 2.  $`\mathrm{}`$

## 4 Estimation of the $`J_W`$ function

Analogously to $`J`$ we want to estimate $`J_W`$ by the ratio of two estimators for the denominator and numerator in Definition 2. The stationarity of $`X`$ and Fubini’s Theorem yield

$$\underset{W}{}𝔽(B(0,r)W_x)𝑑x=𝔼\underset{W}{}\mathrm{𝟙}\{XB(x,r)W\mathrm{}\}𝑑x,$$

and so the denominator of (6) becomes

$$\underset{W}{}[1𝔽(B(0,r)W_x)]𝑑x=|W|\left[1\frac{1}{|W|}𝔼|W\left((XW)B(0,r)\right)|\right],$$
(8)   
where $``$ denotes Minkowski addition.

Applying the Campbell-Mecke formula (Stoyan et al., 1995, (4.4.3)) we find

$`{\displaystyle \underset{W}{}}𝔾(B(0,r)W_x)𝑑x`$ $`=`$ $`{\displaystyle \underset{W}{}}^{0!}(XB(0,r)W_x\mathrm{})𝑑x`$   
$`=`$ $`{\displaystyle \frac{1}{\lambda }}𝔼{\displaystyle \underset{xXW}{}}\mathrm{𝟙}\{XB(x,r)\{x\}W\mathrm{}\}.`$   
Let $`d(x,A)`$ denote the Euclidean distance from a point $`x^d`$ to a set $`A^d`$. The numerator of (6) can then be expressed as

$$\underset{W}{}[1𝔾(B(0,r)W_x)]𝑑x=|W|[1\frac{1}{\lambda |W|}𝔼\underset{xXW}{}\mathrm{𝟙}\{d(x,XW\{x\})r\}].$$
(9)   
The two results (8) and (9) allow uncorrected estimation of $`J_W(r)`$ by

$$\widehat{J}_W(r):=\frac{1\frac{1}{\mathrm{\#}(XW)}\underset{xXW}{}\mathrm{𝟙}\{d(x,XW\{x\})r\}}{1\frac{1}{|W|}|W\left((XW)B(0,r)\right)|}$$
(10)   
which is the uncorrected estimate of the $`J`$ function referred to in the introduction.

Thus the uncorrected estimate of the $`J`$ function, based on the uncorrected (EDF) estimates of $`F`$ and $`G`$, can be thought of as a ratio unbiased estimator of the $`J_W`$ function. As was shown in the previous section, the $`J_W`$ function can be interpreted in the same way as the $`J`$ function. Consequently, the uncorrected estimate of the $`J`$ function, unlike the uncorrected estimates of $`F`$, $`G`$ or $`K`$, can be used directly as an interpretive statistic in classifying deviations from the Poisson process.

## 5 Simulation study

This section reports the results of simulation studies comparing the uncorrected estimator $`\widehat{J}_W`$ with “corrected” estimators $`\widehat{J}`$. In §5.1 we show the results for a single simulated pattern; §5.2 reports the means and variances of the estimators of $`J`$ in a simulation study. These results show that $`\widehat{J}_W`$ typically has smaller variance than the corrected estimators. In §5.3 and §5.4 we consider the power of hypothesis tests based on the uncorrected estimator $`\widehat{J}_W`$. It is not clear, a priori, whether “corrected” estimators $`\widehat{J}`$ or uncorrected estimators $`\widehat{J}_W`$ will yield more powerful tests. There are two competing effects: the variance of $`\widehat{J}_W`$ is smaller than that of $`\widehat{J}`$, but $`J_W`$ is less sensitive than $`J`$ to departures from the Poisson process (say) according to Proposition 2.

For the comparisons which follow, the reduced sample (border method) estimator $`\widehat{J}_{rs}`$ was adopted. However, for the purpose of highlighting the variations which exist between corrected estimates of the $`J`$ function, the Kaplan-Meier estimator $`\widehat{J}_{km}`$ (described in Baddeley and Gill, 1997) is also considered in many cases.

### 5.1 Empirical example

This example highlights the use of $`\widehat{J}_W`$ as a qualitative summary statistic. Consider the point pattern given in Figure 1. This pattern is a realisation of a Matérn Cluster Process, intensity $`\lambda =100`$, cluster radius $`R=0.1`$ and mean number of offspring $`\mu =4`$, observed within an observation window consisting of two rectangular windows, 3.125 by 0.16 units, separated by 0.02 units. (The Matérn Cluster Process is a Neyman-Scott process in which offspring are uniformly distributed in disc of radius $`R`$ about (Poisson) parent points. The number of offspring per parent point is Poisson with mean $`\mu `$ (Stoyan et al., 1995, p. 159)).

Figure 1 also displays the corresponding estimates $`\widehat{J}_W`$ and $`\widehat{J}_{rs}`$ of the given point pattern, together with envelopes of 99 simulations of a binomial process of the same intensity, observed within the same window. The results in this case are marked; the uncorrected estimate suggests strong evidence of clustering in the point pattern, while the corrected estimate appears to suggest no evidence of clustering. Of course, the results are not surprising given the severity of the edge effect introduced by the “censoring” of the middle seventeenth of the window compared to the relatively small bias this introduces. However, it does illustrate the possible benefit of using $`\widehat{J}_W`$ in certain situations.

As the empirical use of $`J_W`$ is the same as the $`J`$ function, readers interested in further examples of the analysis of empirical data are referred to Kerscher (1998), Kerscher et al. (1999), and Kerscher et al. (1998).

### 5.2 Mean and variance

¿From the previous example it is clear that the uncorrected estimate of the J function may be superior in some situations. To examine whether this was true more generally, a number of simulations were conducted to compare the corrected and uncorrected methods across a range of processes. This began with the estimation of the mean and standard deviation of the three estimators ($`\widehat{J}_W`$, $`\widehat{J}_{km}`$, and $`\widehat{J}_{rs}`$) based on $`10,000`$ realisations of a Poisson process with intensity $`100`$, in a unit square window, with the results presented in Figure 2.

With increasing $`r`$, the distributions of the $`\widehat{J}_.`$ (that is, the estimators $`\widehat{J}_W`$, $`\widehat{J}_{km}`$, and $`\widehat{J}_{rs}`$) become skewed to the right; for large $`r`$ there is substantial mass above $`\widehat{J}_.(r)=2`$. As a result all three estimators are positively biased for large values of $`r`$. Empirically it was found that a square root transformation approximately symmetrised the distribution. As expected, the sample standard deviation of the estimates increases with $`r`$, as the denominator of each estimator decreases with $`r`$. However, $`\widehat{J}_W`$ is less biased and has lower variance than $`\widehat{J}_{km}`$ and $`\widehat{J}_{rs}`$.

These simulations were repeated for two processes with more substantial edge effects (namely, a Poisson process of intensity $`\lambda =25`$ in a unit window, and a Poisson process of intensity $`\lambda =10`$ in a $`10`$ by $`1`$ rectangular window). In both cases the results were qualitatively similar to those above.

In addition, further simulations were conducted for point patterns in $`^3`$. Estimates of the means and standard deviations of $`\widehat{J}_{km}`$ and $`\widehat{J}_{rs}`$ based on $`1000`$ realisations in a unit cube were compared for the Poisson process and two alternatives: Matérn hard-core (Stoyan et al., 1995, p. 163) and Matérn cluster processes for a range of parameter values. Some of the key results are presented in Figure 3. As expected, $`\widehat{J}_W`$ is reliable over a wider domain than $`\widehat{J}_{rs}`$. For hard-core processes, the standard deviation of $`\widehat{J}_{rs}`$ is considerably bigger than that of $`\widehat{J}_W`$ and the difference grows with the hard-core radius. For cluster processes the differences are far less apparent, however the overall tendency of $`\widehat{J}_W`$ to have lower variance is also confirmed for this class of processes. Note also that with both alternative processes the mean of $`\widehat{J}_W`$ is bounded by $`1`$ and $`\widehat{J}_{rs}`$, as expected.

It is interesting to note that, unlike the $`J`$ function estimators, the domain of the $`J_W`$ function estimator for a given point process realisation can be easily calculated. The $`J_W`$ estimator is defined for all $`r<r_{F_{max}}`$, where $`r_{F_{max}}`$ is the maximum nearest-point distance (the largest distance, over all points in the window, from a point to the nearest point of the process). Also $`J_W(r)=0`$ for any $`r_{G_{max}}r<r_{F_{max}}`$ if $`r_{G_{max}}<r_{F_{max}}`$, where $`r_{G_{max}}`$ is the maximum nearest-neighbour distance (the largest distance, over all points of the process within the window, from a point of the process to the nearest other point of the process). The value $`r_{F_{max}}`$ is however an upper bound on the domain of both the Reduced Sample and Kaplan-Meier estimators.

### 5.3 The test statistic

We now aim to compare the power of the $`J_W`$ function with the edge corrected estimators of the $`J`$ function in testing the Poisson hypothesis in the two dimensional case. We restricted ourselves to this estimation problem in view of the problems with estimating the range of interaction using the $`J`$ function reported in Kerscher et al. (1999).

The distribution of the following test statistic for each of the three estimators was estimated:

$$\tau =_0^{r_0}\frac{\widehat{J}_.(r)1}{\widehat{\sigma }(r)}𝑑r,$$
(11)   
where $`\widehat{\sigma }`$ denotes the sample standard deviation of $`\widehat{J}_.(r)`$ under the Poisson hypothesis. This form of test statistic was chosen, as opposed to a squared integrand, because of the skewed nature of the distributions of $`\widehat{J}_.`$. The distributions of the test statistics were estimated by a discrete sum and based on $`10,000`$ realisations of a Poisson process. The upper limit of integration $`r_0`$ was chosen to be the 0.9 quantile of the $`F`$ function (for intensity $`100`$, $`r_00.856`$). Having estimated the distribution, the $`0.025`$ and $`0.975`$ quantiles were obtained for use in a two-sided $`5\%`$ significance test for deviation from a Poisson process. One-sided $`5\%`$ significance tests were also constructed to test for clustering or regularity by considering the $`0.05`$ and $`0.95`$ quantiles, respectively.

### 5.4 Power of tests using the various estimators

In order to estimate the power of the hypothesis test described above, realisations from alternative point processes were generated and the proportion of the realisations rejected by the hypothesis test was recorded. The first class of point processes considered was the Matérn hard-core process, with hard-core radius $`R`$. For each of $`22`$ values of $`R`$, $`1000`$ realisations were generated. The proportion of rejections is presented in Figure 4. Note that as $`R0`$ the model approaches the Poisson process, so we expect all power curves to approach $`0.05`$ $`(5\%)`$ as $`R0`$. All three estimators have very similar power curves, with the $`J_W`$ estimator at least as powerful as the two $`J`$ function estimators for all values of $`R`$.

The other class of alternative point processes considered was Matérn’s cluster process. A grid of $`(R,\mu )`$ values was constructed and $`1000`$ realisations were obtained from the corresponding Matérn cluster process. The proportion of rejections for each $`(R,\mu )`$ value is presented in Figure 5. Once again, the curves are very similar, with all three tests performing almost identically. Similar results were obtained for the respective one-sided tests in both cases as well as for the lower intensity $`25`$.

The power tests for the $`10`$ by $`1`$ window support the argument that edge effects are stronger when the boundary is relatively longer. The resulting power function estimates against the Matérn Model II and Matérn cluster process models are also presented in Figures 4 and 5, respectively.

One important observation made while conducting these numerical simulations was that the choice of test statistic had far more impact on the power of the resulting hypothesis test than the choice of $`J`$ function estimator. For a comparison of various test statistics of the $`J`$ function see Thönnes and van Lieshout (1998).