surveillance (version 1.12.1)

R0: Computes reproduction numbers from fitted models

Description

The S3 generic function R0 defined in package surveillance is intended to compute reproduction numbers from fitted epidemic models. The package currently defines a method for the "twinstim" class, which computes expected numbers of infections caused by infected individuals depending on the event type and marks attached to the individual, which contribute to the infection pressure in the epidemic predictor of that class. There is also a method for simulated "epidataCS" (just a wrapper for the "twinstim"-method).

Usage

R0(object, ...)
## S3 method for class 'twinstim':
R0(object, newevents, trimmed = TRUE, newcoef = NULL, ...)
## S3 method for class 'simEpidataCS':
R0(object, trimmed = TRUE, ...)
simpleR0(object, eta = coef(object)[["e.(Intercept)"]],
         eps.s = NULL, eps.t = NULL, newcoef = NULL)

Arguments

object

A fitted epidemic model object for which an R0 method exists.

newevents

an optional data.frame of events for which the reproduction numbers should be calculated. If omitted, it is calculated for the original events from the fit. In this case, if trimmed = TRUE (the default), the result

trimmed

logical indicating if the individual reproduction numbers should be calculated by integrating the epidemic intensities over the observation period and region only (trimmed = TRUE) or over the whole time-space domain R+ x R^2 (

newcoef

the model parameters to use when calculating reproduction numbers. The default (NULL) is to use the MLE coef(object). This argument mainly serves the construction of Monte Carlo confidence intervals by evaluating

...

additional arguments passed to methods. Currently unused for the twinstim method.

eta

a value for the epidemic linear predictor, see details.

eps.s,eps.t

the spatial/temporal radius of interaction. If NULL (the default), the original value from the data is used if this is unique and an error is thrown otherwise.

Value

• For the R0 methods, a numeric vector of estimated reproduction numbers from the fitted model object corresponding to the rows of newevents (if supplied) or the original fitted events including events of the prehistory.

  For simpleR0, a single number (see details).

encoding

latin1

Details

For the "twinstim" class, the individual-specific expected number $\mu_j$ of infections caused by individual (event) $j$ inside its theoretical (untrimmed) spatio-temporal range of interaction given by its eps.t ($\epsilon$) and eps.s ($\delta$) values is defined as follows (cf. Meyer et al, 2012): $$\mu_j = e^{\eta_j} \cdot \int_{b(\bold{0},\delta)} f(\bold{s}) d\bold{s} \cdot \int_0^\epsilon g(t) dt .$$ Here, $b(\bold{0},\delta)$ denotes the disc centred at (0,0)' with radius $\delta$, $\eta_j$ is the epidemic linear predictor, $g(t)$ is the temporal interaction function, and $f(\bold{s})$ is the spatial interaction function. For a type-specific twinstim, there is an additional factor for the number of event types which can be infected by the type of event $j$ and the interaction functions may be type-specific as well.

Alternatively to the equation above, the trimmed (observed) reproduction numbers are obtain by integrating over the observed infectious domains of the individuals, i.e. integrate $f$ over the intersection of the influence region with the observation region W (i.e. over ${ W \cap b(\bold{s}_j,\delta) } - \bold{s}_j$) and $g$ over the intersection of the observed infectious period with the observation period $(t_0;T]$ (i.e. over $(0; \min(T-t_j,\epsilon)]$).

The function simpleR0 computes $$\exp(\eta) \cdot \int_{b(\bold{0},\delta)} f(\bold{s}) d\bold{s} \cdot \int_0^{\epsilon} g(t) dt ,$$ where $\eta$ defaults to $\gamma_0$ disregarding any epidemic effects of types and marks. It is thus only suitable for simple epidemic twinstim models with epidemic = ~1, a diagonal (or secondary diagonal) qmatrix, and type-invariant interaction functions. simpleR0 mainly exists for use by epitest.

(Numerical) Integration is performed exactly as during the fitting of object, for instance object$control.siaf is queried if necessary.

References

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

Examples

## load the 'imdepi' data and a model fit
data("imdepi")
data("imdepifit")

## calculate individual and type-specific reproduction numbers
R0s <- R0(imdepifit)
tapply(R0s, imdepi$events@data[names(R0s), "type"], summary)

## untrimmed R0 for a specific event
R0(imdepifit, newevents=marks(imdepi)[1,], trimmed=FALSE)


### compute a Monte Carlo confidence interval

## use a simpler model with constant 'siaf' for speed
simplefit <- update(imdepifit, epidemic=~type, siaf=NULL, subset = NULL)

## we'd like to compute the mean R0's by event type
meanR0ByType <- function (newcoef) {
    R0events <- R0(simplefit, newcoef=newcoef)
    tapply(R0events, imdepi$events@data[names(R0events),"type"], mean)
}
(meansMLE <- meanR0ByType(newcoef=NULL))

## sample B times from asymptotic multivariate normal of the MLE
B <- 5  # CAVE: toy example! In practice this has to be much larger
set.seed(123)
parsamples <- MASS::mvrnorm(B, mu=coef(simplefit), Sigma=vcov(simplefit))

## for each sample compute the 'meanR0ByType'
meansMC <- apply(parsamples, 1, meanR0ByType)

## get the quantiles and print the result
cisMC <- apply(cbind(meansMLE, meansMC), 1, quantile, probs=c(0.025,0.975))
print(rbind(MLE=meansMLE, cisMC))


### R0 for a simple epidemic model (homogeneous individuals)

mepi1 <- update(simplefit, epidemic = ~1, data = imdepi,
                subset = type == "B", model = TRUE, verbose = FALSE)
## using the default spatial and temporal ranges of interaction
(R0B <- simpleR0(mepi1))  # eps.s=200, eps.t=30
stopifnot(identical(R0B, R0(mepi1, trimmed = FALSE)[[1]]))
## assuming less interaction
simpleR0(mepi1, eps.s = 50, eps.t = 15)