Getting Started with BayesianDEB

Overview

BayesianDEB provides a complete Bayesian framework for Dynamic Energy Budget (DEB) modelling. It wraps pre-written Stan models with a clean R interface for data preparation, model specification, fitting, diagnostics, and visualisation.

The package implements four model types:

Type Stan model Description
individual 2-state (E, V) Single individual growth
growth_repro 3-state (E, V, R) Growth + reproduction
hierarchical 2-state + random effects Multi-individual with partial pooling
debtox 4-state (E, V, R, D) Toxicokinetic-toxicodynamic

Installation

# Install cmdstanr (required backend)
install.packages("cmdstanr",
  repos = c("https://stan-dev.r-universe.dev", getOption("repos")))
cmdstanr::install_cmdstan()

# Install BayesianDEB
# remotes::install_github("sciom/BayesianDEB")

Example: Individual Growth Model

1. Prepare Data

library(BayesianDEB)
data(eisenia_growth)

# Use first individual
df1 <- eisenia_growth[eisenia_growth$id == 1, ]
dat <- bdeb_data(growth = df1, f_food = 1.0)
dat

2. Specify Model

mod <- bdeb_model(
  data   = dat,
  type   = "individual",
  priors = list(
    p_Am    = prior_lognormal(mu = 1.5, sigma = 0.5),
    p_M     = prior_lognormal(mu = -1.0, sigma = 0.5),
    kappa   = prior_beta(a = 3, b = 2),
    sigma_L = prior_halfnormal(sigma = 0.05)
  )
)
mod

Unspecified priors are filled from prior_default("individual").

3. Fit Model

fit <- bdeb_fit(mod, chains = 4, iter_sampling = 1000, seed = 123)
fit

4. Diagnostics

bdeb_diagnose(fit)
plot(fit, type = "trace")
plot(fit, type = "pairs", pars = c("p_Am", "p_M", "kappa"))

5. Posterior Predictive Checks

ppc <- bdeb_ppc(fit, type = "growth")
plot(ppc)

6. Derived Quantities

bdeb_derived(fit, quantities = c("L_inf", "growth_rate"))

7. Trajectory Plot

plot(fit, type = "trajectory", n_draws = 200)

Example: Hierarchical Model (Multiple Individuals)

dat_all <- bdeb_data(growth = eisenia_growth, f_food = 1.0)

mod_hier <- bdeb_model(dat_all, type = "hierarchical")

fit_hier <- bdeb_fit(mod_hier, chains = 4, adapt_delta = 0.9,
                     threads_per_chain = 4)  # within-chain parallelism

bdeb_diagnose(fit_hier)
bdeb_summary(fit_hier, pars = c("mu_log_p_Am", "sigma_log_p_Am", "p_M", "kappa"))

Example: DEBtox Model

data(debtox_growth)

# Concentration mapping
conc_map <- setNames(
  c(0, 20, 80, 200),
  c("1", "2", "3", "4")
)

dat_tox <- bdeb_data(
  growth = debtox_growth,
  concentration = conc_map,
  f_food = 1.0
)

mod_tox <- bdeb_tox(dat_tox, stress = "assimilation")
fit_tox <- bdeb_fit(mod_tox, chains = 4, adapt_delta = 0.95)

# EC50 and NEC
bdeb_ec50(fit_tox)

# Dose-response plot
plot_dose_response(fit_tox)

Prior Specification

BayesianDEB follows the prior recommendations of Kooijman (2010) and the AmP collection (Marques et al., 2018):

# View defaults
prior_default("individual")

# Override specific priors
my_priors <- list(
  p_Am  = prior_lognormal(mu = 2.0, sigma = 0.3),
  kappa = prior_beta(a = 5, b = 2)
)

Observation Models

The default likelihood is Gaussian for growth and negative binomial for reproduction. Switch via observation:

# Robust to outliers
mod <- bdeb_model(dat, type = "individual",
  observation = list(growth = obs_student_t(nu = 5)))

# Multiplicative error
mod <- bdeb_model(dat, type = "individual",
  observation = list(growth = obs_lognormal()))

Further Reading