Longitudinal Mediation

Hi,

I am trying to fit a longitudinal mediation model using brms. Below I attached a reproducible code

n<-50
J<-5
alpha_m<-4
alpha_y<-3
cp<-2
a<-1.75
b<-1.75
b1<-rnorm(n,0,1.5)
b2<-rnorm(n,0,1)
g1<-rnorm(n,0,1.5)
g2<-rnorm(n,0,1)
g3<-rnorm(n,0,1)
X<-matrix(0,n,J)
M<-matrix(0,n,J)
Y<-matrix(0,n,J)
for(i in 1:n)
{
X[i,]<-rnorm(J,0,1)
M[i,]<-alpha_m+b1[i]+(a+b2[i])*X[i,]+rnorm(J,0,1)
Y[i,]<-alpha_y+g1[i]+(cp+g2[i])*X[i,]+(b+g3[i])*M[i,]+rnorm(J,0,1)
}
X.vec1<-X[1,]
M.vec1<-M[1,]
Y.vec1<-Y[1,]
for(i in 2:n)
{
X.vec1<-c(X.vec1,X[i,])
M.vec1<-c(M.vec1,M[i,])
Y.vec1<-c(Y.vec1,Y[i,])
}
data_med<-data.frame(X.vec1,M.vec1,Y.vec1)
data_med$id<-rep(1:n,each=J)
y_mod ← bf(Y.vec1 ~ X.vec1 + M.vec1+(X.vec1+M.vec1|ID|id))
m_mod ← bf(M.vec1 ~ X.vec1+(X.vec1|ID|id))

k_fit_brms<-brm(y_mod+m_mod+set_rescor(FALSE),
                data=data_med,
                cores=2,chains=2)

Question 1.
If I want to get the posterior sample for the parameter Mvec1_X.vec1, I can call the following line of code

a <- posterior_samples(k_fit_brms, pars = "Mvec1_X.vec1")

How can I get the posterior sample for the variance-co-variance parameters, for example, for β€œcor(Yvec1_M.vec1,Mvec1_X.vec1)”

Question 2.
If I have more than one mediator coming from a multivariate normal distribution, how can I jointly fit them allowing for their residual error to be correlated?

  • Operating System: x86_64-pc-linux-gnu (64-bit)
  • brms Version: version 2.3.1

Q1: To get the right names of the parameters, check out parnames(k_fit_brms). However, it is not clear to me what Mvec1_X.vec1 is supposed to be.

Q2: Currently, you can either model all model parts as having residual correlations or none of them. This may change with brms 3.0, though. You may of course extract the stan code via stancode(k_fit_brms) and then amend it to your needs.

Thank you Paul looking forward for brms 3.0. In the mean time I will work on your suggestion.

Hi @paul.buerkner,

Sorry for coming back to my question. In the generated quantities block of the stancode produced by brms , I see
corr_matrix[M_1] Cor_1 = multiply_lower_tri_self_transpose(L_1);

For the purpose of computing direct and indirect effect I was looking for covaraince estimate rather than correlation for the random effect terms in my model above. Is there a way to force brms to let me get these estimates.

Kind regards

You may use method VarCorr().

1 Like