Nothing
data("measlesWeserEms")
## a simple endemic model
measlesFit0 <- hhh4(measlesWeserEms, list(
end = list(f = addSeason2formula(~1), offset = population(measlesWeserEms)),
family = "NegBin1"
))
test_that("endemic-only model has zero-valued Lambda matrix", {
res <- getMaxEV_season(measlesFit0)
expect_equal(res$maxEV.const, 0)
zeromat <- matrix(0, measlesFit0$nUnit, measlesFit0$nUnit)
expect_equal(res$Lambda.const, zeromat)
expect_equal(surveillance:::createLambda(measlesFit0)(2), zeromat)
})
## + AR component
measlesFit1 <- update(measlesFit0, ar = list(f = addSeason2formula(~1)))
test_that("autoregressive model has a diagonal Lambda matrix", {
res <- getMaxEV_season(measlesFit1)
expect_equal(res$Lambda.const, diag(res$maxEV.const, measlesFit1$nUnit))
expect_equal(surveillance:::createLambda(measlesFit1)(2),
diag(res$maxEV.season[2], measlesFit1$nUnit))
})
## + NE component
measlesFit2 <- update(measlesFit1,
ne = list(f = ~1, weights = neighbourhood(measlesWeserEms) == 1)) # symmetric
measlesFit3 <- update(measlesFit2,
ne = list(normalize = TRUE)) # asymmetric
test_that("getMaxEV() and getMaxEV_season() agree", {
expect_equal(getMaxEV_season(measlesFit2)$maxEV.season,
getMaxEV(measlesFit2)[seq_len(measlesWeserEms@freq)])
expect_equal(getMaxEV_season(measlesFit3)$maxEV.season,
getMaxEV(measlesFit3)[seq_len(measlesWeserEms@freq)])
})
## AR within NE + unit-specific epidemic covariate
measlesFit4 <- update(measlesFit0,
ne = list(f = ~pop, weights = (neighbourhood(measlesWeserEms)+1)^-2,
normalize = TRUE),
data = list(pop = population(measlesWeserEms)))
## calculate "nu + Lambda Y_{t-1}" and compare to fitted(object)
check_createLambda <- function (object)
{
mname <- deparse(substitute(object))
model <- terms(object)
means <- meanHHH(object$coefficients, model, subset = seq_len(model$nTime))
expect_equal(means$mean[model$subset,,drop=FALSE], fitted(object),
expected.label = paste0("fitted(", mname, ")"))
Lambda <- surveillance:::createLambda(object)
if (any(object$lags != 1, na.rm = TRUE))
stop("check not implemented for lags != 1")
meansByLambda <- t(vapply(
X = object$control$subset,
FUN = function(t) means$endemic[t,] + Lambda(t) %*% model$response[t-1,],
FUN.VALUE = numeric(object$nUnit), USE.NAMES = FALSE))
expect_equal(meansByLambda, unname(fitted(object)),
expected.label = paste0("fitted(", mname, ")"))
}
test_that("multivariate formulation using Lambda agrees with fitted values", {
check_createLambda(measlesFit0)
check_createLambda(measlesFit1)
check_createLambda(measlesFit2)
check_createLambda(measlesFit3) # failed in surveillance < 1.13.1
check_createLambda(measlesFit4) # failed in surveillance < 1.13.1
})
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.