easystats / bayestestR

:ghost: Utilities for analyzing Bayesian models and posterior distributions
https://easystats.github.io/bayestestR/
GNU General Public License v3.0
562 stars 55 forks source link

Analyzing multiple mediator model using brms and bayestestR #475

Open MichaelWeaver01 opened 2 years ago

MichaelWeaver01 commented 2 years ago

Describe the bug When trying to run a mediation analysis having more than 1 mediation path using brms (2.16.1) and bayestestR (0.11.0), I get the following messages:

Error in .subset2(x, i, exact = exact) : subscript out of bounds
In addition: Warning messages:
1: In response == mediator :
  longer object length is not a multiple of shorter object length
2: In response != mediator :
  longer object length is not a multiple of shorter object length
3: In if (variable %in% colnames(mf)) { :
  the condition has length > 1 and only the first element will be used

Running a simpler model, with only 1 mediator, works fine.

To Reproduce Here is the R code for that produces the brms fit object - the model compiles and runs fine:

y1_model <- bf(CESD_Factor ~ Total_PSS_scaled + PQSIglobal_PSS +
                 Trait_ANX + State_ANX+ PQSI_Global_Score_scaled,
               family = cumulative())

y2_model <- bf(PQSI_Global_Score_scaled ~ Trait_ANX + State_ANX + Total_PSS_scaled ,
               family = student)

y3_model <- bf(Total_PSS_scaled ~ Trait_ANX + State_ANX,
               family = student)

global_pss_cesd_Int_fit <- brm(data = mydata,
                               y1_model + y2_model + y3_model + set_rescor(FALSE),
                               warmup = 1000, iter = 9000, chains = 4,
                               cores = 4)

There are 3 endogenous variables (CESD_Factor, PQSI_Global_Score_scaled, and Total_PSS_scaled) and 2 mediators (PQSI_Global_Score_scaled and Total_PSS_scaled).

I have tried several ways to write the mediation() statement, such as:

pqsi_med <- mediation(global_pss_cesd_Int_fit, response = "CESDFactor")

response <- c("CESD_Factor", "PQSI_Global_Score_scaled")
pqsi_med <- mediation(global_pss_cesd_Int_fit, "PQSI_Global_Score_scaled", response)

but get error messages such as:

Error in .subset2(x, i, exact = exact) : subscript out of bounds
In addition: Warning messages:
1: In response == mediator :
  longer object length is not a multiple of shorter object length
2: In response != mediator :
  longer object length is not a multiple of shorter object length
3: In if (variable %in% colnames(mf)) { :
  the condition has length > 1 and only the first element will be used

** The simpler model (1 mediator) that runs fine:

y_model <- bf(CESD_Factor ~ PQSI_Global_Score_scaled + Total_PSS_scaled + PQSIglobal_PSS,
              family = cumulative())

m_model <- bf(PQSI_Global_Score_scaled ~ Total_PSS_scaled,
              family = student)

global_pss_cesd_Int_fit <- brm(data = mydata,
                           y_model + m_model + set_rescor(FALSE),
                           warmup = 1000, iter = 9000, chains = 4,
                           cores = 4)
mfull <- mediation(global_pss_cesd_Int_fit)

Expected behaviour I expected that analyses of the mediation paths operating through each of the 2 mediators would be produced.

Specifiations (please complete the following information): RStudio 2021.09.0 Build 351 R version 4.1.1 (2021-08-10) -- "Kick Things" brms version 2.16.1 bayestestR version 0.11.0

marioangst commented 2 years ago

I just encountered the same problem, with exactly the same error message as @MichaelWeaver01 encountered.

A slight twist: A brms model composed of three equations and two mediators runs into the error. A brms model composed of two equations (but still two mediators) runs fine.

If needed, I can post a reprex later.

DominiqueMakowski commented 2 years ago

@mattansb any thoughts on that?

mattansb commented 2 years ago

I'm not familiar with the internals here - this was implemented by @strengejacke