Skip to contents

Overview

This vignette covers random effects in numdenom, including random intercepts, random slopes, crossed effects, and nested effects. numdenom uses lme4-style formula syntax with extensions for correlated vs. uncorrelated slopes.

Random Intercepts

Single Grouping Factor

# Random intercept for site
fit <- tratio(
  count | effort ~ depth + (1 | site),
  data = df,
  family = ratiod_poisson_gamma()
)

This estimates a site-specific intercept bsite[i]b_{\text{site}[i]} for each observation ii:

log(μi)=β0+β1depthi+bsite[i] \log(\mu_i) = \beta_0 + \beta_1 \cdot \text{depth}_i + b_{\text{site}[i]}

where bjNormal(0,σsite2)b_j \sim \text{Normal}(0, \sigma_{\text{site}}^2).

Multiple (Crossed) Grouping Factors

When observations are classified by multiple independent grouping factors:

# Crossed random effects: site and year
fit <- tratio(
  count | effort ~ depth + (1 | site) + (1 | year),
  data = df,
  family = ratiod_poisson_gamma()
)

This estimates: - bjsiteNormal(0,σsite2)b^{\text{site}}_j \sim \text{Normal}(0, \sigma_{\text{site}}^2) - bkyearNormal(0,σyear2)b^{\text{year}}_k \sim \text{Normal}(0, \sigma_{\text{year}}^2)

The linear predictor becomes:

log(μi)=β0+β1depthi+bsite[i]site+byear[i]year \log(\mu_i) = \beta_0 + \beta_1 \cdot \text{depth}_i + b^{\text{site}}_{\text{site}[i]} + b^{\text{year}}_{\text{year}[i]}

Nested Random Effects

When groups are hierarchically nested (plots within sites):

# Nested: plots within sites
fit <- tratio(
  count | effort ~ depth + (1 | site/plot),
  data = df,
  family = ratiod_poisson_gamma()
)

This expands to:

# Equivalent explicit form
fit <- tratio(
  count | effort ~ depth + (1 | site) + (1 | site:plot),
  data = df,
  family = ratiod_poisson_gamma()
)

The site:plot interaction creates unique identifiers for each plot-within-site combination.

Deep Nesting

For deeper hierarchies:

# region/site/plot expands to:
# (1 | region) + (1 | region:site) + (1 | region:site:plot)
fit <- tratio(
  count | effort ~ x + (1 | region/site/plot),
  data = df,
  family = ratiod_poisson_gamma()
)

Random Slopes

Random slopes allow the effect of a predictor to vary by group.

Correlated Intercept and Slope

Use single | for correlated random effects:

# Intercept and depth effect vary by site, with correlation
fit <- tratio(
  count | effort ~ depth + (1 + depth | site),
  data = df,
  family = ratiod_poisson_gamma()
)

This estimates:

(b0jb1j)MVN((00),(σ02ρσ0σ1ρσ0σ1σ12)) \begin{pmatrix} b_{0j} \\ b_{1j} \end{pmatrix} \sim \text{MVN}\left(\begin{pmatrix} 0 \\ 0 \end{pmatrix}, \begin{pmatrix} \sigma_0^2 & \rho\sigma_0\sigma_1 \\ \rho\sigma_0\sigma_1 & \sigma_1^2 \end{pmatrix}\right)

The correlation ρ\rho captures whether sites with higher intercepts also tend to have steeper slopes.

Uncorrelated Intercept and Slope

Use double || for independent random effects:

# Intercept and depth effect vary by site, independently
fit <- tratio(
  count | effort ~ depth + (1 + depth || site),
  data = df,
  family = ratiod_poisson_gamma()
)

This is computationally simpler and appropriate when you don’t expect correlation between intercepts and slopes.

Multiple Slopes

# Multiple correlated slopes
fit <- tratio(
  count | effort ~ depth + temp + (1 + depth + temp | site),
  data = df,
  family = ratiod_poisson_gamma()
)

# Multiple uncorrelated slopes
fit <- tratio(
  count | effort ~ depth + temp + (1 + depth + temp || site),
  data = df,
  family = ratiod_poisson_gamma()
)

Formula Expansion

numdenom supports lme4-style formula expansion for random slopes.

Interactions

# x*z expands to x + z + x:z
fit <- tratio(
  count | effort ~ x*z + (1 + x*z | site),
  data = df,
  family = ratiod_poisson_gamma()
)
# Random slopes for x, z, and their interaction

Polynomials

# Quadratic random slope
fit <- tratio(
  count | effort ~ x + I(x^2) + (1 + x + I(x^2) | site),
  data = df,
  family = ratiod_poisson_gamma()
)

# Or using poly()
fit <- tratio(
  count | effort ~ poly(x, 2) + (1 + poly(x, 2) | site),
  data = df,
  family = ratiod_poisson_gamma()
)

Parameter Naming

numdenom cleans parameter names for readability:

Formula Term Parameter Name
x:z x_z
I(x^2) x_pow2
poly(x, 2)1 x_poly1
poly(x, 2)2 x_poly2
log(x) log_x
scale(x) x_scaled

Shared Random Effects

By default, numdenom shares random effects between numerator and denominator when the same grouping factor appears in both. This is the core of the numdenom philosophy (see vignette("philosophy")).

Default: Shared Structure

# Site effect shared between catch and effort processes
fit <- tratio(
  catch | effort ~ depth + (1 | site),
  data = df,
  family = ratiod_poisson_gamma()
)

The site effect enters both linear predictors identically:

log(μcatch)=β0+β1depth+𝐛𝐬𝐢𝐭𝐞 \log(\mu_{\text{catch}}) = \beta_0 + \beta_1 \cdot \text{depth} + \mathbf{b_{\text{site}}} log(μeffort)=α0+𝐛𝐬𝐢𝐭𝐞 \log(\mu_{\text{effort}}) = \alpha_0 + \mathbf{b_{\text{site}}}

Shared effects cancel in the ratio, removing spurious variation from common unmeasured drivers.

Explicit Shared Structure

# Only site is shared, observer is process-specific
fit <- tratio(
  catch | effort ~ (1 | site) + (1 | observer),
  shared = ~ (1 | site),
  data = df,
  family = ratiod_poisson_gamma()
)

Independent Structure (Use with Caution)

# Request independence (triggers warning)
fit <- tratio(
  successes | trials ~ treatment + (1 | block),
  shared = ~ 0,
  data = df,
  family = ratiod_binomial()
)

Extracting Random Effects

# Posterior summaries for random effects
re <- extract_ranef(fit)

# By term
re$site       # Site-level effects
re$site_plot  # Nested plot effects

# For random slopes
re$site$intercept  # Intercept component
re$site$depth      # Slope component

Prior Specification

Default PC Priors

numdenom uses penalized complexity (PC) priors that favor simpler models:

# Default: P(sigma > 1) = 0.01
fit <- tratio(..., priors = ratiod_priors())

Custom Priors

# Expect smaller random effects
fit <- tratio(
  ...,
  priors = ratiod_priors(
    sigma_re_scale = 0.5  # Prior scale for RE standard deviations
  )
)

Latent Factors for Unmeasured Confounders

When you suspect unmeasured confounders are driving both numerator and denominator but don’t have a natural grouping structure, latent factors provide an alternative to random effects.

Single Latent Factor

# Capture shared unmeasured variation
fit <- tratio(
  count | effort ~ depth + (1 | site),
  data = df,
  family = ratiod_poisson_gamma(),
  latent = latent_factor(n_factors = 1)
)

Latent factors are observation-level random effects that enter both linear predictors:

log(μcount)=β0+β1depth+bsite+fiσf \log(\mu_{\text{count}}) = \beta_0 + \beta_1 \cdot \text{depth} + b_{\text{site}} + f_i \cdot \sigma_f log(μeffort)=α0+bsite+fiσf \log(\mu_{\text{effort}}) = \alpha_0 + b_{\text{site}} + f_i \cdot \sigma_f

where fiNormal(0,1)f_i \sim \text{Normal}(0, 1) with sum-to-zero constraint.

Multiple Latent Factors

# Multiple independent unmeasured confounders
fit <- tratio(
  count | effort ~ depth,
  data = df,
  family = ratiod_poisson_gamma(),
  latent = latent_factor(n_factors = 2)
)

Latent Factors vs Random Effects

Feature Random Effects Latent Factors
Grouping Required (site, year, etc.) None needed
Level Group-level Observation-level
Use case Known grouping structure Residual confounding
Dimensionality 1 per group K factors × N observations

Use latent factors when: - You suspect unmeasured confounders but don’t know the grouping structure - Residual variation after random effects is still correlated - You want to capture observation-level shared structure

Prior Specification

# Stronger regularization on latent factors
fit <- tratio(
  count | effort ~ depth,
  data = df,
  family = ratiod_poisson_gamma(),
  latent = latent_factor(
    n_factors = 1,
    prior = prior_pc(U = 0.5, alpha = 0.01)  # P(sigma > 0.5) = 0.01
  )
)

See Also