Type:	Package
Title:	Functions for Spatial Point Pattern Analysis
Version:	1.1-18
Date:	2025-05-01
Depends:	spatstat (≥ 2.0-0)
Imports:	spatstat.geom, spatstat.explore, spatstat.model, spatstat.random
Suggests:	dixon
Description:	Some wrappers, functions and data sets for for spatial point pattern analysis (mainly based on 'spatstat'), used in the book "Introduccion al Analisis Espacial de Datos en Ecologia y Ciencias Ambientales: Metodos y Aplicaciones" and in the papers by De la Cruz et al. (2008) <doi:10.1111/j.0906-7590.2008.05299.x> and Olano et al. (2009) <doi:10.1051/forest:2008074>.
License:	GPL-2 | GPL-3 [expanded from: GPL (≥ 2)]
NeedsCompilation:	yes
Packaged:	2025-05-01 16:04:28 UTC; marcelino
Author:	Marcelino de la Cruz Rot [aut, cre], Philip M. Dixon [ctb], Jose M. Blanco-Moreno [ctb]
Maintainer:	Marcelino de la Cruz Rot <marcelino.delacruz@urjc.es>
Repository:	CRAN
Date/Publication:	2025-05-01 17:00:02 UTC

Functions for spatial point pattern analysis in ecology

Description

This is a summary of the features of ecespa by way of examples of its main functions.

Author(s)

Marcelino de la Cruz Rot, with contributions of Philip M. Dixon and Jose M. Blanco-Moreno and heavily borrowing Baddeley's & Turner's spatstat code.

References

De la Cruz, M. 2006. Introducción al análisis de datos mapeados o algunas de las (muchas) cosas que puedo hacer si tengo coordenadas. Ecosistemas 15 (3): 19-39.

De la Cruz, M. 2008. Métodos para analizar datos puntuales. En: Introducción al Análisis Espacial de Datos en Ecología y Ciencias Ambientales: Métodos y Aplicaciones (eds. Maestre, F. T., Escudero, A. y Bonet, A.), pp 76-127. Asociación Española de Ecología Terrestre, Universidad Rey Juan Carlos y Caja de Ahorros del Mediterráneo, Madrid.

De la Cruz, M., Romao, R.L., Escudero, A. and Maestre, F.T. 2008. Where do seedlings go? A spatio-temporal analysis of early mortality in a semiarid specialist. Ecography, 31(6): 720-730. doi:10.1111/j.0906-7590.2008.05299.x.

Diggle, P. J. 2003. Statistical analysis of spatial point patterns. Arnold, London.

Dixon, P.M. 2002. Nearest-neighbor contingency table analysis of spatial segregation for several species. Ecoscience, 9 (2): 142-151. doi:10.1080/11956860.2002.11682700.

Dixon, P. M. 2002. Ripley's K function. In The encyclopedia of environmetrics (eds. El-Shaarawi, A.H. & Piergorsch, W.W.), pp. 1976-1803. John Wiley & Sons Ltd, NY.

Escudero, A., Romao, R.L., De la Cruz, M. & Maestre, F. 2005. Spatial pattern and neighbour effects on Helianthemum squamatum seedlings in a Mediterranean gypsum community. Journal of Vegetation Science, 16: 383-390. doi:10.1111/j.1654-1103.2005.tb02377.x.

Getis, A. and Franklin, J. 1987. Second-order neighbourhood analysis of mapped point patterns. Ecology 68: 473-477. doi:10.2307/1938452.

Loosmore, N.B. and Ford, E.D. (2006) Statistical inference using the G or K point pattern spatial statistics. Ecology 87, 1925-1931. doi:10.1890/0012-9658(2006)87[1925:SIUTGO]2.0.CO;2.

Lotwick,H.W. & Silverman, B. W. 1982. Methods for analysing spatial processes of several types of points. Journal of the Royal Statistical Society B, 44: 406-413. doi:10.1111/j.2517-6161.1982.tb01221.x.

Olano, J.M., Laskurain, N.A., Escudero, A. and De la Cruz, M. 2009. Why and where adult trees die in a secondary temperate forest? The role of neighbourhood. Annals of Forest Science, 66: 105. doi:10.1051/forest:2008074.

Penttinen, A. 2006. Statistics for Marked Point Patterns. In The Yearbook of the Finnish Statistical Society, pp. 70-91.

Rey-Benayas, J.M., de la Montaña, E., Pérez-Camacho, L., de la Cruz, M., Moreno, D., Parejo, J.L. and Suárez-Seoane, S. 2010. Inter-annual dynamics and spatial congruence of a nocturnal bird assemblage inhabiting a Mediterranean agricultural mosaic. Ardeola,57(2): 303-320.

Syrjala, S. E. 1996. A statistical test for a difference between the spatial distributions of two populations. Ecology 77: 75-80. doi:10.2307/2265656.

Waagepetersen, R. P. 2007. An estimating function approach to inference for inhomogeneous Neymann-Scott processes. Biometrics 63: 252-258. doi:10.1111/j.1541-0420.2006.00667.x.

Examples

## Not run: 

#############################################
### Transfom easily data from a data.frame into the ppp format 
### of spatstat:

data(fig1)

plot(fig1) #typical xyplot

fig1.ppp <- haz.ppp (fig1)

fig1.ppp

plot(fig1.ppp) # point pattern plot of spatstat



#############################################
###  Summarize the joint pattern of points and marks at different scales
###  with the normalized mark-weighted K-function (Penttinen, 2006). 
###  Compare this function in two consecutive cohorts of Helianthemum
###  squamatum seedlings:

 ## Figure 3.10 of De la Cruz (2008):
  
  data(seedlings1)
  
  data(seedlings2)
  
  s1km <- Kmm(seedlings1, r=1:100)
  
  s2km <- Kmm(seedlings2, r=1:100)
  
  plot(s1km, ylime=c(0.6,1.2), lwd=2, maine="", xlabe="r(cm)")

  plot(s2km,  lwd=2, lty=2, add=T )

  abline(h=1, lwd=2, lty=3)
  
  legend(x=60, y=1.2, legend=c("Hs_C1", "Hs_C2", "H0"),
	 lty=c(1, 2, 3), lwd=c(3, 2, 2), bty="n")
 
## A pointwise test of normalized Kmm == 1 for seedlings1:

   s1km.test <- Kmm(seedlings1, r=1:100, nsim=99)

   plot(s1km.test,  xlabe="r(cm)")
   
   
   

#############################################
###  Explore the local relationships between marks and locations (e.g. size 
###  of one cohort of H. squamatum seedlings). Map the marked point pattern 
###  to a random field for visual inspection, with the normalized mark-sum
###  measure (Penttinen, 2006).

data(seedlings1)
   
 seed.m <- marksum(seedlings1, R=25)

 plot(seed.m, what="marksum", sigma = 5)  # raw mark-sum measure; sigma is bandwith for smoothing

 plot(seed.m, what="pointsum", sigma = 5) # point sum measure
   
 plot(seed.m,  what="normalized", dimyx=200, contour=TRUE, sigma = 5) # normalized  mark-sum measure

# the same with added grid

 plot(seed.m,  what="normalized", dimyx=200, contour=TRUE, sigma = 5, grid=TRUE) 



#############################################
###  Test against the null model of "independent labelling",
###  i.e. test asociation/repulsion between  a "fixed" pattern (e.g. adult
###  H. squamatum plants) and a "variable" pattern (e.g. of surviving and 
###  dead seedlings), with 2.5% and 97.5% envelopes of 999 random 
###  labellings (De la Cruz & al. 2008).

data(Helianthemum)


cosa <- K012(Helianthemum, fijo="adultHS", i="deadpl", j="survpl",
             r=seq(0,200,le=201), nsim=999, nrank=25, correction="isotropic")

plot(cosa$k01, sqrt(./pi)-r~r,  col=c(3, 1, 3), lty=c(3, 1, 3), las=1,
         ylab=expression(L[12]), xlim=c(0, 200), 
	 main="adult HS vs. dead seedlings", legend=FALSE)

plot(cosa$k02, sqrt(./pi)-r~r, col=c(3, 1, 3), lty=c(3, 1, 3), las=1, 
         ylab=expression(L[12]), xlim=c(0, 200),
 	 main="adult HS vs. surviving seedlings", legend=FALSE)



#############################################
###  Test differences of agregation and segregation between two patterns, 
###  e.g. surviving and dying H. squamatum seedlings (De la Cruz & al. 2008). 

data(Helianthemum)

cosa12 <- K1K2(Helianthemum, j="deadpl", i="survpl", r=seq(0,200,le=201),
		 nsim=999, nrank=1, correction="isotropic")

plot(cosa12$k1k2, lty=c(2, 1, 2), col=c(2, 1, 2), xlim=c(0, 200),
         main= "survival- death",ylab=expression(K[1]-K[2]), legend=FALSE)

plot(cosa12$k1k12, lty=c(2, 1, 2), col=c(2, 1, 2), xlim=c(0, 200),
	 main="segregation of surviving seedlings",ylab=expression(K[1]-K[12]), legend=FALSE)

plot(cosa12$k2k12, lty=c(2, 1, 2), col=c(2, 1, 2), xlim=c(0, 200),
         main= "segregation of dying seedlings",ylab=expression(K[2]-K[12]), legend=FALSE)


#############################################
###  Test 'univariate' and 'bivariate' point patterns 
###  against non-Poisson (in-)homogeneous models 
###  (De la Cruz and Escudero, submited).

 data(urkiola)

   #####################
   ## univariate example

   # get univariate pp
   I.ppp <- split.ppp(urkiola)$birch

   # estimate inhomogeneous intensity function
   I.lam <- predict (ppm(I.ppp, ~polynom(x,y,2)), type="trend", ngrid=200)

   # Compute and plot envelopes to Kinhom, simulating from an Inhomogeneous
   #  Poisson Process:
   
   I2.env <- envelope( I.ppp,Kinhom, lambda=I.lam, correction="trans", 
                              nsim=99, simulate=expression(rpoispp(I.lam)))
   plot(I2.env, sqrt(./pi)-r~r, xlab="r (metres)", ylab= "L (r)", col=c(1,3,2,2),legend=FALSE) 

   # It seems that there is short scale clustering; let's fit an Inhomogeneous 
   # Poisson Cluster Process: 

   I.ki <- ipc.estK(mippp=I.ppp, lambda=I.lam, correction="trans")

   # Compute and plot envelopes to Kinhom, simulating from the fitted IPCP:

   Ipc.env <- Ki(I.ki, correction="trans", nsim=99, ngrid=200)

   plot (Ipc.env, xlab="r (metres)", ylab= "L (r)")
 
   #####################
   ## bivariate example: test independence between birch and quercus in Urkiola

   J.ppp <- split.ppp(urkiola)$oak
   
   # We want to simulate oak from a homogeneous Poisson model:
   J.ppm <- ppm(J.ppp, trend=~1, interaction=Poisson() )
   
   IJ.env <- Kci (mod1=I.ki, mod2=J.ppm, nsim=99)
   
   plot(IJ.env, type=12)
   
   plot(IJ.env, type=21)




#############################################
###  Simulate envelopes from the fitted values of a logistic model,
###  as in Olano et al. (2009)
   
   
   data(quercusvm)

   # read fitted values from logistic model:
   
   
   probquercus <-c(0.99955463, 0.96563477, 0.97577094, 0.97327199, 0.92437309,
   0.84023396, 0.94926682, 0.89687281, 0.99377915, 0.74157478, 0.95491518,
   0.72366493, 0.66771787, 0.77330148, 0.67569082, 0.9874892, 0.7918891, 
   0.73246803, 0.81614635, 0.66446411, 0.80077908, 0.98290508, 0.54641754,
   0.53546689, 0.73273626, 0.7347013, 0.65559655, 0.89481468, 0.63946334,
   0.62101995, 0.78996371, 0.93179582, 0.80160346, 0.82204428, 0.90050059,
   0.83810669, 0.92153079, 0.47872421, 0.24697004, 0.50680935, 0.6297911, 
   0.46374812, 0.65672284, 0.87951682, 0.35818237, 0.50932432, 0.92293014,
   0.48580241, 0.49692053, 0.52290553, 0.7317549, 0.32445982, 0.30300865,
   0.73599359, 0.6206056, 0.85777043, 0.65758613, 0.50100406, 0.31340849, 
   0.22289286, 0.40002879, 0.29567678, 0.56917817, 0.56866864, 0.27718552,
   0.4910667, 0.47394411, 0.40543788, 0.29571349, 0.30436276, 0.47859015,
   0.31754526, 0.42131675, 0.37468782, 0.73271225, 0.26786274, 0.59506388, 
   0.54801851, 0.38983575, 0.64896835, 0.37282031, 0.67624306, 0.29429766,
   0.29197755, 0.2247629, 0.40697843, 0.17022391, 0.26528042, 0.24373722,
   0.26936163, 0.13052254, 0.19958585, 0.18659692, 0.36686678, 0.47263005,
   0.39557661, 0.68048997, 0.74878567, 0.88352322, 0.93851375)
   
  

   ################################ 
   ## Envelopes for an homogeneous point pattern:
   
   cosap <- Kinhom.log(A=quercusvm, lifemark="0", nsim=99, prob=probquercus)

   plot(cosap)

   
   ################################ 
   ## Envelopes for an inhomogeneous point pattern:
   
   ## First, fit an inhomogeneous Poisson model to alive trees :
   
   quercusalive <- unmark(quercusvm[quercusvm$marks == 0])

    mod2 <- ppm(quercusalive, ~polynom(x,y,2))

    ## Now use mod2 to estimate lambda for K.inhom:
    
    cosapm <- Kinhom.log(A=quercusvm, lifemark="0", prob=probquercus,
                                   nsim=99,  mod=mod2)
    plot(cosapm)




#############################################
###  Test segregation based on the counts in the contingency table
###  of nearest neighbors in a multitype point pattern (Dixon, 2002)

data(swamp)

dixon2002(swamp,nsim=99)




#############################################
###  Fit the Poisson cluster point process to a point pattern with 
###  the method of minimum contrast (Diggle 2003).

data(gypsophylous)

# Estimate K function ("Kobs").

gyps.env <- envelope(gypsophylous, Kest, correction="iso", nsim=99)

plot(gyps.env, sqrt(./pi)-r~r, ylab="L(r)", legend=FALSE)

# Fit Poisson Cluster Process. The limits of integration 
# rmin and rmax are setup to 0 and 60, respectively. 

cosa.pc <- pc.estK(Kobs = gyps.env$obs[gyps.env$r<=60],
		           r = gyps.env$r[gyps.env$r<=60])

# Add fitted Kclust function to the plot.

lines(gyps.env$r,sqrt(Kclust(gyps.env$r, cosa.pc$sigma2,cosa.pc$rho)/pi)-gyps.env$r,
       lty=2, lwd=3, col="purple")

# A kind of pointwise test of the gypsophylous pattern been a realisation
# of the fitted model, simulating with sim.poissonc and using function J (Jest).

gyps.env.sim <- envelope(gypsophylous, Jest, nsim=99,
                    simulate=expression(sim.poissonc(gypsophylous,
		    sigma=sqrt(cosa.pc$sigma2), rho=cosa.pc$rho)))

 plot(gyps.env.sim,  main="",legendpos="bottomleft")






#############################################
###  Compute Syrjala's test for the difference between the spatial 
###  distributions of two populations, as in Rey-Benayas et al. 
### (submited)


 
   data(syr1); data(syr2); data(syr3)
   
   plot(syrjala.test(syr1, syr2, nsim=999)) 
   
   plot(syrjala.test(syr1, syr3, nsim=999)) 
   




## End(Not run)

Spatial point pattern of Helianthemum squamatum adult plants and seedlings

Description

Locations of H. squamatum adult plants and seedlings in a 6 m x 7 m plot over gypsum soil in Chinchón (near Madrid, Spain). These are part of the data collected by Romao (2003) that have been analyzed several times (Escudero et al.2005, De la Cruz 2006, De la Cruz et al. in press. ). The coordinates of the plans are given in cm.

Usage

data(Helianthemum)

Format

An object of class "ppp" of spatstat representing the point pattern of plants locations marked by their type. See ppp.object for details of the format.The dataset has 866 points with the following levels:

adultHS

adult H. squamatum plants

deadp

dying H. squamatum seedlings

survpl

surviving H. squamatum seedlings

Source

Romao, R.L. 2003. Estructura espacial de comunidades de gipsófitos: interacciones bióticas y constricciones abióticas. Tesis Doctoral. Universidad Politécnica de Madrid.

References

De la Cruz, M. 2006. Introducción al análisis de datos mapeados o algunas de las (muchas) cosas que puedo hacer si tengo coordenadas. Ecosistemas 15 (3): 19-39.

De la Cruz, M., Romao, R.L., Escudero, A. and Maestre, F.T. 2008. Where do seedlings go? A spatio-temporal analysis of early mortality in a semiarid specialist. Ecography,31(6): 720-730. doi:10.1111/j.0906-7590.2008.05299.x.

Escudero, A., Romao, R.L., De la Cruz, M. & Maestre, F. 2005. Spatial pattern and neighbour effects on Helianthemum squamatum seedlings in a Mediterranean gypsum community. Journal of Vegetation Science, 16: 383-390. doi:10.1111/j.1654-1103.2005.tb02377.x.

Examples

data(Helianthemum)
plot(Helianthemum)

Tests against 'independent labelling'

Description

Given a "fixed" point pattern and some process that asign labels (I,J) to another "variable" point pattern, K012 estimates the combined bivariate K function between the fixed pattern and every type of the variable pattern, and test that they are independent (i.e. that the labels are randomly assigned, irrespectively of the fixed pattern).

Usage

K012(X, fijo, i, j, nsim = 99, nrank = 1, r = NULL,
	 correction = "isotropic")

Arguments

Details

This test was developped to answer some questions about the spatial pattern of survival and mortality of seedlings and its relationships with adult plants in a plant community (De la Cruz et al. 2008). In order to evaluate the spatial structures of seedlings fates (survive or die), the null hypothesis of random labelling (Cuzick & Edwards 1990, Dixon 2002) would be the appropriate one. This kind of pattern is the result of two hierarchical processes: a first one that generates the pattern of points (seedlings) and other that assign "labels" (i.e. "die", "survive") to the points. On the other hand, to analyze the relationships between the spatial pattern of emerging seedlings and the pattern of adult plants (two patterns that have been generated independently), independence would be the appropriate null hypothesis (Goreaud & Pellisier 2003). However, testing the relationship between the pattern of seedling fates and the pattern of adult plants does not completely fit any of the mentioned hypotheses because, although the pattern of adult plants and the pattern of, e.g., dead seedlings are generated independently, their relationship is conditioned by the dependence of the fate "dead" on the locations of emerging seedlings. This implies that one can not apply the usual technique of toroidal shifting one pattern over the other to test the independence hypothesis. Instead one must permute the label of the focal fate (i.e. survive, die) over the global pattern of seedlings points, keeping the locations and labels of adults fixed. This is the method that K012 uses to build the envelopes. The bivariate K functions are computed with the Lotwick's and Silverman's (1982) combined estimator (Kmulti.ls).

Value

A list with two elements.

Each of the above elements is a fv.object, essentially a data.frame with the following items:

together with the observed corrected estimate of the combined bivariate K function ( iso, trans, border, etc).

Author(s)

Marcelino de la Cruz

References

Cuzick, J. and Edwards, R. 1990. Spatial clustering for inhomogeneous populations (with discussion). Journal of the Royal Statistical Society B 52: 73-104.

De la Cruz, M. 2006. Introducción al análisis de datos mapeados o algunas de las (muchas) cosas que puedo hacer si tengo coordenadas. Ecosistemas 15 (3): 19-39.

De la Cruz, M., Romao, R.L., Escudero, A. and Maestre, F.T. 2008. Where do seedlings go? A spatio-temporal analysis of early mortality in a semiarid specialist. Ecography, 31(6): 720-730. doi:10.1111/j.0906-7590.2008.05299.x.

Dixon, P. M. 2002. Ripley's K function. In The encyclopedia of environmetrics (eds. El-Shaarawi, A.H. & Piergorsch, W.W.), pp. 1976-1803. John Wiley & Sons Ltd, NY.

Goreaud, F. and Pelissier, R. 2003. Avoiding misinterpretation of biotic interactions with the intertype K12-function: population independence vs. random labelling hypotheses. J. Veg. Sci. 14: 681-692.

Lotwick, H. W. & Silverman, B. W. 1982. Methods for analysing spatial processes of several types of points. Journal of the Royal Statistical Society B 44: 406-413.

See Also

dixon2002 for another segregation test, based in the contingency table of counts of nearest neigbors in a marked point pattern.

Examples

data(Helianthemum)


# To get 2.5% and 97.5% envelopes of 999 random labelings.
## as in fig. 10 of De la Cruz (2006), change nsim= 999 and nrank=25
nsim <-19
nrank <-1



## Test asociation/repulsion between the fixed pattern of adult
## H. squamatum plants and the "variable" pattern of surviving and 
## dead seedlings,


cosa <- K012(Helianthemum, fijo="adultHS", i="deadpl", j="survpl",
             r=seq(0,200,le=201), nsim=nsim, nrank=nrank, correction="isotropic")

plot(cosa$k01, sqrt(./pi)-r~r,  col=c(3, 1, 3), lty=c(3, 1, 3), las=1,
         ylab=expression(L[12]), xlim=c(0, 200), 
	 main="adult HS vs. dead seedlings")

plot(cosa$k02, sqrt(./pi)-r~r, col=c(3, 1, 3), lty=c(3, 1, 3), las=1, 
         ylab=expression(L[12]), xlim=c(0, 200),
 	 main="adult HS vs. surviving seedlings")

Differences between univariate and bivariate K-functions

Description

Given two point patterns I and J, K1K2computes the differences between both univariate K-functions (i.e. Ki(r)-Kj(r)) as well as the differences between the univariate and the bivariate K-function (i.e. Ki(r)-Kij(r) and Kj(r)-Kij(r)). It also computes simulation envelopes to test that that the observed differences are within the range expected asuming the random labelling hypothesis.

Usage

K1K2(X, i, j, nsim = 99, nrank = 1, r = NULL,
	 correction = "isotropic")

Arguments

Details

The indiscriminate use of the raw bivariate functions (mainly the K or the L-bivariate functions) in ecological studies for testing the association/ repulsion between different point patterns waste some of the most interesting properties of the K-function. One of them is that under the random labelling hypothesis every individual pattern would be a random thinning of the corresponding bivariate pattern and therefore Ki(r)=Kj(r)= Kij(r)=pi*r^2 (Diggle 2003). Dixon (2002) sugested that some differences of these functions could provide provide interesting ecological information. For example, D(r)= Ki(r)-Kj(r), has an expected value of 0 for all r distances under random labelling and evaluates the differences in the intensity of aggregation of the two point patterns (e.g., in the example bellow, the pattern of drought and herbivory deaths). Other relevant function is D(r) = Ki(r)-Kij(r) and the complementary D(r)= Kj(r)-Kij(r) which evaluate the degree of segregation of every individual pattern, i.e. if every point of the pattern is more -or less- surrounded by other points of the same type than would be expected under the random labelling hypothesis. K1K2 uses K^*ij(r), the combined estimator of Lotwick and Silverman (a weigthed mean of Kij(r) and Kji(r)) as computed by Kmulti.ls.

Value

A list with three elements.

Each of the above elements is a fv.object, essentially a data.frame with the following items:

Author(s)

Marcelino de la Cruz

References

De la Cruz, M. 2006. Introducción al análisis de datos mapeados o algunas de las (muchas) cosas que puedo hacer si tengo coordenadas. Ecosistemas 15 (3): 19-39.

De la Cruz, M., Romao, R.L., Escudero, A. and Maestre, F.T. 2008. Where do seedlings go? A spatio-temporal analysis of early mortality in a semiarid specialist. Ecography,31(6): 720-730. doi:10.1111/j.0906-7590.2008.05299.x.

Diggle, P.J. 2003. Statistical analysis of spatial point patterns. Arnold, London.

Dixon, P. M. 2002. Ripley's K function. In The encyclopedia of environmetrics (eds. El-Shaarawi, A.H. & Piergorsch, W.W.), pp. 1976-1803. John Wiley & Sons Ltd, NY.

Examples

data(Helianthemum)

# set the number of simulations (nsim=199 or larger for real analyses)
nsim<- 19

cosa12 <- K1K2(Helianthemum, j="deadpl", i="survpl", r=seq(0,200,le=201),
		 nsim=nsim, correction="isotropic")

## plots of figure 9 in De la Cruz (2006) (they where made with nsim=999)
plot(cosa12$k1k2, lty=c(2, 1, 2), col=c(2, 1, 2), xlim=c(0, 200),
         main= "survival- death")

plot(cosa12$k1k12, lty=c(2, 1, 2), col=c(2, 1, 2), xlim=c(0, 200),
	 main="segregation of surviving seedlings")

plot(cosa12$k2k12, lty=c(2, 1, 2), col=c(2, 1, 2), xlim=c(0, 200),
         main= "segregation of dying seedlings", legend=FALSE)

Test against non-Poisson (in-)homogeneous models

Description

Functions to automate testing of 'univariate' and 'bivariate' point pattern hypothesis against non-Poisson (in-)homogeneous models.

Usage

Kci(mod1, mod2, correction="trans", nsim=99, ngrid=200, nrep=1e+05,
     r=NULL, simu="both", spctype=1)

Ki(mod1, correction="trans", nsim=99, ngrid=200, nrep=1e+05, r=NULL,
    spctype=1)

## S3 ploth method for objects of class 'ecespa.kci':
## S3 method for class 'ecespa.kci'
plot(x, type=1, q=0.025, kmean=TRUE, add=FALSE, maine=NULL, 
       xlabe=NULL, ylabe=NULL, xlime=NULL, ylime=NULL, 
       lty=c(1,2,3), col=c(1,2,3), lwd=c(1,1,1), ...)

Arguments

Details

These functions are designed to automate the testing of 'univariate' and(/or) 'bivariate' point pattern hypotheses (based on K-functions) against non-Poisson (in-)homogeneous models. These non-Poisson (in-)homogeneous models should have been fitted with pseudolikelihood tools (spatstat ppm function) or with minimum contrast methods (ecespa.minconfit).

Function Ki is designed to test 'univariate' hypotheses. It will compute the (in-)homogeneous K-function (using spatstat Kinhom function) of the point pattern to which the ppm or ecespa.minconfit model has beeen fitted and will compute 'envelopes' simulating from the fitted model. The computed envelopes can be considered as a pointwise test of the point pattern been a realisation of the fitted model.

Function Kci is designed to test 'bivariate' hypotheses. It will compute the (in-)homogeneous cross K-function (using spatstat Kcross.inhom function) and will compute envelopes simulating from the fitted models. As, when dealing with inhomogeneos patterns K12 != K21, Kci will compute both functions. If simu = "both" (default option), Kci will simulate from mod2 to test K12 and from mod1 to test K21. If simu != "both", only mod2 will be simulated. This option may be useful when only K12 is of interest. Function Kci will also compute univariate (in-) homogeneous K-functions and envelopes for each individual point pattern.

The S3 ploth method will plot the selected K-function and envelopes (actually, it will plot the most usual L-function = sqrt[K(r)/pi]-r). The appropriate K function can be selected with the argument type. If type = 1 (default option), it will plot the univariate K function (of the analized model in Ki or of the first model [mod1] in Kci). If type = 2, it will plot the univariate K function of the second model (mod2 in Kci). When type = 12 or type = 21, it will plot respectively K12 or K21. Options type = 112 or type = 221 will graph a kind of 'segregation test' (see K1K2), and will represent de differences K1-K12, K2-K21 and their envelopes.

Value

Both Kci and Ki return an object of class 'ecespa.kci', basically a list with the following items:

Warning

As this implementation involves the use of images as the means of evaluation of the (inhomogeneous) spatial trend, and a mask based on those images will be used as the point pattern window, the "Ripley's" or "isotropic" edge correction can not be employed.

It is usual that during the simulation process some warnings are produced. They are related to some simulated points being rejected as lying outside the specified window.

Note

Even when one of the two point patterns is assumed to be homogeneous Poisson (and, apparently not worth of fitting any model), an homogeneous Poisson model can be easily fitted and passed to Kci with ppm. See the examples.

Author(s)

Marcelino de la Cruz Rot

References

De la Cruz, M. and Escudero, A. 2008. Null models and tools for multivariate heterogeneous point patterns. Submitted.

De Soto, L., Olano, J.M., Rozas, V. and De la Cruz, M. 2009. Release of Juniperus thurifera woodlands from herbivore-mediated arrested succession in Spain. Applied Vegetation Science, 13: 1525. doi:10.1111/j.1654-109X.2009.01045.x

Examples

  
    
    data(urkiola)
    
    # set the number of simulations (nsim=199 or larger for real analyses)
    nsim<- 5

  #####################
   ## univariate case

   # get univariate pp
   I.ppp <- split.ppp(urkiola)$birch

   # estimate inhomogeneous intensity function
   I.lam <- predict (ppm(I.ppp, ~polynom(x,y,2)), type="trend", ngrid=200)

   # Compute and plot envelopes to Kinhom, simulating from an Inhomogeneous
   #  Poisson Process:
   
   I2.env <- envelope( I.ppp,Kinhom, lambda=I.lam, correction="trans", 
                              nsim=nsim, simulate=expression(rpoispp(I.lam)))
   plot(I2.env, sqrt(./pi)-r~r, xlab="r (metres)", ylab= "L (r)", col=c(1,3,2,2),legend=FALSE) 

   # It seems that there is short scale clustering; let's fit an Inhomogeneous 
   # Poisson Cluster Process: 

   I.ki <- ipc.estK(mippp=I.ppp, lambda=I.lam, correction="trans")

   # Compute and plot envelopes to Kinhom, simulating from the fitted IPCP:

   Ipc.env <- Ki(I.ki, correction="trans", nsim=nsim, ngrid=200)

   plot (Ipc.env, xlab="r (metres)", ylab= "L (r)")
 
  ## Not run: 
   #####################
   ## bivariate case: test independence between birch and quercus in Urkiola

   J.ppp <- split.ppp(urkiola)$oak
   
   # We want to simulate oak from a homogeneous Poisson model:
   J.ppm <- ppm(J.ppp, trend=~1, interaction=Poisson() )
   
   IJ.env <- Kci (mod1=I.ki, mod2=J.ppm, nsim=nsim)
   
   plot(IJ.env, type=12)
   
   plot(IJ.env, type=21)
 
## End(Not run)

Simulation envelopes from the fitted values of a logistic model

Description

Computes simulation envelopes for (in-)homogeneous K-function simulating from a vector of probabilitiesn.

Usage

Kinhom.log (A, lambda=NULL, mod=NULL, lifemark="0", prob=NULL,
			r=NULL, nsim=99, correction="trans", ngrid=200)

Arguments

Details

This function is a wrapper to compute the critical envelopes for Monte Carlo test of goodness-of-fit of (in-)homogeneous K functions, simulating from the fittted values of a logistic model (i.e. a binomial GLM with logit link) fitted to the marks ("failure", "success") of a "binomially"-marked point pattern. This is particularly interesting in plant ecology when considering alternatives to the random mortality hypothesis (Kenkel 1988). This hypothesis is usually tested building Monte Carlo envelopes from the "succesful" patterns resulting from a random labelling of a "binomially"-marked point pattern (this is equivalent to a random thinning of the whole pattern irrespective of the marks). As tree mortality is rarely random but instead can be modelled as a function of a certain number of covariates, the most natural alternative to the random mortality hypothesis is the logistic mortality hypothesis, that can be tested thinning the original pattern of trees with retention probabilities defined by the fitted values of a logistic model (Batista and Maguire 1998, Olano et al. 2008).

Kinhom.log will compute the envelopes by thinning the unmarked point pattern A with retention probabilities prob. If no prob vector is provided, all points will be thinned with the same probability ( number of "live" points / number of points ), i.e. Kinhom.log will compute random thinning envelopes.

Kinhom.log will compute envelopes both to homogeneous and inhomogeneous K functions. If no lambda or mode arguments are provided, Kinhom.log assumes that the original pattern is homogeneous and will use a constant lambda to compute the inhomogeneous K (i.e. it will compute the homogeneous K). The most convenient use with inhomogeneous point patterns is to provide the argument mod with an inhomogeneous Poisson model fitted to the original pattern of 'live' points (with spatstat function ppm; see the examples). This model will be used to compute (and to update in the simulations) the inhomogeneous trend (i.e. the "lambda") of the patterns. If the argument lambda is provided but not mod, these lambda will be used as a covariate to fit an inhomogeneous Poisson model that will be used to compute (and to update in the simulations) the inhomogeneous spatial trend.

Kinhom.log will produce an object of class 'ecespa.kci' that can be easily ploted (see the examples). This is accomplished by the S3 ploth method plot.ecespa.kci; it will plot the K-function and its envelopes (actually, it will plot the most usual L-function = sqrt[K(r)/pi]-r).

Value

Kinhom.log returns an object of class ecespa.kci, basically a list with the following items:

Warning

As this implementation involves the use of images as the means of evaluation of the (inhomogeneous) spatial trend, and a mask based on those images will be used as the point pattern window, the "Ripley's" or "isotropic" edge correction can not be employed.

Author(s)

Marcelino de la Cruz Rot

References

Batista, J.L.F. and Maguire, D.A. 1998. Modelling the spatial structure of tropical forests. For. Ecol. Manag., 110: 293-314.

Kenkel, N.C. 1988. Pattern of self-thinning in Jack Pine: testing the random mortality hypothesis. Ecology, 69: 1017-1024.

Olano, J.M., Laskurain, N.A., Escudero, A. and De la Cruz, M. 2009. Why and where adult trees die in a secondary temperate forest? The role of neighbourhood. Annals of Forest Science, 66: 105. doi:10.1051/forest:2008074.

Examples

  
   
   data(quercusvm)
   
   # set the number of simulations (nsim=199 or larger for real analyses)
   nsim<- 19

   # read fitted values from logistic model:
   
   
   probquercus <-c(0.99955463, 0.96563477, 0.97577094, 0.97327199, 0.92437309,
   0.84023396, 0.94926682, 0.89687281, 0.99377915, 0.74157478, 0.95491518,
   0.72366493, 0.66771787, 0.77330148, 0.67569082, 0.9874892, 0.7918891, 
   0.73246803, 0.81614635, 0.66446411, 0.80077908, 0.98290508, 0.54641754,
   0.53546689, 0.73273626, 0.7347013, 0.65559655, 0.89481468, 0.63946334,
   0.62101995, 0.78996371, 0.93179582, 0.80160346, 0.82204428, 0.90050059,
   0.83810669, 0.92153079, 0.47872421, 0.24697004, 0.50680935, 0.6297911, 
   0.46374812, 0.65672284, 0.87951682, 0.35818237, 0.50932432, 0.92293014,
   0.48580241, 0.49692053, 0.52290553, 0.7317549, 0.32445982, 0.30300865,
   0.73599359, 0.6206056, 0.85777043, 0.65758613, 0.50100406, 0.31340849, 
   0.22289286, 0.40002879, 0.29567678, 0.56917817, 0.56866864, 0.27718552,
   0.4910667, 0.47394411, 0.40543788, 0.29571349, 0.30436276, 0.47859015,
   0.31754526, 0.42131675, 0.37468782, 0.73271225, 0.26786274, 0.59506388, 
   0.54801851, 0.38983575, 0.64896835, 0.37282031, 0.67624306, 0.29429766,
   0.29197755, 0.2247629, 0.40697843, 0.17022391, 0.26528042, 0.24373722,
   0.26936163, 0.13052254, 0.19958585, 0.18659692, 0.36686678, 0.47263005,
   0.39557661, 0.68048997, 0.74878567, 0.88352322, 0.93851375)
   
  

   ################################ 
   ## Envelopes for an homogeneous point pattern:
   
   cosap <- Kinhom.log(A=quercusvm, lifemark="0",  prob=probquercus, nsim=nsim)

   plot(cosap)

   
   ################################ 
   ## Envelopes for an inhomogeneous point pattern:
   
   ## First, fit an inhomogeneous Poisson model to alive trees :
   
   quercusalive <- unmark(quercusvm[quercusvm$marks == 0])

    mod2 <- ppm(quercusalive, ~polynom(x,y,2))

    ## Now use mod2 to estimate lambda for K.inhom:
    
    cosapm <- Kinhom.log(A=quercusvm, lifemark="0", prob=probquercus, 
                                   nsim=nsim, mod=mod2)

   
   ################################ 
   ## An example of homogeneous random thinning:
      
   cosa <- Kinhom.log(A=quercusvm, lifemark="0", nsim=nsim)
   
   plot(cosa)
    
    

Mark-weighted K-function

Description

This is a functional data summary for marked point patterns that measures the joint pattern of points and marks at different scales determined by r.

Usage

Kmm(mippp, r = 1:10, nsim=NULL)

## S3 method for ploting objects of class 'ecespa.kmm':
## S3 method for class 'ecespa.kmm'
plot(x, type="Kmm.n", q=0.025, 
            xlime=NULL, ylime=NULL,  maine=NULL, add=FALSE, kmean=TRUE,
            ylabe=NULL, xlabe=NULL, lty=c(1,2,3), col=c(1,2,3), lwd=c(1,1,1),
             ...)

Arguments

Details

Penttinnen (2006) defines Kmm(r), the mark-weighted K-function of a stationary marked point process X, so that

lambda*Kmm(r) = Eo[sum(mo*mn)]/mu^2

where lambda is the intensity of the process, i.e. the expected number of points of X per unit area, Eo[ ] denotes expectation (given that there is a point at the origin); m0 and mn are the marks attached to every two points of the process separated by a distance <= r and mu is the mean mark. It measures the joint pattern of marks and points at the scales determmined by r. If all the marks are set to 1, then lambda*Kmm(r) equals the expected number of additional random points within a distance r of a typical random point of X, i.e. Kmm becomes the conventional Ripley's K-function for unmarked point processes. As the K-function measures clustering or regularity among the points regardless of the marks, one can separate clustering of marks with the normalized weighted K-function

Kmm.normalized(r) = Kmm(r)/K(r)

If the process is independently marked, Kmm(r) equals K(r) so the normalized mark-weighted K-function will equal 1 for all distances r.

If nsim != NULL, Kmm computes 'simulation envelopes' from the simulated point patterns. These are simulated from nsim random permutations of the marks over the points coordinates. This is a kind of pointwise test of Kmm(r) == 1 or normalized Kmm(r) == 1 for a given r.

Value

Kmm returns an object of class 'ecespa.kmm', basically a list with the following items:

Note

This implementation estimates Kmm(r) without any correction of border effects, so it must be used with caution. However, as K(r) is also estimed without correction it migth compensate the border effects on the normalized Kmm-function.

Author(s)

Marcelino de la Cruz Rot

References

De la Cruz, M. 2008. Métodos para analizar datos puntuales. En: Introducción al Análisis Espacial de Datos en Ecología y Ciencias Ambientales: Métodos y Aplicaciones (eds. Maestre, F. T., Escudero, A. y Bonet, A.), pp 76-127. Asociación Española de Ecología Terrestre, Universidad Rey Juan Carlos y Caja de Ahorros del Mediterráneo, Madrid.

Penttinen, A. 2006. Statistics for Marked Point Patterns. In The Yearbook of the Finnish Statistical Society, pp. 70-91.

See Also

markcorr

Examples

  ## Figure 3.10 of De la Cruz (2008):
  # change r to r=1:100
  
   r = seq(1,100, by=5)
  
  data(seedlings1)
  
  data(seedlings2)
  
  s1km <- Kmm(seedlings1, r=r)
  
  s2km <- Kmm(seedlings2, r=r)
  
  plot(s1km, ylime=c(0.6,1.2), lwd=2, maine="", xlabe="r(cm)")

  plot(s2km,  lwd=2, lty=2, add=TRUE )

  abline(h=1, lwd=2, lty=3)
  
  legend(x=60, y=1.2, legend=c("Hs_C1", "Hs_C2", "H0"),
	 lty=c(1, 2, 3), lwd=c(3, 2, 2), bty="n")
## Not run: 
## A pointwise test of normalized Kmm == 1 for seedlings1:

   s1km.test <- Kmm(seedlings1, r=1:100, nsim=99)

   plot(s1km.test,  xlabe="r(cm)")

  
## End(Not run)

Lotwick's and Silverman's combined estimator of the marked K-function

Description

For a multitype point pattern, calculates the combined estimator of the bivariate Kij(r) and Kji(r) functions.

Usage

Kmulti.ls(X, I, J, r = NULL, corre = "isotropic")

Arguments

Details

As a consequence of edge effects, the estimators Kij(r) and Kji(r) of the same bivariate pattern could differ. K^*ij(r) is the combined estimator defined by Lotwick and Silverman (1982) as

nj*Kij(r)+ ni*Kji(r) / (ni + nj) ,

ni and nj being respectively the number of points in I and J.

Value

An object of class "fv" (see fv.object). Essentially a data frame containing numeric columns

.

.

together with a column or columns named "border", "bord.modif", "iso" and/or "trans", according to the selected edge corrections. These columns contain estimates of the function K^*ij(r) obtained by the edge corrections named.

Note

Kmulti.ls is a wrapper for a convenient use of the Kmulti function of spatstat. Please refer to its help page for additional documentation.

Author(s)

Marcelino de la Cruz

References

Lotwick,H.W. & Silverman, B. W. 1982. Methods for analysing spatial processes of several types of points. Journal of the Royal Statistical Society B, 44: 406-413. doi:10.1111/j.2517-6161.1982.tb01221.x.

Examples

data(amacrine)

plot(Kmulti.ls(amacrine, I=amacrine$marks=="on", J=amacrine$marks=="off", 
	 corre="isotropic"), sqrt(./pi)-r~r, main="")

# compare with Kmulti

plot(Kmulti(amacrine, I=amacrine$marks=="on", J=amacrine$marks=="off"),
         sqrt(iso/pi)-r~r, add=TRUE, col=3)

plot(Kmulti(amacrine, J=amacrine$marks=="on", I=amacrine$marks=="off"),
      sqrt(iso/pi)-r~r, add=TRUE, col=4)

Loosmore and Ford Goodness of Fit Test

Description

Performs the Loosmore and Ford (2006) test or the Maximum Absolute Deviation test for a spatial point pattern.

Usage

  LF.gof(X, rmin=NULL, rmax=NULL, na.rm=TRUE)

Arguments

Details

These function perform a tests for goodness-of-fit of a point pattern dataset to a point process model, based on Monte Carlo simulation from the model. The simulations should have been previously computed with the function envelope, applied with the argument savefuns=TRUE in order to save all the simulated functions, required for the computation of the test.

The test, popularized in the ecological field by Loosmore and Ford (2006) is also described in Diggle (2003, page 14), and according to Baddeley and Turner (2005) also in Diggle (1986) and Cressie (1991, page 667, equation (8.5.42)). If the arguments rmin and rmax are set to NULL, the integral of the GoF statistics will be computed over the complete range of r values.

Value

A list with the following components:

u

The GoF statistic, i.e., the value of the integral over the range of r's

p

The p-value of the test

na.count.by.r

Number of NA values for each r. It helps to evaluate the reliability of the computed u's, specially for small r's

Author(s)

Marcelino de la Cruz

References

Cressie, N.A.C. (1991) Statistics for spatial data. John Wiley and Sons, 1991.

Diggle, P. J. (1986). Displaced amacrine cells in the retina of a rabbit : analysis of a bivariate spatial point pattern. J. Neuroscience Methods 18, 115-125.

Diggle, P.J. (2003) Statistical analysis of spatial point patterns, Second edition. Arnold.

Loosmore, N.B. and Ford, E.D. (2006) Statistical inference using the G or K point pattern spatial statistics. Ecology 87, 1925-1931. doi:10.1890/0012-9658(2006)87[1925:SIUTGO]2.0.CO;2.

See Also

dclf.test for an alternative implementation of the test in spatstat.

Dixon (2002) Nearest-neighbor contingency table analysis

Description

dixon2002 is a wrapper to the functions of Dixon (2002) to test spatial segregation for several species by analyzing the counts of the nearest neighbour contingency table for a marked point pattern.

Usage

dixon2002(datos, nsim = 99)

Arguments

Details

A measure of segregation describes the tendency of one species to be associated with itself or with other species. Dixon (2002) proposed a measure of the segregation of species i in a multiespecies spatial pattern as:

S[i] = log{[(N[ii]/(N[i]-N[ii])] / [(N[i]-1)/(N-N[i])]}

where N[i] is the number of individuals of species i, N[ii] is the frequency of species i as neighbor of especies i and N is the total number of locations. Values of S[i] larger than 0 indicate that species i is segregated; the larger the value of S[i], the more extreme the segregation. Values of S[i] less than 0 indicate that species i is is found as neighbor of itself less than expected under random labelling. Values of S[i] close to 0 are consistent with random labelling of the neighbors of species i.

Dixon (2002) also proposed a pairwise segregation index for the off-diagonal elements of the contingency table:

S[ij] = log{[(N[ij]/(N[i]-N[ij])] / [(N[i])/(N-N[j])-1]}

S[ij] is larger than 0 when N[ij], the frequency of neighbors of species j around points of species i, is larger than expected under random labelling and less than 0 when N[ij] is smaller than expected under random labelling.

As a species/neighbor-specific test, Dixon(2002) proposed the statistic

Z[ij] =(N[ij] -EN[ij])/sqrt(Var N[ij])

where j may be the same as i and EN[ij] is the expected count in the contingency table. It has an asymptotic normal distribution with mean 0 and variance 1; its asymptotic p-value can be obtained from the numerical evaluation of the cumulative normal distribution; when the sample size is small, a p-value on the observed counts in each cell (N[ij]) may be obtained by simulation, i.e, by condicting a randomization test.

An overall test of random labelling (i.e. a test that all counts in the k x k nearest-neighbor contingency table are equal to their expected counts) is based on the quadratic form

C = (N-EN)' Sigma^- (N - EN)

where N is the vector of all cell counts in the contingency table, Sigma is the variance-covariance matrix of those counts and Sigma^- is a generalized inverse of Sigma. Under the null hypothesis of random labelling of points, C has a asymptotic Chi-square distribution with k(k-1) degrees of freedom (if the sample sizes are small its distribution should be estimated using Monte-Carlo simulation). P-values are computed from the probability of observing equal or larger values of C. The overall statistic C can be partitioned into k species-specific test statistics C[i]. Each C[i] test if the frequencies of the neighbors of species i are similar to the expected frequencies if the points were randomly labelled. Because the C[i] are not independent Chi-square statistics, they do not sum to the overall C.

Value

A list with the following components:

Warning

The S[i] and S[ij] statistics asume that the spatial nearest-neighbor process is stationary, at least to second order, i.e., have the same sign in every part of the entire plot. A biologically heterogeneous process will violate this asumption.

Author(s)

Philip M. Dixon . Marcelino de la Cruz wrote the wrapper code for the ecespa version.

References

Dixon, P.M. 2002. Nearest-neighbor contingency table analysis of spatial segregation for several species. Ecoscience, 9 (2): 142-151. doi:10.1080/11956860.2002.11682700.

See Also

K012 for another segregation test, based in the differences of univariate and bivariate K-functions. A faster version of this function, with code implemented in FORTRAN it is available in function dixon in dixon.

Examples

data(swamp)

dixon2002(swamp,nsim=99)

Internal ecespa functions.

Description

Internal ecespa functions.

Usage

check(x, v, l1, l2)
ginv(m)
mNNinfo2(n, R, Q)
mNNinfo(xy, label, nnid = NULL)
mNNtest(info, obsN = NULL)
NNid(xy)

Details

These fucntions are called by dixon2002 and are usually not to be called by the user.

Author(s)

Philip M. Dixon

Artificial point data.

Description

The three different point patterns in the figure 3.1 of De la Cruz (2008)

Usage

data(fig1)
data(fig2)
data(fig3)

Format

A data frame with 87 observations on the following 2 variables.

x

x coordinate

y

y coordinate

References

De la Cruz, M. 2008. Métodos para analizar datos puntuales. En: Introducción al Análisis Espacial de Datos en Ecología y Ciencias Ambientales: Métodos y Aplicaciones (eds. Maestre, F. T., Escudero, A. y Bonet, A.), pp 76-127. Asociación Española de Ecología Terrestre, Universidad Rey Juan Carlos y Caja de Ahorros del Mediterráneo, Madrid.

Examples

## Not run: 
data(fig1)

data(fig2)

data(fig3)

# transform to ppp format of spatstat with function haz.ppp:

fig1.ppp <- haz.ppp(fig1)

fig2.ppp <- haz.ppp(fig2)

fig3.ppp <- haz.ppp(fig3)

#Analyses as in Fig.3.2 of De la Cruz (2008). First, compute function K:

cosa1 <- Kest(fig1.ppp)

# Plot different estimators. 
# Fig. 3.2a: 

par("mar"=par("mar")+c(0,1,0,0))

plot(cosa1, col=c(1,0,0,1), lwd=c(2,2,2,2), lty=c(1,1,1,2), 
	 main="")

# Fig. 3.2b:

plot(cosa1, sqrt(./pi)-r~r, col=c(1,0,0,1), lwd=c(2,2,2,2), 
         lty=c(1,1,1,2), main="", ylab="L(r)")

# Fig. 3.2c:

plot(cosa1, .-(pi*r^2)~r, col=c(1,0,0,1), lwd=c(2,2,2,2), 
         lty=c(1,1,1,2), main="", ylab=expression(K(r)-pi*r^2))

# Fig. 3.2d:

plot(cosa1,(./(pi*r^2))-1~r, col=c(1,0,0,1), lwd=c(2,2,2,2), 
	 lty=c(1,1,1,2), main="",
         ylab=expression((K(r)/pi*r^2)-1))

## Analyses as in fig. 3.7 of De la Cruz (2008).
## First, compute function K and pointwise envelopes:

cosa1.env <- envelope(fig1.ppp, Kest)

cosa2.env <- envelope(fig2.ppp, Kest)

cosa3.env <- envelope(fig3.ppp, Kest)

## Plot function L with pointwise envelopes:

plot(cosa1.env,sqrt(./pi)-r~r, lwd=c(1,1,2,2), 
	 lty=c(1,1,3,3), col=c(1,1,1,1), xlab="r", 
         ylab="L(r)", main="", ylim=c(-2,2))

## Add simultaneous envelopes of Ripley (+-1.68 *sqrt(A)/N):

abline(h=1.68*sqrt(area.owin(fig1.ppp$window))/fig1.ppp$n, 
         lty=2, lwd=2)

abline(h=-1.68*sqrt(area.owin(fig1.ppp$window))/fig1.ppp$n, 
         lty=2, lwd=2)

## Plot function L with pointwise envelopes:

plot(cosa2.env,sqrt(./pi)-r~r, lwd=c(1,1,2,2), 
         lty=c(1,1,3,3), col=c(1,1,1,1), xlab="r", 
         ylab="L(r)", main="")

## Add simultaneous envelopes of Ripley:

abline(h=1.68*sqrt(area.owin(fig2.ppp$window))/fig2.ppp$n, 
         lty=2, lwd=2)

abline(h=-1.68*sqrt(area.owin(fig2.ppp$window))/fig2.ppp$n, 
         lty=2, lwd=2)

## Plot function L with pointwise envelopes:

plot(cosa3.env,sqrt(./pi)-r~r, lwd=c(1,1,2,2), 
         lty=c(1,1,3,3), col=c(1,1,1,1), xlab="r", 
         ylab="L(r)", main="")

## Add simultaneous envelopes of Ripley:

abline(h=1.68*sqrt(area.owin(fig3.ppp$window))/fig3.ppp$n, 
         lty=2, lwd=2)

abline(h=-1.68*sqrt(area.owin(fig3.ppp$window))/fig3.ppp$n, 
         lty=2, lwd=2)

## End(Not run)

Neighbourhood density function

Description

Computes and plots the neighbourhood density function, a local version of the K-function defined by Getis and Franklin (1987).

Usage

getis(mippp, nx = 30, ny = 30, R = 10)

## S3 method for class 'ecespa.getis'
plot(x, type="k", dimyx=NULL, xy=NULL, eps=NULL,  color=NULL,
         contour=TRUE, points=TRUE,...)

Arguments

Details

Getis and Franklin (1987) proposed the neigbourhood density function, a local version of Ripley's L- function. Given a spatial point pattern X, the neigbourhood density function associated with the ith point in X is computed by

L[i](r) = sqrt((a/((n-1))*pi))*sum[j]e[i,j])

where the sum is over all points j != i that lie within a distance r of the ith point, a is the area of the observation window, n is the number of points in X, and e[i,j] is the isotropic edge correction term (as described in Kest). The value of L[i](r) can also be interpreted as one of the summands that contributes to the global estimate of the L-function.

The command getis actually computes the local K-function using Kcross. As the main objective of getis is to map the local density function, as sugested by Gestis and Franklin (1987: 476) a grid of points (whose density is controled by nx and ny), is used to accurately estimate the functions in empty or sparse areas. The S3 method plot.ecespa.getis plots the spatial distribution of the local K or L function or other related local statistics, such as n[i](r), the number of neighbor points [= lambda*K[i](r)] or the deviations from the expected value of local L under CSR [= L[i](r) -r]. It some of the arguments dimyx, xy or eps is provided it will use the function interp.im in spatstat package to interpolate the results;otherwise it will plot the estimated values at the origial grid points.

Value

getis gives an object of class ecespa.getis, bassically a list with the following elements:

plot.ecespa.getis plots an interpolated map of the selected local statistics

Note

As plot.ecespa.getis interpolates over rectangular grid of points, it is not apropriate to map irregular windows. In those cases, Smooth.ppp of spatstat can be used to interpolate the local statistics (see examples).

Author(s)

Marcelino de la Cruz Rot

References

Getis, A. and Franklin, J. 1987. Second-order neighbourhood analysis of mapped point patterns. Ecology 68: 473-477. doi:10.2307/1938452.

See Also

localK, a different approach in spatstat.

Examples

  ## Compare with fig. 5b of Getis & Franklin (1987: 476):
  
  data(ponderosa)
  
  #ponderosa12 <- getis(ponderosa, nx = 30, ny = 30, R = 12)
  ponderosa12 <- getis(ponderosa, nx = 20, ny = 20, R = 12)
  
  plot(ponderosa12, type = "l", dimyx=256)
  
## Not run: 
  ## Plot the same, using Smooth.ppp in spatstat
  
  ponderosa.12 <- setmarks(ponderosa, ponderosa12$klocal)
  
  Z <- Smooth(ponderosa.12, sigma=5, dimyx=256)
  
  plot(Z, col=topo.colors(128), main="smoothed neighbourhood density")
  
  contour(Z, add=TRUE)
  
  points(ponderosa, pch=16, cex=0.5) 
  

  ## Example with irregular window:
  
  data(letterR)
  
  X <- rpoispp(50, win=letterR)
  
  X.g <- getis(X, R=0.2)
  
  plot(X.g,dimyx=c(200,100))
  
   ## Plot the same, using Smooth.ppp in spatstat
    X2 <- setmarks(X, X.g$klocal)
  
    Z <- Smooth(X2, sigma=0.05, dimxy=256)
  
    plot(Z, col=topo.colors(128), main="smoothed neighbourhood density")
  
    contour(Z, add=TRUE)
  
    points(X, pch=16, cex=0.5)
  
    
    
## End(Not run)

Spatial point pattern of a plant community

Description

Locations of plants in a gypsophylous plant community in Central Spain. These are part of the data collected by Romao (2003) that have been analyzed several times (Escudero et al. 2005, De la Cruz 2006). The coordinates of the plans are given in cm.

Usage

data(gypsophylous)

Format

An object of class "ppp" of spatstat representing the point pattern of plants locations. See ppp.object for details of the format.

Source

Romao, R.L. 2003. Estructura espacial de comunidades de gipsófitos: interacciones bióticas y constricciones abióticas. Tesis Doctoral. Universidad Politécnica de Madrid.

References

De la Cruz, M. 2006. Introducción al análisis de datos mapeados o algunas de las (muchas) cosas que puedo hacer si tengo coordenadas. Ecosistemas 15 (3): 19-39.

Escudero, A., Romao, R.L., De la Cruz, M. & Maestre, F. 2005. Spatial pattern and neighbour effects on Helianthemum squamatum seedlings in a Mediterranean gypsum community. Journal of Vegetation Science, 16: 383-390. doi:10.1111/j.1654-1103.2005.tb02377.x.

Examples

## Not run: 

data(gypsophylous)

plot(gypsophylous)


## End(Not run)

Easily convert xy data to ppp format

Description

ppp maker for the impatient layman

Usage

haz.ppp(W)

Arguments

Details

This naive function easily transform your xy data to the format required by spatstat (version <2.0). It establishes the window of observation as the rectangle defined by the xy range. It asumes that the first two columns are coordinates x and y, and the third (if any) gives the marks of the points.

Value

A point pattern with the format of spatstat v. < 2.0

Note

This window is an underestimate of the real sampling window of the observed points. For a less biassed estimate, use function ripras in spatstat.

Author(s)

Marcelino de la Cruz

See Also

ppp, as.ppp

Examples

data(fig1)

plot(fig1) #typical xyplot

fig1.ppp <- haz.ppp (fig1)

fig1.ppp

plot(fig1.ppp) # point pattern plot of spatstat

Fit the (In)homogeneous Poisson Cluster Point Process by Minimum Contrast

Description

Fits the (In)homogeneous Poisson Cluster point process to a point pattern dataset by the Method of Minimum Contrast.

Usage

ipc.estK(mippp, lambda = NULL, correction = "iso", r = NULL, sigma2 = NULL,
            rho = NULL, q = 1/4, p = 2)


## S3 method for class 'ecespa.minconfit'
plot(x, type="L", add=FALSE, xlim=NULL, ylim=NULL, lwd=c(1,1),
           lty=c(1,2), col=c(1,2), main=NULL, ...)

Arguments

Details

The algorithm fits the (inhomogeneous) Poisson cluster point process (PCP) to a point pattern, by finding the parameters of the (inhomogeneous) Poisson cluster model which give the closest match between the theoretical K function of the Poisson cluster process and the observed K function. For a concise explanation of the PCP see pc.estK. For a more detailed explanation of the Method of Minimum Contrast, see mincontrast in spatstat or Diggle (2003: 86).

The inhomogeneous PCP can be thought of as a thinned process of an homogeneous PCP, where the spatially varying thinning probability f(s) is related to the spatially varying intensity function lambda(s) as f(s) = lambda(s)/max{lambda(s)} (Waagepetersen, 2007). As the inhomogeneous K function for the IPCP coincides with the (homogeneous) K function for the corresponding homogeneous PCP, the parameters of the underlying homomgeneous PCP can be estimated as those that give the closest match between the theoretical K function for the homogeneous PCP and the empirical inhomogeneous K function for the observed IPCP.

This Poisson cluster process can be simulated with rIPCP.

Value

ipc.estK gives an object of class 'ecespa.minconfit', basically a list with the following components:

Author(s)

Marcelino de la Cruz Rot , inspired by some code of Philip M. Dixon

References

Diggle, P. J. 2003. Statistical analysis of spatial point patterns. Arnold, London.

Waagepetersen, R. P. 2007. An estimating function approach to inference for inhomogeneous Neymann-Scott processes. Biometrics 63: 252-258. doi:10.1111/j.1541-0420.2006.00667.x.

See Also

some functions in spatstat: mincontrast for a more general implementation of the method of mimimum contrast; matclust.estK and lgcp.estK fit other appropriate processes for clustered patterns.

Examples

#####################
## Same example as in pc.estK

data(gypsophylous)

# set the number of simulations (nsim=199 or larger for real analyses)
# and resolution for lambda (ngrid=200  or larger for real analyses)
nsim<- 19
ngrid <- 50

## Estimate K function ("Kobs").

gyps.env <- envelope(gypsophylous, Kest, correction="iso", nsim=nsim)

plot(gyps.env, sqrt(./pi)-r~r, legend=FALSE)

## Fit Poisson Cluster Process. The limits of integration 
## rmin and rmax are setup to 0 and 60, respectively. 

cosa.pc2 <- ipc.estK(gypsophylous, r = gyps.env$r[gyps.env$r<=60])

## Add fitted Kclust function to the plot.


plot(cosa.pc2, add=TRUE, lwd=c(3,3))


## A kind of pointwise test of the gypsophylous pattern been a realisation
## of the fitted model, simulating with rIPCP and using function J (Jest).


gyps.env.sim2 <- envelope(gypsophylous, Jest, nsim=nsim,
                    simulate=expression(rIPCP(cosa.pc2)))


plot(gyps.env.sim2,  main="",legendpos="bottomleft")

## Not run: 

#####################
## Inhomogeneous example

data(urkiola)
   
   #####################
   ## univariate case

   # get univariate pp
   I.ppp <- split.ppp(urkiola)$birch

   # estimate inhomogeneous intensity function
   I.lam <- predict (ppm(I.ppp, ~polynom(x,y,2)), type="trend", ngrid=ngrid)

   # Compute and plot envelopes to Kinhom, simulating from an Inhomogeneous
   #  Poisson Process:
   
   I2.env <- envelope( I.ppp,Kinhom, lambda=I.lam, correction="trans", 
                              nsim=nsim, simulate=expression(rpoispp(I.lam)))
   plot(I2.env, sqrt(./pi)-r~r, xlab="r (metres)", ylab= "L (r)", col=c(1,3,2,2),legend=FALSE) 

   # It seems that there is short scale clustering; let's fit an Inhomogeneous 
   # Poisson Cluster Process: 

   I.ki <- ipc.estK(mippp=I.ppp, lambda=I.lam, correction="trans")

   # Compute and plot envelopes to Kinhom, simulating from the fitted IPCP:

   Ipc.env <- Ki(I.ki, correction="trans", nsim=nsim, ngrid=ngrid)

   plot (Ipc.env, xlab="r (metres)", ylab= "L (r)")
 
   #####################
   ## bivariate case: test independence between birch and quercus in Urkiola

   J.ppp <- split.ppp(urkiola)$oak
   
   # We want to simulate oak from a homogeneous Poisson model:
   J.ppm <- ppm(J.ppp, trend=~1, interaction=Poisson() )
   
   IJ.env <- Kci (mod1=I.ki, mod2=J.ppm, nsim=nsim)
   
   plot(IJ.env, type=12)
   
   plot(IJ.env, type=21)

## End(Not run)

Mark-sum measure

Description

An exploratory data analysis technique for marked point patterns. The marked point pattern is mapped to a random field for visual inspection.

Usage

marksum(mippp, R = 10, nx = 30, ny = 30)

## S3 method for ploting objects of class 'ecespa.marksum':
## S3 method for class 'ecespa.marksum'
plot(x, what="normalized",  contour=FALSE, grid=FALSE,
	ribbon=TRUE,col=NULL ,main=NULL,xlab="",ylab="",...)

Arguments

Details

Penttinen (2006) defines the mark-sum measure as a smoothed summary measuring locally the contribution of points and marks. For any fixed location x within the observational window and a distance R, the mark-sum measure S[R](x) equals the sum of the marks of the points within the circle of radius R with centre in x. The point-sum measure I[R](x) is defined by him as the sum of points within the circle of radius R with centre in x, and describes the contribution of points locally near x. The normalized mark-sum measure describes the contribution of marks near x and is defined (Penttinen, 2006) as

S.normalized[R](x) = S[R](x)/I[R](x)

This implementation of marksum estimates the mark-sum and the point-sum measures in a grid of points whose density is defined by nx and ny.

Value

marksum gives an object of class 'ecespa.marksum'; basically a list with the following elements:

plot.ecespa.marksum plots the selected mark-sum measure.

Author(s)

Marcelino de la Cruz Rot

References

Penttinen, A. 2006. Statistics for Marked Point Patterns. In The Yearbook of the Finnish Statistical Society, pp. 70-91.

See Also

getis, related to the point-sum measure, and markstat for designing different implementations.

Examples

   
 data(seedlings1)
   
 seed.m <- marksum(seedlings1, R=25)

 # raw mark-sum measure; sigma is bandwith for smoothing
 plot(seed.m, what="marksum", sigma = 5)  

 # point sum measure
 plot(seed.m, what="pointsum", sigma = 5) 
   
 # normalized  mark-sum measure
 plot(seed.m,  what="normalized", dimyx=200, contour=TRUE, sigma = 5) 

# the same with added grid and normalized  mark-sum measure
plot(seed.m,  what="normalized", dimyx=200,
      contour=TRUE, sigma = 5, grid=TRUE)

P-value for a discrete distribution on small sample data

Description

Computes the p-value for a two-sided hypothesis test following Dixon's (2002:145) description of the method of Agresti & Min (2001).

Usage

p2colasr(Z)

Arguments

Value

P-value of the two-sided hypothesis test

Note

This function is usually not to be called by the user. It is internally used by dixon2002.

Author(s)

Marcelino de la Cruz Rot

References

Agresti, A. & Min, Y. 2001. On small-sample confidence intervals for parameters in discrete distributions. Biometrics, 57: 963-971.

Dixon, P.M. 2002. Nearest-neighbor contingency table analysis of spatial segregation for several species. Ecoscience, 9(2): 142-151.

Fit the Poisson Cluster Point Process by Minimum Contrast

Description

Fits the Poisson Cluster point process to a point pattern dataset by the Method of Minimum Contrast.

Usage

pc.estK(Kobs, r, sigma2 = NULL, rho = NULL)
Kclust(r, sigma2, rho)

Arguments

Details

The algorithm fits the Poisson cluster point process to a point pattern, by finding the parameters of the Poisson cluster model which give the closest match between the theoretical K function of the Poisson cluster process and the observed K function. For a more detailed explanation of the Method of Minimum Contrast, see mincontrast in spatstat or Diggle (2003: 86).

The Poisson cluster processes are defined by the following postulates (Diggle 2003):

This implementation asumes that the probability distribution p[s] of offspring per parent is a Poisson distribution and that the position of each offspring relative to its parent follows a radially symetric Gaussian distribution with pdf

h(x, y) = [1/(2*pi*sigma^2)]* exp[-(x^2+y^2)/(2*sigma^2)]

The theoretical K-function of this Poisson cluster process is (Diggle, 2003):

pi*r^2 + [1- exp(-r^2/4*sigma^2)]/rho

The command Kclust computes the theoretical K-function of this Poisson cluster process and can be used to find some initial estimates of rho and sigma^2. In any case, the optimization usually finds the correct parameters even without starting values for these parameters.

This Poisson cluster process can be simulated with sim.poissonc.

Value

Note

The exponents p and q of the contrast criterion (see mincontrast) are fixed respectively to p = 2 and q = 1/4. The rmin and rmax limits of integration of the contrast criterion are set up by the sequence of values of r and Kobs passed to pc.estK.

Author(s)

Marcelino de la Cruz Rot, inspired by some code of Philip M. Dixon

References

Diggle, P. J. 2003. Statistical analysis of spatial point patterns. Arnold, London.

See Also

ipc.estK for fitting the inhomogeneous Poisson cluster process; some functions in spatstat ( matclust.estK and lgcp.estK) fit other appropriate processes for clustered patterns; mincontrast performs a more general implementation of the method of mimimum contrast.

Examples

data(gypsophylous)


# set the number of simulations (nsim=199 or larger for real analyses)
nsim<- 19

## Estimate K function ("Kobs").

gyps.env <- envelope(gypsophylous, Kest, correction="iso", nsim=nsim)

plot(gyps.env, sqrt(./pi)-r~r, legend=FALSE)

## Fit Poisson Cluster Process. The limits of integration 
## rmin and rmax are setup to 0 and 60, respectively. 

cosa.pc <- pc.estK(Kobs = gyps.env$obs[gyps.env$r<=60],
		           r = gyps.env$r[gyps.env$r<=60])

## Add fitted Kclust function to the plot.

lines(gyps.env$r,sqrt(Kclust(gyps.env$r, cosa.pc$sigma2,cosa.pc$rho)/pi)-gyps.env$r,
       lty=2, lwd=3, col="purple")

## A kind of pointwise test of the gypsophylous pattern been a realisation
## of the fitted model, simulating with sim.poissonc and using function J (Jest).

gyps.env.sim <- envelope(gypsophylous, Jest, nsim=nsim,
                    simulate=expression(sim.poissonc(gypsophylous,
		    sigma=sqrt(cosa.pc$sigma2), rho=cosa.pc$rho)))

plot(gyps.env.sim,  main="",legendpos="bottomleft")

Alive and dead oak trees

Description

Locations of alive and dead oak trees (Quercus robur) in a secondary wood in Urkiola Natural Park (Basque country, north of Spain). This is part of a more extensive dataset collected and analysed by Laskurain (2008). The coordinates of the trees are given in meters.

Usage

data(quercusvm)

Format

An object of class "ppp" representing the point pattern of tree locations. Entries include

x

Cartesian x-coordinate of tree.

y

Cartesian y-coordinate of tree.

marks

factor with two levels indicating the status of each tree (1 = "alive", 0 = "dead").

See ppp for details of the format of a ppp object.

References

Laskurain, N. A. (2008) Dinámica espacio-temporal de un bosque secundario en el Parque Natural de Urkiola (Bizkaia). Tesis Doctoral. Universidad del País Vasco /Euskal Herriko Unibertsitatea.

Simulate Inhomogeneous Poisson Cluster Process

Description

Generate a random point pattern, a simulated realisation of the Inhomogeneous Poisson Cluster Process.

Usage

rIPCP(x, lambda = NULL, type = 1, lmax = NULL, win = owin(c(0, 1), c(0, 1)), ...)

Arguments

Details

This function simulates the Inhomogeneous Poisson Cluster process from an object of class 'ecespa.minconfit', resulting from fitting an IPCP to some 'original' point pattern using the function ipc.estK. Following the approach of Waagepetersen (2007), the simulation involves a first step in which an homogeneous aggregated pattern is simulated (from the fitted parameters of the 'ecespa.minconfit' object, using function rThomas of spatstat) and a second one in which the homogeneous pattern is thinned with a spatially varying thinning probability f (s) proportional to the spatially varying intensity, i.e. f (s) = lambda(s) / max[lambda(s)]. To obtain a 'final' density similar to that of the original point pattern, a "prethinning" must be performed. There are two alternatives. If the argument 'type' is set equal to '1', the expected number of points per cluster (mu parameter of rThomas is thinned as mu <- mu.0 / mean[f(s)], where mu.0 is the mean number of points per cluster of the original pattern. This alternative produces point patterns most similar to the 'original'. If the argument 'type' is set equal to '2', the fitted intensity of the Poisson process of cluster centres (kappa parameter of rThomas, i.e. the intensity of 'parent' points) is thinned as kappa <- kappa / mean[f(s)]. This alternative produces patterns more uniform than the 'original' and it is provided only for experimental purposes.

Value

A point pattern, with the format of the ppp objects of spatstat.

Author(s)

Marcelino de la Cruz Rot

References

Waagepetersen, R. P. 2007. An estimating function approach to inference for inhomogeneous Neymann-Scott processes. Biometrics 63: 252-258. doi:10.1111/j.1541-0420.2006.00667.x.

See Also

sim.poissonc to simulate homogeneous PCP; rNeymanScott and rThomas in spatstat are the basis of this function

Examples

  
    data(gypsophylous)
  
    plot(gypsophylous) 
    
    ## It 'seems' that the pattern is clustered, so 
    ## fit a Poisson Cluster Process. The limits of integration 
    ## rmin and rmax are setup to 0 and 60, respectively.
    
   cosa.pc2 <- ipc.estK(gypsophylous, r = seq(0, 60, by=0.2))

   ## Create one instance of the fitted PCP:

   pointp <- rIPCP( cosa.pc2)
   
   plot(pointp)
   
   
   
    #####################
    ## Inhomogeneous example

    data(urkiola)

    # get univariate pp
    I.ppp <- split.ppp(urkiola)$birch

    plot(I.ppp)

    #estimate inhomogeneous intensity function
    I.lam <- predict (ppm(I.ppp, ~polynom(x,y,2)), type="trend", ngrid=200)

    # It seems that there is short scale clustering; lets fit an IPCP: 

    I.ki <- ipc.estK(mippp=I.ppp, lambda=I.lam, correction="trans")

    ## Create one instance of the fitted PCP:

    pointpi <- rIPCP( I.ki)
   
    plot(pointpi)
 

Cohorts of Helianthemum squamatum seedlings

Description

Marked point patterns of two consecutive cohorts of seedlings of H. squamatum growing in a gypsophylous plant community in Central Spain. The datasets contains the locations of the seedlings marked with their heigth. Both the coordinates and the heigth of the seedlings are given in cm.

Usage

data(seedlings1)
data(seedlings2)

Format

seedlings1 and seedlings2 are objects of class "ppp" representing the point pattern of seedling locations marked by their heights. See ppp.object for details of the format.

Source

Romao, R.L. 2003. Estructura espacial de comunidades de gipsófitos: interacciones bióticas y constricciones abióticas. Tesis Doctoral. Universidad Politécnica de Madrid.

References

De la Cruz, M. 2006. Introducción al análisis de datos mapeados o algunas de las (muchas) cosas que puedo hacer si tengo coordenadas. Ecosistemas 15 (3): 19-39.

Escudero, A., Romao, R.L., De la Cruz, M. & Maestre, F. 2005. Spatial pattern and neighbour effects on Helianthemum squamatum seedlings in a Mediterranean gypsum community. Journal of Vegetation Science, 16: 383-390. doi:10.1111/j.1654-1103.2005.tb02377.x.

Examples

## Not run: 

 data(seedlings1)
 
 plot(seedlings1)


## End(Not run)

Simulate Poisson Cluster Process

Description

Generate a random point pattern, a simulated realisation of the Poisson Cluster Process

Usage

sim.poissonc(x.ppp, rho, sigma)

Arguments

Details

The Poisson cluster processes are defined by the following postulates (Diggle 2003):

This implementation asumes that the probability distribution p[s] of offspring per parent is a Poisson distribution and that the position of each offspring relative to its parent follows a radially symetric Gaussian distribution with pdf

h(x, y) = [1/(2*pi*sigma^2)]* exp[-(x^2+y^2)/(2*sigma^2)]

Value

The simulated point pattern (an object of class "ppp").

Warning

This implementation simulates only point patterns within rectangular windows. Use ipc.estK to fit and rIPCP (or the spatstat functions) to simulate point patterns within irregular windows.

Note

This function can use the results of pc.estK to simulate point patterns from a fitted model. Be careful as the paramted returned by pc.estK is sigma^2 while sim.poissonc takes its square root, i.e. sigma.

Author(s)

Marcelino de la Cruz Rot

References

Diggle, P.J. 2003. Statistical analysis of spatial point patterns. Arnold, London.

See Also

rIPCP to simulate inhomogeneous PCP; rNeymanScott and rThomas in spatstat

Examples

data(gypsophylous)

# set the number of simulations (nsim=199 or larger for real analyses)
nsim<- 39

## Estimate K function ("Kobs").
gyps.env <- envelope(gypsophylous, Kest, correction="iso", nsim=nsim)

plot(gyps.env, sqrt(./pi)-r~r, legend=FALSE)

## Fit Poisson Cluster Process. The limits of integration 
## rmin and rmax are setup to 0 and 60, respectively. 
cosa.pc <- pc.estK(Kobs = gyps.env$obs[gyps.env$r<=60],
		           r = gyps.env$r[gyps.env$r<=60])

## Add fitted Kclust function to the plot.
lines(gyps.env$r,sqrt(Kclust(gyps.env$r, cosa.pc$sigma2,cosa.pc$rho)/pi)-gyps.env$r,
       lty=2, lwd=3, col="purple")

## A kind of pointwise test of the pattern gypsophilous been a realisation
## of the fitted model, simulating with sim.poissonc and using function J (Jest).

gyps.env.sim <- envelope(gypsophylous, Jest,  nsim=nsim,
                    simulate=expression(sim.poissonc(gypsophylous,
		    sigma=sqrt(cosa.pc$sigma2), rho=cosa.pc$rho)))

plot(gyps.env.sim,  main="")

Tree Species in a Swamp Forest

Description

Locations and botanical classification of trees in a plot in the Savannah River. Locations are given in metres, rounded to the nearest 0.1 metre. The data come from a 1-ha (200 m x 50 m) plot in the Savannah River Site, South Carolina, USA. The 734 mapped stems included 156 Carolina ash (Fraxinus caroliniana), 215 Water tupelo (Nyssa aquatica), 205 Swamp tupelo (Nyssa sylvatica), 98 Bald cypress (Taxodium distichum) and 60 stems of 8 additional species. Although the plots were set up by Bill Good and their spatial patterns described in Good and Whipple(1982), the plots have been maintained and resampled by Rebecca Sharitz and her colleagues of the Savannah River Ecology Laboratory. There are slightly different versions of the Good plot data. Every time the plots are resampled, some errors are corrected. This is mostly a concern for the biologists. The different versions are very similar; they are all very good examples of a marked spatial point pattern.

Usage

data(swamp)

Format

A data frame with 734 observations on the following 3 variables.

x

Cartesian x-coordinate of tree

y

Cartesian y-coordinate of tree

sp

a factor with levels indicating the species of each tree:

References

Dixon, P.M. 2002. Nearest-neighbor contingency table analysis of spatial segregation for several species. Ecoscience, 9 (2): 142-151. doi:10.1080/11956860.2002.11682700.

Good, , B. J. & Whipple, S.A. 1982. Tree spatial patterns: South Carolina bottomland and swamp forest. Bulletin of the Torrey Botanical Club, 109: 529-536.

Jones et al. 1994. Tree population dynamics in seven South Carolina mixed-species forests. Bulletin of the Torrey Botanical Club, 121:360-368.

Examples

data(swamp)
plot(swamp$x,swamp$y, col=as.numeric(swamp$sp),pch=19,
		 xlab="",ylab="",main="Swamp forest")

Syrjala's test for the difference between the spatial distributions of two populations

Description

Computes a two-sample Cramer-von Mises (and Kolmogorov-Smirnov) type test for a difference between the spatial distributions of two populations. It is designed to be sensitive to differences in the way the populations are distributed across the study area but insensitive to differences in abundance between the two populations.

Usage

syrjala0(coords, var1, var2, nsim, R=FALSE)
syrjala(coords = NULL, var1 = NULL, var2 = NULL, nperm = 999)
syrjala.test(ppp1, ppp2, nsim = 999)
## S3 method for class 'syrjala.test'
plot(x, coline=1, ...)
## S3 method for class 'ecespa.syrjala'
plot(x, ...)

Arguments

Details

The null hypothesis of Syrjala's test is that across the study area, the normalized distributions of the two populations are the same (Syrjala, 1996). Population density data are collected at K sampling locations on two populations. Let (xk, yk) denote the coordinates of the kth sampling location ( k= 1,...,K ); let d.i(xk, yk) denote the sample density at the Kth sampling location of the ith population. To construct a test that is independent of the population sizes, the observed density data is first normalized:

gamma.i(xk, yk) = di(xk, yk)/Di,

where Di is the sum of d.i(xk, yk) observations across the K sampling locations. The value of the cumulative distribution function at the location (xk, yk) for the ith population, denoted GAMMA.i(xk, yk), is the sum of all normalized density observations, gamma.i(xk, yk), whose location (x, y) is such that x <= xk and y <= yk. The statistic proposed by Syrjala to test the null hypothesis is the square of the difference between the cumulative distribution functions of the two populations, summed over all sampling locations, that is

psi = sum{GAMMA.1(xk, yk) - GAMMA.2(xk, yk)}^2.

As psi is not invariant with respect to the 'corner' of the rectangle enclosing the study area that is chosen as the origin of the coordinate sytem, psi is computed four times, one with each corner as the origin, and the average psi is employed as the test statistic. The level of significance of the observed psi is determined from its position in the ordered set of test statistic values from all 2^K pairwise permutations (that is approximated from a large number of randomly selected permutations).

Value

Functions syrjala or syrjala0 (with the argument R=FALSE) return an object of class 'syrjala.test'. Functions syrjala.test or syrjala0 (with the argument R=TRUE) return an object of class 'ecespa.syrjala'. In Both cases, the result is a list with the following elements:

Both S3 plot methods plot an histogram with the distribution of the simulated psi's statistics and draws the observed psi as a vertical line.

Warning

The test requires both populations being sampled in exactly the same sampling locations. Althoug this implementation employs ppp's as the supporting data format, this kind of data are not spatial point patterns. They cannot be analysed with the usual tools employed for marked point patterns.

Note

syrjala or syrjala0 (with the argument R=FALSE) implement a Fortran version of Syrjala's test. They run considerably faster than the "whole-R" implementation of syrjala.test or syrjala0 (with the argument R=TRUE). This last implementation is supplied for illustrative purposes and to maintain compability with previous versions of package ecespa. One can use function haz.ppp to easily build the ppp objects from a data.frame with only three columns (x-coordinate, y-coordinate and abundance).

This function has been employed to compute Syrjala's test in Rey-Benayas et al. (2008).

Author(s)

Jose M. Blanco-Moreno jmblanco@ub.edu for the Fortran implementation of Syrjala's original QBasic function, Marcelino de la Cruz Rot for the R version, the wrapping functions and the documentation

References

Rey-Benayas, J.M., de la Montaña, E., Pérez-Camacho, L., de la Cruz, M., Moreno, D., Parejo, J.L. and Suárez-Seoane, S. 2010. Inter-annual dynamics and spatial congruence of a nocturnal bird assemblage inhabiting a Mediterranean agricultural mosaic. Ardeola,57(2): 303-320.

Syrjala, S. E. 1996. A statistical test for a difference between the spatial distributions of two populations. Ecology 77: 75-80. doi:10.2307/2265656.

Examples

   
   data(syr1); data(syr2); data(syr3)
   
   coords <- data.frame(x=syr1$x, y=syr1$y); var1<- syr1$marks; var2<- syr2$marks
   
   (stest<-syrjala(coords, var1, var2, 9999))
   plot(stest)
   
   
   ## Not run: 
   
   plot(syrjala.test(syr1, syr2, nsim=999)) 
   
   plot(syrjala.test(syr1, syr3, nsim=999)) 
   
   syrjala0(coords, var1, var2, 9999)
   
   syrjala0(coords, var1, var2, 999, R=TRUE)
   
   
   coords <- expand.grid(x=1:10,y=1:10)
   var1 <- runif(100)
   var2 <- runif(100)
   syrjala(coords, var1, var2, 9999)
   
   

## End(Not run)

Syrjala test data

Description

Artificial data to exemplify Syrjala's test.

Usage

data(syr1)
data(syr2)
data(syr3)

Format

syr1, syr2 and syr3 are marked point patterns of class "ppp" representing the coordinates of some sampling locations, 'marked' by the value of some parameters (e.g. density of individuals) measured on them. See ppp.object for details of the format. On the other hand, one can use function haz.ppp to easily build ppp objects appropriate for use with syrjala.test.

Examples

data(syr1)