Specify a first-order random walk temporal random effect.
RW1 assumes that the difference between consecutive time points is
normally distributed: phi[t] - phi[t-1] ~ N(0, sigma^2).
This is an intrinsic model (improper prior) that captures smooth temporal trends. It is equivalent to an ICAR model on a 1D chain.
Usage
temporal_rw1(
time_var,
group_var = NULL,
cyclic = FALSE,
shared = TRUE,
parameterization = c("centered", "noncentered")
)Arguments
- time_var
Name of the time variable in data. Can be a formula (e.g.,
~ year) or a character string (e.g.,"year").- group_var
Optional name of grouping variable for panel data. If provided, separate random walks are estimated for each group. Each group's walk sums to zero over time, so it carries the group's trend and not its level; add a group intercept such as
(1 | site)for the levels.- cyclic
Logical; if TRUE, the random walk wraps around (first and last time points are neighbors). Useful for seasonal data.
Logical; if TRUE (default), temporal effect enters both numerator and denominator linear predictors identically.
- parameterization
Coordinate the walk is sampled in.
"centered"(default) samples the effects themselves;"noncentered"samplesz ~ N(0, I)and reaches the effects throughsigma * A^-1 z, whereAstacks the walk's difference operator on1'/sqrt(n), so the variance enters the linear predictor through the transform rather than through a quadratic form. The two describe the same model and differ in the geometry the sampler sees: non-centred helps where the variance is shrunk toward zero and leaves a funnel the centred coordinate cannot descend, and can be worse where the data pin the field. Measured on a 500-observation fixture with no temporal signal, an RW1's bulk ESS fortau_temporalgoes 409 to 805, while an RW2 picks up 23 divergences of 2000 for no gain, which is why this is an argument rather than a default. Withgroup_varthe prior still augments one global constant, so the transform carries the group contrasts as free directions beside the scaled ones and the group levels stay unshrunk, as they are in the centred coordinate.
Details
The RW1 precision matrix Q has the form (for T time points): $$Q[t,t] = 1 \text{ for } t = 1, T$$ $$Q[t,t] = 2 \text{ for } 1 < t < T$$ $$Q[t,t-1] = Q[t-1,t] = -1$$
This is a rank-deficient matrix (rank T-1), so a sum-to-zero constraint is applied for identifiability.
Examples
# Create temporal RW1 specification
temporal_rw1("year")
#> tulpaRatio temporal specification
#> ============================
#>
#> Type: RW1 (First-order Random Walk)
#> Time variable: year
#> Shared: Yes (enters both processes)
temporal_rw1("month", cyclic = TRUE)
#> tulpaRatio temporal specification
#> ============================
#>
#> Type: RW1 (First-order Random Walk)
#> Time variable: month
#> Cyclic: Yes
#> Shared: Yes (enters both processes)
# \donttest{
# Simple temporal trend: a random walk shared by count and effort
set.seed(123)
year <- rep(2010:2021, each = 4)
trend <- cumsum(rnorm(12, sd = 0.25))[year - 2009]
x <- rnorm(48)
df <- data.frame(
year = year,
x = x,
count = rpois(48, exp(3 + 0.3 * x + trend)),
effort = rgamma(48, shape = 5, rate = exp(-trend))
)
fit <- tratio(
count | effort ~ x,
data = df,
family = ratiod_poisson_gamma(),
temporal = temporal_rw1("year"),
mode = "hmc",
control = list(iter = 400, chains = 1, seed = 1)
)
# Panel data: separate trends per site around each site's level
set.seed(124)
site <- rep(paste0("site", 1:5), times = 12)
s <- match(site, paste0("site", 1:5))
t <- rep(1:12, each = 5)
walks <- apply(matrix(rnorm(12 * 5, sd = 0.2), 12, 5), 2, cumsum)
eff <- rnorm(5, sd = 0.5)[s] + walks[cbind(t, s)]
x <- rnorm(60)
df_panel <- data.frame(
year = 2014 + t,
site = site,
x = x,
count = rpois(60, exp(2.5 + 0.3 * x + eff)),
effort = rgamma(60, shape = 3, rate = exp(-eff))
)
# The per-site walks carry the trends, (1 | site) the site levels
fit2 <- tratio(
count | effort ~ x + (1 | site),
data = df_panel,
family = ratiod_poisson_gamma(),
temporal = temporal_rw1("year", group_var = "site"),
mode = "hmc",
control = list(iter = 400, chains = 1, seed = 1)
)
# Cyclic for monthly seasonal data
set.seed(125)
month <- rep(1:12, times = 4)
season <- 0.5 * sin(2 * pi * month / 12)
x <- rnorm(48)
df_monthly <- data.frame(
month = month,
x = x,
count = rpois(48, exp(2.5 + 0.3 * x + season)),
effort = rgamma(48, shape = 2, rate = exp(-season))
)
fit3 <- tratio(
count | effort ~ x,
data = df_monthly,
family = ratiod_poisson_gamma(),
temporal = temporal_rw1("month", cyclic = TRUE),
mode = "hmc",
control = list(iter = 400, chains = 1, seed = 1)
)
# }