plot_moments_by_unit: Plot predictive or stationary moments by unit

In jbracher/hhh4addon: Extensions to endemic-epidemic timeseries modeling from package surveillance

Plot predictive or stationary moments by unit

Description

Plot negative binomial approximation of predictive or stationary distributions. Usually to be used with aggregated predictions (where columns correspond to regions or age groups; no temporal structure kept).

Usage

plot_moments_by_unit(
  mom,
  probs = 1:99/100,
  add_observed = TRUE,
  add_pred_means = TRUE,
  fan.col = colorRampPalette(c("darkgreen", "gray90")),
  pt.col = "red",
  pt.cex = 0.3,
  mean_col = "black",
  mean_lty = "dashed",
  ln = NULL,
  rlab = NULL,
  style = "boxfan",
  space = 0.5,
  add_legend = FALSE,
  probs_legend = c(1, 25, 50, 75, 99)/100,
  ylim = NULL,
  main = NULL,
  xlab = NULL,
  las = NULL,
  axes = TRUE,
  ...
)

Arguments

mom	an object of class predictive_moments_hhh4 or stationary_moments_hhh4, usually an aggregated prediction without temporal structure (aggregated using aggregate_prediction)
probs	probabilities displayed in fanplot (passed to fanplot::fan)
add_observed	should obseved values be added to the plot?
fan.col	color palette for fanplot (passed to fanplot::fan)
pt.col, pt.cex	point color and size for observations
mean_col, mean_lty	line color and type for predictive/stationary means
ln, rlab, style, space	additional arguments passed to fanplot::fan)
add_legend	should a legend with the colour coding be added?
probs_legend	probabilities to be displayed in the legend
ylim, main, xlab, las, axes	usual plotting parameters

Examples

# load data:
data("noroBL")

########
# fit a bivariate model:
controlBL <- list(end = list(f = addSeason2formula(~ -1 + fe(1, unitSpecific = TRUE))),
                  ar = list(f = ~ -1 + fe(1, unitSpecific = TRUE)),
                  ne = list(f = ~ -1 + fe(1, unitSpecific = TRUE)),
                  family = "NegBinM", subset = 2:260) # not a very parsimonious parametrization, but feasible
fitBL <- hhh4(noroBL, control = controlBL)
pred_mom <- predictive_moments(fitBL, t_condition = 260, lgt = 52, return_Sigma = TRUE)
# Sigma is required in order to aggregate predictions.

#########
# perform an aggregation over time: total burden in the two regions
aggr_matr_total_burden <- matrix(rep(c(1, 0, 0, 1), 52), nrow = 2,
                                 dimnames = list(c("Bremen", "Lower Saxony"),
                                                 NULL))
pred_mom_total_burden <- aggregate_moments(pred_mom, aggr_matr_total_burden)
plot_moments_by_unit(pred_mom_total_burden, main = "Total burden 2016", add_legend = TRUE)

jbracher/hhh4addon documentation built on Feb. 16, 2024, 1:45 a.m.