strengejacke / sjstats

Effect size measures and significance tests
https://strengejacke.github.io/sjstats
189 stars 21 forks source link

mediation() with centered mediator #90

Closed franciscowilhelm closed 4 years ago

franciscowilhelm commented 5 years ago

I ran a mediation model in brms with a (person-mean / group-mean) centered predictor. Extracting the mediation model required some tempering with the mediationfunction of sjstats. The problematic thing is that the mediator variable is entered in its centered form in the treatment model, but serves as an uncentered response variable in the mediator model. Hence, the mediator variable has a different name in each of the models.

model_med <- bf(m ~ x_pc + (x_pc|p|id) )
model_out  <- bf(y ~ m_pc + x_pc + (m_pc|p|id) )
fit <- brm(model_med + model_out + set_rescor(FALSE), 
                       data = data, 
                       cores=4)
print(sjstats::mediation(fit, treatment = "x_pc", mediator = "m_pc"), digits=4)
print(sjstats::mediation(fit, treatment = "x_pc", mediator = "m"), digits=4)

Neither of the mediation calls works. I made a rather crude and quick fix to the mediation function in my fork (https://github.com/franciscowilhelm/sjstats/commit/71d20f35130c253cbadcb1db3844d2550f9ae0cb). It would be great if you could implement such a functionality!

Herzlichen Dank Francisco

strengejacke commented 4 years ago

Do you have some sample data to reproduce your example?

franciscowilhelm commented 4 years ago

Using the data from: http://intensivelongitudinal.com/ch9/ch9index.html (Chapter 9 Mplus Files)

df <- read.table("mediation.dat",
                sep = "\t", header = FALSE) %>% 
  select(V1, V16:V18) %>%
  rename(id = V1, x = V16, m = V17, y = V18) %>%
  de_mean(., x, grp = id) %>% de_mean(., m, grp = id)

model_med <- bf(m ~ x_dm + (x_dm|p|id) )
model_out  <- bf(y ~ m_dm + x_dm + (m_dm|p|id) )
fit <- brm(model_med + model_out + set_rescor(FALSE), 
           data = df, 
           cores=4)
print(sjstats::mediation(fit, treatment = "x_dm", mediator = "m_dm"), digits=4)
print(sjstats::mediation(fit, treatment = "x_dm", mediator = "m"), digits=4)
strengejacke commented 4 years ago

Possible in the re-implementation in bayestestR. You need to define responses explicitly, using a named vector for the response-argument.

mediation(fit, response = c(m = "m_dm", y = "y"))
#> # Causal Mediation Analysis for Stan Model
#> 
#>   Treatment: x_dm
#>   Mediator : m_dm
#>   Response : y
#> 
#> Effect                 | Estimate |        89% ETI
#> --------------------------------------------------
#> Direct Effect (ADE)    |    0.100 | [0.066, 0.134]
#> Indirect Effect (ACME) |    0.029 | [0.016, 0.046]
#> Mediator Effect        |    0.156 | [0.106, 0.205]
#> Total Effect           |    0.129 | [0.095, 0.166]
#> 
#> Proportion mediated: 22.60% [11.31%, 33.88%]
franciscowilhelm commented 4 years ago

Great! Thanks for adding this feature.