surveillance (version 1.7-0)

powerlaw: Power-Law Neighbourhood Weights for hhh4 Models

Description

Set up a power-law weight structure for the neighbourhood component of hhh4 models, i.e., without normalization we have $w_{ji} = o_{ji}^{-d}$, where $o_{ji}$ is the order of neighbourhood between regions $i$ and $j$, see nbOrder.

Usage

powerlaw(maxlag, normalize = TRUE, log = FALSE, initial = if (log) 0 else 1)

Arguments

maxlag

single integer specifying up to which order of neighbourhood the weights should be positive; they will be set to 0 for higher orders. If spatial dependence is not to be truncated at some high order, maxlag should be set to the max

normalize

logical indicating if the weights should be normalized such that the rows of the weight matrix sum to 1.

log

logical indicating if the decay parameter should be estimated on the log-scale to ensure positivity.

initial

scalar initial value for the (log-)decay parameter.

Value

• a list which can be passed as a specification of parametric neighbourhood weights in hhh4.

References

Meyer, S. and Held, L. (2013): Modelling power-law spread of infectious diseases. arXiv:1308.5115 arXiv-Link: http://arxiv.org/abs/1308.5115

See Also

nbOrder to determine the matrix of neighbourhood orders from a binary adjacency matrix. siaf.powerlaw for modelling power-law distance decay in twinstim space-time point process models.

Examples

if (requireNamespace("spdep")) {
  data("measles.weser")
  measles <- disProg2sts(measles.weser)

  ## parametric powerlaw weights require neighbourhood orders in the data
  neighbourhood(measles) <- nbOrder(neighbourhood(measles), maxlag=Inf)

  ## a simple hhh4 model with power-law decay of spatial interaction
  m <- list(ar = list(f = ~ 1),
            ne = list(f = ~ 1, weights = powerlaw(maxlag=4, log=FALSE)),
            end = list(f = addSeason2formula(~-1 + ri(), S=1, period=52),
                       offset = population(measles)),
            family = "NegBin1", verbose=TRUE)

  ## fit the model
  set.seed(1)  # random intercepts are initialized randomly
  fit <- hhh4(measles, m)
  summary(fit)  # "neweights.d" is the decay parameter
}