Predicting response correlations in multivariate models

[mvuorre]

Currently, brms syntax allows predicting the residual SD sigma:

bf(y ~ x + (x|id), sigma ~ x)

And multivariate response models:

bf(cbind(y1, y2) ~ x + (x|id))

But multiple response models can't include predictors for sigma, nor is it possible to model the response correlation on x (rescor(y1, y2)).

Something like this would be easy to use (if it's easy / possible to implement, that's another issue!):

bf(cbind(y1, y2) ~ x + (x|id), sigma ~ x)

where sigma could include residuals for y1, y2, and the response correlation rescor(y1, y2). Alternatively, the syntax could require / allow users to specify predictors on these three independently (perhaps that isn't too important?).

Thoughts? This is probably related to issue #3.

[paul-buerkner]

I agree that this is desireble and it will certainly be possible once multivariate models are generally implemented.

Implementing this separately from #3 would likely require too much special case coding, which I try to avoid in order to keep the package maintainable.

[paul-buerkner]

In the new dev version, sigma can now be predicted in multivariate models as well. However, this is not (yet) possible with rescor. The reason is that rescor is a correlation matrix with more than one unique entry in models with more than two response variables. It remains unclear to me, how to predict a correlation matrix or parts of it. Any suggestions are welcome!

[jeremy-koster]

In the new dev version, sigma can now be predicted in multivariate models as well.

Is there a vignette that includes the syntax for predicting sigma in multivariate models? Thanks.

[paul-buerkner]

vignette("brms_distreg") explains the syntax. See also https://vuorre.netlify.com/post/2017/how-to-compare-two-groups-with-robust-bayesian-estimation-using-r-stan-and-brms/

[jeremy-koster]

So to confirm, code for a bivariate response model in which residual variances are expected to vary by gender would look something like the following:

	bf_y1 <- bf (y1 ~ gender) + lf(sigma ~ gender)
	bf_y2 <- bf (y2 ~ gender) + lf(sigma ~ gender)
	m1 <- brm ( bf_y1 + bf_y2, data = dat, chains = 2, cores = 2)

Is that right?

And then all reported sigmas are on the log scale, correct?

[paul-buerkner]

looks correct to me.

suldauna [email protected] schrieb am Di., 17. Apr. 2018, 20:51:

So to confirm, code for a bivariate response model in which residual variances are expected to vary by gender would look something like the following:

bf_y1 <- bf (y1 ~ gender) + lf(sigma ~ gender) bf_y2 <- bf (y2 ~ gender) + lf(sigma ~ gender) m1 <- brm ( bf_y1 + bf_y2, data = dat, chains = 2, cores = 2)

Is that right?

And then all reported sigmas are on the log scale, correct?

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/paul-buerkner/brms/issues/254#issuecomment-382102283, or mute the thread https://github.com/notifications/unsubscribe-auth/AMVtAGtR41dOZGNA4PIM1zJ28xQl6FkWks5tpjnIgaJpZM4O6ePs .

[jeremy-koster]

Regarding the above question about predicting a correlation matrix, perhaps that would be less problematic in a bivariate probit model in which the respective sigmas are fixed at 1 and only the covariance is allowed to vary freely.

[paul-buerkner]

It would be easier in any bivariate model, where we only have one correlation. Still, I am looking for a general solution in brms and I have yet to find a good way of doing it.

[stephensrmmartin]

@paul-buerkner This isn't for "predicting" correlations in a matrix, per se, but if partial pooling of a correlation matrix is desired across groups - Ben goodrich suggested a method to me that worked fairly well. See thread here: http://discourse.mc-stan.org/t/hierarchical-prior-for-partial-pooling-on-correlation-matrices/4852

[paul-buerkner]

Thanks Stephen, I will take a look!

[jackobailey]

@paul-buerkner, this is just an idle thought, but could you model the correlation matrix by assuming that it comes from an LKJ distribution (like we use as a prior anyway) and then modelling the LKJ distribution's scalar parameter using a log link (to keep it positive)?

[paul-buerkner]

Yes, but that would prevent us from having a separate model for each correlation and it would only control the scale not the location of the prior.

[jackobailey]

I might be missing something, but in that case why not simply create a vector of correlation parameters the same length as the number of comparisons between response variables and use them to fill in the correlation matrix? Then you could let users model, say, cor_a_b or cor_b_c, etc.

[paul-buerkner]

Then we don't ensure that the created correlation matrix is actually postive definite, a prerequisite for it being an actually valid correlation matrix.

[jeremy-koster]

Some co-authors and I currently are working on implementing unconventional correlation matrices that fulfill the positive definite requirement, but it's a special case written in Stan that I'm not sure would generalize easily to the brm() function.

[paul-buerkner]

Great! Once you have something you can share, I am happy to take a look and see if we can generalize it to the needs of brms (or other regression frameworks).

[bwiernik]

Here is an approach to modeling correlation matrices that comes from the meta-analysis literature https://onlinelibrary.wiley.com/doi/full/10.1002/jrsm.1214

[shannonajw]

@paul-buerkner This might be useful for a way to model correlation matrices: https://gitlab.com/dschrage/rcovreg