hhh4_simulate_plot: Plot Simulations from `"hhh4"` Models

Description

Arrays of simulated counts from simulate.hhh4 can be visualized as final size boxplots, individual or average time series, or fan charts (using the fanplot package). An aggregate-method is also available.

Usage

# S3 method for hhh4sims
plot(x, ...)
# S3 method for hhh4sims
aggregate(x, units = TRUE, time = FALSE, ..., drop = FALSE)
as.hhh4simslist(x, ...)
# S3 method for hhh4simslist
plot(x, type = c("size", "time", "fan"), ...,
     groups = NULL, par.settings = list())
# S3 method for hhh4simslist
aggregate(x, units = TRUE, time = FALSE, ..., drop = FALSE)
plotHHH4sims_size(x, horizontal = TRUE, trafo = NULL, observed = TRUE,
                  names = base::names(x), ...)
plotHHH4sims_time(x, average = mean, individual = length(x) == 1,
    conf.level = if (individual) 0.95 else NULL,
    matplot.args = list(), initial.args = list(), legend = length(x) > 1,
    xlim = NULL, ylim = NULL, add = FALSE, ...)
plotHHH4sims_fan(x, which = 1,
    fan.args = list(), observed.args = list(), initial.args = list(),
    means.args = NULL, key.args = NULL, xlim = NULL, ylim = NULL,
    add = FALSE, xaxis = list(), ...)

Arguments

x: an object of class "hhh4sims" (as resulting from the simulate-method for "hhh4" models if simplify = TRUE was set), or an "hhh4simslist", i.e., a list of such simulations potentially obtained from different model fits (using the same simulation period).
type: a character string indicating the summary plot to produce.
...: further arguments passed to methods.
groups: an optional factor to produce stratified plots by groups of units. The special setting groups = TRUE is a convenient shortcut for one plot by unit.
par.settings: a list of graphical parameters for par. Sensible defaults for mfrow, mar and las will be applied unless overridden or !is.list(par.settings).
horizontal: a logical indicating if the boxplots of the final size distributions should be horizontal (the default).
trafo: an optional transformation function from the scales package, e.g., sqrt_trans.
observed: a logical indicating if a line and axis value for the observed size of the epidemic should be added to the plot. Alternatively, a list with graphical parameters can be specified to modify the default values.
names: a character vector of names for x.
average: scalar-valued function to apply to the simulated counts at each time point.
individual: a logical indicating if the individual simulations should be shown as well.
conf.level: a scalar in (0,1), which determines the level of the pointwise quantiles obtained from the simulated counts at each time point. A value of NULL disables the confidence interval.
matplot.args: a list of graphical parameters for matlines.
initial.args: if a list (of graphical parameters for lines), a bar for the initial number of cases is added to the plot.
legend: a logical, a character vector (providing names for x), or a list of parameters for legend.
xlim,ylim: vectors of length 2 determining the axis limits.
add: a logical indicating if the (mean) simulated time series or the fan chart, respectively, should be added to an existing plot.
which: a single integer or a character string selecting the model in x for which to produce the fan chart. This is only relevant if x is a "hhh4simslist" of simulations from multiple models. Defaults to the first model.
fan.args: a list of graphical parameters for the fan, e.g., to employ a different colorRampPalette as fan.col, or to enable contour lines via ln.
observed.args: if a list (of graphical parameters for lines), the originally observed counts are added to the plot.
means.args: if a list (of graphical parameters for lines), the point forecasts are added to the plot (by default as a white line within the fan).
key.args: if a list, a color key (in fan's "boxfan"-style) is added to the fan chart. The list may include positioning parameters start (the x-position) and ylim (the y-range of the color key), space to modify the width of the boxfan, and rlab to modify the labels. The color key is disabled by default. An alternative way of labeling the quantiles is via the argument ln in fan.args, see the Examples.
xaxis: if a list of arguments for addFormattedXAxis, that function is used to draw the time axis, otherwise a default x-axis is drawn.
units: a logical indicating aggregation over units. Can also be a factor (or something convertible to a factor using as.factor) to aggregate groups of units.
time: a logical indicating if the counts should be summed over the whole simulation period.
drop: a logical indicating if the unit dimension and the "hhh4sims" (or "hhh4simslist") class should be dropped after aggregating over (groups of) units.

Author

Sebastian Meyer

Examples

Run this code

### univariate example
data("salmAllOnset")

## fit a hhh4 model to the first 13 years
salmModel <- list(end = list(f = addSeason2formula(~1 + t)),
                  ar = list(f = ~1), family = "NegBin1", subset = 2:678)
salmFit <- hhh4(salmAllOnset, salmModel)

## simulate the next 20 weeks ahead
salmSims <- simulate(salmFit, nsim = 300, seed = 3, subset = 678 + seq_len(20),
                     y.start = observed(salmAllOnset)[678,])

## compare final size distribution to observed value
summary(aggregate(salmSims, time = TRUE))  # summary of simulated values
plot(salmSims, type = "size")

## individual and average simulated time series with a confidence interval
plot(salmSims, type = "time", main = "20-weeks-ahead simulation")

## fan chart based on the quantiles of the simulated counts at each time point
## point forecasts are represented by a white line within the fan
if (requireNamespace("fanplot")) {
    plot(salmSims, type = "fan", main = "20-weeks-ahead simulation",
         fan.args = list(ln = 1:9/10), means.args = list())
}


### multivariate example
data("measlesWeserEms")

## fit a hhh4 model to the first year
measlesModel <- list(
    end = list(f = addSeason2formula(~1), offset = population(measlesWeserEms)),
    ar = list(f = ~1),
    ne = list(f = ~1 + log(pop),
        weights = W_powerlaw(maxlag = 5, normalize = TRUE)),
    family = "NegBin1", subset = 2:52,
    data = list(pop = population(measlesWeserEms)))
measlesFit1 <- hhh4(measlesWeserEms, control = measlesModel)

## use a Poisson distribution instead (just for comparison)
measlesFit2 <- update(measlesFit1, family = "Poisson")

## simulate realizations from these models during the second year
measlesSims <- lapply(X = list(NegBin = measlesFit1, Poisson = measlesFit2),
                      FUN = simulate, nsim = 50, seed = 1, subset = 53:104,
                      y.start = observed(measlesWeserEms)[52,])

## final size of the first model
plot(measlesSims[[1]])

## stratified by groups of districts
mygroups <- factor(substr(colnames(measlesWeserEms), 4, 4))
apply(aggregate(measlesSims[[1]], time = TRUE, units = mygroups), 1, summary)
plot(measlesSims[[1]], groups = mygroups)

## a class and plot-method for a list of simulations from different models
measlesSims <- as.hhh4simslist(measlesSims)
plot(measlesSims)

## simulated time series
plot(measlesSims, type = "time", individual = TRUE, ylim = c(0, 80))

## fan charts
if (requireNamespace("fanplot")) {
    opar <- par(mfrow = c(2,1))
    plot(measlesSims, type = "fan", which = 1, ylim = c(0, 80), main = "NegBin",
         key.args = list())
    plot(measlesSims, type = "fan", which = 2, ylim = c(0, 80), main = "Poisson")
    par(opar)
}

Run the code above in your browser using DataLab