Package 'hhh4ZI'

object of class "sts" containing the (multivariate) count data time series.

a list containing the model specification and control arguments, the parts relating to hhh4 model is the same as in surveillance::hhh4:

• ar Model for the autoregressive component given as list with the following components:

  • f = ~ -1 a formula specifying $\log(\lambda_{it})$.

  • offset = 1 optional multiplicative offset, either 1 or a matrix of the same dimension as observed(stsObj).

  • lag = 1 a positive integer meaning autoregression on $y_{i,t-lag}$.

• ne Model for the neighbour-driven component given as list with the following components:

  • f = ~ -1 a formula specifying $\log(\phi_{it})$.

  • offset = 1 optional multiplicative offset, either 1 or a matrix of the same dimension as observed(stsObj).

  • lag = 1 a non-negative integer meaning dependency on $y_{j,t-lag}$.

  • weights = neighbourhood(stsObj) == 1 neighbourhood weights $w_{ji}$. The default corresponds to the original formulation by Held et al (2005), i.e., the spatio-temporal component incorporates an unweighted sum over the lagged cases of the first-order neighbours. See Paul et al (2008) and Meyer and Held (2014) for alternative specifications, e.g., W_powerlaw. Time-varying weights are possible by specifying an array of dim() c(nUnits, nUnits, nTime), where nUnits=ncol(stsObj) and nTime=nrow(stsObj).

  • scale = NULL optional matrix of the same dimensions as weights (or a vector of length ncol(stsObj)) to scale the weights to scale * weights.

  • normalize = FALSE logical indicating if the (scaled) weights should be normalized such that each row sums to 1.

• end Model for the endemic component given as list with the following components

  • f = ~ 1 a formula specifying $\log(\nu_{it})$.

  • offset = 1 optional multiplicative offset $e_{it}$, either 1 or a matrix of the same dimension as observed(stsObj).

• zi Model for the zero inflation component given as list with the following components:

  • f = ~ -1 a formula specifying $logit(\gamma_{it})$.

  • lag = 1 a vector of positive integers meaning autoregression on the respective $y_{i,t-lag}$ terms. Use NULL or integer(0) to exclude AR terms.

  • lag.unitSpecific logical indicating if the autoregressive parameter in the zero inflation part is unit specific.

• family Distributional family – the Negative Binomial distribution. The overdispersion parameter can be assumed to be the same for all units ("NegBin1"), to vary freely over all units ("NegBinM"), or to be shared by some units (specified by a factor of length ncol(stsObj) such that its number of levels determines the number of overdispersion parameters). Note that "NegBinM" is equivalent to factor(colnames(stsObj), levels = colnames(stsObj)).

• subset Typically 2:nrow(obs) if model contains autoregression.

• optimizer a list of three lists of control arguments.

  The "stop" list specifies two criteria for the outer optimization of regression and variance parameters: the relative tolerance for parameter change using the criterion max(abs(x[i+1]-x[i])) / max(abs(x[i])), and the maximum number niter of outer iterations.

  Control arguments for the single optimizers are specified in the lists named "regression" and "variance". method="nlminb" is the default optimizer for regression update, however, the methods from optim may also be specified (as well as "nlm" but that one is not recommended here). For the variance updates, only Nelder-Mead optimization (method="Nelder-Mead") is provided. All other elements of these two lists are passed as control arguments to the chosen method, e.g., if method="nlminb" adding iter.max=50 increases the maximum number of inner iterations from 20 (default) to 50.

• verbose non-negative integer (usually in the range 0:3) specifying the amount of tracing information to be output during optimization.

• start a list of initial parameter values replacing initial values set via fe and ri. Since surveillance 1.8-2, named vectors are matched against the coefficient names in the model (where unmatched start values are silently ignored), and need not be complete, e.g., start = list(fixed = c("-log(overdisp)" = 0.5)) (default: 2) for a family = "NegBin1" model. In contrast, an unnamed start vector must specify the full set of parameters as used by the model.

• data a named list of covariates that are to be included as fixed effects (see fe) in any of the 3 component formulae. By default, the time variable t is available and used for seasonal effects created by addSeason2formula. In general, covariates in this list can be either vectors of length nrow(stsObj) interpreted as time-varying but common across all units, or matrices of the same dimension as the disease counts observed(stsObj).

• keep.terms logical indicating if the terms object used in the fit is to be kept as part of the returned object. This is usually not necessary, since the terms object is reconstructed by the terms-method for class "hhh4ZI" if necessary (based on stsObj and control, which are both part of the returned "hhh4ZI" object).

logical (or a subset of c("numDeriv", "maxLik")), indicating if (how) the implemented analytical score vector and Fisher information matrix should be checked against numerical derivatives at the parameter starting values, using the packages numDeriv and/or maxLik. If activated, hhh4 will return a list containing the analytical and numerical derivatives for comparison (no ML estimation will be performed). This is mainly intended for internal use by the package developers.

fitted hhh4 model (class "hhh4").

numeric vector of length 2 specifying the time range in which to compute one-step-ahead predictions (for the time points tp[1]+1, ..., tp[2]+1). If a single time index is specified, it is interpreted as tp[1], and tp[2] is set to the penultimate time point of result$control$subset.

unused (argument of the generic).

The default "rolling" procedure sequentially refits the model up to each time point in tp and computes the one-step-ahead predictions for the respective next time point. The alternative types are no true one-step-ahead predictions but much faster: "first" will refit the model for the first time point tp[1] only and use this specific fit to calculate all subsequent predictions, whereas "final" will just use result to calculate these. The latter case thus gives nothing else than a subset of result$fitted.values if the tp's are part of the fitted subset result$control$subset.

Which initial parameter values should be used when successively refitting the model to subsets of the data (up to time point tp[1], up to tp[1]+1, ...) if type="rolling"? Default ("current") is to use the parameter estimates from the previous time point, and "final" means to always use the estimates from result as initial values. Alternatively, which.start can be a list of start values as expected by hhh4, which then replace the corresponding estimates from result as initial values. This argument is ignored for “non-rolling” types.

logical indicating if parameter estimates and log-likelihoods from the successive fits should be returned.

non-negative integer (usually in the range 0:3) specifying the amount of tracing information to output. During hhh4 model updates, the following verbosity is used: 0 if cores > 1, otherwise verbose-1 if there is more than one time point to predict, otherwise verbose.

the number of cores to use when computing the predictions for the set of time points tp in parallel (with mclapply). Note that parallelization is not possible in the default setting type="rolling" and which.start="current" (use which.start="final" for this to work).

an object of class "oneStepAhead" or "hhh4".

numeric vector of probabilities with values in [0,1].

an object of class "oneStepAhead".

unused (argument of the generic).

required confidence level of the prediction interval.

single integer or character selecting a unit for which to produce the plot.

x-coordinate of the first prediction. If start=NULL (default), this is derived from x.

an object of class "oneStepAhead_hhh4ZI" or "hhh4ZI", respectively.

integer or character vector indexing the units for which to compute the PIT histogram. By default, all units are considered.

further arguments passed to pit.default.

subset of time points for which to produce the PIT histogram. Defaults to the subset used for fitting the model.

a fitted hhh4ZI object.

type of plot: either "fitted" component means of selected units along time along with the observed counts, or "season"ality plots of the model components and the epidemic dominant eigenvalue (which may also be plotted along overall time by type="maxEV", especially if the model contains time-varying neighbourhood weights or unit-specific epidemic effects), or "maps" of the fitted mean components averaged over time, or a map of estimated region-specific random intercepts ("ri") of a specific model component. The latter two require x$stsObj to contain a map.

For plotHHH4_season and plotHHH4_maxEV, one or more hhh4-fits, or a single list of these. Otherwise further arguments passed on to other functions.
For the plot-method these go to the specific plot type function.
plotHHH4_fitted passes them to plotHHH4_fitted1, which is called sequentially for every unit in units.
plotHHH4_maps and plotHHH4_ri pass additional arguments to spplot, and plotHHH4_neweights to the plotter.

integer or character vector specifying a single unit or possibly multiple units to plot. It indexes colnames(x$stsObj).
In plotHHH4_fitted, units=NULL plots all units.
In the seasonality plot, selection of a unit is only relevant if the model contains unit-specific intercepts or seasonality terms.

main title(s) for the selected unit(s) / components. If NULL (default), plotHHH4_fitted1 will use the appropriate element of colnames(x$stsObj), whereas plotHHH4_season uses default titles.

length 3 vectors specifying the fill and border colors for the endemic, autoregressive, and spatio-temporal component polygons (in this order).

style specifications for the dots drawn to represent the observed counts. pch=NA can be used to disable these dots.

list of graphical parameters for par. Sensible defaults for mfrow, mar and las will be applied unless overridden or !is.list(par.settings).

Integer vector specifying in which of the length(units) frames the legend should be drawn. If a logical vector is supplied, which(legend) determines the frame selection, i.e., the default is to drawn the legend in the first (upper left) frame only, and legend=FALSE results in no legend being drawn.

list of arguments for legend, e.g., to modify the default positioning list(x="topright", inset=0.02).

logical indicating if the legend should contain a line for the dots corresponding to observed counts.

if TRUE or (a permutation of) colnames(x$stsObj), the fitted mean will be decomposed into the contributions from each single unit and the endemic part instead of the default endemic + AR + neighbours decomposition.

logical indicating if the fitted components should be summed over all units to be compared with the total observed counts at each time point. If total=TRUE, the units/unit argument is ignored.

(internal) use different component means than those estimated and available from x.

time range to plot specified by vectors of length two in the form c(year,number), see "sts".

if this is a list (of arguments for addFormattedXAxis), the time axis is nicely labelled similar to stsplot_time. Note that in this case or if xaxis = NA, the basic time indexes 1:nrow(x$stsObj) will be used as x coordinates, which is different from the long-standing default (xaxis = NULL) with a real time scale.

numeric vector of length 2 specifying the x-axis range. The default (NULL) is to plot the complete time range (type="fitted") or period (type="season"), respectively.

y-axis range. For type="fitted", this defaults to c(0,max(observed(x$stsObj)[,unit])). For type="season", ylim must be a list of length length(components) specifying the range for every component plot, or a named list to customize only a subset of these. If only one ylim is specified, it will be recycled for all components plots.

axis labels. For plotHHH4_season, ylab specifies the y-axis labels for all components in a list (similar to ylim). If NULL or incomplete, default mathematical expressions are used. If a single name is supplied such as the default ylab="" (to omit y-axis labels), it is used for all components.

logical indicating if dots for zero observed counts should be omitted. Especially useful if there are too many.

a character vector specifying the components of the mean for which to produce maps. By default, the overall mean and all three components are shown.

a logical indicating whether the component maps should display proportions of the total mean instead of absolute numbers.

a numeric vector of length length(which) (recycled as necessary) specifying upper limits for the color keys of the maps, using a lower limit of 0. A missing element (NA) means to use a map-specific color key only covering the range of the values in that map (can be useful for prop = TRUE). The default zmax = NULL means to use the same scale for the component maps and a separate scale for the map showing the overall mean.

a vector of colors used to encode the fitted component means (see levelplot). For plotHHH4_maps, the length of this color vector also determines the number of levels, using 10 heat colors by default.

determines if and how regions are labeled, see layout.labels.

optional list of additional layout items, see spplot.

an object inheriting from "SpatialPolygons" with row.names covering colnames(x).

list of line style specifications passed to matplot, e.g., lty, lwd, col.

list of line style specifications (e.g., lty or col) passed to abline when drawing the reference line (h=1) in plots of seasonal effects (if intercept=FALSE) and of the dominant eigenvalue. The reference line is omitted if refline.args is not a list.

component for which to plot the estimated region-specific random intercepts. Must partially match one of colnames(ranef(x, tomatrix=TRUE)).

logical indicating whether to exp-transform the color-key axis labels to show the multiplicative effect of the region-specific random intercept on the respective component. Axis labels are then computed using log_breaks from package scales (if that is available) or axisTicks (as a fallback) respecting the colorkey$tick.number setting (default: 7). The default is FALSE.

a numeric vector of breaks for the color levels (see levelplot), or a list specifying the number of breaks n (default: 10) and their range (default: range of the random effects, extended to be symmetric around 0). In the latter case, breaks are equally spaced (on the original, non-exp scale of the random intercepts). If exp=TRUE, custom breaks (or range) need to be given on the exp-scale.

a Boolean indicating whether to draw the color key. Alternatively, a list specifying how to draw it, see levelplot.

list of graphical parameters for sp.polygons, applied to regions with missing random intercepts, i.e., not included in the model. Such extra regions won't be plotted if !is.list(gpar.missing).

character vector of component names, i.e., a subset of c("ar", "ne", "end"), for which to plot the estimated seasonality. If NULL (the default), only components which appear in any of the models in ... are plotted.
A seasonality plot of the epidemic dominant eigenvalue is also available by including "maxEV" in components, but it only supports models without epidemic covariates/offsets.

logical indicating whether to include the global intercept. For plotHHH4_season, the default (FALSE) means to plot seasonality as a multiplicative effect on the respective component. Multiplication by the intercept only makes sense if there are no further (non-centered) covariates/offsets in the component.

a numeric value giving the (longest) period of the harmonic terms in the model. This usually coincides with the freq of the data (the default), but needs to be adjusted if the model contains harmonics with a longer periodicity.

the (name of a) function used to produce the plot of weights (a numeric vector) as a function of neighbourhood order (a factor variable). It is called as plotter(Weight ~ Distance, ...) and defaults to boxplot. A useful alternative is, e.g., stripplot from package lattice.

vector of neighbourhood orders to be excluded from plotting (passed to factor). By default, the neighbourhood weight for order 0 is not shown, which is usually zero anyway.

maximum order of neighbourhood to be assumed when computing the nbOrder matrix. This additional step is necessary iff neighbourhood(x$stsObj) only specifies a binary adjacency matrix.

an object of class "oneStepAhead_hhh4ZI" or "hhh4ZI", respectively.

character vector determining which scores to compute. The package surveillance implements the following proper scoring rules: logarithmic score ("logs"), ranked probability score ("rps"), Dawid-Sebastiani score ("dss"), and squared error score ("ses"). The normalized SES ("nses") is also available but it is improper and hence not computed by default.
It is possible to name own scoring rules in which. These must be functions of (x, mu, size), vectorized in all arguments (time x unit matrices) except that size is NULL in case of a Poisson model. See the available scoring rules for guidance, e.g., dss.

integer or character vector indexing the units for which to compute the scores. By default, all units are considered.

logical indicating if the function should also return sign(x-mu), i.e., the sign of the difference between the observed counts and corresponding predictions. This does not really make sense when averaging over multiple units with individual=FALSE.

logical indicating if the individual scores of the units should be returned. By default (FALSE), the individual scores are averaged over all units.

unused (argument of the generic).

subset of time points for which to compute the scores.

of class "hhh4ZI".

number of time series to simulate. Defaults to 1.

an object specifying how the random number generator should be initialized for simulation (via set.seed). The initial state will also be stored as an attribute "seed" of the result. The original state of the .Random.seed will be restored at the end of the simulation. By default (NULL), neither initialization nor recovery will be done. This behaviour is copied from the simulate.lm method.

vector or matrix (with ncol(object$stsObj) columns) with starting counts for the epidemic components. If NULL, the observed means in the respective units of the data in object during subset are used.

time period in which to simulate data. Defaults to (and cannot exceed) the whole period defined by the underlying "sts" object.

coefficients used for simulation from the model in object. Default is to use the fitted parameters. Note that the coefs-vector must be in the same order and scaling as coef(object), which especially means reparamPsi = TRUE (as per default when using the coef-method to extract the parameters). The overdispersion parameter in coefs is the inverse of the dispersion parameter size in rnbinom.

character vector indicating which components of the fitted model object should be active during simulation. For instance, a simulation with components="end" is solely based on the fitted endemic mean.

logical indicating if only the simulated counts (TRUE) or the full "sts" object (FALSE) should be returned for every replicate. By default a full "sts" object is returned iff nsim=1.

unused (argument of the generic).

a fitted "hhh4ZI" model.

elements modifying the original control list.

a named list of numeric vectors (with names "ar", "end", etc) to adjust the number of harmonics in the model components via addSeason2formula, or NULL (meaning no modification of seasonal terms).

refit on a subset of the data up to that time point (used by oneStepAhead).

logical indicating if coef(object) should be used as starting values for the new fit.

logical indicating if the updated hhh4ZI call should be evaluated or returned.

a "hhh4" fit.

either a list of specifications for the zi component of hhh4ZI or a hhh4 control list with added zi specification.

further (non-zi) control elements.