twimstimThe 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").
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), ...)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).
object of class "twinstim" containing the fitted model.
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)\).
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.
character string determining if/how the the interaction function should be scaled. Possible choices are:
multiplication by the epidemic intercept.
division by the value at 0 distance such that the function starts at 1.
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.
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.
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").
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.
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.
number of samples for the "MC" (Monte Carlo)
confidence interval.
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).
vector of colours to use for the function point estimates of the different types.
vector of colours to use for the confidence intervals of the different types.
alpha transparency value (as relative opacity) used for the conf.B
sampled interaction functions in case conf.level = NA
numeric vectors of length two specifying the line width and type of point estimates (first element) and confidence limits (second element), respectively.
graphical settings for step function
kernels. These can be logical (as in plot.stepfun) or
lists of graphical parameters.
add to an existing plot?
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.
labels for the axes with NULL providing sensible defaults.
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.
Sebastian Meyer
Meyer, S., Elias, J. and Höhle, M. (2012): A space-time conditional intensity model for invasive meningococcal disease occurrence. Biometrics, 68, 607-616. tools:::Rd_expr_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. tools:::Rd_expr_doi("10.18637/jss.v077.i11")
plot.twinstim, which calls this function.
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
Run the code above in your browser using DataLab