surveillance (version 1.5-4)

twinstim_intensityplot: Plotting Intensities of Infection over Time or Space

Description

intensityplot method to plot the evolution of the total infection intensity, its epidemic proportion or its endemic proportion over time or space (integrated over the other dimension) of fitted twinstim models (or simEpidataCS). The "simEpidataCS"-method is just a wrapper around intensityplot.twinstim by making the "simEpidataCS" object "twinstim"-compatible, i.e. enriching it by the required model components and environment.

The intensity.twinstim auxiliary function returns functions which calculate the endemic or epidemic intensity at a specific time point or location (integrated over the other dimension).

Usage

## S3 method for class 'twinstim':
intensityplot(x,
    which = c("epidemic proportion", "endemic proportion", "total intensity"),
    aggregate = c("time", "space"), types = 1:nrow(x$qmatrix),
    tiles, tiles.idcol = NULL, plot = TRUE, add = FALSE,
    tgrid = 101, rug.opts = list(),
    sgrid = 128, polygons.args = list(), points.args = list(),
    cex.fun = sqrt, ...)
## S3 method for class 'simEpidataCS':
intensityplot(x, ...)
intensity.twinstim(x,
    aggregate = c("time", "space"), types = 1:nrow(x$qmatrix), 
    tiles, tiles.idcol = NULL)

Arguments

x

an object of class "twinstim" or "simEpidataCS", respectively.

which

"epidemic proportion", "endemic proportion", or "total intensity". Partial matching is applied. Determines whether to plot the path of the total intensity or its epidemic or endemic proportions over ti

aggregate

One of "time" or "space". The former results in a plot of the evolution of which as a function of time (integrated over the observation region $\bold{W}$), whereas the latter produces a spplot

types

event types to aggregate. By default, all types of events are aggregated, but one could also be interested in only one specific type or a subset of event types.

tiles

object of class SpatialPolygons representing the decomposition of $\bold{W}$ into different regions (as used in the corresponding stgrid of the "epidataCS

tiles.idcol

either a column index for tiles@data, or NULL (default), which refers to the "ID" slot of the polygons. The resulting ID's must correspond to the ID's used in the tile column of the stg

plot

logical indicating if a plot is desired, which defaults to TRUE. Otherwise, a function will be returned, which takes a vector of time points (if aggregate = "time") or a matrix of coordinates (if aggregate = "sp

add

logical. If TRUE and aggregate = "time", paths are added to the current plot, using lines. This does not work for aggregate = "space".

tgrid

either a numeric vector of time points when to evaluate which, or a scalar representing the desired number of evaluation points in the observation interval $[t_0, T]$. This argument is unused for aggregate = "space".

rug.opts

if a list, its elements are passed as arguments to the function rug, which will mark the time points of the events if aggregate = "time" (it is unused in the spatial case); otherwise (e

sgrid

either an object of class "SpatialPixels" (or coercible to that class) representing the locations where to evaluate which, or a scalar representing the total number of points of a grid const

polygons.args

if a list, its elements are passed as arguments to sp.polygons, which will add tiles to the plot if aggregate = "space" (it is unused for the temporal plot). By def

points.args

if a list, its elements are passed as arguments to sp.points, which will add the event locations to the plot if aggregate = "space" (it is unused for the temporal plot). By defaul

cex.fun

function which takes a vector of counts of events at each unique location and returns a (vector of) cex value(s) for the sizes of the points at the event locations used in points.args. Defaults to the sqrt()

...

further arguments passed to plot or lines (if aggregate = "time"), or to spplot (if aggregate = "space"). For intensityplot.simEpidataCS

Value

• If plot = FALSE or aggregate = "time", a function is returned, which takes a vector of time points (if aggregate = "time") or a matrix of coordinates (if aggregate = "space"), and returns which on this grid.

  If plot = TRUE and aggregate = "space", the trellis.object containg the spatial plot is returned.

See Also

plot.twinstim, which calls this function.

Examples

data(imdepi)
data(imdepifit)

# for the intensityplot we need the model environment, which can be
# easily added by the intelligent update method (no need to refit the model)
imdepifit <- update(imdepifit, model=TRUE)

## path of the total intensity
opar <- par(mfrow=c(2,1))
intensityplot(imdepifit, which="total intensity",
              aggregate="time", tgrid=500)
plot(imdepi, breaks=100)
par(opar)

## path of the proportion of the epidemic component by event
intensityplot(imdepifit, which="epidemic proportion",
              aggregate="time", tgrid=500, types=1)
intensityplot(imdepifit, which="epidemic proportion",
              aggregate="time", tgrid=500, types=2, add=TRUE, col=2)
legend("topright", legend=levels(imdepi$events$type), lty=1, col=1:2,
       title = "event type")

if (surveillance.options("allExamples") && require("maptools"))
{
  ## spatial shape of the intensity (aggregated over time)
  # 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"))

  # total intensity (using a rather sparse 'sgrid' for speed)
  intensityplot(imdepifit, which="total intensity",
                aggregate="space", tiles=districtsD, sgrid=500)

  # epidemic proportion by type
  maps_epiprop <- lapply(1:2, function (type) {
      intensityplot(imdepifit, which="epidemic", aggregate="space",
                    types=type, tiles=districtsD, sgrid=1000,
                    at=seq(0,1,by=0.1), col.regions=rev(heat.colors(20)))
  })
  plot(maps_epiprop[[1]], split=c(1,1,2,1), more=TRUE)
  plot(maps_epiprop[[2]], split=c(2,1,2,1))
}