2
$\begingroup$

For this example, I am using the data "appendix_example1_wide.SUPP.FINAL.csv" posted here.

In the paper, the authors use the to MCMCglmm package fit a multivariate multilevel model. Of particular interest, they use the following code to estimate covariances between random effects (e.g., random intercepts from each outcome):

prior1 <- list(R = list(V = .2*diag(2), nu = 3),
               G = list(G1 = list(V = diag(c(.2,.2,.1,.1)), nu = 5)))

mult.mcmc <- MCMCglmm(
  cbind(y1, y2) ~ -1 + trait + trait:tx + trait:time + 
  trait:time:tx,
  data = data1,
  random = ~ us(-1 + trait + trait:time):pid,
  rcov = ~ us(trait):units,
  family = c("gaussian","gaussian"),
  nitt = 25000, burnin = 5000, 
  prior = prior1
)
summary(mult.mcmc)

Is it possible to specify the same model using the brms package? Here's my attempt (that did not even closely reproduce the estimates from the paper):

mod1 <- brm(cbind(y1,y2) ~ 1 + time*tx + (1 + time|pid),data=data)
$\endgroup$
1
$\begingroup$

Update: I was missing a small piece of code. This code does the trick by specifying the covariances between the random effects:

 mod1 <- brm(cbind(y1,y2) ~ 1 + time*tx + (1 + time|p|pid), 
        data=data)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.