STIKhat {stpp}

Estimation of the Space-Time Inhomogeneous K-function

Description

Compute an estimate of the Space-Time Inhomogeneous K-function.

Usage

STIKhat(xyt, s.region, t.region, dist, times, lambda, 
correction="isotropic", infectious=FALSE) 

Arguments

xyt

Coordinates and times (x,y,t) of the point pattern.

s.region

Two-column matrix specifying polygonal region containing all data locations. If s.region is missing, the bounding box of xyt[,1:2] is considered.

t.region

Vector containing the minimum and maximum values of the time interval. If t.region is missing, the range of xyt[,3] is considered.

dist

Vector of distances u at which K(u,v) is computed. If missing, the maximum of dist is given by \min(S_x,S_y)/4, where S_x and S_y represent the maximum width and height of the bounding box of s.region.

times

Vector of times v at which K(u,v) is computed. If missing, the maximum of times is given by (T_{\max} - T_{\min})/4, where T_{\min} and T_{\max} are the minimum and maximum of the time interval T.

lambda

Vector of values of the space-time intensity function evaluated at the points (x,y,t) in S\times T. If lambda is missing, the estimate of the space-time K-function is computed as for the homogeneous case (Diggle et al., 1995), i.e. considering n/|S \times T| as an estimate of the space-time intensity.

correction

A character vector specifying the edge correction(s) to be applied among "isotropic", "border", "modified.border", "translate" and "none" (see Details). The default is "isotropic".

infectious

Logical value. If TRUE, only future events are considered and the isotropic edge correction method is used. See Details.

Details

Gabriel (2014) proposes the following unbiased estimator for the STIK-function, based on data giving the locations of events x_i: i=1,\ldots,n on a spatio-temporal region S\times T, where S is an arbitrary polygon and T is a time interval:

\widehat{K}(u,v)=\sum_{i=1}^{n}\sum_{j\neq i}\frac{1}{w_{ij}}\frac{1}{\lambda(x_i)\lambda(x_j)}\mathbf{1}_{\lbrace \|s_i - s_j\| \leq u \ ; \ |t_i - t_j| \leq v \rbrace},

where \lambda(x_i) is the intensity at x_i = (s_i,t_i) and w_{ij} is an edge correction factor to deal with spatial-temporal edge effects. The edge correction methods implemented are:

isotropic: w_{ij} = |S \times T| w_{ij}^{(t)} w_{ij}^{(s)}, where the temporal edge correction factor w_{ij}^{(t)} = 1 if both ends of the interval of length 2 |t_i - t_j| centred at t_i lie within T and w_{ij}^{(t)}=1/2 otherwise and w_{ij}^{(s)} is the proportion of the circumference of a circle centred at the location s_i with radius \|s_i -s_j\| lying in S (also called Ripley's edge correction factor).

border: w_{ij}=\frac{\sum_{j=1}^{n}\mathbf{1}\lbrace d(s_j,S)>u \ ; \ d(t_j,T) >v\rbrace/\lambda(x_j)}{\mathbf{1}_{\lbrace d(s_i,S) > u \ ; \ d(t_i,T) >v \rbrace}}, where d(s_i,S) denotes the distance between s_i and the boundary of S and d(t_i,T) the distance between t_i and the boundary of T.

modified.border: w_{ij} = \frac{|S_{\ominus u}|\times|T_{\ominus v}|}{\mathbf{1}_{\lbrace d(s_i,S) > u \ ; \ d(t_i,T) >v \rbrace}}, where S_{\ominus u} and T_{\ominus v} are the eroded spatial and temporal region respectively, obtained by trimming off a margin of width u and v from the border of the original region.

translate: w_{ij} =|S \cap S_{s_i-s_j}| \times |T \cap T_{t_i-t_j}|, where S_{s_i-s_j} and T_{t_i-t_j} are the translated spatial and temporal regions.

none: No edge correction is performed and w_{ij}=|S \times T|.

If parameter infectious = TRUE, ony future events are considered and the estimator is, using an isotropic edge correction factor (Gabriel and Diggle, 2009):

\widehat{K}(u,v)=\frac{1}{|S\times T|}\frac{n}{n_v}\sum_{i=1}^{n_v}\sum_{j=1; j > i}^{n_v} \frac{1}{w_{ij}} \frac{1}{\lambda(x_i) \lambda(x_j)}\mathbf{1}_{\left\lbrace u_{ij} \leq u\right\rbrace}\mathbf{1}_{\left\lbrace t_j - t_i \leq v \right\rbrace}.

In this equation, the points x_i=(s_i, t_i) are ordered so that t_i < t_{i+1}, with ties due to round-off error broken by randomly unrounding if necessary. To deal with temporal edge-effects, for each v, n_v denotes the number of events for which t_i \leq T_1 -v, with T=[T_0,T_1]. To deal with spatial edge-effects, we use Ripley's method.

If lambda is missing in argument, STIKhat computes an estimate of the space-time (homogeneous) K-function:

\widehat{K}(u,v)=\frac{|S\times T|}{n_v(n-1)} \sum_{i=1}^{n_v}\sum_{j=1;j>i}^{n_v}\frac{1}{w_{ij}}\mathbf{1}_{\lbrace u_{ij}\leq u \rbrace}\mathbf{1}_{\lbrace t_j - t_i \leq v \rbrace}

Value

A list containing:

Khat

ndist x ntimes matrix containing values of \hat{K}_{ST}(u,v)..

Ktheo

ndist x ntimes matrix containing theoretical values for a Poisson process; \pi u^2 v for K and 2 \pi u^2 v for K^*.

dist, times, infectious

Parameters passed in argument.

correction

The name(s) of the edge correction method(s) passed in argument.

Author(s)

Edith Gabriel <[e-mail]>

References

Baddeley A., Moller J. and Waagepetersen R. (2000). Non- and semi-parametric estimation of interaction in inhomogeneous point patterns. Statistica Neerlandica, 54, 329–350.

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

Diggle P. , Chedwynd A., Haggkvist R. and Morris S. (1995). Second-order analysis of space-time clustering. Statistical Methods in Medical Research, 4, 124–136.

Gabriel E., Diggle P. (2009). Second-order analysis of inhomogeneous spatio-temporal point process data. Statistica Neerlandica, 63, 43–51.

Gabriel E., Rowlingson B., Diggle P. (2013). stpp: an R package for plotting, simulating and analyzing Spatio-Temporal Point Patterns. Journal of Statistical Software, 53(2), 1–29.

Gabriel E. (2014). Estimating second-order characteristics of inhomogeneous spatio-temporal point processes: influence of edge correction methods and intensity estimates. Methodology and computing in Applied Probabillity, 16(2), 411–431.

Examples

# First example

data(fmd)
data(northcumbria)
FMD<-as.3dpoints(fmd[,1]/1000,fmd[,2]/1000,fmd[,3])
Northcumbria=northcumbria/1000

# estimation of the temporal intensity
Mt<-density(FMD[,3],n=1000)
mut<-Mt$y[findInterval(FMD[,3],Mt$x)]*dim(FMD)[1]

# estimation of the spatial intensity
h<-mse2d(as.points(FMD[,1:2]), Northcumbria, nsmse=50, range=4)
h<-h$h[which.min(h$mse)]
Ms<-kernel2d(as.points(FMD[,1:2]), Northcumbria, h, nx=5000, ny=5000)
atx<-findInterval(x=FMD[,1],vec=Ms$x)
aty<-findInterval(x=FMD[,2],vec=Ms$y)
mhat<-NULL
for(i in 1:length(atx)) mhat<-c(mhat,Ms$z[atx[i],aty[i]])

# estimation of the STIK function
u <- seq(0,10,by=1)
v <- seq(0,15,by=1)
stik1 <- STIKhat(xyt=FMD, s.region=northcumbria/1000,t.region=c(1,200), 
lambda=mhat*mut/dim(FMD)[1], dist=u, times=v, infectious=TRUE)

# plotting the estimation
plotK(stik1)
plotK(stik1,type="persp",theta=-65,phi=35)
 
# Second example

xyt=rpp(lambda=200)
stik2=STIKhat(xyt$xyt,dist=seq(0,0.16,by=0.02),
times=seq(0,0.16,by=0.02),correction=c("border","translate"))
plotK(stik2,type="contour",legend=TRUE,which="translate")