egf {epigrowthfit} | R Documentation |
Fit Nonlinear Mixed Effects Models of Epidemic Growth
Description
Fits nonlinear mixed effects models of epidemic growth to collections of one or more disease incidence time series.
Usage
egf(model, ...)
## S3 method for class 'egf_model'
egf(model,
formula_ts,
formula_windows,
formula_parameters = list(),
formula_priors = list(),
data_ts,
data_windows,
subset_ts = NULL,
subset_windows = NULL,
select_windows = NULL,
na_action_ts = c("fail", "pass"),
na_action_windows = c("fail", "omit"),
control = egf_control(),
init = list(),
map = list(),
fit = TRUE,
se = FALSE,
...)
Arguments
model |
an R object specifying a top level nonlinear model,
typically of class |
formula_ts |
a formula of the form |
formula_windows |
a formula of the form |
formula_parameters |
a list of formulae of the form |
formula_priors |
a list of formulae of the form |
data_ts , data_windows |
data frames, lists, or environments to be searched for variables
named in the corresponding |
subset_ts , subset_windows |
expressions to be evaluated in the corresponding |
select_windows |
an expression indicating additional variables in |
na_action_ts |
a character string affecting the handling of |
na_action_windows |
a character string affecting the handling of |
control |
an |
init |
a named list of numeric vectors with possible elements |
map |
a named list of factors with possible elements |
fit |
a logical. If |
se |
a logical. If |
... |
additional arguments passed from or to other methods. |
Details
Users attempting to set arguments formula_priors
, init
,
and map
should know the structure of the bottom level parameter
vector. It is described under topic egf-class
.
If
formula_ts = cbind(time, x) ~ ts1 formula_windows = cbind(start, end) ~ ts2
then it is expected that time
, start
, and end
(after coercion to numeric) measure time on the same scale.
To be precise, numeric times should have a common unit of measure
and, at least within time series, represent displacements from a
common reference time.
These conditions will always hold if time
, start
, and
end
all evaluate to Date
, POSIXct
,
or POSIXlt
vectors.
When day of week effects are estimated, numeric times are interpreted
as numbers of days since midnight on January 1, 1970, so that time
points can be mapped unambiguously to days of week.
Furthermore, in this case, time
(after coercion to numeric) is
required to be integer-valued with one day spacing in all time series.
This means that
isTRUE(all.equal(time, round(time))) && all(range(diff(round(time))) == 1)
must be TRUE
in each level of ts1
.
These conditions ensure that intervals between successive time points
represent exactly one day of week.
Value
A list inheriting from class egf
.
See topic egf-class
for class documentation.
See Also
The many methods for class egf
,
listed by methods(class = "egf")
.
Examples
## Simulate 'N' incidence time series exhibiting exponential growth
set.seed(180149L)
N <- 10L
f <- function(time, r, c0) {
lambda <- diff(exp(log(c0) + r * time))
c(NA, rpois(lambda, lambda))
}
time <- seq.int(0, 40, 1)
r <- rlnorm(N, -3.2, 0.2)
c0 <- rlnorm(N, 6, 0.2)
data_ts <-
data.frame(country = gl(N, length(time), labels = LETTERS[1:N]),
time = rep.int(time, N),
x = unlist(Map(f, time = list(time), r = r, c0 = c0)))
rm(f, time)
## Define fitting windows (here, two per time series)
data_windows <-
data.frame(country = gl(N, 1L, 2L * N, labels = LETTERS[1:N]),
wave = gl(2L, 10L),
start = c(sample(seq.int(0, 5, 1), N, TRUE),
sample(seq.int(20, 25, 1), N, TRUE)),
end = c(sample(seq.int(15, 20, 1), N, TRUE),
sample(seq.int(35, 40, 1), N, TRUE)))
## Estimate the generative model
m1 <-
egf(model = egf_model(curve = "exponential", family = "pois"),
formula_ts = cbind(time, x) ~ country,
formula_windows = cbind(start, end) ~ country,
formula_parameters = ~(1 | country:wave),
data_ts = data_ts,
data_windows = data_windows,
se = TRUE)
## Re-estimate the generative model with:
## * Gaussian prior on beta[1L]
## * LKJ prior on all random effect covariance matrices
## (here there happens to be just one)
## * initial value of 'theta' set explicitly
## * theta[3L] fixed at initial value
m2 <-
update(m1,
formula_priors = list(beta[1L] ~ Normal(mu = -3, sigma = 1),
Sigma ~ LKJ(eta = 2)),
init = list(theta = c(log(0.5), log(0.5), 0)),
map = list(theta = 3L))