A toolbox for fitting complex spatial point process models using integrated nested Laplace approximation (INLA)

A toolbox for fitting complex spatial point process models using integrated nested Laplace approximation (INLA)

[ [    [    [ University of St Andrews, University of Tromsø and Norwegian University of Science and Technology J. B. Illian
Centre for Research into Ecological
 and Environmental Modelling
The Observatory
University of St Andrews
St Andrews KY16 9LZ
S. H. Sørbye
Department of Mathematics
 and Statistics
University of Tromsø
9037 Tromsø
H. Rue
Department of Mathematical Sciences
Norwegian University of Science
 and Technology
7491 Trondheim
\smonth10 \syear2011\smonth11 \syear2011
\smonth10 \syear2011\smonth11 \syear2011
\smonth10 \syear2011\smonth11 \syear2011

This paper develops methodology that provides a toolbox for routinely fitting complex models to realistic spatial point pattern data. We consider models that are based on log-Gaussian Cox processes and include local interaction in these by considering constructed covariates. This enables us to use integrated nested Laplace approximation and to considerably speed up the inferential task. In addition, methods for model comparison and model assessment facilitate the modelling process. The performance of the approach is assessed in a simulation study. To demonstrate the versatility of the approach, models are fitted to two rather different examples, a large rainforest data set with covariates and a point pattern with multiple marks.


10.1214/11-AOAS530 \volume6 \issue4 2012 \firstpage1499 \lastpage1530


Fitting complex spatial point process models with INLA


A]\fnmsJanine B. \snmIllian\correflabel=e1]janine@mcs.st-and.ac.uk, B]\fnmsSigrunn H. \snmSørbye and C]\fnmsHåvard \snmRue

Cox processes \kwdmarked point patterns \kwdmodel assessment \kwdmodel comparison.

1 Introduction

1.1 Complex point process models

These days a large variety of complex statistical models can be fitted routinely to complex data sets as a result of widely accessible high-level statistical software, such as R [R Development Core Team (2009)] or winbugs [Lunn et al. (2000)]. For instance, the nonspecialist user can estimate parameters in generalized linear mixed models or run a Gibbs sampler to fit a model in a Bayesian setting, and expert programming skills are no longer required. Researchers from many different disciplines are now able to analyze their data with sufficiently complex methods rather than resorting to simpler yet nonappropriate methods. In addition, methods for the assessment of a model’s fit as well as for the comparison of different models are widely used in practical applications.

The routine fitting of spatial point process models to complex data sets, however, is still in its infancy. This is despite a rapidly improving technology that facilitates data collection, and a growing awareness of the importance and relevance of small-scale spatial information. Spatially explicit data sets have become increasingly available in many areas of science, including plant ecology [Burslem, Garwood and Thomas (2001); Law et al. (2001)], animal ecology [Forchhammer and Boomsma (1995, 1998)], geosciences [Naylor et al. (2009); Ogata (1999)], molecular genetics [Hardy and Vekemans (2002)], evolution [Johnson and Boerlijst (2002)] and game theory [Killingback and Doebeli (1996)], with the aim of answering a similarly broad range of scientific questions. Currently, these data sets are often analyzed with methods that do not make full use of the available spatially explicit information. Hence, there is a need for making existing point process methodology available to applied scientists by facilitating the fitting of suitable models.

In addition, real data sets are often more complex than the classical data sets that have been analyzed with point process methodology in the past. They often consist of the exact spatial locations of the objects or events of interest, and of further information on these objects, that is, potentially dependent qualitative as well as quantitative marks or spatial covariates [Burslem, Garwood and Thomas (2001); Moore et al. (2010)]. There is an interest in fitting complex joint models to the marks (or the covariates) as well as to the point pattern. So far, the statistical literature has discussed few examples of complex point process models of this type.

There have been previous advances in facilitating routine model fitting for spatial point processes, in particular, for Gibbs processes. Most markedly, the work by Baddeley and Turner (2000) has facilitated the routine fitting of Gibbs point processes based on an approximation of the pseudolikelihood to avoid the issue of intractable normalizing constants [Berman and Turner (1992); Lawson (1992)] as well as the approximate likelihood approach by Huang and Ogata (1999). Work by Baddeley et al. (2005) and Stoyan and Grabarnik (1991) has provided methods for model assessment for some Gibbs processes. Many of these have been made readily available through the library spatstat for R [Baddeley and Turner (2005)].

However, most Gibbs process models considered in the literature are relatively simple in comparison to models that are commonly used in the context of other types of data. In an attempt to generalize the approach in Baddeley and Turner (2005), Illian and Hendrichsen (2010) include random effects in Gibbs point processes but more complex models, such as hierarchical models or models including quantitative marks, currently cannot be fitted in this framework. Similarly, methods for model comparison or assessment considered in Baddeley et al. (2005) and Stoyan and Grabarnik (1991) are restricted to relatively simple models. Furthermore, both estimation based on maximum likelihood and that based on pseudolikelihood are approximate so that inference is not straightforward. The approximations become less reliable with increasing interaction strength [Baddeley and Turner (2000)].

Cox processes are another, flexible, class of spatial point process models [Møller and Waagepetersen (2007)], assuming a stochastic spatial trend makes them particularly realistic and relevant in applications. Even though many theoretical results have been discussed in the literature for these [Møller and Waagepetersen (2004)], the practical fitting of Cox point process models to point pattern data remains difficult due to intractable likelihoods. Fitting a Cox process to data is often based on Markov chain Monte Carlo (MCMC) methods. These require expert programming skills and can be very time-consuming both to tune and to run [Møller and Waagepetersen (2004)] so that fitting complex models can easily become computationally prohibitive. For simple models, fast minimum contrast approaches to parameter estimation have been discussed [Møller and Waagepetersen (2007)].

However, approaches to routinely fitting Cox process models have been discussed very little in the literature; similarly, methods for model comparison or assessment for Cox processes have rarely been discussed in the literature [Illian and Rue (2010); Illian et al. (2012)]. To the authors’ knowledge, Cox processes have not been used outside the statistical literature to answer concrete scientific questions. Within the statistical literature Cox process models have focused on the analysis of relatively small spatial patterns in terms of the locations of individual species. Very few attempts have been made at fitting models to both the pattern and the marks [Ho and Stoyan (2008); Myllymäki and Penttinen (2009)], in particular, not to patterns with multiple dependent continuous marks, and joint models of covariates and patterns have not been considered.

This paper addresses two issues. It develops complex joint models and, at the same time, provides methods facilitating the routine fitting of these models. This provides a toolbox that allows applied researchers to appropriately analyze realistic point pattern data sets. We consider joint models of both the spatial pattern and associated marks as well as of the spatial pattern and covariates. Using a Bayesian approach, we provide modern model fitting methodology for complex spatial point pattern data similar to what is common in other areas of statistics and has become a standard in many areas of application, including methods for model comparison and validation. The approach is based on integrated nested Laplace approximation (INLA) [Rue, Martino and Chopin (2009)], which speeds up parameter estimation substantially so that Cox processes can be fitted within feasible time. In order to make the methods accessible to nonspecialists, an R package that may be used to run INLA is available and contains generic functions for fitting spatial point process models; see http://www.r-inla.org/.

1.2 Cox processes with local spatial structure

Applied researchers are aware that spatial behavior tends to vary at a number of spatial scales as a result of different underlying mechanisms that drive the pattern [Wiegand et al. (2007); Latimer et al. (2009)]. Local spatial behavior is often of specific interest but the spatial structure also varies on a larger spatial scale due to the influence of observed or unobserved spatial covariates. Cox processes model spatial patterns relative to observed or unobserved spatial trends and would be ideal models for these data sets.

However, Cox processes typically do not consider spatial structures at different spatial scales within the same model. More specifically, a specific strength of spatial point process models is their ability to take into account detailed information at very small spatial scales contained in spatial point pattern data, in terms of the local structure formed by an individual and its neighbors. So far, Cox processes have often been used to relate the locations of individuals to environmental variation, phenomena that typically operate on larger spatial scales. However, different mechanisms operate at a smaller spatial scale. Spatial point data sets are often collected with a specific interest in the local behavior of individuals, such as spatial interaction or local clustering [Law et al. (2001); Latimer et al. (2009)].

We consider an approach to fitting Cox process models that reflects both the local spatial structure and spatial behavior at a larger spatial scale by using a constructed covariate together with spatial effects that account for spatial behavior at different spatial scales. This approach is assessed in a simulation study and we also discuss issues specific to this approach that arise when several spatial scales are accounted for in a model.

This paper is structured as follows. The general methodology is introduced in Section 2. In Section 3 we investigate the idea of mimicking local spatial behavior by using constructed covariates in a simulation study in the context of (artificial) data with known spatial structures and inspect patterns resulting from the fitted models. Section 4 discusses a joint model of a large point pattern and two empirical covariates along with a constructed covariate and fits this to a rainforest data set. A hierarchical approach is considered in Section 5, where both (multiple) marks and the underlying pattern are included in a joint model and fitted to a data set of eucalyptus trees and koalas foraging on these trees.

2 Methods

2.1 Spatial point process models

Spatial point processes have been discussed in detail in the literature; see Stoyan, Kendall and Mecke (1995), van Lieshout (2000), Diggle (2003), Møller and Waagepetersen (2004, 2007) and Illian et al. (2008). Here we aim at modeling a spatial point pattern , regarding it as a realization from a spatial point process . For simplicity we consider only point processes in , but the approaches can be generalized to point patterns in higher dimensions.

We refer the reader to the literature for information on different (classes of) spatial point process models such as the simple Poisson process, the standard null model of complete spatial randomness, as well as the rich class of Gibbs (or Markov) processes [van Lieshout (2000)]. Here, we discuss the class of Cox processes, in particular, log-Gaussian Cox processes. Cox processes lend themselves well to modeling spatial point pattern data with spatially varying environmental conditions [Møller and Waagepetersen (2007)], as they model spatial patterns based on an underlying (or latent) random field that describes the random intensity, assuming independence given this field. In other words, given the random field, the point pattern forms a Poisson process. Log-Gaussian Cox processes as considered, for example, in Møller, Syversveen and Waagepetersen (1998) and Møller and Waagepetersen (2004, 2007), are a particularly flexible class, where has the form and is a Gaussian random field, . Other examples of Cox processes include shot-noise Cox processes [Møller and Waagepetersen (2004)].

Here, we consider a general class of complex spatial point process models based on log-Gaussian Cox processes that allows the joint modeling of spatial patterns along with marks and covariates. We include both small and larger scale spatial behavior, using a constructed covariate and additional spatial effects. The resulting models can be regarded as latent Gaussian models and, hence, INLA can be used for parameter estimation and model fitting.

2.2 Integrated nested Laplace approximation (INLA)

Cox processes are a special case of the very general class of latent Gaussian models, models of an outcome variable that assume independence conditional on some underlying latent field and hyperparameters . Rue, Martino and Chopin (2009) show that if has a sparse precision matrix and the number of hyperparameters is small (i.e., 7), inference based on INLA is fast.

The main aim of the INLA approach is to approximate the posteriors of interest, that is, the marginal posteriors for the latent field and the marginal posteriors for the hyperparameters , and use these to calculate posterior means, variances, etc. These posteriors can be written as


The nested formulation is used to compute by approximating and , and then to use numerical integration to integrate out . This is feasible, since the dimension of is small. Similarly, is calculated by approximating and integrating out .

The marginal posterior in equations (1) and (2) can be calculated using the Laplace approximation

where is the Gaussian approximation to the full conditional of , and is the mode of the full conditional for , for a given . This makes sense, since the full conditional of a zero mean Gauss Markov random field can often be well approximated by a Gaussian distribution by matching the mode and the curvature at the mode [Rue and Held (2005)]. Further details are given in Rue, Martino and Chopin (2009) who show that the nested approach yields a very accurate approximation if applied to latent Gaussian models. As a result, the time required for fitting these models is substantially reduced.

2.3 Fitting log-Gaussian Cox processes with INLA

The class of latent Gaussian models comprises log-Gaussian Cox processes and, hence, the INLA approach may be applied to fit these. Specifically, the observation window is discretized into grid cells , each with area , . The points in the pattern can then be described by with , where denotes the observed number of points in grid cell . We condition on the point pattern and, conditionally on , we have


see Rue, Martino and Chopin (2009).

We model as


where the functions are spatially structured effects that reflect large scale spatial variation in the pattern. These effects are modeled using a second-order random walk on a lattice, using vague gamma priors for the hyperparameter and constrained to sum to zero [Rue and Held (2005)]. In the models that we discuss below, the spatially structured effects relate to observed and unobserved spatial covariates as discussed in the examples in Sections 4 and 5. Including spatial covariates directly in the model as fixed effects in addition to the random effects is straightforward. For simplicity, we omit these in equation (4) since this is not relevant in the specific data sets and models discussed below. denotes a spatially unstructured zero-mean Gaussian i.i.d. error term, using a gamma prior for the precision.

Further, denotes a constructed covariate. Constructed covariates are summary characteristics defined for any location in the observation window reflecting inter-individual spatial behavior such as local interaction or competition. We assume that this behavior operates at a smaller spatial scale than spatial aggregation due to (observed or unobserved) spatial covariates, and hence the spatially structured effects. The use of constructed covariates yields models with local spatial interaction within the flexible class of log-Gaussian Cox process models. It avoids issues with intractable normalizing constants that are common in the context of Gibbs processes [Møller and Waagepetersen (2004)], since the covariates operate directly on the intensity of the pattern rather than on the density or the conditional intensity [Schoenberg (2005)].

The functional relationship between the outcome variable and the constructed covariate is typically not obvious and might often not be linear. We thus estimate this relationship explicitly by a smooth function and inspect this estimate to gain further information on the form of the spatial dependence. This function will be modeled as a first-order random walk, also constrained to sum to zero.

The constructed covariate considered in this paper is based on the nearest point distance, which is simple and fast to compute. Specifically, for each center point of the grid cells we find the distance to the nearest point in the pattern outside this grid cell as


where denotes the center point of cell and denotes the Euclidean distance. Defined this way, the constructed covariate can be used both to model local repulsion and local clustering.

During the modeling process, methods for model comparison based on the deviance information criterion (DIC) [Spiegelhalter et al. (2002)] may be used to compare different models with different levels of complexity. Furthermore, both the (estimated) spatially structured field and the error field in (4) may be used to assess the model fit. The spatially structured effect may be used to reveal remaining spatial structure that is unexplained by the current model and the unstructured effects may be interpreted as a spatial residual. This provides a method for model assessment akin to residuals in, for example, linear models.

This approach yields a toolbox for fitting, comparing and assessing realistically complex models of spatial point pattern data. We show that different types of flexible models can be fitted to point pattern data with complex structures using the INLA approach within reasonable computation time. This includes joint models of large point patterns and covariates operating on a large spatial scale and local clustering (Section 4) as well as of a pattern with several dependent marks which also depend on the pattern (Section 5).

2.4 Issues of spatial scale

In the natural world, different mechanisms operate at different spatial scales [Steffan-Dewenter et al. (2002)], and hence are reflected in a spatial pattern at these scales. It is crucial to bear this in mind during the analysis of spatial data derived from nature, including spatial point pattern data. Some mechanisms, such as seed dispersal in plants or territorial behavior in animals, may operate at a local spatial scale, while others, such as aggregation resulting from an association with certain environmental covariates, operate on the scale of the variation in these covariates, and hence often on a larger spatial scale. In addition, a spatial scale that is relevant in one application may not be relevant for a different data set. Hence, the analysis of a spatial point pattern always involves a consideration of the appropriate spatial scales at which mechanisms of interest may operate, regardless of the concrete analysis methods. Even as early as at the outset of a study, when an appropriately sized observation window has to be chosen, relevant spatial scales operating in the system of interest have to be taken into consideration.

During the analysis the researcher has to carefully decide if variation at a specific scale constitutes noise or whether it reflects a true signal. It is hence crucial to be aware of which mechanisms operate at which spatial scales prior to any spatial data analysis. This may be done based on either background knowledge (such as existing data on dispersal distances in plants or the sizes of home ranges in territorial animals) or common sense.

In the models we discuss here, we explicitly take mechanisms operating at several different scales into account and have to choose these sensibly, based on knowledge of the systems. The spatially structured effect reflects spatial autocorrelation at a large spatial scale, whereas the constructed covariate is used to describe small scale inter-individual behavior. In addition, since we grid the data in this approach, the number of grid cells clearly determines the spatial resolution, especially at a small scale, and is clearly linked to computational costs and the extent to which information is lost through gridding the data. In the following, we discuss issues related to each of these three parts of the models where spatial scale is relevant.

A spatially structured effect is typically included in a spatial model as a spatially structured error term, that is, in order to account for any spatial autocorrelation unexplained by covariates in the model. INLA currently supports the 2nd order random walk on a lattice as a model for this, with a gamma prior for the variance of the spatially structured effect. The choice of this prior determines the smoothness of the spatial effect and through this, the spatial scale at which it operates. This prior has to be chosen carefully to avoid overfitting. This is particularly crucial in the context of spatial point patterns with relatively small numbers of points, where the gridded data are typically rather sparse [Illian et al. (2012)]. If the spatial effect is chosen to be too coarse, it explains the spatial variation at too small a scale, resulting in a coarse estimate of the spatially structured effect. This estimate would perfectly explain every single data point, resulting in overfitting rather than in a model of a generally interpretable trend. Given the role of the spatially structured effect, it appears plausible to choose the prior so that the spatial effect operates at a similar spatial scale as the covariate. Problems can occur when the spatially structured effect operates at a smaller scale than the covariate, as it is then likely to explain the data better than the covariates, rendering the model rather useless. In the absence of covariate data, background knowledge on spatial scales may aid in choosing the prior.

Small scale inter-individual spatial behavior is modeled by the constructed covariate. As mentioned, this is done to account for local spatial behavior if this is of specific interest in the application. Again, there is a danger of overfitting, especially since the constructed covariate is estimated directly from the data. We discuss the practicality of using a spatial constructed covariate in detail in Section 3 and only point out here that it has to be carefully chosen, if possible with appropriate knowledge of the specific system the data have been derived from.

The choice of prior for the spatially structured effect is strongly related to the choice of grid size. However, in our experience the overall results often do not change substantially when the grid size was varied within reason. In applications, the locations of the modeled objects as well as spatial covariates are sometimes given on a grid with a fixed resolution. We recommend using a grid that is not finer than that given by the data in the analysis.

3 Using a constructed covariate to account for local spatial structure—a simulation study

In Section 4 we use a constructed covariate primarily to incorporate local spatial structure into a model, while accounting for spatial variation at a larger spatial scale. To illustrate the use of the given constructed covariate and to assess the performance of the resulting models, we simulate point patterns from various classical point-process models. Note, however, that we do not aim at explicitly estimating the parameters of these models but at assessing (i) whether known spatial structures may be detected through the use of the constructed covariate, as suggested here, and (ii) whether simulations from the fitted models generate patterns with similar characteristics. In the applications we have in mind, such as those discussed in the example in Section 4, the data structure is typically more complicated.

For the purpose of this simulation study we consider three different situations: patterns with local repulsion (Section 3.1), patterns with local clustering (Section 3.2) and patterns with local clustering in the presence of a larger-scale spatial trend (Section 3.3). We generate example patterns from different point process models with these properties on the unit square. For all simulation results this observation window has been discretized into a grid.

In Sections 3.1 and 3.2 we initially assume that there is no large-scale spatial variation, with the aim of inspecting only the constructed covariate, and we consider


using the notation in Section 2.3. In Section 3.3 we consider both small- and large-scale spatial structures by including a spatially structured effect in addition to the constructed covariate and


To evaluate a fitted model, we apply the Metropolis algorithm [Metropolis et al. (1953)] to simulate patterns from these models and then compare characteristics of the simulated patterns with the generated example patterns. More specifically, for and , denote the joint distribution of given the latent field , by

where the mean . For a given example pattern, we first apply INLA to find the estimate of the latent field for all grid cells. To evaluate the estimated function of the constructed covariate for all arguments, we apply the splinefun command in R to perform cubic spline interpolation of the original data points. Using the Metropolis algorithm, we assume an initial pattern , which is randomly scattered in the unit square, having the same number of points as the original pattern. The th step of the algorithm is performed by randomly selecting one point of the pattern and proposing to move this point to a new position drawn uniformly in the unit square. The proposal is accepted with probability

where denotes the resulting grid cell counts for . The simulated patterns in Sections 3.13.3 each result from iterations of the algorithm.

(a) (b)
(c) (d)
(e) (f)
Figure 1: Simulated Strauss process with medium repulsion (, , ) (a), the associated constructed covariate for this pattern (b), the estimated functional relationship between the outcome and the constructed covariate (c), a pattern simulated from the fitted model after 100,000 iterations (d), the estimated -function for the original pattern (solid line) and for the simulated pattern (dashed line) (e) and simulation envelopes for the -function for 50 simulated patterns (f).

3.1 Modeling repulsion

To inspect the performance of the constructed covariate for repulsion, we generate patterns from a homogeneous Strauss process [Strauss (1975)] on the unit square, with medium repulsion (intensity parameter), (interaction parameter) and interaction radius [see Figure 1(a) for an example]. We then fit a model to the pattern as in equation (6) using the constructed covariate in (5) [Figure 1(b)]. The shape of the estimated functional relationship between the constructed covariate and the outcome variable is shown in Figure 1(c). This function illustrates that the intensity in a grid cell is influenced by the calculated distance in (5), as higher distances will give higher intensities. Thus, the intensity is positively related to the value of the constructed covariate, clearly reflecting repulsion. At larger distances (0.05) the function levels out distinctly, indicating that beyond these distances the covariate and the intensity are unrelated, that is, the spatial pattern shows random behavior. In other words, the functional relationship not only characterizes the pattern as regular but also correctly identifies the interaction distance as 0.05.

The pattern resulting from the Metropolis–Hastings algorithm [Figure 1(d)] shows very similar characteristics to those in the original pattern. This indicates that the model based on the nearest point constructed covariate in equation (5) captures adequately the spatial information contained in the repulsive pattern.

The estimated -function [Besag (1977)] for the simulated pattern and the original pattern confirm this impression, as they look very similar [Figure 1(e)]. Additionally, we have calculated simulation envelopes for the -function of Strauss processes with the given parameter values, using 50 simulated patterns and 100,000 iterations of the Metropolis algorithm for each pattern [Figure 1(f)]. We notice that the estimated -functions of the original patterns are well within the simulation envelopes for all distances.

3.2 Modeling clustering

In order to assess the performance of the model in (6) in the context of clustered patterns, we generate patterns from a homogeneous Thomas process [Neyman and Scott (1952)] in the unit square, with parameters (the intensity of the Poisson process of cluster centers), (the standard deviation of the distance of a process point from the cluster center) and (the expected number of points per cluster) [see Figure 2(a) for an example]. We fit the model in equation (6) using the constructed covariate in (5) [Figure 2(b)]. The shape of the estimated functional relationship between the constructed covariate and the outcome variable [Figure 2(c)] now indicates that the intensity is negatively related to the value of the constructed covariate as the intensities increase for smaller distances, reflecting local clustering. At larger distances (0.1) the function levels out, indicating that at these distances the covariate and the intensity are unrelated.

(a) (b)
(c) (d)
(e) (f)
Figure 2: Simulated Thomas process with parameters , and (a), the associated constructed covariate for this pattern (b), the estimated functional relationship between the outcome and the constructed covariate (c), a pattern simulated from the fitted model after 100,000 iterations (d), the estimated L-function for the original pattern (solid line) and the simulated pattern (dashed line) (e) and simulation envelopes for the L-function for 50 simulated patterns (f).

The pattern simulated from the fitted model [Figure 2(d)] shows that the constructed covariate introduces some clustering in the model. However, the resulting pattern shows fewer and less distinct clusters than the original pattern. Similarly, the estimated -function for the pattern simulated from the fitted model shows a weaker local clustering effect than the original pattern [Figure 2(e)]. This is also illustrated by the simulation envelopes for 50 patterns of the fitted model which do not include the true -function [Figure 2(f)].

3.3 Modeling small scale clustering in the presence of large-scale inhomogeneity

So far, we have considered constructed covariates only for patterns with local interaction to illustrate their use. In applications, however, different mechanisms operate at different spatial scales. Patterns may be locally clustered, for example, due to dispersal mechanisms, but may also show aggregation at a larger spatial scale, for example, due to dependence on underlying observed or unobserved covariates. Hence, the main reason for using constructed covariates in the data example in Section 4 is to distinguish behavior at different spatial resolutions, in order to provide information on mechanisms operating at different spatial scales.

We illustrate the use of constructed covariates in this context by generating an inhomogeneous, locally clustered pattern mimicking a situation where different mechanisms have caused local clustering and large scale inhomogeneity. In applications, the inhomogeneity may be modeled using suitable spatially varying covariates or assuming an unobserved spatial variation or both. We generate patterns from an inhomogeneous Thomas process with parameters and and a simple trend function for the intensity of parent points given by . Each pattern is then superimposed with a pattern generated from an inhomogeneous Poisson process with trend function [Figure 3(a)].

(a) (b)
(c) (d)
(e) (f)
Figure 3: Realization of an inhomogeneous Thomas process with parameters , and trend function superimposed on an inhomogeneous Poisson process with trend function (a), the associated constructed covariate for this pattern (b), the estimated functional relationship between the outcome and the constructed covariate (c), the estimated spatially structured effect (d), a pattern simulated from the fitted model after 100,000 iterations (e) and the inhomogeneous L-function for the original pattern (solid line) and the simulated pattern (dashed line) (f).

We again use the constructed covariate in (5) [see Figure 3(b)] and fit the model in (7). The inspection of the functional relationship between the constructed covariate and the outcome [Figure 3(c)] shows that at small values of the covariate the intensity is negatively related to the constructed covariate, reflecting clustering at smaller distances. The estimated spatially structured effect picks up the larger-scale spatial behavior [Figure 3(d)]. Patterns simulated from the fitted model look quite similar to the original pattern [Figure 3(e)]. However, local clustering is slightly stronger in the original pattern than in the simulated pattern [Figure 3(f)].

This is again confirmed by the simulation envelopes for the simulated patterns from the fitted model, as shown in Figure 4. The mean estimated L-function for the generated patterns is very close to the upper edge of the simulation envelopes and partly outside, indicating that the fitted model does not reflect the strength of clustering sufficiently well.

Figure 4: Inhomogeneous Thomas process. Simulation envelopes for 50 patterns generated from the fitted model using 100,000 iterations, the inhomogeneous -function for a Poisson process (bold solid line), the mean of the inhomogeneous -function for the generated (solid) and simulated (dashed) patterns.

3.4 Discussion on constructed covariates

With the aim of assessing the performance of models with constructed covariates reflecting small scale inter-individual spatial behavior, we consider a number of simulated point patterns for three different scenarios: repulsion, clustering and small-scale clustering in the presence of large scale inhomogeneity. In all cases, the local spatial structure can be clearly identified. The constructed covariate does not only take account of local spatial structures but also characterizes the spatial behavior. The functional form of the dependence of the intensity on the constructed covariate clearly reflects the character of the local behavior.

This section presents only a small part of an extensive simulation study; the results shown here are typical examples. We have run simulations from the same models as above with different sets of parameters and have obtained essentially the same results. Further, fitting the model in equation (6) to patterns simulated from a homogeneous Poisson process resulted in a nonsignificant functional relationship, that is, the modeling approach does not pick up spurious clustering or regularity.

The approach allows us to fit models that take into account small-scale spatial behavior, regularity as well as clustering, in the context of log-Gaussian Cox processes, that is, as latent Gaussian models. Since these can be fitted using the INLA approach, fitting is fast and exact. In addition, we avoid some of the typical problems that arise with Gibbs process models, that is, we do not face issues of intractable normalizing constants, and regular as well as clustered patterns may be modeled.

However, the simulations also show that the approach of using constructed covariates works clearly better with repulsive patterns than with clustered patterns. This is akin to similar issues with Gibbs processes, where repulsive patterns are less problematic to model than clustered patterns. Certainly, this is related to the fact that it is difficult to tell apart clustering from inhomogeneity [Diggle (2003)]. When working with constructed covariates the issues highlighted, that is, that local clustering may have been underestimated, have to be taken into account, especially in the interpretation of results.

Certainly, the constructed covariate in equation (5) that we consider here is not the only possible choice. A covariate based on distance to the nearest point is likely to be rather “short-sighted,” so that other constructed covariates might be more suitable for detecting specific spatial structures. In particular, taking into account these limitations, it is not surprising that patterns simulated from models show less clustering than the original data. More general covariates such as the distance to the th nearest point may be considered. Other covariates, such as the local intensity or the number of points within a fixed interaction radius from a location , are certainly also suitable. A nice property of the given constructed covariate based on nearest-point distance is that it is parameter-free. For this reason, it is not necessary to choose explicitly the resolution of the local spatial behavior, for example, as an interaction radius. Also, note that since the distance to the nearest point in point pattern for a location may be interpreted as a graph associated with , other constructed covariates based on different types of graphs [Rajala and Illian (2012)] may also be used as constructed covariates. Similarly, an approach based on morphological functions may be used for this purpose. Note that one could also consider constructed marks based on first or second order summary characteristics [Illian et al. (2008)] that are defined only for the points in the pattern and include these in the model.

Distinguishing spatial behavior at different spatial scales is clearly an ill-posed problem, since the behavior at one spatial scale is not independent of that at different spatial scales [Diggle (2003)]. The approach we take here will not always be able to distinguish clustering at different scales. However, different mechanisms that operate at very similar spatial scales are likely to be nonidentifiable by any method, irrespective of the choice of model or the constructed covariate. Constructed covariates hence only provide useful results when the processes they are meant to describe operate at a spatial scale that is distinctly smaller than the larger scale processes in the same model.

Admittedly, the use of constructed covariates is of a rather subjective and ad hoc nature. Clearly, in applications the covariates have to be constructed carefully, depending on the questions of interest; different types of constructed covariates may be suitable in different contexts. However, similarly subjective decisions are usually made when a model is fitted that is purely based on empirical covariates, as these have been specifically chosen as potentially influencing the outcome variable, based on background knowledge. In addition, due to the apparent danger of overfitting, constructed covariates should only be used if there is an interest in the local spatial behavior in a specific data set and if there is reason to believe that small- and large-scale spatial behavior are operating at scales that are different enough to make them identifiable.

4 Joint model of a point pattern and environmental covariates

4.1 Modeling approach

In this example we consider a point pattern , where the number of points is potentially very large and several spatial covariates have been measured. The point pattern is assumed to depend on one or several (observed or unobserved) environmental covariates for which data exist. In the application that we have in mind the values of these have been observed in a few locations that are typically different from locations of the objects that form the pattern. In previous modeling attempts the values of the covariates in the locations of the objects are then either interpolated or modeled separately so that (estimated) values are used for locations were the covariates have not been observed. However, these covariates are likely to have been collected with both sampling and measurement error. In the specific case we consider here (see Section 4.2) they concern soil properties, which are measured much less reliably than the topography covariates in models such as those in Waagepetersen (2007), Waagepetersen and Guan (2009). In addition, it is less clear for soil variables than for topography covariates if these influence the presence of trees, or whether the presence of trees impacts on the soil variables. Whereas models in which the soil variables are considered fixed and not modeled alongside the pattern, the model we deal with here does not make any assumption on the direction of this influence.

As a result, we suggest a joint model of the covariates along with the pattern that uses the original (noninterpolated) data on the covariates and accounts for measurement error. That is, we fit the model in equation (4) to and jointly fit a model to the covariates. The pattern and the covariates are linked by joint spatial fields. An additional spatially structured effect is used to detect any remaining spatial structures in the pattern that cannot be explained by the joint fields with the covariates.

In the case of we fit the following model, where the pattern is modeled as


and the covariates as




where and are the observed covariates in grid cells where the covariates have been measured and missing where they have not been measured. represents the function of the constructed covariate (5). and are spatially structured effects, that is, reflect a random field for each of the covariates and reflects spatial autocorrelation in the pattern unexplained by the covariates; and are spatially unstructured fields used to account for measurement or sampling error.

In addition to the spatial effect reflecting the empirical covariates, which are likely to have an impact on the larger scale spatial behavior, we use the constructed covariate to account for local clustering. In the application we have in mind (see Section 4.2) this clustering is a result of seed-dispersal mechanisms operating on a much smaller spatial scale than that of the aggregation of individuals due to an association with environmental covariates.

4.2 Application to example data set

4.2.1 The rainforest data

Some extraordinarily detailed multi-species maps are being collected in tropical forests as part of an international effort to gain greater understanding of these ecosystems [Condit (1998); Hubbell et al. (1999); Burslem, Garwood and Thomas (2001); Hubbell, Condit and Foster (2005)]. These data comprise the locations of all trees with diameters at breast height (dbh) 1 cm or greater, a measure of the size of the trees (dbh), and the species identity of the trees. The data usually amount to several hundred thousand trees in large (25 ha or 50 ha) plots that have not been subject to any sustained disturbance such as logging. The spatial distribution of these trees is likely to be determined by both spatially varying environmental conditions and local dispersal.

Recently, spatial point process methodology has been applied to analyze some of these data sets [Law et al. (2009); Wiegand et al. (2007)] using nonparametric descriptive methods as well as explicit models [Waagepetersen (2007); Guan (2008); Waagepetersen and Guan (2009); Yue and Loh (2011)]. Rue, Martino and Chopin (2009) model the spatial pattern formed by a tropical rain forest tree species on the underlying environmental conditions and use the INLA approach to fit the model.

We analyze a data set that is similar to those discussed in the above references. Since the spatial structure in a forest reflects dispersal mechanisms as well as association with environmental conditions, we include a constructed covariate to account for local clustering. The model is fitted to a data set from a 50 ha forest dynamics plot at Pasoh Forest Reserve, Peninsular Malaysia. This study focuses on the species Aporusa microstachya consisting of 7416 individuals [Figure 5(a)]. The environmental covariates have been observed in 83 locations that are distinct from the locations of the trees [Figure 5(b)]. The plot lies in a forest that has never been logged with very narrow streams on almost flat land. The data collected in 1995 are used here when the plot contained 320,903 stems from 817 species. The species is the most common small tree on the plot. It is of interest if this species, as an aluminium accumulator, covaries with magnesium availability, as aluminium uptake might constrain its capacity to take up nutrient cations such as magnesium. In addition, its covariation with phosphorus is considered here as the element thought to be the nutrient primarily limiting forest productivity and individual tree growth in tropical forests [Burslem, personal communication (February 2011)].

(a) (b)
Figure 5: Spatial pattern of the species Aporusa microstachya in Pasoh Forest Reserve, Peninsular Malaysia and locations where soil variables have been measured.
(a) (b)
(c) (d)
(e) (f)
Figure 6: Rainforest data. Top panels: The estimated spatially structured effect for the covariates phosphorus (a) and magnesium (b). Middle panels: The calculated constructed covariate (c) and the estimated function of the constructed covariate (d). Bottom panels: The estimated spatially structured effect for the pattern with (e) and without the constructed covariate term in the model (f).

4.2.2 Results

We run the full model as described in equations (8) to (10), in which the observation area is discretized into grid cells. The spatial effect of the two empirical covariates, phosphorus and magnesium , are displayed in Figure 6(a) and (b). We notice that these effects are very smooth, but we have to remember that the covariate information is sparse and only available in grid cells. In terms of DIC, the empirical covariate terms explain some spatial structure of the pattern as DIC increases from 15,379 to 15,440 if these two terms are not included. High phosphorus seems to coincide with low tree density and a similar, but less clear, pattern emerges for magnesium. Currently, the ecological literature cannot explain these results, but they could be related to resource partitioning along axes of soil nutrient availability [Burslem, personal communication (September 2011), John et al. (2007)]. In addition, it is currently also unclear if the soil properties cause an aggregation of trees, as they provide suitable growing conditions, or whether a high tree intensity leads to low levels of magnesium or phosphorus resulting from the chemical composition of the leaf litter.

The plot of the constructed covariate in Figure 6(c) illustrates the resolution of the local clustering represented by it. The resulting estimated function of the constructed covariate is shown in Figure 6(d), which indicates that it accounts for clustering of up to a distance of 15 metres. The estimated spatial effect for the pattern is given in Figure 6(e), while Figure 6(f) displays the estimated spatially structured effect if the constructed covariate is left out of the full model. This last figure shows clear local structure in the spatial effect and might give a model which is overfitted to the actual pattern. Including the constructed covariate, the local structure of the spatial effect is removed, making the spatial effect smoother. This indicates that spatial behavior at a local scale has been picked up by the constructed covariate. In this way the model can account for spatial structures at different scales. The two unstructured spatial fields in equations (9) and (10) do not show any particular pattern (results not shown). Fitting this model took 55 minutes to run (2.66 GHz Intel Core i7 processor).

4.3 Discussion on rainforest data

In this section we consider a log Gaussian Cox process model and fit it jointly to a point pattern data set with a large number of points and two covariates that have been observed at a relatively small number of points within the plot. Waagepetersen (2007) and Waagepetersen and Guan (2009) model the patterns formed by rainforest tree species with this data structure, using Thomas processes to include local clustering resulting from seed dispersal. This approximate approach is based on the minimum contrast method for parameter estimation. Rue, Martino and Chopin (2009) consider the same data in the context of Cox processes to demonstrate that log-Gaussian Cox processes can be fitted conveniently to a large spatial point pattern using INLA relative to environmental covariates which are assumed to be known everywhere and fixed. In many typical applications, however, the values of spatial covariates in the location of the points forming the point process are not known. Similarly, the direction of the relationship between soil properties and tree presence may be not clear. We generalize the approach in Rue, Martino and Chopin (2009) here and fit a joint model of the pattern and the covariates. This approach distinguishes between locations where the values of the covariates are available but potentially subject to measurement error and those where they are not. In addition, it does not assume that the soil variables impact on the pattern but not vice versa. We also consider a constructed covariate that reflects local clustering as a result of local seed dispersal, as discussed above.

The given approach accommodates model comparison and model assessment, both of which are of practical value in many applications. An inspection of the estimated spatially structured effect in Figure 6(e) indicates that some spatial structure still remains in the point pattern which cannot be explained by the current model, that is, the current model can still be improved on. Hence, judging by Figure 6(e), it might be possible to improve the model by including further covariates and the structure of the estimated spatial effect might be used to suggest a suitable covariate. Previous approaches to fitting a model to these data [Waagepetersen (2007); Waagepetersen and Guan (2009)] neither have been able to reveal the shortcomings of the models nor to provide mechanisms that help identify covariates that might improve the model.

The function of a constructed covariate [Figure 6(d)], which reflects local clustering up to a distance of 15 metres, may be interpreted as a seed dispersal kernel. Biological research has shown that this species is likely to be dispersed primarily by small understorey birds that feed in the canopy and mostly drop the seeds beneath the parent tree. Since trees of the species Aporusa microstachyathese are relatively small, 15 m reflect the maximum radius of the tree crown [Burslem, Garwood and Thomas (2001)].

The approach discussed here can be extended easily to allow more complex models to be fitted, such as a model of both the spatial pattern and associated marks, along the lines of the model discussed in Section 5. For instance, this may include a model of both the spatial pattern and the size and the growth of the trees. Here, both size and growth might depend on the spatial pattern and growth might also depend on size.

5 Modeling marks and pattern in a marked point pattern with multiple marks

Modeling the behavior of individuals in space based simply on the individuals’ locations and ignoring their properties is certainly a gross over-simplification for many systems. In practice, researchers hence often collect data on the locations of the individuals along with data on additional properties, that is, marks. In this section we discuss a marked point pattern with several dependent marks, which also depend on the spatial pattern, and consider a joint model of the marks and the pattern. Models where marks depend on the point pattern have recently been considered in the literature [Menezes (2005); Ho and Stoyan (2008); Myllymäki and Penttinen (2009)]. Also note the work by Diggle, Menezes and Su (2010), where a point process with intensity dependent marks is used in the context of preferential sampling in geostatistics. The model we fit here is more general than these related models, since we model multiple dependent marks jointly with the pattern.

5.1 Data structure and modeling approach

We analyze a spatial point pattern together with several types of nonindependent associated marks. We consider only two marks and here, but the approach can be generalized in a straightforward way to include more than two marks. The are assumed to follow an exponential family distribution with parameter vector and to depend on the intensity of the point pattern, while the are assumed to follow a (different) exponential family distribution with parameter vector and to depend both on the intensity of the point pattern and on the marks . Without loss of generality, the parameters and are the location parameters of the distributions and , respectively.

We discretize the observation window as discussed in Section 2.3, and for the spatial pattern we assume the model


using the same notation as in (4). For the marks, we construct a model where the marks depend on the pattern by assuming that they depend on the same spatially structured effect . Specifically, we assume that with


where is another error term. The marks are assumed to depend both on the spatial pattern through and on the marks . We thus have that with


where denotes another error term.

Figure 7: Spatial pattern formed by the locations of the eucalyptus trees in the koala data set; the diameters of the circles reflect the value of the leaf marks (a) and the frequency marks (b), respectively.

5.2 Application to example data set

5.2.1 Koala data

Koalas are arboreal marsupial herbivores native to Australia with a very low metabolic rate. They rest motionless for about 18 to 20 hours a day, sleeping most of that time. They feed selectively and live almost entirely on eucalyptus leaves. Whereas these leaves are poisonous to most other species, the koala gut has adapted to digest them. It is likely that the animals preferentially forage leaves that are high in nutrients and low in toxins as an extreme example of evolutionary adaptation. An understanding of the koala-eucalyptus interaction is crucial for conservation efforts [Moore et al. (2010)].

The data have been collected in a study conducted at the Koala Conservation Centre on Phillip Island, near Melbourne, Australia. For each of 915 trees within a reserve enclosed by a koala-proof fence (Figure 7), information on the leaf chemistry and on the frequency of koala visits has been collected. The leaf chemistry is summarized in a measure of the palatability of the leaves (“leaf mark” ). Palatability is assumed to depend on the intensity of the point pattern. In addition, “frequency marks” describe for each tree the diurnal tree use by individual koalas collected at monthly intervals between 1993 and March 2004. The are assumed to depend on the intensity of the point pattern as well as on the leaf marks.

There are no additional covariate data available for the given data set. Hence, for the locations of the trees we use the model in (11) with notation as above. For the leaf and frequency marks we use the models in equations (12) and (13), respectively. The leaf marks are assumed to follow a normal distribution and the frequency marks a Poisson distribution, that is, and .

5.2.2 Results

With these distributional assumptions for the marks, we fit a joint model as given in equations (11)–(13) to the data set. The results are based on an observation window discretized into 1571 grid cells. In order to fit spatial effects, we embed this area within a rectangular area. For the constructed covariate, we perform a simple edge correction for the distances in (5), assuming missing values in grid cells in which the distance from the center point to the border is shorter than the nearest-point distance.

Model Terms DIC Time (s)
1. Only error terms
2. Add intercepts
3. Add fixed covariate ()
4. Add spatial effect
– Only for pattern
– For pattern and leaf marks
– For pattern and frequency marks
– For pattern and both marks (final model)
5. Add constructed covariate
Table 1: DIC values and computation time for different fitted models for the koala data

When fitting complex models it can be useful to apply a stepwise procedure to study the impact of each term in the model. Table 1 illustrates DIC-values and computation time (in seconds) of models of increasing complexity. In the first three steps we initially run a model with only error terms and then add intercepts and the fixed covariate for the frequency marks. Step 4 illustrates the effect of adding the spatial effect in modeling the pattern together with one or both of the two marks, in which DIC decreases to 6943. Inclusion of the constructed covariate in (11) does not improve the model fitting for this data set. This is not surprising, as the original pattern does not seem to exhibit any strong local clustering effect and as a result the estimated function of the constructed covariate is not significantly different from 0.

(a) (b)
(c) (d)
Figure 8: Plots of the estimated common spatial effect (a) and the three unstructured effects , and (b)–(d) for the koala data set.


Table 2: Posterior means and credible intervals for parameters in the koala model
Parameter Mean credible interval

The estimated common spatial effect [Figure 8(a)] represents spatial autocorrelation present in the pattern and the marks which might be the result of related environmental processes such as nutrient levels in the soil. The estimated parameter value for and have opposite signs (Table 2). The negative sign for indicates that palatability is low where the trees are aggregated, which might have been caused by competition for soil nutrients in these areas. The positive sign for reflects that the koalas are more likely to be present in areas with higher intensity. Recalling that the data have been accumulated over time, this might be due to the koalas being more likely to change from one tree to a neighboring tree where the trees are aggregated. The mean of the posterior density for the parameter in the final model is 1.38, indicating a significant positive influence of palatability on the frequency of koala visits to the trees. The three unstructured terms are given in Figures 8(b)–(d). A slight trend in the residuals for the leaf marks may be observed in Figure 8(c), with lower values toward the bottom left probably reflecting an inhomogeneity that cannot be accounted for by the joint spatial effect .

5.3 Discussion of koala data

The example considered in this section is a marked Cox process model, that is, a model of both the spatial pattern and two types of dependent marks, providing information on the spatial pattern at the same time as on the marks and their dependence. In cases where the marks are of primary scientific interest, one could view this approach as a model of the marks which implicitly takes the spatial dependence into account by modeling it alongside the marks. The model we use here is similar to approaches taken in Menezes (2005), Ho and Stoyan (2008), Myllymäki and Penttinen (2009). Since our approach is very flexible, it can easily be generalized to allow for separate spatially structured effects for the pattern and the marks and to include additional empirical covariates; these have not been available here. Hence, using the approach considered here, we are able to fit easily a complex spatial point process model to a marked point pattern and to assess its suitability for a specific data set.

Marked point pattern data sets where data on marks are likely to depend on an underlying spatial pattern are not uncommon. Within ecology, for instance, metapopulation data [Hanski and Gilpin (1997)] typically consist of the locations of subpopulations and their properties, and have a similar structure to the data set considered here. These data sets may be modeled using a similar approach and it is straightforward to fit related but more complex models, including empirical covariates or temporal replicates. Similarly, marks are available for the rainforest data discussed in Section 4. As mentioned there, a model that includes the marks of the trees may also be fitted using the approach discussed here.

6 Discussion

Researchers outside the statistical community have become familiar with fitting a large range of different models to complex data sets using software available in R. This paper provides a very flexible framework for routinely fitting models to complex spatial point pattern data with little computational effort using models that account for both local and global spatial behavior. We consider complex data examples and demonstrate how marks as well as covariates can be included in a joint model. That is, we consider a situation where the marks and the covariates can be modeled along with the pattern and show that it is computationally feasible to do so. We can take account of local spatial structure by using a constructed covariate, which we discuss in detail in Section 3.

The two models discussed here indicate that our approach can be applied in a wide range of situations and is flexible enough to facilitate the fitting of other even more complex models. It is feasible to fit several related models to realistically complex data sets if necessary, and to use the DIC to aid the choice of covariates. The posterior distributions of the estimated parameters can be used to assess the significance of the influence of different covariates in the models. Through the use of a structured spatial effect and an unstructured spatial effect it is possible to assess the quality of the model fit. Specifically, the structured spatial effect can be used to reveal spatial correlations in the data that have not been explained with the covariates and may help researchers identify suitable covariates to incorporate into the model. Spatially unstructured effects may be used to account for and identify extreme observations such as locations where covariate values have been collected with a particularly strong measurement error.

There is an extensive literature on descriptive and nonparametric approaches to the analysis of spatial point patterns, specifically on (functional) summary characteristics describing first and second order spatial behavior, in particular, on Ripley’s -function [Ripley (1976)] and the pair correlation function [Stoyan, Kendall and Mecke (1995)]. In both the statistical and the applied literature these have been discussed far more frequently than likelihood based modeling approaches, and provide an elegant means for characterizing the properties of spatial patterns [Illian et al. (2008)]. A thorough analysis of a spatial point pattern typically includes an extensive exploratory analysis and in many cases it may even seem unnecessary to continue the analysis and fit a spatial point process model to a pattern. An exploratory analysis based on functional summary characteristics, such as Ripley’s -function or the pair-correlation function, considers spatial behavior at a multitude of spatial scales, making this approach particularly appealing. However, with increasing complexity of the data, it becomes less obvious how suitable summary characteristics should be defined for these, and a point process model may be a suitable alternative. For example, it is not obvious how one would jointly analyze the two different marks together with the pattern in the koala data set based on summary characteristics. However, as discussed in Section 5, it is straightforward to do this with a hierarchical model. In addition, most exploratory analysis tools assume the process to be first-order stationary or at least second-order reweighted stationary [Baddeley et al. (2000)]—a situation that is both rare and difficult to assess in applications, in particular, in the context of realistic and complex data sets. The approach discussed here does not make any assumptions about stationarity but explicitly includes spatial trends into the model.

In the literature, local spatial behavior has often been modelled by a Gibbs process. Large-scale spatial behavior may be incorporated into a Gibbs process model as a parametric or nonparametric, yet deterministic, trend, while it is treated as a stochastic process in itself here. Modeling the spatial trend in a Gibbs process hence often assumes that an explicit and deterministic model of the trend as a function of location (and spatial covariates) is known [Baddeley and Turner (2005)]. Even in the nonparametric situation, the estimated values of the underlying spatial trend are considered fixed values, which are subject neither to stochastic variation nor to measurement error. Since it is based on a latent random field, the approach discussed here differs substantially from the Gibbs process approach and assumes a hierarchical, doubly stochastic structure. This very flexible class of point processes provides models of local spatial behavior relative to an underlying large-scale spatial trend. In realistic applications this spatial trend is not known. Values of the covariates that are continuous in space are typically not known everywhere and have been interpolated. It is likely that spatial trends exist in the data that cannot be accounted for by the covariates. The spatial trend is hence not regarded as deterministic but assumed to be a random field. This approach allows to jointly model the covariate and the spatial pattern as in the model used for the rainforest example data set. Clearly, unlike Gibbs processes, log Gaussian Cox processes do not allow second order inter-individual interactions to be included in a model. In a situation where these are of primary interest, Cox processes are certainly not suitable.

In order to make model fitting feasible, the continuous Gaussian random field is approximated here by a discrete Gauss Markov random field. While this is computationally elegant, one might argue that this approximation is not justified and is too coarse, resulting in an unnecessary loss of information. Clearly, since any model only has a finite representation in a computer, model fitting approaches often work with some degree of discretization. However, and more importantly, Lindgren, Rue and Lindström (2011) show that there is an explicit link between a large class of covariance functions (and hence the Gaussian random field based on these) and Gauss Markov random fields, clearly pointing out that the approximation is indeed justified. In addition, based on the results discussed in Lindgren, Rue and Lindström (2011), the approaches taken in this paper may be extended to avoid the computationally wasteful need of having to use a regular grid [Illian and Simpson (2011)]. Illian et al. (2012) also mention the issue of complex boundaries structures that are particularly relevant for point process data sets where the observation window has been chosen to align with natural boundaries that may impact on pattern. While this is clearly not an issue for the rainforest data set since the boundaries have been chosen arbitrarily, the koala data set, however, has been observed in an observation window surrounded by a koala proof fence. This fence does probably not impact on the locations of the trees nor the leaf chemistry but might increase the frequency of koala visits near the fence. The approach in Lindgren, Rue and Lindström (2011) may be used to define varying boundary conditions for different parts of the data set, and hence allow for more realistic modeling for data sets with complicated boundary structures.

In summary, the methodology discussed here, together with the R library R-INLA (http://www.r-inla.org/), makes complex spatial point process models accessible to scientists outside the statistical sciences and provides them with a toolbox for routinely fitting and assessing the fit of suitable and realistic point process models to complex spatial point pattern data.


Some of the ideas relevant to the rainforest data were developed during a working group on “Spatial analysis of tropical forest biodiversity” funded by the Natural Environment Research Council and English Nature through the NERC Centre for Population Biology and UK Population Biology Network. The data were collected by the Center for Tropical Forest Science and the Forest Institute of Malaysia funded by the U.S. National Science Foundation, the Smithsonian Tropical Research Institution and the National Institute of Environmental Studies (Japan).

We would like to thank David Burslem, University of Aberdeen, and Richard Law, University of York, for introducing the rainforest data into the statistical community and for many in-depth discussions over the last few years. We also thank Colin Beale, University of York, and Ben Moore, James Hutton Institute, Aberdeen, for extended discussions on the koala data.

The authors also gratefully acknowledge the financial support of Research Councils UK for Illian.


  • Baddeley, Møller and Waagepetersen (2000) {barticle}[mr] \bauthor\bsnmBaddeley, \bfnmA.\binitsA., \bauthor\bsnmMøller, \bfnmJ.\binitsJ. \AND\bauthor\bsnmWaagepetersen, \bfnmR.\binitsR. (\byear2000). \btitleNon- and semi-parametric estimation of interaction in inhomogeneous point patterns. \bjournalStat. Neerl. \bvolume54 \bpages329–350. \bptokimsref \endbibitem
  • Baddeley and Turner (2000) {barticle}[mr] \bauthor\bsnmBaddeley, \bfnmAdrian\binitsA. \AND\bauthor\bsnmTurner, \bfnmRolf\binitsR. (\byear2000). \btitlePractical maximum pseudolikelihood for spatial point patterns (with discussion). \bjournalAust. N. Z. J. Stat. \bvolume42 \bpages283–322. \biddoi=10.1111/1467-842X.00128, issn=1369-1473, mr=1794056 \bptnotecheck related\bptokimsref \endbibitem
  • Baddeley and Turner (2005) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmBaddeley, \bfnmA. J.\binitsA. J. \AND\bauthor\bsnmTurner, \bfnmR.\binitsR. (\byear2005). \btitleSpatstat: An R package for analyzing spatial point patterns. \bjournalJournal of Statistical Software \bvolume12 \bpages1–42. \bptokimsref \endbibitem
  • Baddeley et al. (2005) {barticle}[mr] \bauthor\bsnmBaddeley, \bfnmA.\binitsA., \bauthor\bsnmTurner, \bfnmR.\binitsR., \bauthor\bsnmMøller, \bfnmJ.\binitsJ. \AND\bauthor\bsnmHazelton, \bfnmM.\binitsM. (\byear2005). \btitleResidual analysis for spatial point processes. \bjournalJ. R. Stat. Soc. Ser. B Stat. Methodol. \bvolume67 \bpages617–666. \biddoi=10.1111/j.1467-9868.2005.00519.x, issn=1369-7412, mr=2210685 \bptnotecheck related\bptokimsref \endbibitem
  • Berman and Turner (1992) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmBerman, \bfnmM.\binitsM. \AND\bauthor\bsnmTurner, \bfnmR.\binitsR. (\byear1992). \btitleApproximating point process likelihoods with GLIM. \bjournalApplied Statistics \bvolume41 \bpages31–38. \bptokimsref \endbibitem
  • Besag (1977) {barticle}[mr] \bauthor\bsnmBesag, \bfnmJ. E.\binitsJ. E. (\byear1977). \btitleContribution to the discussion of Dr. Ripley’s paper. \bjournalJ. Roy. Statist. Soc. Ser. B \bvolume39 \bpages193–195. \bptokimsref \endbibitem
  • Burslem, Garwood and Thomas (2001) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmBurslem, \bfnmD. F. R. P.\binitsD. F. R. P., \bauthor\bsnmGarwood, \bfnmN. C.\binitsN. C. \AND\bauthor\bsnmThomas, \bfnmS. C.\binitsS. C. (\byear2001). \btitleTropical forest diversity—the plot thickens. \bjournalScience \bvolume291 \bpages606–607. \bptokimsref \endbibitem
  • Condit (1998) {bmisc}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmCondit, \bfnmR.\binitsR. (\byear1998). \bhowpublishedTropical Forest Census Plots. Springer and R. G. Landes Company, Berlin and Georgetown, TX. \bptokimsref \endbibitem
  • Diggle (2003) {bbook}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmDiggle, \bfnmPeter J.\binitsP. J. (\byear2003). \btitleStatistical Analysis of Spatial Point Patterns, \bedition2nd ed. \bpublisherHodder Arnold, \baddressLondon. \bptokimsref \endbibitem
  • Diggle, Menezes and Su (2010) {barticle}[mr] \bauthor\bsnmDiggle, \bfnmPeter J.\binitsP. J., \bauthor\bsnmMenezes, \bfnmRaquel\binitsR. \AND\bauthor\bsnmSu, \bfnmTing-li\binitsT.-l. (\byear2010). \btitleGeostatistical inference under preferential sampling. \bjournalJ. R. Stat. Soc. Ser. C. Appl. Stat. \bvolume59 \bpages191–232. \biddoi=10.1111/j.1467-9876.2009.00701.x, issn=0035-9254, mr=2744471 \bptnotecheck related\bptokimsref \endbibitem
  • Forchhammer and Boomsma (1995) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmForchhammer, \bfnmM. C.\binitsM. C. \AND\bauthor\bsnmBoomsma, \bfnmJ.\binitsJ. (\byear1995). \btitleForaging strategies and seasonal diet optimization of muskoxen in West Greenland. \bjournalOecologia \bvolume104 \bpages169–180. \bptokimsref \endbibitem
  • Forchhammer and Boomsma (1998) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmForchhammer, \bfnmM. C.\binitsM. C. \AND\bauthor\bsnmBoomsma, \bfnmJ.\binitsJ. (\byear1998). \btitleOptimal mating strategies in nonterritorial ungulates: A general model tested on muskoxen. \bjournalBehavioural Ecology \bvolume9 \bpages136–143. \bptokimsref \endbibitem
  • Guan (2008) {barticle}[mr] \bauthor\bsnmGuan, \bfnmYongtao\binitsY. (\byear2008). \btitleOn consistent nonparametric intensity estimation for inhomogeneous spatial point processes. \bjournalJ. Amer. Statist. Assoc. \bvolume103 \bpages1238–1247. \biddoi=10.1198/016214508000000526, issn=0162-1459, mr=2528839 \bptokimsref \endbibitem
  • Hanski and Gilpin (1997) {bbook}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmHanski, \bfnmI. A.\binitsI. A. \AND\bauthor\bsnmGilpin, \bfnmM. E.\binitsM. E. (\byear1997). \btitleMetapopulation Biology: Ecology, Genetics and Evolution. \bpublisherAcademic Press, \baddressSan Diego. \bptokimsref \endbibitem
  • Hardy and Vekemans (2002) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmHardy, \bfnmO. J.\binitsO. J. \AND\bauthor\bsnmVekemans, \bfnmX.\binitsX. (\byear2002). \btitleSPAGEDi: A versatile computer program to analyse spatial genetic structure at the individual or population levels. \bjournalMolecular Ecology Notes \bvolume2 \bpages618–620. \bptokimsref \endbibitem
  • Ho and Stoyan (2008) {barticle}[mr] \bauthor\bsnmHo, \bfnmLai Ping\binitsL. P. \AND\bauthor\bsnmStoyan, \bfnmD.\binitsD. (\byear2008). \btitleModelling marked point patterns by intensity-marked Cox processes. \bjournalStatist. Probab. Lett. \bvolume78 \bpages1194–1199. \biddoi=10.1016/j.spl.2007.11.013, issn=0167-7152, mr=2441462 \bptokimsref \endbibitem
  • Huang and Ogata (1999) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmHuang, \bfnmF.\binitsF. \AND\bauthor\bsnmOgata, \bfnmY.\binitsY. (\byear1999). \btitleImprovements of the maximum pseudo-likelihood estimators in various spatial statistical models. \bjournalJ. Comput. Graph. Statist. \bvolume8 \bpages519–530. \bptokimsref \endbibitem
  • Hubbell, Condit and Foster (2005) {bmisc}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmHubbell, \bfnmS. P.\binitsS. P., \bauthor\bsnmCondit, \bfnmR.\binitsR. \AND\bauthor\bsnmFoster, \bfnmR. B.\binitsR. B. (\byear2005). \bhowpublishedBarro Colorado forest census plot data. Available at http://www.ctfs.si.edu/datasets/bci. \bptokimsref \endbibitem
  • Hubbell et al. (1999) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmHubbell, \bfnmS. P.\binitsS. P., \bauthor\bsnmFoster, \bfnmR. B.\binitsR. B., \bauthor\bsnmO’Brien, \bfnmS. T.\binitsS. T., \bauthor\bsnmHarms, \bfnmK. E.\binitsK. E., \bauthor\bsnmCondit, \bfnmR.\binitsR., \bauthor\bsnmWechsler, \bfnmB.\binitsB., \bauthor\bsnmWright, \bfnmS. J.\binitsS. J. \AND\bauthor\bparticleLoo de \bsnmLao, \bfnmS.\binitsS. (\byear1999). \btitleLight gap disturbances, recruitment limitation, and tree diversity in a neotropical forest. \bjournalScience \bvolume283 \bpages554–557. \bptokimsref \endbibitem
  • Illian and Hendrichsen (2010) {barticle}[mr] \bauthor\bsnmIllian, \bfnmJanine B.\binitsJ. B. \AND\bauthor\bsnmHendrichsen, \bfnmDitte K.\binitsD. K. (\byear2010). \btitleGibbs point process models with mixed effects. \bjournalEnvironmetrics \bvolume21 \bpages341–353. \biddoi=10.1002/env.1008, issn=1180-4009, mr=2842247 \bptokimsref \endbibitem
  • Illian and Rue (2010) {bmisc}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmIllian, \bfnmJ. B.\binitsJ. B. \AND\bauthor\bsnmRue, \bfnmH.\binitsH. (\byear2010). \bhowpublishedA toolbox for fitting complex spatial point process models using integrated Laplace transformation (INLA). Technical report, Trondheim Univ. \bptokimsref \endbibitem
  • Illian and Simpson (2011) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmIllian, \bfnmJ. B.\binitsJ. B. \AND\bauthor\bsnmSimpson, \bfnmD.\binitsD. (\byear2011). \btitleComment on Lindgren et al., an explicit link between Gaussian fields and Gaussian Markov random fields: The SPDE approach. \bjournalJ. Roy. Statist. Soc. Ser. B \bvolume73 \bpages423–498. \bptokimsref \endbibitem
  • Illian et al. (2008) {bbook}[mr] \bauthor\bsnmIllian, \bfnmJanine\binitsJ., \bauthor\bsnmPenttinen, \bfnmAntti\binitsA., \bauthor\bsnmStoyan, \bfnmHelga\binitsH. \AND\bauthor\bsnmStoyan, \bfnmDietrich\binitsD. (\byear2008). \btitleStatistical Analysis and Modelling of Spatial Point Patterns. \bpublisherWiley, \baddressChichester. \bidmr=2384630 \bptokimsref \endbibitem
  • Illian et al. (2012) {bmisc}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmIllian, \bfnmJ. B.\binitsJ. B., \bauthor\bsnmSørbye, \bfnmS. H.\binitsS. H., \bauthor\bsnmRue, \bfnmH.\binitsH. \AND\bauthor\bsnmHendrichsen, \bfnmD. K.\binitsD. K. (\byear2012). \bhowpublishedFitting a log Gaussian Cox process with temporally varying effects—a case study. Journal of Environmental Statistics. To appear. \bptokimsref \endbibitem
  • John et al. (2007) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmJohn, \bfnmR. C.\binitsR. C., \bauthor\bsnmDalling, \bfnmJ. W.\binitsJ. W., \bauthor\bsnmHarms, \bfnmK. E.\binitsK. E., \bauthor\bsnmYavitt, \bfnmJ. B.\binitsJ. B., \bauthor\bsnmStallard, \bfnmR. F.\binitsR. F., \bauthor\bsnmMirabello, \bfnmM.\binitsM., \bauthor\bsnmHubbell, \bfnmS. P.\binitsS. P., \bauthor\bsnmValencia, \bfnmR.\binitsR., \bauthor\bsnmNavarrete, \bfnmH.\binitsH., \bauthor\bsnmVallejo, \bfnmM.\binitsM. \AND\bauthor\bsnmFoster, \bfnmR. B.\binitsR. B. (\byear2007). \btitleSoil nutrients influence spatial distributions of tropical tree species. \bjournalProc. Natl. Acad. Sci. USA \bvolume104 \bpages864–869. \bptokimsref \endbibitem
  • Johnson and Boerlijst (2002) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmJohnson, \bfnmC. R.\binitsC. R. \AND\bauthor\bsnmBoerlijst, \bfnmM. C.\binitsM. C. (\byear2002). \btitleSelection at the level of the community: The importance of spatial structure. \bjournalTrends in Ecology & Evolution \bvolume17 \bpages83–90. \bptokimsref \endbibitem
  • Killingback and Doebeli (1996) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmKillingback, \bfnmT.\binitsT. \AND\bauthor\bsnmDoebeli, \bfnmM.\binitsM. (\byear1996). \btitleSpatial evolutionary game theory: Hawks and doves revisited. \bjournalProc. R. Soc. Lond. Ser. B \bvolume263 \bpages1135–1144. \bptokimsref \endbibitem
  • Latimer et al. (2009) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmLatimer, \bfnmA. M.\binitsA. M., \bauthor\bsnmBanerjee, \bfnmS.\binitsS., \bauthor\bsnmSang, \bfnmS.\binitsS., \bauthor\bsnmMosher, \bfnmE. S.\binitsE. S. \AND\bauthor\bsnmSilanderJr., \bfnmJ. A.\binitsJ. A. (\byear2009). \btitleHierarchical models facilitate spatial analysis of large data sets: A case study on invasive plant species in the northeastern United States. \bjournalEcology Letters \bvolume12 \bpages144–154. \bptokimsref \endbibitem
  • Law et al. (2001) {bincollection}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmLaw, \bfnmR.\binitsR., \bauthor\bsnmPurves, \bfnmD. W.\binitsD. W., \bauthor\bsnmMurrell, \bfnmD. J.\binitsD. J. \AND\bauthor\bsnmDieckmann, \bfnmU.\binitsU. (\byear2001). \btitleCauses and effects of small scale spatial structure in plant populations. In \bbooktitleIntegrating Ecology and Evolution in a Spatial Context (\beditor\bfnmJ.\binitsJ. \bsnmSilvertown \AND\beditor\bfnmJ.\binitsJ. \bsnmAntonovics, eds.) \bpages21–44. \bpublisherBlackwell, \baddressOxford. \bptokimsref \endbibitem
  • Law et al. (2009) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmLaw, \bfnmR.\binitsR., \bauthor\bsnmIllian, \bfnmJ. B.\binitsJ. B., \bauthor\bsnmBurslem, \bfnmD. F. R. P.\binitsD. F. R. P., \bauthor\bsnmGratzer, \bfnmG.\binitsG., \bauthor\bsnmGunatilleke, \bfnmC. V. S.\binitsC. V. S. \AND\bauthor\bsnmGunatilleke, \bfnmI. A. U. N.\binitsI. A. U. N. (\byear2009). \btitleEcological information from spatial patterns of plants: Insights from point process theory. \bjournalJournal of Ecology \bvolume97 \bpages616–628. \bptokimsref \endbibitem
  • Lawson (1992) {bincollection}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmLawson, \bfnmA.\binitsA. (\byear1992). \btitleOn fitting non-stationary Markov point process models on GLIM. In \bbooktitleComputational Statistics, Volume 1 (\beditor\bfnmY.\binitsY. \bsnmDodge \AND\beditor\bfnmJ.\binitsJ. \bsnmWhittaker, eds.) \bpages35–40. \bpublisherPhysica Verlag, \baddressHeidelberg. \bptokimsref \endbibitem
  • Lindgren, Rue and Lindström (2011) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmLindgren, \bfnmF.\binitsF., \bauthor\bsnmRue, \bfnmH.\binitsH. \AND\bauthor\bsnmLindström, \bfnmJ.\binitsJ. (\byear2011). \btitleAn explicit link between Gaussian fields and Gaussian Markov random fields: The SPDE approach. \bjournalJ. Roy. Statist. Soc. Ser. B \bvolume73 \bpages423–498. \bptokimsref \endbibitem
  • Lunn et al. (2000) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmLunn, \bfnmD. J.\binitsD. J., \bauthor\bsnmThomas, \bfnmA.\binitsA., \bauthor\bsnmBest, \bfnmN.\binitsN. \AND\bauthor\bsnmSpiegelhalter, \bfnmD.\binitsD. (\byear2000). \btitleWinBUGS—a Bayesian modelling framework: Concepts, structure, and extensibility. \bjournalStat. Comput. \bvolume10 \bpages325–337. \bptokimsref \endbibitem
  • Menezes (2005) {bmisc}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmMenezes, \bfnmR.\binitsR. (\byear2005). \bhowpublishedAssessing spatial dependency under non-standard sampling. Ph.D. thesis, Universidad de Santiago de Compostela, Santiago de Compostela, Spain. \bptokimsref \endbibitem
  • Metropolis et al. (1953) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmMetropolis, \bfnmN.\binitsN., \bauthor\bsnmRosenbluth, \bfnmA. W.\binitsA. W., \bauthor\bsnmRosenbluth, \bfnmM. N.\binitsM. N., \bauthor\bsnmTeller, \bfnmA. H.\binitsA. H. \AND\bauthor\bsnmTeller, \bfnmE.\binitsE. (\byear1953). \btitleEquations of state calculations by fast computing machines. \bjournalJournal of Chemical Physics \bvolume6 \bpages1087–1092. \bptokimsref \endbibitem
  • Møller, Syversveen and Waagepetersen (1998) {barticle}[mr] \bauthor\bsnmMøller, \bfnmJesper\binitsJ., \bauthor\bsnmSyversveen, \bfnmAnne Randi\binitsA. R. \AND\bauthor\bsnmWaagepetersen, \bfnmRasmus Plenge\binitsR. P. (\byear1998). \btitleLog Gaussian Cox processes. \bjournalScand. J. Stat. \bvolume25 \bpages451–482. \biddoi=10.1111/1467-9469.00115, issn=0303-6898, mr=1650019 \bptokimsref \endbibitem 
  • Møller and Waagepetersen (2004) {bbook}[mr] \bauthor\bsnmMøller, \bfnmJesper\binitsJ. \AND\bauthor\bsnmWaagepetersen, \bfnmRasmus Plenge\binitsR. P. (\byear2004). \btitleStatistical Inference and Simulation for Spatial Point Processes. \bseriesMonographs on Statistics and Applied Probability \bvolume100. \bpublisherChapman & Hall/CRC, \baddressBoca Raton, FL. \bidmr=2004226 \bptokimsref \endbibitem
  • Møller and Waagepetersen (2007) {barticle}[mr] \bauthor\bsnmMøller, \bfnmJesper\binitsJ. \AND\bauthor\bsnmWaagepetersen, \bfnmRasmus P.\binitsR. P. (\byear2007). \btitleModern statistics for spatial point processes. \bjournalScand. J. Stat. \bvolume34 \bpages643–711. \biddoi=10.1111/j.1467-9469.2007.00569.x, issn=0303-6898, mr=2392447 \bptnotecheck related\bptokimsref \endbibitem
  • Moore et al. (2010) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmMoore, \bfnmB. D.\binitsB. D., \bauthor\bsnmLawler, \bfnmI. R.\binitsI. R., \bauthor\bsnmWallis, \bfnmI. R.\binitsI. R., \bauthor\bsnmBeale, \bfnmC. M.\binitsC. M. \AND\bauthor\bsnmFoley, \bfnmW. J.\binitsW. J. (\byear2010). \btitlePalatability mapping: A koala’s eye view of spatial variation in habitat quality. \bjournalEcology \bvolume91 \bpages3165–3176. \bptokimsref \endbibitem
  • Myllymäki and Penttinen (2009) {barticle}[mr] \bauthor\bsnmMyllymäki, \bfnmMari\binitsM. \AND\bauthor\bsnmPenttinen, \bfnmAntti\binitsA. (\byear2009). \btitleConditionally heteroscedastic intensity-dependent marking of log Gaussian Cox processes. \bjournalStat. Neerl. \bvolume63 \bpages450–473. \biddoi=10.1111/j.1467-9574.2009.00433.x, issn=0039-0402, mr=2598981 \bptokimsref \endbibitem
  • Naylor et al. (2009) {bmisc}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmNaylor, \bfnmM.\binitsM., \bauthor\bsnmGreenhough, \bfnmJ.\binitsJ., \bauthor\bsnmMcCloskey, \bfnmJ.\binitsJ., \bauthor\bsnmBell, \bfnmA. F.\binitsA. F. \AND\bauthor\bsnmMain, \bfnmI. G.\binitsI. G. (\byear2009). \bhowpublishedStatistical evaluation of characteristic earthquakes in the frequency-magnitude distributions of sumatra and other subduction zone regions. Geophysical Research Letters 36 DOI:\doiurl10.1029/2009GL040460. \bptokimsref \endbibitem
  • Neyman and Scott (1952) {barticle}[mr] \bauthor\bsnmNeyman, \bfnmJ.\binitsJ. \AND\bauthor\bsnmScott, \bfnmE. L.\binitsE. L. (\byear1952). \btitleA theory of the spatial distribution of galaxies. \bjournalAstrophys. J. \bvolume116 \bpages144–163. \bidissn=0004-637X, mr=0053640 \bptokimsref \endbibitem
  • Ogata (1999) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmOgata, \bfnmY.\binitsY. (\byear1999). \btitleSeismicity analysis through point-process modeling: A review. \bjournalPure and Applied Geophysics \bvolume155 \bpages471–507. \bptokimsref \endbibitem
  • R Development Core Team (2009) {bmisc}[auto:STB—2012/03/21—07:41:58] \borganizationR Development Core Team. (\byear2009). \bhowpublishedR: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. Available at http://www.R-project.org. ISBN 3-900051-07-0. \bptokimsref \endbibitem
  • Rajala and Illian (2012) {bmisc}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmRajala, \bfnmT. A.\binitsT. A. \AND\bauthor\bsnmIllian, \bfnmJ. B.\binitsJ. B. (\byear2012). \bhowpublishedA family of spatial biodiversity measures based on graphs. Environ. Ecol. Stat. To appear. \bptokimsref \endbibitem
  • Ripley (1976) {barticle}[mr] \bauthor\bsnmRipley, \bfnmB. D.\binitsB. D. (\byear1976). \btitleThe second-order analysis of stationary point processes. \bjournalJ. Appl. Probab. \bvolume13 \bpages255–266. \bidissn=0021-9002, mr=0402918 \bptokimsref \endbibitem
  • Rue and Held (2005) {bbook}[mr] \bauthor\bsnmRue, \bfnmHåvard\binitsH. \AND\bauthor\bsnmHeld, \bfnmLeonhard\binitsL. (\byear2005). \btitleGaussian Markov Random Fields: Theory and Applications. \bseriesMonographs on Statistics and Applied Probability \bvolume104. \bpublisherChapman & Hall/CRC, \baddressBoca Raton, FL. \biddoi=10.1201/9780203492024, mr=2130347 \bptokimsref \endbibitem
  • Rue, Martino and Chopin (2009) {barticle}[mr] \bauthor\bsnmRue, \bfnmHåvard\binitsH., \bauthor\bsnmMartino, \bfnmSara\binitsS. \AND\bauthor\bsnmChopin, \bfnmNicolas\binitsN. (\byear2009). \btitleApproximate Bayesian inference for latent Gaussian models by using integrated nested Laplace approximations. \bjournalJ. R. Stat. Soc. Ser. B Stat. Methodol. \bvolume71 \bpages319–392. \biddoi=10.1111/j.1467-9868.2008.00700.x, issn=1369-7412, mr=2649602 \bptnotecheck related\bptokimsref \endbibitem
  • Schoenberg (2005) {barticle}[mr] \bauthor\bsnmSchoenberg, \bfnmFrederic Paik\binitsF. P. (\byear2005). \btitleConsistent parametric estimation of the intensity of a spatial-temporal point process. \bjournalJ. Statist. Plann. Inference \bvolume128 \bpages79–93. \biddoi=10.1016/j.jspi.2003.09.027, issn=0378-3758, mr=2110179 \bptokimsref \endbibitem
  • Spiegelhalter et al. (2002) {barticle}[mr] \bauthor\bsnmSpiegelhalter, \bfnmDavid J.\binitsD. J., \bauthor\bsnmBest, \bfnmNicola G.\binitsN. G., \bauthor\bsnmCarlin, \bfnmBradley P.\binitsB. P. \AND\bauthor\bparticlevan der \bsnmLinde, \bfnmAngelika\binitsA. (\byear2002). \btitleBayesian measures of model complexity and fit. \bjournalJ. R. Stat. Soc. Ser. B Stat. Methodol. \bvolume64 \bpages583–639. \biddoi=10.1111/1467-9868.00353, issn=1369-7412, mr=1979380 \bptnotecheck related\bptokimsref \endbibitem
  • Steffan-Dewenter et al. (2002) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmSteffan-Dewenter, \bfnmI.\binitsI., \bauthor\bsnmMünzenberg, \bfnmU.\binitsU., \bauthor\bsnmThies, \bfnmC.\binitsC. \AND\bauthor\bsnmTscharntke, \bfnmT.\binitsT. (\byear2002). \btitleScale dependent effects of landscape context on three pollinator guilds. \bjournalEcology \bvolume83 \bpages1421–1432. \bptokimsref \endbibitem
  • Stoyan and Grabarnik (1991) {barticle}[mr] \bauthor\bsnmStoyan, \bfnmDietrich\binitsD. \AND\bauthor\bsnmGrabarnik, \bfnmPavel\binitsP. (\byear1991). \btitleSecond-order characteristics for stochastic structures connected with Gibbs point processes. \bjournalMath. Nachr. \bvolume151 \bpages95–100. \biddoi=10.1002/mana.19911510108, issn=0025-584X, mr=1121200 \bptokimsref \endbibitem
  • Stoyan, Kendall and Mecke (1995) {bbook}[mr] \bauthor\bsnmStoyan, \bfnmD.\binitsD., \bauthor\bsnmKendall, \bfnmW. S.\binitsW. S. \AND\bauthor\bsnmMecke, \bfnmJ.\binitsJ. (\byear1995). \btitleStochastic Geometry and Its Applications, \bedition2nd ed. \bpublisherWiley, \baddressLondon. \bptokimsref \endbibitem
  • Strauss (1975) {barticle}[mr] \bauthor\bsnmStrauss, \bfnmDavid J.\binitsD. J. (\byear1975). \btitleA model for clustering. \bjournalBiometrika \bvolume62 \bpages467–475. \bidissn=0006-3444, mr=0383493 \bptokimsref \endbibitem
  • van Lieshout (2000) {bbook}[mr] \bauthor\bparticlevan \bsnmLieshout, \bfnmM. N. M.\binitsM. N. M. (\byear2000). \btitleMarkov Point Processes and Their Applications. \bpublisherImperial College Press, \baddressLondon. \biddoi=10.1142/9781860949760, mr=1789230 \bptokimsref \endbibitem
  • Waagepetersen (2007) {barticle}[mr] \bauthor\bsnmWaagepetersen, \bfnmRasmus Plenge\binitsR. P. (\byear2007). \btitleAn estimating function approach to inference for inhomogeneous Neyman–Scott processes. \bjournalBiometrics \bvolume63 \bpages252–258, 315. \biddoi=10.1111/j.1541-0420.2006.00667.x, issn=0006-341X, mr=2345595 \bptokimsref \endbibitem
  • Waagepetersen and Guan (2009) {barticle}[mr] \bauthor\bsnmWaagepetersen, \bfnmRasmus\binitsR. \AND\bauthor\bsnmGuan, \bfnmYongtao\binitsY. (\byear2009). \btitleTwo-step estimation for inhomogeneous spatial point processes. \bjournalJ. R. Stat. Soc. Ser. B Stat. Methodol. \bvolume71 \bpages685–702. \biddoi=10.1111/j.1467-9868.2008.00702.x, issn=1369-7412, mr=2749914 \bptokimsref \endbibitem
  • Wiegand et al. (2007) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmWiegand, \bfnmT.\binitsT., \bauthor\bsnmGunatilleke, \bfnmS.\binitsS., \bauthor\bsnmGunatilleke, \bfnmN.\binitsN. \AND\bauthor\bsnmOkuda, \bfnmT.\binitsT. (\byear2007). \btitleAnalysing the spatial structure of a Sri Lankan tree species with multiple scales of clustering. \bjournalEcology \bvolume88 \bpages3088–3012. \bptokimsref \endbibitem 
  • Yue and Loh (2011) {barticle}[auto:STB—2012/03/21—07:41:58] \bauthor\bsnmYue, \bfnmY.\binitsY. \AND\bauthor\bsnmLoh, \bfnmJ. M.\binitsJ. M. (\byear2011). \btitleBayesian semiparametric intensity estimation for inhomogeneous spatial point processes. \bjournalBiometrics \bvolume67 \bpages937–946. \bptokimsref \endbibitem
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description