surveillance (version 1.5-4)

twinstim_simulation: Simulation of a Self-Exciting Spatio-Temporal Point Process

Description

The function simEpidataCS simulates events of a self-exiciting spatio-temporal point process of the "twinstim" class. Simulation works via Ogata's modified thinning of the conditional intensity as described in Meyer et al. (2012).

The simulate method for objects of class "twinstim" simulates new epidemic data using the model and the parameter estimates of the fitted object.

Usage

simEpidataCS(endemic, epidemic, siaf, tiaf, qmatrix, rmarks,
    events, stgrid, tiles, beta0, beta, gamma, siafpars, tiafpars,
    t0 = stgrid$start[1], T = tail(stgrid$stop,1), nEvents = 1e5,
    control.siaf = list(F=list(), Deriv=list()),
    W = NULL, trace = 5, nCircle2Poly = 32, gmax = NULL, .allocate = 500,
    .skipChecks = FALSE, .onlyEvents = FALSE)
## S3 method for class 'twinstim':
simulate(object, nsim = 1, seed = NULL, data, tiles,
    rmarks = NULL, t0 = NULL, T = NULL, nEvents = 1e5,
    control.siaf = object$control.siaf,
    W = NULL, trace = FALSE, nCircle2Poly = NULL, gmax = NULL,
    .allocate = 500, simplify = TRUE, ...)

Arguments

endemic

see twinstim. Note that type-specific endemic intercepts are specified by beta0 here, not by the term (1|type).

epidemic

see twinstim. Marks appearing in this formula must be returned by the generating function rmarks.

siaf

see twinstim. In addition to what is required for fitting with twinstim, the siaf specification must also contain the element simulate, a function which d

tiaf

e.g. what is returned by the generating function tiaf.constant or tiaf.exponential. See also tw

qmatrix

see epidataCS. Note that this square matrix and its dimnames determine the number and names of the different event types. In the simplest case, there is only a single type of even

rmarks

function of single time (1st arg) and location (2nd arg) returning a one-row data.frame of marks (named according to the variables in epidemic) for an event at this point. This must include the columns eps.s

events

NULL or missing (default) in case of an empty prehistory, or a SpatialPointsDataFrame containing events of the prehistory (-Inf;t0] of the process (required

stgrid

see as.epidataCS. Simulation only works inside the spatial and temporal range of stgrid.

tiles

object inheriting from SpatialPolygons with row.names equal to the tile names in stgrid and having the same proj4string as events and W. This is necessary

beta0,beta,gamma,siafpars,tiafpars

these are the parameter subvectors of the twinstim. beta and gamma must be given in the same order as they appear in endemic and epidemic, respectively. beta0 is ei

t0

events having occured during (-Inf;t0] are regarded as part of the prehistory $H_0$ of the process. The time point t0 must be an element of stgrid$start. For simEpidataCS, by default, and als

T, nEvents

simulate a maximum of nEvents events up to time T, then stop. For simEpidataCS, by default, and also if T=NULL, T equals the last stop time in stgrid (it cannot be g

W

see as.epidataCS. Must have the same proj4string as events and tiles. If not specified (NULL), W is generated automatically

trace

logical (or integer) indicating if (or how often) the current simulation status should be cated. For the simulate.twinstim method, trace currently only applies to the first of the nsim simul

.allocate

number of rows (events) to initially allocate for the event history; defaults to 500. Each time the simulated epidemic exceeds the allocated space, the event data.frame will be enlarged by .allocate rows.

.skipChecks,.onlyEvents

these logical arguments are not meant to be set by the user. They are used by the simulate-method for twinstim objects.

object

an object of class "twinstim".

nsim

number of epidemics (i.e. spatio-temporal point patterns inheriting from class "epidataCS") to simulate. Defaults to 1 when the result is a simple object inheriting from class "simEpidataCS" (as if simEpidataCS

seed

an object specifying if and how the random number generator should be initialized (seeded). By default (NULL), the state of the random generator is not changed and the value of .Random.seed prior to simulation is stor

data

an object of class "epidataCS", usually the one to which the "twinstim" object was fitted.

simplify

logical. It is strongly recommended to set simplify = TRUE (default) if nsim is large. This saves space and computation time, because for each simulated epidemic only the events component is saved. All ot

control.siaf

see twinstim.

nCircle2Poly

see as.epidataCS. For simulate.twinstim, NULL means to use the same value as for data.

gmax

maximum value the temporal interaction function tiaf$g can attain. If NULL, then it is assumed as the maximum value of the type-specific values at 0, i.e. max(tiaf$g(rep.int(0,nTypes), tiafpars, 1:nTypes))

...

unused (arguments of the generic).

Value

• The function simEpidataCS returns a simulated epidemic of class "simEpidataCS", which enhances the class "epidataCS" by the following additional components known from objects of class "twinstim": timeRange, formula, coefficients, npars, call, runtime.

  The simulate.twinstim method has some additional attributes set on its result: call, seed, simplified, and runtime with their obvious meanings. Furthermore, if nsim > 1, it returns an object of class "simEpidataCSlist", the form of which depends on the value of simplify: if simplify = FALSE, then the return value is just a list of sequential simulations, each of class "simEpidataCS". However, if simplify = TRUE, then the sequential simulations share all components but the simulated events, i.e. the result is a list with the same components as a single object of class "simEpidataCS", but with events replaced by an eventsList containing the events returned by each of the simulations.

  The stgrid component of the returned "simEpidataCS" will be truncated to the actual end of the simulation, which might be $nEvents is reached during simulation. CAVE: Currently, simplify=TRUE in simulate.twinstim ignores that multiple simulated epidemics (nsim > 1) may have different stgrid time ranges. In a "simEpidataCSlist", the stgrid shared by all of the simulated epidemics is just the stgrid returned by the first simulation.

encoding

latin1

References

Douglas, D. H. and Peucker, T. K. (1973): Algorithms for the reduction of the number of points required to represent a digitized line or its caricature. Cartographica: The International Journal for Geographic Information and Geovisualization, 10, 112-122

Harrower, M. and Bloch, M. (2006): MapShaper.org: A Map Generalization Web Service. IEEE Computer Graphics and Applications, 26(4), 22-27. DOI-Link: http://dx.doi.org/10.1109/MCG.2006.85

Meyer, S., Elias, J. and H�{oe}hle, M. (2012): A space-time conditional intensity model for invasive meningococcal disease occurrence. Biometrics, 68, 607-616. DOI-Link: http://dx.doi.org/10.1111/j.1541-0420.2011.01684.x

Meyer, S. (2010): Spatio-Temporal Infectious Disease Epidemiology based on Point Processes. Master's Thesis, Ludwig-Maximilians-Universit�{ae}t M�{ue}nchen. Available as http://epub.ub.uni-muenchen.de/11703/

See Also

The plot.epidataCS and animate.epidataCS methods for plotting and animating continuous-space epidemic data, respectively, which also work for simulated epidemics (by inheritance).

Function twinstim for fitting spatio-temporal conditional intensity models to epidemic data.

Examples

data(imdepi)
data(imdepifit)

## load borders of Germany's districts (obtained from the
## Bundesamt f�r Kartographie und Geod�sie, Frankfurt am Main, Germany,
## www.geodatenzentrum.de), simplified by the "special" Visvalingam
## algorithm (level=60\%) using www.MapShaper.org:
load(system.file("shapes", "districtsD.RData", package="surveillance"))

## map of Germany's districts
plot(districtsD)
plot(stateD, add=TRUE, border=2, lwd=2)
# 'stateD' is the union of the polygons in 'districtsD'
# (created by maptools::unionSpatialPolygons)

if (require("gpclib")) {
    oopt <- surveillance.options(gpclib=TRUE)

    ## simulate 2 realizations (during a VERY short period -- for speed)
    ## considering events from data(imdepi) before t=31 as pre-history
    mysims <- simulate(imdepifit, nsim=2, seed=1, data=imdepi,
	               tiles=districtsD, W=stateD, t0=31, T=61, trace=FALSE,
                       nCircle2Poly=16, .allocate=100, simplify=TRUE)

    ## extract the first realization -> object of class simEpidataCS
    mysim1 <- mysims[[2]]
    summary(mysim1)
    plot(mysim1, aggregate="space")

    ## plot both epidemics using the plot-method for simEpidataCSlist's
    plot(mysims, aggregate="time", subset=type=="B")

    surveillance.options(oopt)
}