Multivariate analysis (original) (raw)

Hi everybody,

I am new to brms, and bayesian statistics. I would like to dig a bit deeper though, and I tried applying it to my own data.

I am running a multivariate model. I have 3 metrics, which I follow over Time (factor), and which I want to compare between 2 Groups (factor).
The metrics are heavily correlated, and the distributions are definitely not normal (rather skewed). It is not an option to apply a transformation on the data.

I found in the online material how to run a multivariate model in brms (here: Estimating Multivariate Models with brms). However, this approach doesn’t support non normal underlying distributions.

Would it be an option to construct the formula like this?

Form1 ← bf(Response1 ~ 1 + Group * Time + (1|Participant)) + skew_normal()
Form2 ← bf(Response2 ~ 1 + Group * Time + (1|Participant)) + skew_normal()
Form3 ← bf(Response3 ~ 1 + Group * Time + (1|Participant)) + skew_normal()

ModelFormula ← Form1 + Form2 + Form3

But how do I account for the fact that the data are correlated? And how can I correct for multiple comparisons?

Thanks for your reply!
A

wpetry March 21, 2024, 4:54pm 2

This will fit a multivariate model, but it won’t include residual correlations between the responses (rescor = FALSE) because brms doesn’t have a multivariate skew-normal.

Depending on your data and goals, you might be able to exploit correlations in the effects of Participant by specifying (1 | p | Participant). The letter p doesn’t mean anything special, it’s just used as a grouping marker that let’s brms know which group-level effects should be correlated. I think the argument here might be that participants have consistent, correlated responses to 1-3 that mean that they differ in the set of intercepts their specific time series will have.

I’m not sure where the multiple comparisons issue would arise with your problem. Andrew Gelman has a blog post that summarizes thinking on multiple comparisons in Bayesian analyses.

Solomon March 22, 2024, 2:27pm 3

As a half-measure, you can allow correlations in your upper-level deviation terms with syntax like this:

Form1 <- bf(Response1 ~ 1 + Group * Time + (1 |p| Participant)) + skew_normal()
Form2 <- bf(Response2 ~ 1 + Group * Time + (1 |p| Participant)) + skew_normal()
Form3 <- bf(Response3 ~ 1 + Group * Time + (1 |p| Participant)) + skew_normal()

Notice the double || with some key in between. p does not need to be in your data, and you can use other characters like X or whatever. The trick is that by placing that p in all 3 of your formulas, you will allow the random intercepts for Response1, Response2, and Response3 to correlate. You can learn more about this in the brmsformula section of the brms reference manual (link).

Hi, @Solomon and @wpetry, thanks for your helpful advice!

I have a question in a similar vein, so replying here instead of opening a new thread.
I’d like to jointly estimate the survival and the rate of recurrent events (hospitalizations) for participants (needed for a while-alive estimand). They are likely negative correlated.

formula_hosp_negbinomial <- bf(
    n_re | rate(survival_time) ~ trt + ecog_fstcnt,
    family = negbinomial()
    )

    formula_death_weibull <- bf(
    survival_time | cens(censoring) ~ trt + ecog_fstcnt,
    family = weibull()
    )

Could I use brms to somehow account for this correlation beyond what is captured by the covariates? I don’t think a random intercept at the patient level can work, given that I only have 1 observation per patient.

Solomon April 30, 2025, 4:13pm 6

I’m not aware that brms can combine different likelihoods like this.

@Solomon thanks for the swift reply!

jroon May 1, 2025, 6:33am 8

This sounds like the kind of problem a multi-state survival model would be suited to. Alas while I know they exist, I’ve no experience with them myself, and I don’t know if they are implemented in Stan. There’s an explanation here: flexsurv: A Platform for Parametric Survival Modeling in R - PMC

Interesting. I know that Frank Harrell has developed packages for multi-state models which IIRC are based on stan, so I will look into that.

I also found this publication in which the authors address the correlation between survival and the recurrent event rate with a shared frailty term. The authors also provide their rstan code.

lldb14 May 2, 2025, 8:08pm 10

Who is the likelihood of brms?
What is the marginal correlation between responses for brms?
Who is the posterior of brms?
How does brms write when it has a random effect a posteriori?
How does brms induce correlation between responses?