`superthin`

takes a space-time point pattern and conditional intensity model and calculates a set of super-thinned residuals for further analysis.

1 |

`X` |
A “ |

`cifunction` |
A function returning the value of the conditional intensity at all points in |

`theta` |
Optional: A vector of parameters to be passed to |

`k` |
The super-thinning rate. |

`lambda` |
Optional: A vector of conditional intensities at each point in |

Super-thinned residuals (Clements et. al. (2012)) is a type of transformation based residuals for space-time point processes based on both thinned residuals (see Schoenberg (2003)) and superposed residuals (see Bremaud (1981)). The residuals consist of a set of points that should be homogeneous Poisson, with rate `k`

, if the model for the conditional intensity is correct. Any patterns or inter-point interaction in the residuals indicates a lack of fit of the model. To test for homogeneity, a commonly used tool is Ripley's K-function, a version of which can be found in the `spatstat`

package.

Super-thinned residuals are found as follows:

1. The super-thinning rate `k`

is specified. This rate determines the amount of thinning and superposition conducted, and also determines the final rate of the super-thinned residual point process.

2. All observed points in `X`

where *lambda_hat <* `k`

are automatically kept.

3. All points in `X`

where *lambda_hat >=* `k`

are kept with probability `k`

*/lambda_hat*.

4. In all space-time locations where *λ <* `k`

, points are simulated with rate `k`

*- lambda_hat*.

The result should be a homogeneous Poisson process with rate `k`

if the model is correct.

The conditional intensity function, `cifunction`

, should take `X`

as the first argument, and an optional `theta`

as the second argument, and return a vector of conditional intensity estimates with length equal to the number of points in `X`

, i.e. the length of `X$x`

. `cifunction`

is required, while `lambda`

is optional. `lambda`

eliminates the need for `superthin`

to calculate the conditional intensity at each observed point in `X`

.

If `k`

is not specified, the default is the mean of *lambda_hat* estimated by the total number of points divided by the volume of the space-time window.

Outputs an object of class “`superthin`

”, which is a list of

`X ` |
An object of class “ |

`k ` |
The super-thinning rate. |

`residuals ` |
A data frame consisting of the x, y, and t coordinates of the super-thinned residuals. |

`super ` |
A data frame consisting of the x, y, and t coordinates of the superposed points. |

`keep1 ` |
A data frame consisting of the x, y, and t coordinates of the automatically retained points. |

`keep2 ` |
A data frame consisting of the x, y, and t coordinates of the points remaining after the thinning has taken place. |

`deleted ` |
A data frame consisting of the x, y, and t coordinates of the points removed during the thinning process. |

Robert Clements

Bremaud, P. *Point Processes and Queues: Martingale Dynamics.* SpringerVerlag, New York, 1981.

Clements, R.A., Schoenberg, F.P., and Veen, A. (2012) Evaluation of space-time point process models using super-thinning. *Environmetrics*, to appear.

Schoenberg, F.P. (2003) Multi-dimensional residuals analysis of point process models for earthquake occurrences. *Journal of the American Statistical Association,* **98**, 789–795.

Clements, R.A., Schoenberg, F.P., and Schorlemmer, D. (2011) Residual analysis methods for space-time point processes with applications to earthquake forecast models in California. *Annals of Applied Statistics,* **5**, Number 4, 2549–2571.

`stpp`

, `thinresid`

, `supresid`

, `spatstat`

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 | ```
#===> load simulated data <===#
data(simdata)
X <- stpp(simdata$x, simdata$y, simdata$t)
#===> define conditional intensity function <===#
ci1 <- function(X, theta){theta[1]*exp(-theta[2]*X$x -
theta[3]*X$y - theta[4]*X$t)} #correct model
stresiduals1 <- superthin(X, ci1, theta = c(3000, 2, 2, 2), k = 250)
stresiduals2 <- superthin(X, ci1, theta = c(2500, 5, 5, 10), k = 250)
#===> plot results <===#
par(mfrow = c(1,2))
plot(stresiduals1)
plot(stresiduals2)
summary(stresiduals1)
summary(stresiduals2)
``` |

Questions? Problems? Suggestions? Tweet to @rdrrHQ or email at ian@mutexlabs.com.

All documentation is copyright its authors; we didn't write any of that.