Fits a distributed lag nonlinear model (DLNM) from an epidemic-level design matrix and stores the standardized EpiExposure metadata required by downstream prediction, effect summarization, scenario simulation, and model comparison functions.
Usage
fit_epidlnm(
data,
model_engine,
family,
random_effect = NULL,
random_effect_prior = NULL,
epiexposure_spec = NULL,
basis_objects = NULL,
spatial_effect = NULL,
spatial_structure = "matern",
spatial_group = NULL,
...
)Arguments
- data
A data.frame containing `y_model` and the fitted DLNM design columns. For all engines except `bdlnm`, cross-basis columns must follow the package convention `cb_<variable>_<index>`. For `bdlnm`, epidemic-level, matrix-form cross-basis objects are supplied through `basis_objects` instead. Data prepared by `prepare_response()` may carry the `response_family_name` attribute; when present, it must agree with `family`.
- model_engine
Character scalar identifying the modeling engine. One of `"glm"`, `"glmmTMB"`, `"gam"`, `"gamm"`, `"gls"`, `"spamm"`, `"brms"`, `"inla"`, or `"bdlnm"`.
- family
Response distribution supplied as a supported canonical name or as a family object from which a supported canonical family and link can be identified. The supported EpiExposure v1 families are `"beta"`, `"binomial"`, `"poisson"`, `"gamma"`, `"gaussian"`, and `"negative_binomial"`. The canonical negative-binomial model is the quadratic-variance NB2 parameterization. NB1 and ordinal models are not supported in this version and produce explicit errors.
- random_effect
Optional character scalar naming one grouping column. In EpiExposure v1, this argument represents a conventional random intercept only. Random slopes, nested or crossed random-effect specifications, and arbitrary engine-specific random-effect expressions are not supported.
For the Bayesian `brms` engine, the random intercept is fitted as `(1 | group)`. For the INLA-backed `inla` and `bdlnm` engines, it is fitted as `INLA::f(group, model = "iid")`. EpiExposure v1 does not expose other INLA latent structures through this argument, including `"rw1"`, `"rw2"`, `"ar1"`, `"besag"`, `"bym"`, `"bym2"`, or SPDE-based effects.
INLA-backed fits internally re-index arbitrary grouping labels to consecutive integers while retaining the original column name in EpiExposure metadata. `glm` and `gls` do not support this argument. `gamm` requires it in the current EpiExposure interface because its generated fixed formula contains no other smooth or random term.
- random_effect_prior
Optional named list defining the hyperprior for the precision of the random intercept when `model_engine = "inla"` or `"bdlnm"`. The list is passed unchanged to the `hyper` argument of the internally generated `INLA::f(..., model = "iid")` term.
For example, `list(prec = list(prior = "pc.prec", param = c(1, 0.01)))`. `NULL` retains the engine's default prior. A non-`NULL` value requires `random_effect` and is currently supported only for the INLA-backed engines.
This argument customizes the hyperprior of the supported IID random intercept; it does not select or configure other INLA latent models. Structures such as `"rw1"`, `"rw2"`, `"ar1"`, `"besag"`, `"bym"`, `"bym2"`, and SPDE effects are not supported by the EpiExposure v1 random-effect interface.
- epiexposure_spec
Named list describing the cross-basis construction for every fitted exposure. This metadata is required so downstream functions reconstruct exactly the fitted exposure-lag basis rather than re-inferring `df`, knots, or basis functions. Each element must contain `max_lag`, `argvar`, and `arglag`.
Under the EpiExposure exact-profile contract, each `max_lag` must be one non-negative integer and all fitted exposures must use the same `max_lag`. Consequently, every original exposure profile used to build the epidemic-level design must contain exactly `max_lag + 1` observations for every fitted exposure. `fit_epidlnm()` receives the already collapsed epidemic-level design and therefore validates the common fitted `max_lag`; exact original profile length is validated upstream by the exposure/design construction functions.
- basis_objects
Optional named list of `dlnm::crossbasis()` objects. For `model_engine = "bdlnm"`, these must be the epidemic-level, matrix-form cross-basis objects containing exactly one row per model row. They are normally recovered from the `epiexposure_bdlnm_basis_objects` attribute produced by `build_design()`, but may also be supplied explicitly. For other engines, when supplied, their names must match the fitted exposure variables.
- spatial_effect
NULL (default) or two distinct names of numeric, finite coordinate columns in the epidemic-level design. Only `model_engine = "spamm"` supports spatial autocorrelation. Coordinates may repeat across different epidemics; at least two distinct pairs are needed.
- spatial_structure
Spatial correlation structure; currently only `"matern"` (case-insensitive) when `spatial_effect` is supplied. This argument does not change non-spatial fits.
- spatial_group
NULL for one shared spatial field, or the name of one grouping column (e.g., `"year"`) for independent Matérn realizations with shared correlation parameters. Requires `spatial_effect`; the source column may be factor, character, or integer and is converted in a local fit copy. `random_effect` remains a separate, conventional random intercept and both effects can be fitted together. Coordinates and spatial groups are not used as epidemic-profile identifiers.
- ...
Named additional arguments passed to the selected engine. Core arguments managed by EpiExposure (`formula`/`model`, `data`, `family`, and the engine-specific random-effect argument) cannot be supplied again in `...`. For INLA-backed fits (`inla` and `bdlnm`), `control.compute$config = TRUE` is required by the EpiExposure uncertainty contract and is added automatically when absent. If a conflicting value is supplied, the function stops.
Value
A fitted model object with a strict metadata contract. The following attributes are attached for downstream functions: `epiexposure_family_name`, `epiexposure_link`, `epiexposure_engine`, `epiexposure_cb_cols`, `epiexposure_vars`, `epiexposure_data_template`, `epiexposure_id_col`, `epiexposure_spec`, and `epiexposure_basis_objects`. Additional attributes record the standardized family parameterization, random-intercept structure, the INLA latent random-effect model in `epiexposure_random_effect_model`, the optional `epiexposure_random_effect_prior`, and the spatial specification (`epiexposure_spatial_effect`, `epiexposure_spatial_structure`, `epiexposure_spatial_group`, and `epiexposure_spatial_term`), common fitted `max_lag`, expected profile length (`max_lag + 1`), the exact-profile contract, the validated temporal step in `epiexposure_time_step`, and the EpiExposure v1 prediction contract.
Details
## Engine harmonization
`fit_epidlnm()` standardizes the *statistical target* across engines without requiring every engine to use the same computational implementation. The model is fitted using the selected engine, while downstream EpiExposure v1 predictions are defined as population-level expected responses. Random effects may therefore affect estimation during model fitting, but downstream predictions exclude group-specific random-effect contributions.
The random-intercept translation used at fitting is:
* `glmmTMB`, `brms`, and `spaMM`: `(1 | group)`; * `gam`: `s(group, bs = "re")`; * `gamm`: `random = list(group = ~1)`; * `INLA` and `bdlnm`: `f(group, model = "iid")`, optionally extended to `f(group, model = "iid", hyper = random_effect_prior)` when a custom random-effect hyperprior is supplied.
## Random-effect scope in EpiExposure v1
The EpiExposure v1 random-effect interface is intentionally restricted to one conventional random intercept. For `brms`, this corresponds to `(1 | group)`; for `inla` and `bdlnm`, it corresponds to `f(group, model = "iid")`.
The interface does not currently construct structured Bayesian latent effects such as first- or second-order random walks (`"rw1"` or `"rw2"`), autoregressive effects (`"ar1"`), areal spatial effects (`"besag"`, `"bym"`, or `"bym2"`), or continuous spatial effects based on SPDE models. These structures may be available in the underlying INLA framework but are outside the harmonized fitting, prediction, validation, and uncertainty contract of EpiExposure v1.
The `spatial_effect` interface described below is separate from these Bayesian latent structures and currently supports Matérn spatial covariance through `spaMM` only.
`glm` and `gls` are fixed-effect engines in this interface and reject a non-`NULL` `random_effect`. The `mgcv::gamm()` implementation uses the engine's native `random` argument rather than lme4-style syntax. Because `gamm()` uses PQL for non-Gaussian responses and is specifically known to be problematic for binary data, a warning is issued for the binomial family.
## Spatial covariance (spaMM only)
Spatial terms are explicitly distinct from `random_effect`: the former model correlated Matérn fields, whereas the latter is an ordinary intercept. `spatial_effect = c("x_coord", "y_coord")` adds `Matern(1 | x_coord + y_coord)`; adding `spatial_group = "year"` instead uses `Matern(1 | x_coord + y_coord leaves every pre-existing engine and random-intercept fit unchanged. Coordinates and spatial grouping affect estimation, but not EpiExposure's default fixed-component prediction, which sets all non-fixed contributions to zero rather than integrating them over their distributions.
## Families and links
Canonical families are translated to engine-native representations. The canonical `negative_binomial` family is standardized to an NB2-type model (`Var(Y) = mu + mu^2 / shape`) where supported. An explicit NB1 request is rejected instead of being silently converted.
Automatic family availability is engine-specific:
* `glm`: Gaussian, Binomial, Poisson, Gamma; * `glmmTMB`: all six EpiExposure v1 families; * `gam`: all six EpiExposure v1 families; * `gamm`: Gaussian, Binomial, Poisson, Gamma; * `gls`: Gaussian with identity link only; * `spaMM`: all six EpiExposure v1 families; * `brms`: all six EpiExposure v1 families; * `INLA`: all six EpiExposure v1 families, subject to the installed INLA likelihood implementation; * `bdlnm`: all six EpiExposure v1 families through its INLA backend.
When a family object supplies a non-default link, EpiExposure records that exact link in `epiexposure_link` and reconstructs the corresponding engine-native family. For INLA-based engines the same link is forwarded via `control.family$control.link$model`. Unsupported family-link combinations fail explicitly in EpiExposure or in the selected engine; no Gaussian or identity-link fallback is used.
To specify a non-default link, supply a supported family object, such as family = stats::Gamma(link = "inverse"); character family names use the EpiExposure default link.
For `gam`, Beta and negative-binomial models use `mgcv` extended families. When needed, `method = "REML"` is supplied by default to satisfy the supported fitting route; an explicitly supplied incompatible method is rejected. A GAM random intercept also defaults to REML unless the user explicitly selects another method.
For `model_engine = "gamm"`, the returned object contains the native `gam` and `lme` components and inherits from `"epiexposure_gamm"`. Calling `summary()` returns the fixed/population GAM summary by default; use `component = "lme"` or `"both"` to inspect the mixed-model component.
## Response-scale validation
To avoid engine-specific silent coercion, `y_model` is checked against the canonical family before fitting: Beta responses must lie strictly inside `(0, 1)`; Binomial responses must be coded `0/1`; Poisson and negative- binomial responses must be non-negative integer counts; Gamma responses must be strictly positive; Gaussian responses need only be finite numeric values.
## Common lag and exact-profile contract
EpiExposure models use one common retrospective lag window for every fitted exposure. If the common maximum lag is `L`, all original exposure profiles must contain exactly
$$ L + 1 $$
equally spaced observations per epidemiological unit and per exposure. Exposures with different fitted `max_lag` values are rejected here.
Because `fit_epidlnm()` is called after `build_design()` has collapsed each complete profile to one epidemic-level model row, the original long-format row count is no longer available at this stage. Exact `L + 1` temporal coverage must therefore be enforced by `define_exposures()` and `build_design()`. This function records the common lag/profile contract in the fitted model metadata so downstream functions can enforce it.
## Metadata and downstream prediction contract
The fitted exposure specification is part of the model contract, not an optional hint. Downstream functions must use the stored family, link, engine, cross-basis mapping, exposure specification, and data template. They should not infer a link from the family, rebuild a cross-basis from user-supplied degrees of freedom, align coefficients by position, or substitute a Gaussian/identity model when metadata are missing.
EpiExposure v1 prediction functions are intended to return the expected outcome under the fitted model, not a newly simulated observation. They use population-level predictions: random effects are included in model fitting when requested here, but group-specific random effects are not included in downstream predictions. This function itself only fits the model and does not create outcome predictions or uncertainty summaries.
For `glm` and `glmmTMB`, models are fitted with the fully evaluated EpiExposure design, but their stored calls are compacted after fitting. The original `data` expression is retained for display instead of the complete evaluated data frame and its EpiExposure attributes. For `glm`, the stored function and family expressions are also compacted to prevent `summary()` from printing the evaluated function definition and complete family object. These changes affect only the stored calls used for display and re-evaluation; they do not alter fitted coefficients, likelihoods, covariance matrices, predictions, or EpiExposure metadata.
Bayesian engines may generate posterior coefficient draws as part of their native fitting procedure (notably `bdlnm`). Those stored coefficient draws are model output and are available for later draw-by-draw uncertainty propagation; they are not used here as a single deterministic prediction.
Ordinal outcomes are deliberately excluded from EpiExposure v1 so fitting, prediction, uncertainty, performance metrics, and ensemble behavior remain harmonized across the supported response families.
Examples
# Minimal epidemic-level design with one exposure basis column
dat <- data.frame(
y_model = c(0.12, 0.18, 0.25, 0.31, 0.38, 0.44, 0.51, 0.57),
cb_rain_1 = c(-1.5, -1.0, -0.5, 0, 0.5, 1.0, 1.5, 2.0)
)
# Cross-basis specification corresponding to the fitted exposure
spec <- list(
rain = list(
max_lag = 0L,
argvar = list(fun = "lin"),
arglag = list(fun = "integer")
)
)
fit <- fit_epidlnm(
data = dat,
model_engine = "glm",
family = "gaussian",
epiexposure_spec = spec
)
stats::coef(fit)
#> (Intercept) cb_rain_1
#> 0.3126190 0.1295238
attr(fit, "epiexposure_engine")
#> [1] "glm"
attr(fit, "epiexposure_family_name")
#> [1] "gaussian"
attr(fit, "epiexposure_vars")
#> [1] "rain"