Applying BayestestR to longitudinal mediation models fitted using brms package

Modeling
, ,
1

Greetings to Stan community!

Thanks so much to the help from the Stan community, I managed to successfully (?) run a basic longitudinal mediation model.

And I tried to extract the direct and indirect effects, failing in replicating this suggestion from the change lab.
(unfortunately, ‘med_post’ function did not give an output product that starts with ‘cor-’ like in ‘cor_id__fwkdiscw_fwkstrcw__freldiscw_fwkdiscw’, as suggested in the link)

I however found a fascinating package called ‘BayestestR’ which automatically gives direct, and indirect effects that seem to be based on the causal mediation effect framework.

And using this package, I managed to obtain direct and indirect effects!
And although the package is very convenient in obtaining effects out of brms based models, I’m seeking your advice on whether it is okay to plug longitudinal mediation models from brms into mediation function of the BayestestR package. (I’m just very grateful it worked but now I’m having concerns whether its use can be justified in my dissertation!)

I would really appreciate your help and advice!!

(please find below the replicable code in R :)!!)

library(brms)
library(bmlm)
#use dataset from 'bmlm' package for a reproducible example
df <- BLch9

#create a binary mediator and covariate
df$mediator_dum <- rbinom(nrow(df), 1, 0.5)
df$covariate_gender <- rbinom(nrow(df), 1, 0.5)

#implement longitudinal mediation with a covariate, that addresses autoregressive error (AR(1))

#Dependent variable (time varying, continuous) = fwkstrs
#Independent variable (time varying, continuous) = fwkdis
#Mediator (time varying, dichotomous) = mediator_dum
#Covariate (time invariant, dichotomous) = covariate_gender

a<-bf(fwkstrs ~ time + mediator_dum + fwkdis + covariate_gender +
      (1 + covariate_gender|id), family = "gaussian")

b<-bf(mediator_dum ~ time + fwkdis + covariate_gender +
      (1 + covariate_gender|id), family = "bernoulli")

#temporal autocorrelation included ('AR(1)')
model <- brm(a + b + set_rescor(FALSE), data = df, ar(time = time, gr = id, p=1))

library(bayestestR)
mediation(model = model, treatment = "fwkdis", 
          mediator = "mediator_dum", response = NULL, centrality = "mean", 
          ci = 0.95)