Skip to contents

Plot the Spatial or Temporal Interaction Function of a twimstim

twinstim_iafplot.Rd

The function plots the fitted temporal or (isotropic) spatial interaction function of a twinstim object. The implementation is illustrated in Meyer et al. (2017, Section 3), see vignette("twinstim").

Usage

iafplot(object, which = c("siaf", "tiaf"), types = NULL,
        scaled = c("intercept", "standardized", "no"), truncated = FALSE,
        log = "", conf.type = if (length(pars) > 1) "MC" else "parbounds",
        conf.level = 0.95, conf.B = 999, xgrid = 101,
        col.estimate = rainbow(length(types)), col.conf = col.estimate,
        alpha.B = 0.15, lwd = c(3,1), lty = c(1,2),
        verticals = FALSE, do.points = FALSE,
        add = FALSE, xlim = NULL, ylim = NULL, xlab = NULL, ylab = NULL,
        legend = !add && (length(types) > 1), ...)

Arguments

object

object of class "twinstim" containing the fitted model.

which

argument indicating which of the two interaction functions to plot. Possible values are "siaf" (default) for the spatial interaction \(f(x)\) as a function of the distance \(x\), and "tiaf" for the temporal interaction function \(g(t)\).

types

integer vector indicating for which event types the interaction function should be plotted in case of a marked "twinstim". The default types=NULL checks if the interaction function is type-specific: if so, types=1:nrow(object$qmatrix) is used, otherwise types=1.

scaled

character string determining if/how the the interaction function should be scaled. Possible choices are:

"intercept":

multiplication by the epidemic intercept.

"standardized":

division by the value at 0 distance such that the function starts at 1.

"no":

no scaling.

The first one is the default and required for the comparison of estimated interaction functions from different models. For backward compatibility, scaled can also be a boolean, where TRUE refers to "intercept" scaling and FALSE to "no" scaling.

truncated

logical indicating if the plotted interaction function should take the maximum range of interaction (eps.t/eps.s) into account, i.e., drop to zero at that point (if it is finite after all). If there is no common range of interaction, a rug indicating the various ranges will be added to the plot if truncated=TRUE. If truncated is a scalar, this value is used as the point eps where the function drops to 0.

log

a character string passed to plot.default indicating which axes should be logarithmic. If add=TRUE, log is set according to par("xlog") and par("ylog").

conf.type

type of confidence interval to produce.
If conf.type="MC" (or "bootstrap"), conf.B parameter vectors are sampled from the asymptotic (multivariate) normal distribution of the ML estimate of the interaction function parameters; the interaction function is then evaluated on the xgrid (i.e. temporal or spatial distances from the host) for each parameter realization to obtain a conf.level confidence interval at each point of the xgrid (or to plot the interaction functions of all Monte-Carlo samples if conf.level=NA). Note that the resulting plot is .Random.seed-dependent for the Monte-Carlo type of confidence interval.
If conf.type="parbounds", the conf.level Wald confidence intervals for the interaction function parameters are calculated and the interaction function is evaluated on the xgrid (distances from the host) for all combinations of the bounds of the parameters and the point-wise extremes of those functions are plotted. This type of confidence interval is only valid in case of a single parameter, i.e. scaled + nsiafpars == 1, but could also be used as a rough indication if the Monte-Carlo approach takes too long. A warning is thrown if the "parbounds" type is used for multiple parameters.
If conf.type="none" or NA or NULL, no confidence interval will be calculated.

conf.level

the confidence level required. For conf.type = "MC" it may also be specified as NA, in which case all conf.B sampled functions will be plotted with transparency value given by alpha.B.

conf.B

number of samples for the "MC" (Monte Carlo) confidence interval.

xgrid

either a numeric vector of x-values (distances from the host) where to evaluate which, or a scalar representing the desired number of evaluation points in the interval c(0,xlim[2]).
If the interaction function is a step function (siaf.step or tiaf.step), xgrid is ignored and internally set to c(0, knots).

col.estimate

vector of colours to use for the function point estimates of the different types.

col.conf

vector of colours to use for the confidence intervals of the different types.

alpha.B

alpha transparency value (as relative opacity) used for the conf.B sampled interaction functions in case conf.level = NA

lwd, lty

numeric vectors of length two specifying the line width and type of point estimates (first element) and confidence limits (second element), respectively.

verticals,do.points

graphical settings for step function kernels. These can be logical (as in plot.stepfun) or lists of graphical parameters.

add

add to an existing plot?

xlim, ylim

vectors of length two containing the x- and y-axis limit of the plot. The default y-axis range (ylim=NULL) is from 0 to the value of the (scaled) interaction function at \(x = 0\). The default x-axis (xlim=NULL) starts at 0, and the upper limit is determined as follows (in decreasing order of precedence):

  • If xgrid is a vector of evaluation points, xlim[2] is set to max(xgrid).

  • eps.t/eps.s if it is unique and finite.

  • If the interaction function is a step function with maxRange<Inf, i.e. it drops to 0 at maxRange, xlim[2] is set to maxRange.

  • Otherwise, it is set to the length of the observation period (which="tiaf") or the diagonal length of the bounding box of the observation region (which="siaf"), respectively.

xlab, ylab

labels for the axes with NULL providing sensible defaults.

legend

logical indicating if a legend for the types should be added. It can also be a list of arguments passed to legend to tweak the default settings.

...

additional arguments passed to the default plot method.

Value

A plot is created – see e.g. Figure 3(b) in Meyer et al. (2012).

The function invisibly returns a matrix of the plotted values of the interaction function (evaluated on xgrid, by type). The first column of the matrix contains the distance \(x\), and the remaining length(types) columns contain the (scaled) function values for each type.

The pointwise confidence intervals of the interaction functions are returned in similar matrices as attributes: if length(types)==1, there is a single attribute "CI", whereas for multiple types, the attributes are named paste0("CI.",typeNames) (where the typeNames are retrieved from object$qmatrix).

References

Meyer, S., Elias, J. and Höhle, M. (2012): A space-time conditional intensity model for invasive meningococcal disease occurrence. Biometrics, 68, 607-616. doi:10.1111/j.1541-0420.2011.01684.x

Meyer, S., Held, L. and Höhle, M. (2017): Spatio-temporal analysis of epidemic phenomena using the R package surveillance. Journal of Statistical Software, 77 (11), 1-55. doi:10.18637/jss.v077.i11

Author

Sebastian Meyer

See also

plot.twinstim, which calls this function.

Examples

data("imdepifit")

iafplot(imdepifit, "tiaf", scaled=FALSE)   # tiaf.constant(), not very exciting
iafplot(imdepifit, "siaf", scaled=FALSE)

# scaled version uses a Monte-Carlo-CI
set.seed(1)  # result depends on .Random.seed
iafplot(imdepifit, "siaf", scaled=TRUE, conf.type="MC", conf.B=199,
        col.conf=gray(0.4), conf.level=NA)  # show MC samples