Ensuring equivalent results between univariate vs. multivariate brms models

Modeling
1

Hi,

I’m comparing two approaches in brms: Fitting separate univariate models versus one multivariate model—without any residual or random-effect correlations (using set_rescor(FALSE)). My goal is to confirm that, aside from seed-related stochastic variability, all parameters estimates (e.g., in terms of location and scale) are equivalent across the two methods. This is crucial because I plan to use the multivariate approach to minimize runtime when dealing with many outcomes.

Below is a simple illustration:

library(brms)

set.seed(123)

x <- rnorm(100)
dt <- data.frame(x = x,
                 y1 = 1 + 2*x + rnorm(100),
                 y2 = 1 + 2*x + rnorm(100))

# Univariate models
fit_uni1 <- brm(y1 ~ x, data = DT)
fit_uni2 <- brm(y2 ~ x, data = DT)

# Multivariate model with no residual or random-effect correlations
fit_multi <- brm(bf(mvbind(y1, y2) ~ x) + set_rescor(FALSE), data = dt)

Any feedback or confirmation that these two approaches are equivalent would be much appreciated. Thanks!

2

Looks right, to me. Just consider setting the seed argument in brm() to make the results more reproducible. Otherwise, since your goal is to reduce run time, set cores = 4 to run all HMC chains in parallel, rather than sequentially, as is the default.

3

It looks like this is your first post. Welcome to the Stan forums, @Mat!