I have a dataset is ~ 600 participants each measured on four dichotomous outcomes. These are meaningfully clinically (e.g. death), and it is expected that outcomes will be correlated. The purpose of the model is ultimately to inform a patient level simulation for a health economic model, so capturing this correlation will be important. From reading this forum and a quick lit scan (e.g. https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6290916/), I was thinking a multivariate GLM linked by patient level random effects might be the solution modeled for ex as
o1 <- bf(out1 ~ (1 | p | id))
o2 <- bf(out2 ~ (1 | p | id))
m1 <- brm(o1 + o2...)
This works in the sense that I end up with correlated predictions, but I am finding standard errors to be about 2-10 times as large as the univariate models. Using posterior_predict and allow_new_levels on a new data set the correlated model gives me average rates across the sample with 95% percent intervals from 0.04 to 0.72. The same interval for the univariate model is 0.17 - 0.27.
I feel like the issue here is that most of these mGLMs are applied to longitudinal scenarios and that it may not translate to one observation per patient but thought Iād give a quick post here before giving up.
To clarify, by univariate model, you mean brm(o1, data = your_data) ?
Does that happen already for two outcomes? Or do the wide posteriors appear when you use all the four outcomes? One possible thing is that the correlation is not actually common across all the outcomes and hence the model has trouble fitting your data, resulting in wide posteriors.
Finally, I have doubts whether the y ~ 1 | id effect is at least theoretically identifiable for logistic regression. But I would refer to @paul.buerkner for confirmation.
Do you have a reason why having a more compact model, such as outcome_value ~ outcome_type + (1 || id) would be problematic?
1 Like
If Id is an identifier of individual observations then this will probably not work.
1 Like
This is what I was afraid of. Is there a solution to link these distributions in this case?
There may be some in the literature but I dont know them out of mind right now.
I was able to find a SAS example here (SAS/STAT(R) 14.1 User's Guide) which shows how patient level random effect is used to model dependence between a binary and poisson variable each measured once per patient. This is what I was hoping my original code would mimic. Would it be possible to fit this style of model in brms? I have similar binary/poisson dataset for another analysis.
Since my original analysis is just multiple binomial outcomes I have converted to long format and taken this approach. Thank you for the suggestion!
1 Like