Package 'spatstat.model'

Description

The spatstat.model package belongs to the spatstat family of packages. It contains the core functionality for parametric statistical modelling of spatial data.

Details

spatstat is a family of R packages for the statistical analysis of spatial data. Its main focus is the analysis of spatial patterns of points in two-dimensional space.

The original spatstat package has now been split into several sub-packages.

This sub-package spatstat.model contains all the main user-level functions that perform parametric statistical modelling of spatial data.

(The main exception is that functions for linear networks are in the separate sub-package spatstat.linnet.)

Structure of the spatstat family

The orginal spatstat package grew to be very large. It has now been divided into several sub-packages:

• spatstat.utils containing basic utilities

• spatstat.sparse containing linear algebra utilities

• spatstat.data containing datasets

• spatstat.univar containing functions for estimating probability distributions of random variables

• spatstat.geom containing geometrical objects and geometrical operations

• spatstat.explore containing the functionality for exploratory analysis and nonparametric modelling of spatial data

• spatstat.model containing the main functionality for parametric modelling, analysis and inference for spatial data

• spatstat.linnet containing functions for spatial data on a linear network

• spatstat, which simply loads the other sub-packages listed above, and provides documentation.

When you install spatstat, these sub-packages are also installed. Then if you load the spatstat package by typing library(spatstat), the other sub-packages listed above will automatically be loaded or imported.

For an overview of all the functions available in the sub-packages of spatstat, see the help file for "spatstat-package" in the spatstat package.

Additionally there are several extension packages:

• spatstat.gui for interactive graphics

• spatstat.local for local likelihood (including geographically weighted regression)

• spatstat.Knet for additional, computationally efficient code for linear networks

• spatstat.sphere (under development) for spatial data on a sphere, including spatial data on the earth's surface

The extension packages must be installed separately and loaded explicitly if needed. They also have separate documentation.

Overview of Functionality in spatstat.model

The spatstat family of packages is designed to support a complete statistical analysis of spatial data. It supports

• creation, manipulation and plotting of point patterns;

• exploratory data analysis;

• spatial random sampling;

• simulation of point process models;

• parametric model-fitting;

• non-parametric smoothing and regression;

• formal inference (hypothesis tests, confidence intervals);

• model diagnostics.

For an overview, see the help file for "spatstat-package" in the spatstat package.

Following is a list of the functionality provided in the spatstat.model package only.

To simulate a random point pattern:

Functions for generating random point patterns are now contained in the spatstat.random package.

Exploratory analysis

Exploratory graphics, smoothing, and exploratory analysis of spatial data are now provided in the spatstat.explore package.

Model fitting (Cox and cluster models)

Cluster process models (with homogeneous or inhomogeneous intensity) and Cox processes can be fitted by the function kppm. Its result is an object of class "kppm". The fitted model can be printed, plotted, predicted, simulated and updated.

For model selection, you can also use the generic functions step, drop1 and AIC on fitted point process models. For variable selection, see sdr.

The theoretical models can also be simulated, for any choice of parameter values, using rThomas, rMatClust, rCauchy, rVarGamma, and rLGCP.

Lower-level fitting functions include:

Model fitting (Poisson and Gibbs models)

Poisson point processes are the simplest models for point patterns. A Poisson model assumes that the points are stochastically independent. It may allow the points to have a non-uniform spatial density. The special case of a Poisson process with a uniform spatial density is often called Complete Spatial Randomness.

Poisson point processes are included in the more general class of Gibbs point process models. In a Gibbs model, there is interaction or dependence between points. Many different types of interaction can be specified.

For a detailed explanation of how to fit Poisson or Gibbs point process models to point pattern data using spatstat, see Baddeley and Turner (2005b) or Baddeley (2008).

To fit a Poisson or Gibbs point process model:

Model fitting in spatstat is performed mainly by the function ppm. Its result is an object of class "ppm".

Here are some examples, where X is a point pattern (class "ppp"):

It is also possible to fit models that depend on other covariates.

Manipulating the fitted model:

For model selection, you can also use the generic functions step, drop1 and AIC on fitted point process models. For variable selection, see sdr.

See spatstat.options to control plotting of fitted model.

To specify a point process model:

The first order “trend” of the model is determined by an R language formula. The formula specifies the form of the logarithm of the trend.

The higher order (“interaction”) components are described by an object of class "interact". Such objects are created by:

Note that it is also possible to combine several such interactions using Hybrid.

Simulation and goodness-of-fit for fitted models:

Model fitting (determinantal point process models)

Code for fitting determinantal point process models has recently been added to spatstat.

For information, see the help file for dppm.

Model fitting (spatial logistic regression)

Pixel-based spatial logistic regression is an alternative technique for analysing spatial point patterns that is widely used in Geographical Information Systems. It is approximately equivalent to fitting a Poisson point process model.

In pixel-based logistic regression, the spatial domain is divided into small pixels, the presence or absence of a data point in each pixel is recorded, and logistic regression is used to model the presence/absence indicators as a function of any covariates.

Facilities for performing spatial logistic regression are provided in spatstat for comparison purposes.

Fitting a spatial logistic regression

Spatial logistic regression is performed by the function slrm. Its result is an object of class "slrm". There are many methods for this class, including methods for print, fitted, predict, simulate, anova, coef, logLik, terms, update, formula and vcov.

For example, if X is a point pattern (class "ppp"):

Manipulating a fitted spatial logistic regression

There are many other undocumented methods for this class, including methods for print, update, formula and terms. Stepwise model selection is possible using step or stepAIC. For variable selection, see sdr.

Simulation

There are many ways to generate a random point pattern, line segment pattern, pixel image or tessellation in spatstat.

Random point patterns: Functions for random generation are now contained in the spatstat.random package.

See also varblock for estimating the variance of a summary statistic by block resampling, and lohboot for another bootstrap technique.

Fitted point process models:

If you have fitted a point process model to a point pattern dataset, the fitted model can be simulated.

Cluster process models are fitted by the function kppm yielding an object of class "kppm". To generate one or more simulated realisations of this fitted model, use simulate.kppm.

Gibbs point process models are fitted by the function ppm yielding an object of class "ppm". To generate a simulated realisation of this fitted model, use rmh.ppm. To generate one or more simulated realisations of the fitted model, use simulate.ppm.

Other random patterns: Functions for random generation are now contained in the spatstat.random package.

Simulation-based inference

Simulation-based inference including simulation envelopes and hypothesis tests is now supported by the package spatstat.explore.

Sensitivity diagnostics:

Classical measures of model sensitivity such as leverage and influence have been adapted to point process models.

Diagnostics for covariate effect:

Classical diagnostics for covariate effects have been adapted to point process models.

Residual diagnostics:

Residuals for a fitted point process model, and diagnostic plots based on the residuals, were introduced in Baddeley et al (2005) and Baddeley, Rubak and Moller (2011).

Type demo(diagnose) for a demonstration of the diagnostics features.

Resampling and randomisation procedures

You can build your own tests based on randomisation and resampling using the following capabilities:

Licence

This library and its documentation are usable under the terms of the "GNU General Public License", a copy of which is distributed with the package.

Acknowledgements

Kasper Klitgaard Berthelsen, Ottmar Cronie, Tilman Davies, Julian Gilbey, Yongtao Guan, Ute Hahn, Kassel Hingee, Abdollah Jalilian, Marie-Colette van Lieshout, Greg McSwiggan, Tuomas Rajala, Suman Rakshit, Dominic Schuhmacher, Rasmus Waagepetersen and Hangsheng Wang made substantial contributions of code.

For comments, corrections, bug alerts and suggestions, we thank Monsuru Adepeju, Corey Anderson, Ang Qi Wei, Ryan Arellano, Jens Astrom, Robert Aue, Marcel Austenfeld, Sandro Azaele, Malissa Baddeley, Guy Bayegnak, Colin Beale, Melanie Bell, Thomas Bendtsen, Ricardo Bernhardt, Andrew Bevan, Brad Biggerstaff, Anders Bilgrau, Leanne Bischof, Christophe Biscio, Roger Bivand, Jose M. Blanco Moreno, Florent Bonneu, Jordan Brown, Ian Buller, Julian Burgos, Simon Byers, Ya-Mei Chang, Jianbao Chen, Igor Chernayavsky, Y.C. Chin, Bjarke Christensen, Lucia Cobo Sanchez, Jean-Francois Coeurjolly, Kim Colyvas, Hadrien Commenges, Rochelle Constantine, Robin Corria Ainslie, Richard Cotton, Marcelino de la Cruz, Peter Dalgaard, Mario D'Antuono, Sourav Das, Peter Diggle, Patrick Donnelly, Ian Dryden, Stephen Eglen, Ahmed El-Gabbas, Belarmain Fandohan, Olivier Flores, David Ford, Peter Forbes, Shane Frank, Janet Franklin, Funwi-Gabga Neba, Oscar Garcia, Agnes Gault, Jonas Geldmann, Marc Genton, Shaaban Ghalandarayeshi, Jason Goldstick, Pavel Grabarnik, C. Graf, Ute Hahn, Andrew Hardegen, Martin Bogsted Hansen, Martin Hazelton, Juha Heikkinen, Mandy Hering, Markus Herrmann, Maximilian Hesselbarth, Paul Hewson, Hamidreza Heydarian, Kurt Hornik, Philipp Hunziker, Jack Hywood, Ross Ihaka, Cenk Icos, Aruna Jammalamadaka, Robert John-Chandran, Devin Johnson, Mahdieh Khanmohammadi, Bob Klaver, Lily Kozmian-Ledward, Peter Kovesi, Mike Kuhn, Jeff Laake, Robert Lamb, Frederic Lavancier, Tom Lawrence, Tomas Lazauskas, Jonathan Lee, George Leser, Angela Li, Li Haitao, George Limitsios, Andrew Lister, Nestor Luambua, Ben Madin, Martin Maechler, Kiran Marchikanti, Jeff Marcus, Robert Mark, Peter McCullagh, Monia Mahling, Jorge Mateu Mahiques, Ulf Mehlig, Frederico Mestre, Sebastian Wastl Meyer, Mi Xiangcheng, Lore De Middeleer, Robin Milne, Enrique Miranda, Jesper Moller, Annie Mollie, Ines Moncada, Mehdi Moradi, Virginia Morera Pujol, Erika Mudrak, Gopalan Nair, Nader Najari, Nicoletta Nava, Linda Stougaard Nielsen, Felipe Nunes, Jens Randel Nyengaard, Jens Oehlschlaegel, Thierry Onkelinx, Sean O'Riordan, Evgeni Parilov, Jeff Picka, Nicolas Picard, Tim Pollington, Mike Porter, Sergiy Protsiv, Adrian Raftery, Ben Ramage, Pablo Ramon, Xavier Raynaud, Nicholas Read, Matt Reiter, Ian Renner, Tom Richardson, Brian Ripley, Ted Rosenbaum, Barry Rowlingson, Jason Rudokas, Tyler Rudolph, John Rudge, Christopher Ryan, Farzaneh Safavimanesh, Aila Sarkka, Cody Schank, Katja Schladitz, Sebastian Schutte, Bryan Scott, Olivia Semboli, Francois Semecurbe, Vadim Shcherbakov, Shen Guochun, Shi Peijian, Harold-Jeffrey Ship, Tammy L Silva, Ida-Maria Sintorn, Yong Song, Malte Spiess, Mark Stevenson, Kaspar Stucki, Jan Sulavik, Michael Sumner, P. Surovy, Ben Taylor, Thordis Linda Thorarinsdottir, Leigh Torres, Berwin Turlach, Torben Tvedebrink, Kevin Ummer, Medha Uppala, Andrew van Burgel, Tobias Verbeke, Mikko Vihtakari, Alexendre Villers, Fabrice Vinatier, Maximilian Vogtland, Sasha Voss, Sven Wagner, Hao Wang, H. Wendrock, Jan Wild, Carl G. Witthoft, Selene Wong, Maxime Woringer, Luke Yates, Mike Zamboni and Achim Zeileis.

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

Description

Computes the coordinates for an Added Variable Plot for a fitted point process model.

Usage

addvar(model, covariate, ...,
                   subregion=NULL,
                   bw="nrd0", adjust=1,
                   from=NULL, to=NULL, n=512,
                   bw.input = c("points", "quad"),
                   bw.restrict = FALSE,
                   covname, crosscheck=FALSE)

Arguments

Details

This command generates the plot coordinates for an Added Variable Plot for a spatial point process model.

Added Variable Plots (Cox, 1958, sec 4.5; Wang, 1985) are commonly used in linear models and generalized linear models, to decide whether a model with response $y$ and predictors $x$ would be improved by including another predictor $z$.

In a (generalised) linear model with response $y$ and predictors $x$, the Added Variable Plot for a new covariate $z$ is a plot of the smoothed Pearson residuals from the original model against the scaled residuals from a weighted linear regression of $z$ on $x$. If this plot has nonzero slope, then the new covariate $z$ is needed. For general advice see Cook and Weisberg(1999); Harrell (2001).

Essentially the same technique can be used for a spatial point process model (Baddeley et al, 2012).

The argument model should be a fitted spatial point process model (object of class "ppm").

The argument covariate identifies the covariate that is to be considered for addition to the model. It should be either a pixel image (object of class "im") or a function(x,y) giving the values of the covariate at any spatial location. Alternatively covariate may be a character string, giving the name of a covariate that was supplied (in the covariates argument to ppm) when the model was fitted, but was not used in the model.

The result of addvar(model, covariate) is an object belonging to the classes "addvar" and "fv". Plot this object to generate the added variable plot.

Note that the plot method shows the pointwise significance bands for a test of the null model, i.e. the null hypothesis that the new covariate has no effect.

The smoothing bandwidth is controlled by the arguments bw, adjust, bw.input and bw.restrict. If bw is a numeric value, then the bandwidth is taken to be adjust * bw. If bw is a string representing a bandwidth selection rule (recognised by density.default) then the bandwidth is selected by this rule.

The data used for automatic bandwidth selection are specified by bw.input and bw.restrict. If bw.input="points" (the default) then bandwidth selection is based on the covariate values at the points of the original point pattern dataset to which the model was fitted. If bw.input="quad" then bandwidth selection is based on the covariate values at every quadrature point used to fit the model. If bw.restrict=TRUE then the bandwidth selection is performed using only data from inside the subregion.

Value

An object of class "addvar" containing the coordinates for the added variable plot. There is a plot method.

Slow computation

In a large dataset, computation can be very slow if the default settings are used, because the smoothing bandwidth is selected automatically. To avoid this, specify a numerical value for the bandwidth bw. One strategy is to use a coarser subset of the data to select bw automatically. The selected bandwidth can be read off the print output for addvar.

Internal data

The return value has an attribute "spatial" which contains the internal data: the computed values of the residuals, and of all relevant covariates, at each quadrature point of the model. It is an object of class "ppp" with a data frame of marks.

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected], Ya-Mei Chang and Yong Song.

References

Baddeley, A., Chang, Y.-M., Song, Y. and Turner, R. (2013) Residual diagnostics for covariate effects in spatial point process models. Journal of Computational and Graphical Statistics, 22, 886–905.

Cook, R.D. and Weisberg, S. (1999) Applied regression, including computing and graphics. New York: Wiley.

Cox, D.R. (1958) Planning of Experiments. New York: Wiley.

Harrell, F. (2001) Regression Modeling Strategies. New York: Springer.

Wang, P. (1985) Adding a variable in generalized linear models. Technometrics 27, 273–276.

See Also

parres, rhohat, rho2hat.

Examples

X <-  rpoispp(function(x,y){exp(3+3*x)})
  model <- ppm(X, ~y)
  adv <- addvar(model, "x")
  plot(adv)
  adv <- addvar(model, "x", subregion=square(0.5))

Description

Performs analysis of deviance for one or more point process models fitted to replicated point pattern data.

Usage

## S3 method for class 'mppm'
anova(object, ...,
                  test=NULL, adjust=TRUE,
                  fine=FALSE, warn=TRUE)

Arguments

Details

This is a method for anova for comparing several fitted point process models of class "mppm", usually generated by the model-fitting function mppm).

If the fitted models are all Poisson point processes, then this function performs an Analysis of Deviance of the fitted models. The output shows the deviance differences (i.e. 2 times log likelihood ratio), the difference in degrees of freedom, and (if test="Chi") the two-sided p-values for the chi-squared tests. Their interpretation is very similar to that in anova.glm.

If some of the fitted models are not Poisson point processes, the ‘deviance’ differences in this table are 'pseudo-deviances' equal to 2 times the differences in the maximised values of the log pseudolikelihood (see ppm). It is not valid to compare these values to the chi-squared distribution. In this case, if adjust=TRUE (the default), the pseudo-deviances will be adjusted using the method of Pace et al (2011) and Baddeley, Turner and Rubak (2015) so that the chi-squared test is valid. It is strongly advisable to perform this adjustment.

The argument test determines which hypothesis test, if any, will be performed to compare the models. The argument test should be a character string, partially matching one of "Chisq", "F" or "Cp", or NULL. The first option "Chisq" gives the likelihood ratio test based on the asymptotic chi-squared distribution of the deviance difference. The meaning of the other options is explained in anova.glm.

Value

An object of class "anova", or NULL.

Random effects models are currently not supported

For models with random effects (i.e. where the call to mppm included the argument random), analysis of deviance is currently not supported, due to changes in the nlme package. We will try to find a solution.

Error messages

An error message that reports system is computationally singular indicates that the determinant of the Fisher information matrix of one of the models was either too large or too small for reliable numerical calculation. See vcov.ppm for suggestions on how to handle this.

Author(s)

Adrian Baddeley, Ida-Maria Sintorn and Leanne Bischoff. Implemented by Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

References

Baddeley, A., Rubak, E. and Turner, R. (2015) Spatial Point Patterns: Methodology and Applications with R. Chapman and Hall/CRC Press.

Baddeley, A., Turner, R. and Rubak, E. (2015) Adjusted composite likelihood ratio test for Gibbs point processes. Journal of Statistical Computation and Simulation 86 (5) 922–941. DOI: 10.1080/00949655.2015.1044530.

Pace, L., Salvan, A. and Sartori, N. (2011) Adjusting composite likelihood ratio statistics. Statistica Sinica 21, 129–148.

See Also

mppm

Examples

H <- hyperframe(X=waterstriders)
 #' test for loglinear trend in x coordinate
 mod0 <- mppm(X~1, data=H, Poisson())
 modx <- mppm(X~x, data=H, Poisson())
 anova(mod0, modx, test="Chi")
 # not significant
 anova(modx, test="Chi")
 # not significant

 #' test for inhibition
 mod0S <- mppm(X~1, data=H, Strauss(2))
 anova(mod0, mod0S, test="Chi")
 # significant! 

 #' test for trend after accounting for inhibition
 modxS <- mppm(X~x, data=H, Strauss(2))
 anova(mod0S, modxS, test="Chi")
 # not significant

Description

Performs analysis of deviance for one or more fitted point process models.

Usage

## S3 method for class 'ppm'
anova(object, ..., test=NULL,
                      adjust=TRUE, warn=TRUE, fine=FALSE)

Arguments

Details

This is a method for anova for fitted point process models (objects of class "ppm", usually generated by the model-fitting function ppm).

If the fitted models are all Poisson point processes, then by default, this function performs an Analysis of Deviance of the fitted models. The output shows the deviance differences (i.e. 2 times log likelihood ratio), the difference in degrees of freedom, and (if test="Chi" or test="LRT") the two-sided p-values for the chi-squared tests. Their interpretation is very similar to that in anova.glm. If test="Rao" or test="score", the score test (Rao, 1948) is performed instead.

If some of the fitted models are not Poisson point processes, the ‘deviance’ differences in this table are 'pseudo-deviances' equal to 2 times the differences in the maximised values of the log pseudolikelihood (see ppm). It is not valid to compare these values to the chi-squared distribution. In this case, if adjust=TRUE (the default), the pseudo-deviances will be adjusted using the method of Pace et al (2011) and Baddeley et al (2015) so that the chi-squared test is valid. It is strongly advisable to perform this adjustment.

Value

An object of class "anova", or NULL.

Errors and warnings

models not nested:

There may be an error message that the models are not “nested”. For an Analysis of Deviance the models must be nested, i.e. one model must be a special case of the other. For example the point process model with formula ~x is a special case of the model with formula ~x+y, so these models are nested. However the two point process models with formulae ~x and ~y are not nested.

If you get this error message and you believe that the models should be nested, the problem may be the inability of R to recognise that the two formulae are nested. Try modifying the formulae to make their relationship more obvious.

different sizes of dataset:

There may be an error message from anova.glmlist that “models were not all fitted to the same size of dataset”. This implies that the models were fitted using different quadrature schemes (see quadscheme) and/or with different edge corrections or different values of the border edge correction distance rbord.

To ensure that models are comparable, check the following:

• the models must all have been fitted to the same point pattern dataset, in the same window.

• all models must have been fitted by the same fitting method as specified by the argument method in ppm.

• If some of the models depend on covariates, then they should all have been fitted using the same list of covariates, and using allcovar=TRUE to ensure that the same quadrature scheme is used.

• all models must have been fitted using the same edge correction as specified by the arguments correction and rbord. If you did not specify the value of rbord, then it may have taken a different value for different models. The default value of rbord is equal to zero for a Poisson model, and otherwise equals the reach (interaction distance) of the interaction term (see reach). To ensure that the models are comparable, set rbord to equal the maximum reach of the interactions that you are fitting.

Error messages

An error message that reports system is computationally singular indicates that the determinant of the Fisher information matrix of one of the models was either too large or too small for reliable numerical calculation. See vcov.ppm for suggestions on how to handle this.

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

References

Baddeley, A., Turner, R. and Rubak, E. (2015) Adjusted composite likelihood ratio test for Gibbs point processes. Journal of Statistical Computation and Simulation 86 (5) 922–941. DOI: 10.1080/00949655.2015.1044530.

Pace, L., Salvan, A. and Sartori, N. (2011) Adjusting composite likelihood ratio statistics. Statistica Sinica 21, 129–148.

Rao, C.R. (1948) Large sample tests of statistical hypotheses concerning several parameters with applications to problems of estimation. Proceedings of the Cambridge Philosophical Society 44, 50–57.

See Also

ppm, vcov.ppm

Examples

mod0 <- ppm(swedishpines ~1)
 modx <- ppm(swedishpines ~x)
 # Likelihood ratio test
 anova(mod0, modx, test="Chi")
 # Score test
 anova(mod0, modx, test="Rao")

 # Single argument
 modxy <- ppm(swedishpines ~x + y)
 anova(modxy, test="Chi")

 # Adjusted composite likelihood ratio test
 modP <- ppm(swedishpines ~1, rbord=9)
 modS <- ppm(swedishpines ~1, Strauss(9))
 anova(modP, modS, test="Chi")

Description

Performs Analysis of Deviance for two or more fitted Spatial Logistic Regression models.

Usage

## S3 method for class 'slrm'
anova(object, ..., test = NULL)

Arguments

Details

This is a method for anova for fitted spatial logistic regression models (objects of class "slrm", usually obtained from the function slrm).

The output shows the deviance differences (i.e. 2 times log likelihood ratio), the difference in degrees of freedom, and (if test="Chi") the two-sided $p$-values for the chi-squared tests. Their interpretation is very similar to that in anova.glm.

Value

An object of class "anova", inheriting from class "data.frame", representing the analysis of deviance table.

Author(s)

Adrian Baddeley [email protected] and Rolf Turner [email protected]

See Also

slrm

Examples

X <- rpoispp(42)
  fit0 <- slrm(X ~ 1)
  fit1 <- slrm(X ~ x+y)
  anova(fit0, fit1, test="Chi")

Description

Creates an instance of the Area Interaction point process model (Widom-Rowlinson penetrable spheres model) which can then be fitted to point pattern data.

Usage

AreaInter(r)

Arguments

Details

This function defines the interpoint interaction structure of a point process called the Widom-Rowlinson penetrable sphere model or area-interaction process. It can be used to fit this model to point pattern data.

The function ppm(), which fits point process models to point pattern data, requires an argument of class "interact" describing the interpoint interaction structure of the model to be fitted. The appropriate description of the area interaction structure is yielded by the function AreaInter(). See the examples below.

In standard form, the area-interaction process (Widom and Rowlinson, 1970; Baddeley and Van Lieshout, 1995) with disc radius $r$, intensity parameter $\kappa$ and interaction parameter $\gamma$ is a point process with probability density

$f(x_1,\ldots,x_n) = \alpha \kappa^{n(x)} \gamma^{-A(x)}$

for a point pattern $x$, where $x_1,\ldots,x_n$ represent the points of the pattern, $n(x)$ is the number of points in the pattern, and $A(x)$ is the area of the region formed by the union of discs of radius $r$ centred at the points $x_1,\ldots,x_n$. Here $\alpha$ is a normalising constant.

The interaction parameter $\gamma$ can be any positive number. If $\gamma = 1$ then the model reduces to a Poisson process with intensity $\kappa$. If $\gamma < 1$ then the process is regular, while if $\gamma > 1$ the process is clustered. Thus, an area interaction process can be used to model either clustered or regular point patterns. Two points interact if the distance between them is less than $2r$.

The standard form of the model, shown above, is a little complicated to interpret in practical applications. For example, each isolated point of the pattern $x$ contributes a factor $\kappa \gamma^{-\pi r^2}$ to the probability density.

In spatstat, the model is parametrised in a different form, which is easier to interpret. In canonical scale-free form, the probability density is rewritten as

$f(x_1,\ldots,x_n) = \alpha \beta^{n(x)} \eta^{-C(x)}$

where $\beta$ is the new intensity parameter, $\eta$ is the new interaction parameter, and $C(x) = B(x) - n(x)$ is the interaction potential. Here

$B(x) = \frac{A(x)}{\pi r^2}$

is the normalised area (so that the discs have unit area). In this formulation, each isolated point of the pattern contributes a factor $\beta$ to the probability density (so the first order trend is $\beta$). The quantity $C(x)$ is a true interaction potential, in the sense that $C(x) = 0$ if the point pattern $x$ does not contain any points that lie close together (closer than $2r$ units apart).

When a new point $u$ is added to an existing point pattern $x$, the rescaled potential $-C(x)$ increases by a value between 0 and 1. The increase is zero if $u$ is not close to any point of $x$. The increase is 1 if the disc of radius $r$ centred at $u$ is completely contained in the union of discs of radius $r$ centred at the data points $x_i$. Thus, the increase in potential is a measure of how close the new point $u$ is to the existing pattern $x$. Addition of the point $u$ contributes a factor $\beta \eta^\delta$ to the probability density, where $\delta$ is the increase in potential.

The old parameters $\kappa,\gamma$ of the standard form are related to the new parameters $\beta,\eta$ of the canonical scale-free form, by

$\beta = \kappa \gamma^{-\pi r^2} = \kappa /\eta$

and

$\eta = \gamma^{\pi r^2}$

provided $\gamma$ and $\kappa$ are positive and finite.

In the canonical scale-free form, the parameter $\eta$ can take any nonnegative value. The value $\eta = 1$ again corresponds to a Poisson process, with intensity $\beta$. If $\eta < 1$ then the process is regular, while if $\eta > 1$ the process is clustered. The value $\eta = 0$ corresponds to a hard core process with hard core radius $r$ (interaction distance $2r$).

The nonstationary area interaction process is similar except that the contribution of each individual point $x_i$ is a function $\beta(x_i)$ of location, rather than a constant beta.

Note the only argument of AreaInter() is the disc radius r. When r is fixed, the model becomes an exponential family. The canonical parameters $\log(\beta)$ and $\log(\eta)$ are estimated by ppm(), not fixed in AreaInter().

Value

An object of class "interact" describing the interpoint interaction structure of the area-interaction process with disc radius $r$.

Warnings

The interaction distance of this process is equal to 2 * r. Two discs of radius r overlap if their centres are closer than 2 * r units apart.

The estimate of the interaction parameter $\eta$ is unreliable if the interaction radius r is too small or too large. In these situations the model is approximately Poisson so that $\eta$ is unidentifiable. As a rule of thumb, one can inspect the empty space function of the data, computed by Fest. The value $F(r)$ of the empty space function at the interaction radius r should be between 0.2 and 0.8.

Author(s)

Adrian Baddeley [email protected] and Rolf Turner [email protected]

References

Baddeley, A.J. and Van Lieshout, M.N.M. (1995). Area-interaction point processes. Annals of the Institute of Statistical Mathematics 47 (1995) 601–619.

Widom, B. and Rowlinson, J.S. (1970). New model for the study of liquid-vapor phase transitions. The Journal of Chemical Physics 52 (1970) 1670–1684.

See Also

ppm, pairwise.family, ppm.object

ragsAreaInter and rmh for simulation of area-interaction models.

Examples

# prints a sensible description of itself
   AreaInter(r=0.1)

   # Note the reach is twice the radius
   reach(AreaInter(r=1))

   # Fit the stationary area interaction process to Swedish Pines data
   ppm(swedishpines ~1, AreaInter(r=7))

   # Fit the stationary area interaction process to `cells'
   ppm(cells ~1, AreaInter(r=0.06))
   # eta=0 indicates hard core process.

   # Fit a nonstationary area interaction with log-cubic polynomial trend
   
     ppm(swedishpines ~polynom(x/10,y/10,3), AreaInter(r=7))

Description

Converts an object of class "leverage.ppm" to a function of the $x$ and $y$ coordinates.

Usage

## S3 method for class 'leverage.ppm'
as.function(x, ...)

Arguments

Details

An object of class "leverage.ppm" represents the leverage function of a fitted point process model. This command converts the object to a function(x,y) where the arguments x and y are (vectors of) spatial coordinates. This function returns the leverage values at the specified locations (calculated by referring to the nearest location where the leverage has been computed).

Value

A function in the R language, also belonging to the class "funxy".

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

See Also

as.im.leverage.ppm

Examples

X <- rpoispp(function(x,y) { exp(3+3*x) })
  fit <- ppm(X ~x+y)
  lev <- leverage(fit)
  f <- as.function(lev)
  
  f(0.2, 0.3)  # evaluate at (x,y) coordinates
  y <- f(X)    # evaluate at a point pattern

Description

Converts fitted model into a function table (an object of class "fv").

Usage

## S3 method for class 'kppm'
as.fv(x)

  ## S3 method for class 'dppm'
as.fv(x)

  ## S3 method for class 'minconfit'
as.fv(x)

Arguments

Details

The generic command as.fv converts data x, that could be interpreted as the values of a function, into a function value table (object of the class "fv" as described in fv.object). This object can then be plotted easily using plot.fv.

Objects of class "kppm" (and related classes) represent a model that has been fitted to a dataset by computing a summary function of the dataset and matching it to the corresponding summary function of the model. The methods for as.fv for classes "kppm", "dppm" and "minconfit" extract this information: the result is a function table containing the observed summary function and the best fit summary function.

Value

An object of class "fv" (see fv.object).

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected]

Examples

as.fv(kppm(redwood))

Description

Extracts the interpoint interaction structure from a point pattern model.

Usage

as.interact(object)
## S3 method for class 'fii'
as.interact(object)
## S3 method for class 'interact'
as.interact(object)
## S3 method for class 'ppm'
as.interact(object)

Arguments

Details

The function as.interact extracts the interpoint interaction structure from a suitable object.

An object of class "interact" describes an interpoint interaction structure, before it has been fitted to point pattern data. The irregular parameters of the interaction (such as the interaction range) are fixed, but the regular parameters (such as interaction strength) are undetermined. Objects of this class are created by the functions Poisson, Strauss and so on. The main use of such objects is in a call to ppm.

The function as.interact is generic, with methods for the classes "ppm", "fii" and "interact". The result is an object of class "interact" which can be printed.

Value

An object of class "interact" representing the interpoint interaction. This object can be printed and plotted.

Note on parameters

This function does not extract the fitted coefficients of the interaction. To extract the fitted interaction including the fitted coefficients, use fitin.

Author(s)

Adrian Baddeley [email protected] and Rolf Turner [email protected]

See Also

fitin, ppm.

Examples

model <- ppm(cells ~1, Strauss(0.07))
   f <- as.interact(model)
   f

Description

Converts a measure into a layered object.

Usage

## S3 method for class 'msr'
as.layered(X)

Arguments

Details

This function converts the object X into an object of class "layered".

It is a method for the generic as.layered for the class of measures.

If X is a vector-valued measure, then as.layered(X) consists of several layers, each containing a scalar-valued measure.

Value

An object of class "layered" (see layered).

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

See Also

as.layered, msr.

Examples

P <- rpoispp(100)
   fit <- ppm(P ~ x+y)
   rs <- residuals(fit, type="score")
   as.layered(rs)

Description

Converts data specifying an observation window in any of several formats, into an object of class "owin".

Usage

## S3 method for class 'ppm'
as.owin(W, ..., from=c("points", "covariates"), fatal=TRUE)

 ## S3 method for class 'kppm'
as.owin(W, ..., from=c("points", "covariates"), fatal=TRUE)

 ## S3 method for class 'dppm'
as.owin(W, ..., from=c("points", "covariates"), fatal=TRUE)

 ## S3 method for class 'slrm'
as.owin(W, ..., from=c("points", "covariates"))

 ## S3 method for class 'msr'
as.owin(W, ..., fatal=TRUE)

Arguments

Details

The class "owin" is a way of specifying the observation window for a point pattern. See owin.object for an overview.

The generic function as.owin converts data in any of several formats into an object of class "owin" for use by the spatstat package. The function as.owin is generic, with methods for different classes of objects, and a default method.

The argument W may be

• an object of class "owin"

• a structure with entries xrange, yrange specifying the $x$ and $y$ dimensions of a rectangle

• a structure with entries named xmin, xmax, ymin, ymax (in any order) specifying the $x$ and $y$ dimensions of a rectangle. This will accept objects of class bbox in the sf package.

• a numeric vector of length 4 (interpreted as (xmin, xmax, ymin, ymax) in that order) specifying the $x$ and $y$ dimensions of a rectangle

• a structure with entries named xl, xu, yl, yu (in any order) specifying the $x$ and $y$ dimensions of a rectangle as (xmin, xmax) = (xl, xu) and (ymin, ymax) = (yl, yu). This will accept objects of class spp used in the Venables and Ripley spatial package.

• an object of class "ppp" representing a point pattern. In this case, the object's window structure will be extracted.

• an object of class "psp" representing a line segment pattern. In this case, the object's window structure will be extracted.

• an object of class "tess" representing a tessellation. In this case, the object's window structure will be extracted.

• an object of class "quad" representing a quadrature scheme. In this case, the window of the data component will be extracted.

• an object of class "im" representing a pixel image. In this case, a window of type "mask" will be returned, with the same pixel raster coordinates as the image. An image pixel value of NA, signifying that the pixel lies outside the window, is transformed into the logical value FALSE, which is the corresponding convention for window masks.

• an object of class "ppm", "kppm", "slrm" or "dppm" representing a fitted point process model. In this case, if from="data" (the default), as.owin extracts the original point pattern data to which the model was fitted, and returns the observation window of this point pattern. If from="covariates" then as.owin extracts the covariate images to which the model was fitted, and returns a binary mask window that specifies the pixel locations.

• an object of class "lpp" representing a point pattern on a linear network. In this case, as.owin extracts the linear network and returns a window containing this network.

• an object of class "lppm" representing a fitted point process model on a linear network. In this case, as.owin extracts the linear network and returns a window containing this network.

• A data.frame with exactly three columns. Each row of the data frame corresponds to one pixel. Each row contains the $x$ and $y$ coordinates of a pixel, and a logical value indicating whether the pixel lies inside the window.

• A data.frame with exactly two columns. Each row of the data frame contains the $x$ and $y$ coordinates of a pixel that lies inside the window.

• an object of class "distfun", "nnfun" or "funxy" representing a function of spatial location, defined on a spatial domain. The spatial domain of the function will be extracted.

• an object of class "rmhmodel" representing a point process model that can be simulated using rmh. The window (spatial domain) of the model will be extracted. The window may be NULL in some circumstances (indicating that the simulation window has not yet been determined). This is not treated as an error, because the argument fatal defaults to FALSE for this method.

• an object of class "layered" representing a list of spatial objects. See layered. In this case, as.owin will be applied to each of the objects in the list, and the union of these windows will be returned.

• an object of some other suitable class from another package. For full details, see vignette('shapefiles').

If the argument W is not in one of these formats and cannot be converted to a window, then an error will be generated (if fatal=TRUE) or a value of NULL will be returned (if fatal=FALSE).

When W is a data frame, the argument step can be used to specify the pixel grid spacing; otherwise, the spacing will be guessed from the data.

Value

An object of class "owin" (see owin.object) specifying an observation window.

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

See Also

as.owin, as.owin.rmhmodel, as.owin.lpp.

owin.object, owin.

Additional methods for as.owin may be provided by other packages outside the spatstat family.

Examples

fit <- ppm(cells ~ 1)
  as.owin(fit)

Description

Extracts the fitted point process model from some kind of fitted model.

Usage

as.ppm(object)

## S3 method for class 'ppm'
as.ppm(object)

## S3 method for class 'profilepl'
as.ppm(object)

## S3 method for class 'kppm'
as.ppm(object)

## S3 method for class 'dppm'
as.ppm(object)

## S3 method for class 'rppm'
as.ppm(object)

Arguments

Details

The function as.ppm extracts the fitted point process model (of class "ppm") from a suitable object.

The function as.ppm is generic, with methods for the classes "ppm", "profilepl", "kppm", "dppm" and "rppm", and possibly for other classes.

For the class "profilepl" of models fitted by maximum profile pseudolikelihood, the method as.ppm.profilepl extracts the fitted point process model (with the optimal values of the irregular parameters).

For the class "kppm" of models fitted by minimum contrast (or Palm or composite likelihood) using Waagepetersen's two-step estimation procedure (see kppm), the method as.ppm.kppm extracts the Poisson point process model that is fitted in the first stage of the procedure.

The behaviour for the class "dppm" is analogous to the "kppm" case above.

For the class "rppm" of models fitted by recursive partitioning (regression trees), the method as.ppm.rppm extracts the corresponding loglinear model that is fitted in the first stage of the procedure (whose purpose is merely to identify and evaluate the explanatory variables).

Value

An object of class "ppm".

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

See Also

ppm, profilepl.

Examples

# fit a model by profile maximum pseudolikelihood
   rvals <- data.frame(r=(1:10)/100)
   pfit <- profilepl(rvals, Strauss, cells, ~1)
   # extract the fitted model
   fit <- as.ppm(pfit)

Description

Compute the AUC (area under the Receiver Operating Characteristic curve) for a fitted point process model.

Usage

## S3 method for class 'ppm'
auc(X, ...)

## S3 method for class 'kppm'
auc(X, ...)

## S3 method for class 'slrm'
auc(X, ...)

Arguments

Details

This command computes the AUC, the area under the Receiver Operating Characteristic curve. The ROC itself is computed by roc.

For a fitted point process model X, the AUC measures the ability of the fitted model intensity to separate the spatial domain into areas of high and low density of points. Suppose $\lambda(u)$ is the intensity function of the model. The AUC is the probability that $\lambda(x_i) > \lambda(U)$. That is, AUC is the probability that a randomly-selected data point has higher predicted intensity than does a randomly-selected spatial location. The AUC is not a measure of the goodness-of-fit of the model (Lobo et al, 2007).

(For spatial logistic regression models (class "slrm") replace “intensity” by “probability of presence” in the text above.)

Value

Numeric. For auc.ppm, auc.kppm and auc.lppm, the result is a numeric vector of length 2 giving the AUC value and the theoretically expected AUC value for this model.

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

References

Lobo, J.M., Jimenez-Valverde, A. and Real, R. (2007) AUC: a misleading measure of the performance of predictive distribution models. Global Ecology and Biogeography 17(2) 145–151.

Nam, B.-H. and D'Agostino, R. (2002) Discrimination index, the area under the ROC curve. Pages 267–279 in Huber-Carol, C., Balakrishnan, N., Nikulin, M.S. and Mesbah, M., Goodness-of-fit tests and model validity, Birkhauser, Basel.

See Also

roc

Examples

fit <- ppm(swedishpines ~ x+y)
  auc(fit)

Description

Creates an instance of the Baddeley-Geyer point process model, defined as a hybrid of several Geyer interactions. The model can then be fitted to point pattern data.

Usage

BadGey(r, sat)

Arguments

Details

This is Baddeley's generalisation of the Geyer saturation point process model, described in Geyer, to a process with multiple interaction distances.

The BadGey point process with interaction radii $r_1,\ldots,r_k$, saturation thresholds $s_1,\ldots,s_k$, intensity parameter $\beta$ and interaction parameters $\gamma_1,\ldots,gamma_k$, is the point process in which each point $x_i$ in the pattern $X$ contributes a factor

$\beta \gamma_1^{v_1(x_i, X)} \ldots gamma_k^{v_k(x_i,X)}$

to the probability density of the point pattern, where

$v_j(x_i, X) = \min( s_j, t_j(x_i,X) )$

where $t_j(x_i, X)$ denotes the number of points in the pattern $X$ which lie within a distance $r_j$ from the point $x_i$.

BadGey is used to fit this model to data. The function ppm(), which fits point process models to point pattern data, requires an argument of class "interact" describing the interpoint interaction structure of the model to be fitted. The appropriate description of the piecewise constant Saturated pairwise interaction is yielded by the function BadGey(). See the examples below.

The argument r specifies the vector of interaction distances. The entries of r must be strictly increasing, positive numbers.

The argument sat specifies the vector of saturation parameters that are applied to the point counts $t_j(x_i, X)$. It should be a vector of the same length as r, and its entries should be nonnegative numbers. Thus sat[1] is applied to the count of points within a distance r[1], and sat[2] to the count of points within a distance r[2], etc. Alternatively sat may be a single number, and this saturation value will be applied to every count.

Infinite values of the saturation parameters are also permitted; in this case $v_j(x_i,X) = t_j(x_i,X)$ and there is effectively no ‘saturation’ for the distance range in question. If all the saturation parameters are set to Inf then the model is effectively a pairwise interaction process, equivalent to PairPiece (however the interaction parameters $\gamma$ obtained from BadGey have a complicated relationship to the interaction parameters $\gamma$ obtained from PairPiece).

If r is a single number, this model is virtually equivalent to the Geyer process, see Geyer.

Value

An object of class "interact" describing the interpoint interaction structure of a point process.

Hybrids

A ‘hybrid’ interaction is one which is built by combining several different interactions (Baddeley et al, 2013). The BadGey interaction can be described as a hybrid of several Geyer interactions.

The Hybrid command can be used to build hybrids of any interactions. If the Hybrid operator is applied to several Geyer models, the result is equivalent to a BadGey model. This can be useful for incremental model selection.

Author(s)

Adrian Baddeley [email protected] and Rolf Turner [email protected] in collaboration with Hao Wang and Jeff Picka

References

Baddeley, A., Turner, R., Mateu, J. and Bevan, A. (2013) Hybrids of Gibbs point process models and their implementation. Journal of Statistical Software 55:11, 1–43. DOI: 10.18637/jss.v055.i11

See Also

ppm, pairsat.family, Geyer, PairPiece, SatPiece, Hybrid

Examples

BadGey(c(0.1,0.2), c(1,1))
   # prints a sensible description of itself
   BadGey(c(0.1,0.2), 1)

   # fit a stationary Baddeley-Geyer model
   ppm(cells ~1, BadGey(c(0.07, 0.1, 0.13), 2))

   # nonstationary process with log-cubic polynomial trend
   
     ppm(cells ~polynom(x,y,3), BadGey(c(0.07, 0.1, 0.13), 2))

Description

Applies a first-order bias correction to a fitted model.

Usage

bc(fit, ...)

  ## S3 method for class 'ppm'
bc(fit, ..., nfine = 256)

Arguments

Details

This command applies the first order Newton-Raphson bias correction method of Baddeley and Turner (2014, sec 4.2) to a fitted model. The function bc is generic, with a method for fitted point process models of class "ppm".

A fine grid of locations, of dimensions nfine * nfine or nfine[2] * nfine[1], is created over the original window of the data, and the intensity or conditional intensity of the fitted model is calculated on this grid. The result is used to update the fitted model parameters once by a Newton-Raphson update.

This is only useful if the quadrature points used to fit the original model fit are coarser than the grid of points specified by nfine.

Value

A numeric vector, of the same length as coef(fit), giving updated values for the fitted model coefficients.

Author(s)

Adrian Baddeley [email protected] and Rolf Turner [email protected].

References

Baddeley, A. and Turner, R. (2014) Bias correction for parameter estimates of spatial point process models. Journal of Statistical Computation and Simulation 84, 1621–1643. DOI: 10.1080/00949655.2012.755976

See Also

rex

Examples

fit <- ppm(cells ~ x, Strauss(0.07))
  coef(fit)
  if(!interactive()) {
    bc(fit, nfine=64)
  } else {
    bc(fit)
  }

Description

Tests the goodness-of-fit of a Poisson point process model using methods of Berman (1986).

Usage

## S3 method for class 'ppm'
berman.test(model, covariate,
                         which = c("Z1", "Z2"),
               alternative = c("two.sided", "less", "greater"), ...)

Arguments

Details

These functions perform a goodness-of-fit test of a Poisson point process model fitted to point pattern data. The observed distribution of the values of a spatial covariate at the data points, and the predicted distribution of the same values under the model, are compared using either of two test statistics $Z_1$ and $Z_2$ proposed by Berman (1986). The $Z_1$ test is also known as the Lawson-Waller test.

The function berman.test is generic, with methods for point patterns ("ppp" or "lpp") and point process models ("ppm" or "lppm").

• If X is a point pattern dataset (object of class "ppp" or "lpp"), then berman.test(X, ...) performs a goodness-of-fit test of the uniform Poisson point process (Complete Spatial Randomness, CSR) for this dataset.

• If model is a fitted point process model (object of class "ppm" or "lppm") then berman.test(model, ...) performs a test of goodness-of-fit for this fitted model. In this case, model should be a Poisson point process.

The test is performed by comparing the observed distribution of the values of a spatial covariate at the data points, and the predicted distribution of the same covariate under the model. Thus, you must nominate a spatial covariate for this test.

The argument covariate should be either a function(x,y) or a pixel image (object of class "im" containing the values of a spatial function. If covariate is an image, it should have numeric values, and its domain should cover the observation window of the model. If covariate is a function, it should expect two arguments x and y which are vectors of coordinates, and it should return a numeric vector of the same length as x and y.

First the original data point pattern is extracted from model. The values of the covariate at these data points are collected.

Next the values of the covariate at all locations in the observation window are evaluated. The point process intensity of the fitted model is also evaluated at all locations in the window.

• If which="Z1", the test statistic $Z_1$ is computed as follows. The sum $S$ of the covariate values at all data points is evaluated. The predicted mean $\mu$ and variance $\sigma^2$ of $S$ are computed from the values of the covariate at all locations in the window. Then we compute $Z_1 = (S-\mu)/\sigma$. Closely-related tests were proposed independently by Waller et al (1993) and Lawson (1993) so this test is often termed the Lawson-Waller test in epidemiological literature.

• If which="Z2", the test statistic $Z_2$ is computed as follows. The values of the covariate at all locations in the observation window, weighted by the point process intensity, are compiled into a cumulative distribution function $F$. The probability integral transformation is then applied: the values of the covariate at the original data points are transformed by the predicted cumulative distribution function $F$ into numbers between 0 and 1. If the model is correct, these numbers are i.i.d. uniform random numbers. The standardised sample mean of these numbers is the statistic $Z_2$.

In both cases the null distribution of the test statistic is the standard normal distribution, approximately.

The return value is an object of class "htest" containing the results of the hypothesis test. The print method for this class gives an informative summary of the test outcome.

Value

An object of class "htest" (hypothesis test) and also of class "bermantest", containing the results of the test. The return value can be plotted (by plot.bermantest) or printed to give an informative summary of the test.

Warning

The meaning of a one-sided test must be carefully scrutinised: see the printed output.

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

References

Berman, M. (1986) Testing for spatial association between a point process and another stochastic process. Applied Statistics 35, 54–62.

Lawson, A.B. (1993) On the analysis of mortality events around a prespecified fixed point. Journal of the Royal Statistical Society, Series A 156 (3) 363–377.

Waller, L., Turnbull, B., Clark, L.C. and Nasca, P. (1992) Chronic Disease Surveillance and testing of clustering of disease and exposure: Application to leukaemia incidence and TCE-contaminated dumpsites in upstate New York. Environmetrics 3, 281–300.

See Also

cdf.test, quadrat.test, ppm

Examples

# Berman's data
   X <- copper$SouthPoints
   L <- copper$SouthLines
   D <- distmap(L, eps=1)
   # test of fitted model
   fit <- ppm(X ~ x+y)
   berman.test(fit, D)

Description

Fits the Neyman-Scott Cluster point process with Cauchy kernel to a point pattern dataset by the Method of Minimum Contrast.

Usage

cauchy.estK(X, startpar=c(kappa=1,scale=1), lambda=NULL,
            q = 1/4, p = 2, rmin = NULL, rmax = NULL, ...)

Arguments

Details

This algorithm fits the Neyman-Scott cluster point process model with Cauchy kernel to a point pattern dataset by the Method of Minimum Contrast, using the $K$ function.

The argument X can be either

a point pattern:

An object of class "ppp" representing a point pattern dataset. The $K$ function of the point pattern will be computed using Kest, and the method of minimum contrast will be applied to this.

a summary statistic:

An object of class "fv" containing the values of a summary statistic, computed for a point pattern dataset. The summary statistic should be the $K$ function, and this object should have been obtained by a call to Kest or one of its relatives.

The algorithm fits the Neyman-Scott cluster point process with Cauchy kernel to X, by finding the parameters of the Matern Cluster model which give the closest match between the theoretical $K$ function of the Matern Cluster process and the observed $K$ function. For a more detailed explanation of the Method of Minimum Contrast, see mincontrast.

The model is described in Jalilian et al (2013). It is a cluster process formed by taking a pattern of parent points, generated according to a Poisson process with intensity $\kappa$, and around each parent point, generating a random number of offspring points, such that the number of offspring of each parent is a Poisson random variable with mean $\mu$, and the locations of the offspring points of one parent follow a common distribution described in Jalilian et al (2013).

If the argument lambda is provided, then this is used as the value of the point process intensity $\lambda$. Otherwise, if X is a point pattern, then $\lambda$ will be estimated from X. If X is a summary statistic and lambda is missing, then the intensity $\lambda$ cannot be estimated, and the parameter $\mu$ will be returned as NA.

The remaining arguments rmin,rmax,q,p control the method of minimum contrast; see mincontrast.

The corresponding model can be simulated using rCauchy.

For computational reasons, the optimisation procedure uses the parameter eta2, which is equivalent to 4 * scale^2 where scale is the scale parameter for the model as used in rCauchy.

Homogeneous or inhomogeneous Neyman-Scott/Cauchy models can also be fitted using the function kppm and the fitted models can be simulated using simulate.kppm.

The optimisation algorithm can be controlled through the additional arguments "..." which are passed to the optimisation function optim. For example, to constrain the parameter values to a certain range, use the argument method="L-BFGS-B" to select an optimisation algorithm that respects box constraints, and use the arguments lower and upper to specify (vectors of) minimum and maximum values for each parameter.

Value

An object of class "minconfit". There are methods for printing and plotting this object. It contains the following main components:

Author(s)

Abdollah Jalilian and Rasmus Waagepetersen. Adapted for spatstat by Adrian Baddeley [email protected]

References

Ghorbani, M. (2013) Cauchy cluster process. Metrika 76, 697–706.

Jalilian, A., Guan, Y. and Waagepetersen, R. (2013) Decomposition of variance for spatial Cox processes. Scandinavian Journal of Statistics 40, 119-137.

Waagepetersen, R. (2007) An estimating function approach to inference for inhomogeneous Neyman-Scott processes. Biometrics 63, 252–258.

See Also

kppm, cauchy.estpcf, lgcp.estK, thomas.estK, vargamma.estK, mincontrast, Kest, Kmodel.

rCauchy to simulate the model.

Examples

u <- cauchy.estK(redwood)
    u
    plot(u)

Description

Fits the Neyman-Scott Cluster point process with Cauchy kernel to a point pattern dataset by the Method of Minimum Contrast, using the pair correlation function.

Usage

cauchy.estpcf(X, startpar=c(kappa=1,scale=1), lambda=NULL,
            q = 1/4, p = 2, rmin = NULL, rmax = NULL, ...,
            pcfargs = list())

Arguments

Details

This algorithm fits the Neyman-Scott cluster point process model with Cauchy kernel to a point pattern dataset by the Method of Minimum Contrast, using the pair correlation function.

The argument X can be either

a point pattern:

An object of class "ppp" representing a point pattern dataset. The pair correlation function of the point pattern will be computed using pcf, and the method of minimum contrast will be applied to this.

a summary statistic:

An object of class "fv" containing the values of a summary statistic, computed for a point pattern dataset. The summary statistic should be the pair correlation function, and this object should have been obtained by a call to pcf or one of its relatives.

The algorithm fits the Neyman-Scott cluster point process with Cauchy kernel to X, by finding the parameters of the Matern Cluster model which give the closest match between the theoretical pair correlation function of the Matern Cluster process and the observed pair correlation function. For a more detailed explanation of the Method of Minimum Contrast, see mincontrast.

The model is described in Jalilian et al (2013). It is a cluster process formed by taking a pattern of parent points, generated according to a Poisson process with intensity $\kappa$, and around each parent point, generating a random number of offspring points, such that the number of offspring of each parent is a Poisson random variable with mean $\mu$, and the locations of the offspring points of one parent follow a common distribution described in Jalilian et al (2013).

If the argument lambda is provided, then this is used as the value of the point process intensity $\lambda$. Otherwise, if X is a point pattern, then $\lambda$ will be estimated from X. If X is a summary statistic and lambda is missing, then the intensity $\lambda$ cannot be estimated, and the parameter $\mu$ will be returned as NA.

The remaining arguments rmin,rmax,q,p control the method of minimum contrast; see mincontrast.

The corresponding model can be simulated using rCauchy.

For computational reasons, the optimisation procedure internally uses the parameter eta2, which is equivalent to 4 * scale^2 where scale is the scale parameter for the model as used in rCauchy.

Homogeneous or inhomogeneous Neyman-Scott/Cauchy models can also be fitted using the function kppm and the fitted models can be simulated using simulate.kppm.

The optimisation algorithm can be controlled through the additional arguments "..." which are passed to the optimisation function optim. For example, to constrain the parameter values to a certain range, use the argument method="L-BFGS-B" to select an optimisation algorithm that respects box constraints, and use the arguments lower and upper to specify (vectors of) minimum and maximum values for each parameter.

Value

An object of class "minconfit". There are methods for printing and plotting this object. It contains the following main components:

Author(s)

Abdollah Jalilian and Rasmus Waagepetersen. Adapted for spatstat by Adrian Baddeley [email protected]

References

Ghorbani, M. (2013) Cauchy cluster process. Metrika 76, 697–706.

Jalilian, A., Guan, Y. and Waagepetersen, R. (2013) Decomposition of variance for spatial Cox processes. Scandinavian Journal of Statistics 40, 119-137.

Waagepetersen, R. (2007) An estimating function approach to inference for inhomogeneous Neyman-Scott processes. Biometrics 63, 252–258.

See Also

kppm, cauchy.estK, lgcp.estpcf, thomas.estpcf, vargamma.estpcf, mincontrast, pcf, pcfmodel.

rCauchy to simulate the model.

Examples

u <- cauchy.estpcf(redwood)
    u
    plot(u, legendpos="topright")

Description

Performs a spatial distribution test of a point process model fitted to multiple spatial point patterns. The test compares the observed and predicted distributions of the values of a spatial covariate, using either the Kolmogorov-Smirnov, Cramer-von Mises or Anderson-Darling test of goodness-of-fit.

Usage

## S3 method for class 'mppm'
cdf.test(model, covariate, test=c("ks", "cvm", "ad"), ...,
            nsim=19, verbose=TRUE, interpolate=FALSE, fast=TRUE, jitter=TRUE)

Arguments

Details

This function is a method for the generic function cdf.test for the class mppm.

This function performs a goodness-of-fit test of a point process model that has been fitted to multiple point patterns. The observed distribution of the values of a spatial covariate at the data points, and the predicted distribution of the same values under the model, are compared using the Kolmogorov-Smirnov, Cramer-von Mises or Anderson-Darling test of goodness-of-fit. These are exact tests if the model is Poisson; otherwise, for a Gibbs model, a Monte Carlo p-value is computed by generating simulated realisations of the model and applying the selected goodness-of-fit test to each simulation.

The argument model should be a fitted point process model fitted to multiple point patterns (object of class "mppm").

The argument covariate contains the values of a spatial function. It can be

• a function(x,y)

• a pixel image (object of class "im"

• a list of function(x,y), one for each point pattern

• a list of pixel images, one for each point pattern

• a hyperframe (see hyperframe) of which the first column will be taken as containing the covariate

• a character string giving the name of one of the covariates in model

• one of the character strings "x" or "y", indicating the spatial coordinates.

If covariate is an image, it should have numeric values, and its domain should cover the observation window of the model. If covariate is a function, it should expect two arguments x and y which are vectors of coordinates, and it should return a numeric vector of the same length as x and y.

First the original data point pattern is extracted from model. The values of the covariate at these data points are collected.

The predicted distribution of the values of the covariate under the fitted model is computed as follows. The values of the covariate at all locations in the observation window are evaluated, weighted according to the point process intensity of the fitted model, and compiled into a cumulative distribution function $F$ using ewcdf.

The probability integral transformation is then applied: the values of the covariate at the original data points are transformed by the predicted cumulative distribution function $F$ into numbers between 0 and 1. If the model is correct, these numbers are i.i.d. uniform random numbers. A goodness-of-fit test of the uniform distribution is applied to these numbers using ks.test, cvm.test or ad.test.

The argument interpolate determines how pixel values will be handled when covariate is a pixel image. The value of the covariate at a data point is obtained by looking up the value of the nearest pixel if interpolate=FALSE, or by linearly interpolating between the values of the four nearest pixels if interpolate=TRUE. Linear interpolation is slower, but is sometimes necessary to avoid tied values of the covariate arising when the pixel grid is coarse.

If model is a Poisson point process, then the Kolmogorov-Smirnov, Cramer-von Mises and Anderson-Darling tests are theoretically exact. This test was apparently first described (in the context of spatial data, and for Kolmogorov-Smirnov) by Berman (1986). See also Baddeley et al (2005).

If model is not a Poisson point process, then the Kolmogorov-Smirnov, Cramer-von Mises and Anderson-Darling tests are biased. Instead they are used as the basis of a Monte Carlo test. First nsim simulated realisations of the model will be generated. Each simulated realisation consists of a list of simulated point patterns, one for each of the original data patterns. This can take a very long time. The model is then re-fitted to each simulation, and the refitted model is subjected to the goodness-of-fit test described above. A Monte Carlo p-value is then computed by comparing the p-value of the original test with the p-values obtained from the simulations.

Value

An object of class "cdftest" and "htest" containing the results of the test. See cdf.test for details.

Author(s)

Adrian Baddeley [email protected], Ida-Maria Sintorn and Leanne Bischoff. Implemented by Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

References

Baddeley, A., Rubak, E. and Turner, R. (2015) Spatial Point Patterns: Methodology and Applications with R. Chapman and Hall/CRC Press.

Baddeley, A., Turner, R., Moller, J. and Hazelton, M. (2005) Residual analysis for spatial point processes. Journal of the Royal Statistical Society, Series B 67, 617–666.

Berman, M. (1986) Testing for spatial association between a point process and another stochastic process. Applied Statistics 35, 54–62.

See Also

cdf.test, quadrat.test, mppm

Examples

# three i.i.d. realisations of nonuniform Poisson process
   lambda <- as.im(function(x,y) { 200 * exp(x) }, square(1))
   dat <- hyperframe(X=list(rpoispp(lambda), rpoispp(lambda), rpoispp(lambda)))

   # fit uniform Poisson process
   fit0 <- mppm(X~1, dat)
   # fit correct nonuniform Poisson process
   fit1 <- mppm(X~x, dat)

   # test wrong model
   cdf.test(fit0, "x")
   # test right model
   cdf.test(fit1, "x")

   # Gibbs model
   fitGibbs <- update(fit0, interaction=Strauss(0.05))
   ns <- if(interactive()) 19 else 2
   cdf.test(fitGibbs, "x", nsim=ns)

Description

Performs a test of goodness-of-fit of a point process model. The observed and predicted distributions of the values of a spatial covariate are compared using either the Kolmogorov-Smirnov test, Cramer-von Mises test or Anderson-Darling test. For non-Poisson models, a Monte Carlo test is used.

Usage

## S3 method for class 'ppm'
cdf.test(model, covariate,  test=c("ks", "cvm", "ad"), ...,
          interpolate=TRUE, jitter=TRUE, nsim=99, verbose=TRUE)

## S3 method for class 'slrm'
cdf.test(model, covariate,  test=c("ks", "cvm", "ad"), ...,
          modelname=NULL, covname=NULL)

Arguments

Details

These functions perform a goodness-of-fit test of a Poisson or Gibbs point process model fitted to point pattern data. The observed distribution of the values of a spatial covariate at the data points, and the predicted distribution of the same values under the model, are compared using the Kolmogorov-Smirnov test, the Cramer-von Mises test or the Anderson-Darling test. For Gibbs models, a Monte Carlo test is performed using these test statistics.

The function cdf.test is generic, with methods for point patterns ("ppp" or "lpp"), point process models ("ppm" or "lppm") and spatial logistic regression models ("slrm").

• If X is a point pattern dataset (object of class "ppp"), then cdf.test(X, ...) performs a goodness-of-fit test of the uniform Poisson point process (Complete Spatial Randomness, CSR) for this dataset. For a multitype point pattern, the uniform intensity is assumed to depend on the type of point (sometimes called Complete Spatial Randomness and Independence, CSRI).

• If model is a fitted point process model (object of class "ppm" or "lppm") then cdf.test(model, ...) performs a test of goodness-of-fit for this fitted model.

• If model is a fitted spatial logistic regression (object of class "slrm") then cdf.test(model, ...) performs a test of goodness-of-fit for this fitted model.

The test is performed by comparing the observed distribution of the values of a spatial covariate at the data points, and the predicted distribution of the same covariate under the model, using a classical goodness-of-fit test. Thus, you must nominate a spatial covariate for this test.

If X is a point pattern that does not have marks, the argument covariate should be either a function(x,y) or a pixel image (object of class "im" containing the values of a spatial function, or one of the characters "x" or "y" indicating the Cartesian coordinates. If covariate is an image, it should have numeric values, and its domain should cover the observation window of the model. If covariate is a function, it should expect two arguments x and y which are vectors of coordinates, and it should return a numeric vector of the same length as x and y.

If X is a multitype point pattern, the argument covariate can be either a function(x,y,marks), or a pixel image, or a list of pixel images corresponding to each possible mark value, or one of the characters "x" or "y" indicating the Cartesian coordinates.

First the original data point pattern is extracted from model. The values of the covariate at these data points are collected.

The predicted distribution of the values of the covariate under the fitted model is computed as follows. The values of the covariate at all locations in the observation window are evaluated, weighted according to the point process intensity of the fitted model, and compiled into a cumulative distribution function $F$ using ewcdf.

The probability integral transformation is then applied: the values of the covariate at the original data points are transformed by the predicted cumulative distribution function $F$ into numbers between 0 and 1. If the model is correct, these numbers are i.i.d. uniform random numbers. The A goodness-of-fit test of the uniform distribution is applied to these numbers using stats::ks.test, goftest::cvm.test or goftest::ad.test.

This test was apparently first described (in the context of spatial data, and using Kolmogorov-Smirnov) by Berman (1986). See also Baddeley et al (2005).

If model is not a Poisson process, then a Monte Carlo test is performed, by generating nsim point patterns which are simulated realisations of the model, re-fitting the model to each simulated point pattern, and calculating the test statistic for each fitted model. The Monte Carlo $p$ value is determined by comparing the simulated values of the test statistic with the value for the original data.

The return value is an object of class "htest" containing the results of the hypothesis test. The print method for this class gives an informative summary of the test outcome.

The return value also belongs to the class "cdftest" for which there is a plot method plot.cdftest. The plot method displays the empirical cumulative distribution function of the covariate at the data points, and the predicted cumulative distribution function of the covariate under the model, plotted against the value of the covariate.

The argument jitter controls whether covariate values are randomly perturbed, in order to avoid ties. If the original data contains any ties in the covariate (i.e. points with equal values of the covariate), and if jitter=FALSE, then the Kolmogorov-Smirnov test implemented in ks.test will issue a warning that it cannot calculate the exact $p$-value. To avoid this, if jitter=TRUE each value of the covariate will be perturbed by adding a small random value. The perturbations are normally distributed with standard deviation equal to one hundredth of the range of values of the covariate. This prevents ties, and the $p$-value is still correct. There is a very slight loss of power.

Value

An object of class "htest" containing the results of the test. See ks.test for details. The return value can be printed to give an informative summary of the test.

The value also belongs to the class "cdftest" for which there is a plot method.

Warning

The outcome of the test involves a small amount of random variability, because (by default) the coordinates are randomly perturbed to avoid tied values. Hence, if cdf.test is executed twice, the $p$-values will not be exactly the same. To avoid this behaviour, set jitter=FALSE.

Author(s)

Adrian Baddeley [email protected] and Rolf Turner [email protected]

References

Baddeley, A., Turner, R., Moller, J. and Hazelton, M. (2005) Residual analysis for spatial point processes. Journal of the Royal Statistical Society, Series B 67, 617–666.

Berman, M. (1986) Testing for spatial association between a point process and another stochastic process. Applied Statistics 35, 54–62.

See Also

plot.cdftest, quadrat.test, berman.test, ks.test, cvm.test, ad.test, ppm

Examples

op <- options(useFancyQuotes=FALSE)


   # fit inhomogeneous Poisson model and test
   model <- ppm(nztrees ~x)
   cdf.test(model, "x")

   if(interactive()) {
     # synthetic data: nonuniform Poisson process
     X <- rpoispp(function(x,y) { 100 * exp(x) }, win=square(1))

     # fit uniform Poisson process
     fit0 <- ppm(X ~1)
     # fit correct nonuniform Poisson process
     fit1 <- ppm(X ~x)

     # test wrong model
     cdf.test(fit0, "x")
     # test right model
     cdf.test(fit1, "x")
   }

   # multitype point pattern
   yimage <- as.im(function(x,y){y}, W=Window(amacrine))
   cdf.test(ppm(amacrine ~marks+y), yimage)

   options(op)

Description

Low-level functions to count the number of close pairs of points.

Usage

closepaircounts(X, r)

crosspaircounts(X, Y, r)

Arguments

Details

These are the efficient low-level functions used by spatstat to count close pairs of points in a point pattern or between two point patterns.

closepaircounts(X,r) counts the number of neighbours for each point in the pattern X. That is, for each point X[i], it counts the number of other points X[j] with j != i such that d(X[i],X[j]) <= r where d denotes Euclidean distance. The result is an integer vector v such that v[i] is the number of neighbours of X[i].

crosspaircounts(X,Y,r) counts, for each point in the pattern X, the number of neighbours in the pattern Y. That is, for each point X[i], it counts the number of points Y[j] such that d(X[i],X[j]) <= r. The result is an integer vector v such that v[i] is the number of neighbours of X[i] in the pattern Y.

Value

An integer vector of length equal to the number of points in X.

Warning about accuracy

The results of these functions may not agree exactly with the correct answer (as calculated by a human) and may not be consistent between different computers and different installations of R. The discrepancies arise in marginal cases where the interpoint distance is equal to, or very close to, the threshold rmax.

Floating-point numbers in a computer are not mathematical Real Numbers: they are approximations using finite-precision binary arithmetic. The approximation is accurate to a tolerance of about .Machine$double.eps.

If the true interpoint distance $d$ and the threshold rmax are equal, or if their difference is no more than .Machine$double.eps, the result may be incorrect.

Author(s)

Adrian Baddeley [email protected] and Rolf Turner [email protected]

See Also

closepairs to identify all close pairs of points.

Examples

a <- closepaircounts(cells, 0.1)
   sum(a)
   Y <- split(amacrine)
   b <- crosspaircounts(Y$on, Y$off, 0.1)

Description

Calculate the superposition of cluster kernels at the location of a point pattern.

Usage

## S3 method for class 'kppm'
clusterfield(model, locations = NULL, ...)

Arguments

Details

The function clusterfield is generic, with a method for "kppm" (described here) and methods for "character" and "function".

The method clusterfield.kppm extracts the relevant information from the fitted model and calls clusterfield.function.

The calculations are performed by density.ppp and ... arguments are passed thereto for control over the pixel resolution etc. (These arguments are then passed on to pixellate.ppp and as.mask.)

Value

A pixel image (object of class "im").

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

See Also

kppm,

clusterfield.

Examples

fit <- kppm(redwood~1, "Thomas")
  Z <- clusterfield(fit, eps = 0.01)

Description

Fit a homogeneous or inhomogeneous cluster process or Cox point process model to a point pattern by the Method of Minimum Contrast.

Usage

clusterfit(X, clusters, lambda = NULL, startpar = NULL, ...,
           q = 1/4, p = 2, rmin = NULL, rmax = NULL,
           ctrl=list(q=q, p=p, rmin=rmin, rmax=rmax),
           statistic = NULL, statargs = NULL, algorithm="Nelder-Mead",
           verbose=FALSE, pspace=NULL)

Arguments

Details

This function fits the clustering parameters of a cluster or Cox point process model by the Method of Minimum Contrast, that is, by matching the theoretical $K$-function of the model to the empirical $K$-function of the data, as explained in mincontrast.

If statistic="pcf" (or X appears to be an estimated pair correlation function) then instead of using the $K$-function, the algorithm will use the pair correlation function.

If X is a point pattern of class "ppp" an estimate of the summary statistic specfied by statistic (defaults to "K") is first computed before minimum contrast estimation is carried out as described above. In this case the argument statargs can be used for controlling the summary statistic estimation. The precise algorithm for computing the summary statistic depends on whether the intensity specification (lambda) is:

homogeneous:

If lambda is NUll or a single numeric the pattern is considered homogeneous and either Kest or pcf is invoked. In this case lambda is not used for anything when estimating the summary statistic.

inhomogeneous:

If lambda is a pixel image (object of class "im"), a fitted point process model (object of class "ppm" or "kppm") or a function(x,y) the pattern is considered inhomogeneous. In this case either Kinhom or pcfinhom is invoked with lambda as an argument.

After the clustering parameters of the model have been estimated by minimum contrast lambda (if non-null) is used to compute the additional model parameter $\mu$.

The algorithm parameters q,p,rmax,rmin are described in the help for mincontrast. They may be provided either as individually-named arguments, or as entries in the list ctrl. The individually-named arguments q,p,rmax,rmin override the entries in the list ctrl.

Value

An object of class "minconfit". There are methods for printing and plotting this object. See mincontrast.

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

References

Diggle, P.J. and Gratton, R.J. (1984) Monte Carlo methods of inference for implicit statistical models. Journal of the Royal Statistical Society, series B 46, 193 – 212.

Moller, J. and Waagepetersen, R. (2003). Statistical Inference and Simulation for Spatial Point Processes. Chapman and Hall/CRC, Boca Raton.

Waagepetersen, R. (2007). An estimating function approach to inference for inhomogeneous Neyman-Scott processes. Biometrics 63 (2007) 252–258.

See Also

kppm

Examples

fit <- clusterfit(redwood, "Thomas")
  fit
  if(interactive()){
    plot(fit)
  }
  K <- Kest(redwood)
  fit2 <- clusterfit(K, "MatClust")

Description

Given a fitted cluster point process model, this command returns the probability density of the cluster offspring.

Usage

## S3 method for class 'kppm'
clusterkernel(model, ...)

Arguments

Details

Given a cluster point process model, this command returns a function(x,y) giving the two-dimensional probability density of the cluster offspring points assuming a cluster parent located at the origin.

The function clusterkernel is generic, with methods for class "kppm" (described here) and "character" (described in clusterkernel.character).

Value

A function in the R language with arguments x,y,....

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

See Also

clusterkernel.character, clusterfield, kppm

Examples

fit <- kppm(redwood ~ x, "MatClust")
  f <- clusterkernel(fit)
  f(0.05, 0.02)

Description

Given a cluster point process model, this command returns a value beyond which the the probability density of the cluster offspring is neglible.

Usage

## S3 method for class 'kppm'
clusterradius(model, ..., thresh = NULL, precision = FALSE)

Arguments

Details

Given a cluster model this function by default returns the effective range of the model with the given parameters as used in spatstat. For the Matern cluster model (see e.g. rMatClust) this is simply the finite radius of the offsring density given by the paramter scale irrespective of other options given to this function. The remaining models in spatstat have infinite theoretical range, and an effective finite value is given as follows: For the Thomas model (see e.g. rThomas the default is 4*scale where scale is the scale or standard deviation parameter of the model. If thresh is given the value is instead found as described for the other models below.

For the Cauchy model (see e.g. rCauchy) and the Variance Gamma (Bessel) model (see e.g. rVarGamma) the value of thresh defaults to 0.001, and then this is used to compute the range numerically as follows. If $k(x,y)=k_0(r)$ with $r=\sqrt(x^2+y^2)$ denotes the isotropic cluster kernel then $f(r) = 2 \pi r k_0(r)$ is the density function of the offspring distance from the parent. The range is determined as the value of $r$ where $f(r)$ falls below thresh times $k_0(r)$.

If precision=TRUE the precision related to the chosen range is returned as an attribute. Here the precision is defined as the polar integral of the kernel from distance 0 to the calculated range. Ideally this should be close to the value 1 which would be obtained for the true theretical infinite range.

Value

A positive numeric.

Additionally, the precision related to this range value is returned as an attribute "prec", if precision=TRUE.

Author(s)

Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

See Also

clusterkernel, kppm, rMatClust, rThomas, rCauchy, rVarGamma, rNeymanScott.

Examples

fit <- kppm(redwood ~ x, "MatClust")
  clusterradius(fit)

Description

Given a point process model fitted to a list of point patterns, extract the coefficients of the fitted model. A method for coef.

Usage

## S3 method for class 'mppm'
coef(object, ...)

Arguments

Details

This function is a method for the generic function coef.

The argument object must be a fitted point process model (object of class "mppm") produced by the fitting algorithm mppm). This represents a point process model that has been fitted to a list of several point pattern datasets. See mppm for information.

This function extracts the vector of coefficients of the fitted model. This is the estimate of the parameter vector $\theta$ such that the conditional intensity of the model is of the form

$\lambda(u,x) = \exp(\theta S(u,x))$

where $S(u,x)$ is a (vector-valued) statistic.

For example, if the model object is the uniform Poisson process, then coef(object) will yield a single value (named "(Intercept)") which is the logarithm of the fitted intensity of the Poisson process.

If the fitted model includes random effects (i.e. if the argument random was specified in the call to mppm), then the fitted coefficients are different for each point pattern in the original data, so coef(object) is a data frame with one row for each point pattern, and one column for each parameter. Use fixef.mppm to extract the vector of fixed effect coefficients, and ranef.mppm to extract the random effect coefficients at each level.

Use print.mppm to print a more useful description of the fitted model.

Value

Either a vector containing the fitted coefficients, or a data frame containing the fitted coefficients for each point pattern.

Author(s)

Adrian Baddeley, Ida-Maria Sintorn and Leanne Bischoff. Implemented in spatstat by Adrian Baddeley [email protected], Rolf Turner [email protected] and Ege Rubak [email protected].

References

Baddeley, A., Rubak, E. and Turner, R. (2015) Spatial Point Patterns: Methodology and Applications with R. Chapman and Hall/CRC Press.

See Also

fixef.mppm and ranef.mppm for the fixed and random effect coefficients in a model that includes random effects.

print.mppm, mppm

Examples

H <- hyperframe(X=waterstriders)

    fit.Poisson <- mppm(X ~ 1, H)
    coef(fit.Poisson)

    # The single entry "(Intercept)" 
    # is the log of the fitted intensity of the Poisson process

    fit.Strauss <- mppm(X~1, H, Strauss(7))
    coef(fit.Strauss)

    # The two entries "(Intercept)" and "Interaction"
    # are respectively log(beta) and log(gamma)
    # in the usual notation for Strauss(beta, gamma, r)

    # Tweak data to exaggerate differences
    H$X[[1]] <- rthin(H$X[[1]], 0.3)
    # Model with random effects
    fitran <- mppm(X ~ 1, H, random=~1|id)
    coef(fitran)

Description

Given a point process model fitted to a point pattern, extract the coefficients of the fitted model. A method for coef.

Usage

## S3 method for class 'ppm'
coef(object, ...)

Arguments

Details

This function is a method for the generic function coef.

The argument object must be a fitted point process model (object of class "ppm"). Such objects are produced by the maximum pseudolikelihood fitting algorithm ppm).

This function extracts the vector of coefficients of the fitted model. This is the estimate of the parameter vector $\theta$ such that the conditional intensity of the model is of the form

$\lambda(u,x) = \exp(\theta S(u,x))$

where $S(u,x)$ is a (vector-valued) statistic.

For example, if the model object is the uniform Poisson process, then coef(object) will yield a single value (named "(Intercept)") which is the logarithm of the fitted intensity of the Poisson process.

Use print.ppm to print a more useful description of the fitted model.

Value

A vector containing the fitted coefficients.

Author(s)

Adrian Baddeley [email protected] and Rolf Turner [email protected]

See Also

print.ppm, ppm.object, ppm

Examples

poi <- ppm(cells, ~1, Poisson())
    coef(poi)
    # This is the log of the fitted intensity of the Poisson process

    stra <- ppm(cells, ~1, Strauss(r=0.07))
    coef(stra)

    # The two entries "(Intercept)" and "Interaction"
    # are respectively log(beta) and log(gamma)
    # in the usual notation for Strauss(beta, gamma, r)

Description

Extracts the coefficients (parameters) from a fitted Spatial Logistic Regression model.

Usage

## S3 method for class 'slrm'
coef(object, ...)

Arguments

Details

This is a method for coef for fitted spatial logistic regression models (objects of class "slrm", usually obtained from the function slrm).

It extracts the fitted canonical parameters, i.e.\ the coefficients in the linear predictor of the spatial logistic regression.

Value

Numeric vector of coefficients.

Author(s)

Adrian Baddeley [email protected] and Rolf Turner [email protected]

See Also

slrm

Examples

X <- rpoispp(42)
  fit <- slrm(X ~ x+y)
  coef(fit)

Description

Compares several fitted point process models using the same residual diagnostic.

Usage

compareFit(object, Fun, r = NULL, breaks = NULL, ...,
         trend = ~1, interaction = Poisson(), rbord = NULL,
         modelnames = NULL, same = NULL, different = NULL)

Arguments

Details

This is a convenient way to collect diagnostic information for several different point process models fitted to the same point pattern dataset, or for point process models of the same form fitted to several different datasets, etc.

The first argument, object, is usually a list of fitted point process models (objects of class "ppm"), obtained from the model-fitting function ppm.

For convenience, object can also be a list of point patterns (objects of class "ppp"). In that case, point process models will be fitted to each of the point pattern datasets, by calling ppm using the arguments trend (for the first order trend), interaction (for the interpoint interaction) and rbord (for the erosion distance in the border correction for the pseudolikelihood). See ppm for details of these arguments.

Alternatively object can be a single point pattern (object of class "ppp") and one or more of the arguments trend, interaction or rbord can be a list. In this case, point process models will be fitted to the same point pattern dataset, using each of the model specifications listed.

The diagnostic function Fun will be applied to each of the point process models. The results will be collected into a single function value table. The modelnames are used to label the results from each fitted model.

Value

Function value table (object of class "fv").

Author(s)

Ege Rubak [email protected], Adrian Baddeley [email protected] and Jesper Moller.

See Also

ppm, Kcom, Kres, Gcom, Gres, psst, psstA, psstG, collapse.fv

Examples

nd <- 40
   
   ilist <- list(Poisson(), Geyer(7, 2), Strauss(7))
   iname <- c("Poisson", "Geyer", "Strauss")
   
   K <- compareFit(swedishpines, Kcom, interaction=ilist, rbord=9,
            correction="translate",
            same="trans", different="tcom", modelnames=iname, nd=nd)
   K

Description

Creates an instance of the Connected Component point process model which can then be fitted to point pattern data.

Usage

Concom(r)

Arguments

Details

This function defines the interpoint interaction structure of a point process called the connected component process. It can be used to fit this model to point pattern data.

The function ppm(), which fits point process models to point pattern data, requires an argument of class "interact" describing the interpoint interaction structure of the model to be fitted. The appropriate description of the connected component interaction is yielded by the function Concom(). See the examples below.

In standard form, the connected component process (Baddeley and Moller, 1989) with disc radius $r$, intensity parameter $\kappa$ and interaction parameter $\gamma$ is a point process with probability density

$f(x_1,\ldots,x_n) = \alpha \kappa^{n(x)} \gamma^{-C(x)}$

for a point pattern $x$, where $x_1,\ldots,x_n$ represent the points of the pattern, $n(x)$ is the number of points in the pattern, and $C(x)$ is defined below. Here $\alpha$ is a normalising constant.

To define the term C(x), suppose that we construct a planar graph by drawing an edge between each pair of points $x_i,x_j$ which are less than $r$ units apart. Two points belong to the same connected component of this graph if they are joined by a path in the graph. Then $C(x)$ is the number of connected components of the graph.