rlgcp: Generate log-Gaussian Cox point patterns

In stpp-GitHub/stpp: Space-Time Point Pattern Simulation, Visualisation and Analysis

Generate log-Gaussian Cox point patterns

Description

Generate one (or several) realisation(s) of the log-Gaussian cox process in a region S\times T.

Usage

 rlgcp(s.region, t.region, replace=TRUE, npoints=NULL, nsim=1, 
 nx=100, ny=100, nt=100,separable=TRUE,model="exponential",
 param=c(1,1,1,1,1,1), scale=c(1,1),var.grf=1,mean.grf=0,
 lmax=NULL,discrete.time=FALSE,exact=FALSE,anisotropy=FALSE,ani.pars=NULL)

Arguments

s.region	Two-column matrix specifying polygonal region containing all data locations. If s.region is missing, the unit square is considered.
t.region	Vector containing the minimum and maximum values of the time interval. If t.region is missing, the interval [0,1] is considered.
npoints	Number of points to simulate. If NULL, the number of points is from a Poisson distribution with mean the double integral of \lambda(s,t) over s.region and t.region.
nsim	number of simulations to generate. Default is 1.
separable	Logical. If TRUE, the covariance function of the Gaussian random field is separable.
model	Vector of length 1 or 2 specifying the model(s) of covariance of the Gaussian random field. If separable=TRUE and model is of length 2, then the elements of model define the spatial and temporal covariances respectively. If separable=TRUE and model is of length 1, then the spatial and temporal covariances belongs to the same class of covariances, among "matern", "exponential", "stable", "cauchy" and "wave" (see Details). If separable=FALSE, model must be of length 1 and is either "gneiting" or "cesare" (see Details).
param	(\alpha_1,\alpha_2,\alpha_3,\alpha_4,\alpha_5,\alpha_6). Vector of parameters of the covariance function (see Details).
scale	Vector of length 2 defining the spatial and temporal scale.
var.grf	Variance of the Gaussian random field.
mean.grf	Mean of the Gaussian random field.
replace	Logical allowing times repeat.
nx, ny, nt	Define the size of the 3-D grid on which the intensity is evaluated.
lmax	Upper bound for the value of \lambda(x,y,t).
discrete.time	If TRUE, times belong to N, otherwise belong to R^+.
exact	logical allowing exact simulation of Gaussian random fields (see manual for details).
anisotropy	If TRUE, simulate an anisotropic point pattern. Currently only implemented for separable covariance functions.
ani.pars	Vector of length 2, the anisotropy angle and the anisotropy ratio, respectively (see details).

Details

We implemented stationary, isotropic spatio-temporal covariance functions.

Separable covariance functions

c(h,t) = c_s(\| h \|) \, c_t(|t|) , h \in S, t \in T

The purely spatial and purely temporal covariance functions can be:

• Exponential: c(r) = \exp(-r),

• Stable: c(r) = \exp(-r^\alpha), \alpha \in [0,2],

• Cauchy: c(r) = (1+r^2)^{-\alpha}, \alpha >0,

• Wave: c(r) = \frac{\sin(r)}{r} if r>0, c(0)=1,

• Matern: c(r) = \frac{(\alpha r)^\nu}{2^{\nu-1}\Gamma(\nu)} {\cal K}_{\nu}(\alpha r), \nu > 0 and \alpha > 0.

  {\cal K}_{\nu} is the modified Bessel function of second kind:

  {\cal K}_{\nu}(x) = \frac{\pi}{2} \frac{I_{-\nu}(x) - I_{\nu}(x)}{\sin(\pi \nu)},

  with I_{\nu}(x) = \left( \frac{x}{2} \right)^{\nu} \sum_{k=0}^\infty \frac{1}{k! \Gamma(\nu+k+1)} \left( \frac{x}{2} \right)^{2k}.

  The parameters \alpha_1 and \alpha_2 correspond to the parameters of the spatial and temporal covariance respectively. For the Matern model, the parameters \alpha_1, \alpha_3 and \alpha_2, \alpha_4 correspond to the parameters \nu, \alpha of the spatial and temporal covariance.

Non-separable covariance functions

The spatio-temporal covariance function can be:

• Gneiting: c(h,t) = \psi (t^2/\beta_2)^{-\alpha_6} \phi \left(\frac{h^2/ \beta_1}{\psi (t^2/\beta_2)} \right), \beta_1, \beta_2 >0 are spatial and temporal scales respectively,

  • If \alpha_2=1, \phi(r) is the Stable model.

  • if \alpha_2=2, \phi(r) is the Cauchy model.

  • If \alpha_2=3, \phi(r) is the Matern model.

  • If \alpha_5=1, \psi(r) = (r^{\alpha_3}+ 1)^{\alpha_4},

  • If \alpha_5=2, \psi(r) = (\alpha_4^{-1} r^{\alpha_3} +1)/(r^{\alpha_3}+1),

  • If \alpha_5=3, \psi(r) = \log(r^{\alpha_3} + \alpha_4)/\log {\alpha_4},

  The parameter \alpha_1 is the respective parameter for the model of \phi(\cdot), \alpha_3 \in (0,1], \alpha_4 \in (0,1] and \alpha_6 \geq 2.

• cesare: c(h,t) = \left( 1 + (h/\beta_1)^{\alpha_1} + (t/\beta_2)^{\alpha_2} \right)^{-\alpha_3}, \beta_1, \beta_2 >0, \alpha_1, \alpha_2 \in [1,2] and \alpha_3 \geq 3/2.

We also implemented anisotropic Log-Gaussian Cox processes. We considered geometric spatial anisotropy (see Moller and Toftaker, 2014). In this case the covariance function is elliptical and anisotropy is characterized by two parameters: the anisotropy angle 0 \leq \theta < \pi and the anisotropy ratio 0 < \delta \leq 1 of the minor axis 2 \omega \delta and the major axis 2 \omega.

C(h,t)=C_0\left( \sqrt{h \Sigma^{-1} h'},t \right), \ h \in R^2.

Value

A list containing:

xyt	Matrix (or list of matrices if nsim>1) containing the points (x,y,t) of the simulated point pattern. xyt (or any element of the list if nsim>1) is an object of the class stpp.
s.region, t.region	parameters passed in argument.
Lambda	nx \times ny \times nt array (or list of array if nsim>1) of the intensity.

Author(s)

Edith Gabriel <[e-mail]>, Peter J Diggle.

References

Chan, G. and Wood A. (1997). An algorithm for simulating stationary Gaussian random fields. Applied Statistics, Algorithm Section, 46, 171–181.

Chan, G. and Wood A. (1999). Simulation of stationary Gaussian vector fields. Statistics and Computing, 9, 265–268.

Gneiting T. (2002). Nonseparable, stationary covariance functions for space-time data. Journal of the American Statistical Association, 97, 590–600.

Moller J. and Toftaker H. (2014). Geometric anisotropic spatial point pattern analysis and Cox processes. Scandinavian Journal of Statistics, 41, 414–435.

See Also

plot.stpp, animation and stan for plotting space-time point patterns.

Examples

# non separable covariance function: 
lgcp1 <- rlgcp(npoints=200, nx=50, ny=50, nt=50, separable=FALSE, 
model="gneiting", param=c(1,1,1,1,1,2), var.grf=1, mean.grf=0)
N <- lgcp1$Lambda[,,1];for(j in 2:(dim(lgcp1$Lambda)[3])){N <-
N+lgcp1$Lambda[,,j]}
image(N,col=grey((1000:1)/1000));box()
animation(lgcp1$xyt, cex=0.8, runtime=10, add=TRUE, prevalent="orange")

# separable covariance function: 
lgcp2 <- rlgcp(npoints=200, nx=50, ny=50, nt=50, separable=TRUE, 
model="exponential", param=c(1,1,1,1,1,2), var.grf=2, mean.grf=-0.5*2)
N <- lgcp2$Lambda[,,1];for(j in 2:(dim(lgcp2$Lambda)[3])){N <-
N+lgcp2$Lambda[,,j]}
image(N,col=grey((1000:1)/1000));box()
animation(lgcp2$xyt, cex=0.8, pch=20, runtime=10, add=TRUE,
prevalent="orange")

# anisotropic

sigma2=0.5
simlgcp <- rlgcp(npoints=500,nx=250, ny=200, nt=50,separable=TRUE,
s.region=matrix(c(0,2,2,0,0,0,0.5,0.5),byrow=FALSE,ncol=2), model="exponential", 
param=c(1,1,1,1,1,2), var.grf=sigma2, mean.grf=-0.5*sigma2,anisotropy = TRUE,
ani.pars = c(pi/4,0.1))

N <- simlgcp$Lambda[,,1];for(j in 2:dim(simlgcp$Lambda)[3]){N <- N+simlgcp$Lambda[,,j]}
image(x=simlgcp$grid[[1]]$x,y=simlgcp$grid[[1]]$y,z=N,col=grey((1000:1)/1000));box()
points(simlgcp$xyt[,1:2],pch=19,cex=0.25,col=2)

stpp-GitHub/stpp documentation built on April 19, 2024, 4:21 p.m.