STIKhat: 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)

Value

A list containing:

Khat: ndist x ntimes matrix containing values of ${\hat{K}}_{S T} (u, v) .$ .
Ktheo: ndist x ntimes matrix containing theoretical values for a Poisson process; $π u^{2} v$ for $K$ and $2 π 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.

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.

Author

Edith Gabriel <edith.gabriel@inrae.fr>

Details

Gabriel (2014) proposes the following unbiased estimator for the STIK-function, based on data giving the locations of events $x_{i} : i = 1, \dots, n$ on a spatio-temporal region $S \times T$ , where $S$ is an arbitrary polygon and $T$ is a time interval: $\hat{K} (u, v) = \sum_{i = 1}^{n} \sum_{j \neq i} \frac{1}{w_{i j}} \frac{1}{λ (x_{i}) λ (x_{j})} 1_{{‖ s_{i} - s_{j} ‖ \leq u; | t_{i} - t_{j} | \leq v}},$ where $λ (x_{i})$ is the intensity at $x_{i} = (s_{i}, t_{i})$ and $w_{i j}$ is an edge correction factor to deal with spatial-temporal edge effects. The edge correction methods implemented are:

isotropic: $w_{i j} = | S \times T | w_{i j}^{(t)} w_{i j}^{(s)}$ , where the temporal edge correction factor $w_{i j}^{(t)} = 1$ if both ends of the interval of length $2 | t_{i} - t_{j} |$ centred at $t_{i}$ lie within $T$ and $w_{i j}^{(t)} = 1 / 2$ otherwise and $w_{i j}^{(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_{i j} = \frac{\sum_{j = 1}^{n} 1 {d (s_{j}, S) > u; d (t_{j}, T) > v} / λ (x_{j})}{1_{{d (s_{i}, S) > u; d (t_{i}, T) > v}}}$ , 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_{i j} = \frac{| S_{⊖ u} | \times | T_{⊖ v} |}{1_{{d (s_{i}, S) > u; d (t_{i}, T) > v}}}$ , where $S_{⊖ u}$ and $T_{⊖ 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_{i j} = | 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_{i j} = | 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): $\hat{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_{i j}} \frac{1}{λ (x_{i}) λ (x_{j})} 1_{{u_{i j} \leq u}} 1_{{t_{j} - t_{i} \leq v}} .$

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: $\hat{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_{i j}} 1_{{u_{i j} \leq u}} 1_{{t_{j} - t_{i} \leq v}}$

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

Run this code

# \donttest{
# 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")

Run the code above in your browser using DataLab

Join us for RADAR: AI Edition