Multilevel path anaylysis in lavaan: Testing for mediation at both Level 1 and Level 2
Overview
For this example, we are testing a model that incorporates mediation at both Level 1 and Level 2. We will include compositional variables (i.g., group means of Level 1) variables as predictors of the mediator at Level 2 and (indirectly) the outcome variable (also at Level 2). This discussion picks up from previous demonstrations on multilevel path analysis with lavaan found at https://mikesquanthub.blogspot.com/2024/04/multilevel-path-analysis-in-lavaan.html and https://mikesquanthub.blogspot.com/2024/04/multilevel-path-model-using-lavaan.html.
The example is adapted from a study by Xie et al. (2015), “Linking Colleague Support to Employees’ Promotive Voice: A Moderated Mediation Model” published in PLOS One. Although we will not be testing the authors’ model, we will rely on the dataset they published in the Supplementary Materials on the journal site for their article (https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0132123).
I have constructed a smaller dataset from the original data and made some changes in naming, etc. You can download a copy of the .csv file containing the data here. The .csv file is named prvoiceV5. Once you have downloaded the data to your computer, you will need to import it into RStudio for use.
For our analysis, we are modeling data based on \(n_1\)=162 employees (level 1 units) nested within \(n_2\)=51 work teams (level 2 units). At the employee level, we are predicting promotive voice as a function of felt obligation, psychological safety, and team identification. Felt obligation is treated as a mediator of the effects of psychological safety and team identification on promotive voice.
At level 2, we are modeling the effects of the compositional effects of psychological safety and team identification on felt obligation and promotive voice at the team-level.
Here is a figure of the model we will be specifying and testing using lavaan.
Let’s take a look at the names of the
variables in our dataset.
names(prvoiceV5)## [1] "X" "teamcode"
## [3] "colleaguesupport" "feltobl"
## [5] "gender" "age"
## [7] "eduation" "tenure"
## [9] "teamsize" "teamID"
## [11] "subgroupformationindividual" "psysafe"
## [13] "prvoice"The variables we will be using when running our analysis include:
- prvoice (level 1 outcome; promotive voice)
- feltobl (level 1 mediator; felt obligation)
- teamID (level 1 antecedent; team identification)
- psysafe (level 1 antecedent; psychological safety)
All variables are treated as continuous in the analysis. [Keep in mind too that lavaan currently does not allow for categorical endogenous variables]
Before carrying out our main analyses, we will want to generate cluster means (i.e., team-level means) for our exogenous variables to serve as predictors at Level 2.
library(mitml)## *** This is beta software. Please report any bugs!
## *** See the NEWS file for recent changes.prvoiceV5$teamID_M <- clusterMeans(prvoiceV5$teamID, cluster = prvoiceV5$teamcode)
prvoiceV5$psysafe_M <- clusterMeans(prvoiceV5$psysafe, cluster = prvoiceV5$teamcode)
names(prvoiceV5) ## [1] "X" "teamcode"
## [3] "colleaguesupport" "feltobl"
## [5] "gender" "age"
## [7] "eduation" "tenure"
## [9] "teamsize" "teamID"
## [11] "subgroupformationindividual" "psysafe"
## [13] "prvoice" "teamID_M"
## [15] "psysafe_M"head(prvoiceV5)## X teamcode colleaguesupport feltobl gender age eduation tenure teamsize
## 1 1 2 4.50 5.00 1 2 4 28 10
## 2 2 2 4.50 4.20 1 2 4 12 10
## 3 3 2 3.50 4.45 1 3 3 23 10
## 4 4 4 3.75 4.20 0 2 2 3 4
## 5 5 4 5.00 4.20 0 2 3 18 4
## 6 6 4 2.25 3.00 1 2 3 18 4
## teamID subgroupformationindividual psysafe prvoice teamID_M psysafe_M
## 1 6.50 2.00 6.57 4.2 6.553333 6.046667
## 2 6.33 3.50 6.00 3.8 6.553333 6.046667
## 3 6.83 3.00 5.57 5.4 6.553333 6.046667
## 4 5.50 5.50 5.71 4.0 5.113333 4.236667
## 5 4.67 4.50 4.14 5.2 5.113333 4.236667
## 6 5.17 2.75 2.86 6.0 5.113333 4.236667We see here that both compositional variables have been added to the data frame.
The main analyses
library(lavaan)model1a <- '
level: 1 # specification of level 1 model
prvoice ~ feltobl
feltobl ~ psysafe + teamID
level: 2 # specification of level 2 model
prvoice ~ feltobl
feltobl ~ psysafe_M + teamID_M
'
The following code will fit the model and generate our output.
Make sure you include the cluster argument (see below) and set it equal
to the name of the clustering variable (enclosed in quotes). [The output
will contain warning messages about lack of variance within clusters for
certain variables. This is simply a warning. But I do not reproduce it
here.]
fit1 <- sem(model1a, data = prvoiceV5, cluster = "teamcode")
summary(fit1, fit.measures=T, rsquare=T)## lavaan 0.6.13 ended normally after 51 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 12
##
## Number of observations 162
## Number of clusters [teamcode] 51
##
## Model Test User Model:
##
## Test statistic 2.502
## Degrees of freedom 4
## P-value (Chi-square) 0.644
##
## Model Test Baseline Model:
##
## Test statistic 127.043
## Degrees of freedom 10
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.032
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -343.478
## Loglikelihood unrestricted model (H1) -342.227
##
## Akaike (AIC) 710.956
## Bayesian (BIC) 748.008
## Sample-size adjusted Bayesian (SABIC) 710.018
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent confidence interval - lower 0.000
## 90 Percent confidence interval - upper 0.095
## P-value H_0: RMSEA <= 0.050 0.788
## P-value H_0: RMSEA >= 0.080 0.091
##
## Standardized Root Mean Square Residual (corr metric):
##
## SRMR (within covariance matrix) 0.026
## SRMR (between covariance matrix) 0.075
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Observed
## Observed information based on Hessian
##
##
## Level 1 [within]:
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## prvoice ~
## feltobl 0.430 0.109 3.926 0.000
## feltobl ~
## psysafe 0.365 0.051 7.179 0.000
## teamID 0.285 0.051 5.559 0.000
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .prvoice 0.000
## .feltobl 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .prvoice 0.625 0.084 7.434 0.000
## .feltobl 0.229 0.029 7.933 0.000
##
## R-Square:
## Estimate
## prvoice 0.123
## feltobl 0.518
##
##
## Level 2 [teamcode]:
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## prvoice ~
## feltobl 3.945 2.372 1.663 0.096
## feltobl ~
## psysafe_M -0.041 0.071 -0.575 0.565
## teamID_M -0.077 0.092 -0.838 0.402
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .prvoice 0.067 2.699 0.025 0.980
## .feltobl 1.457 0.452 3.223 0.001
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .prvoice 0.328 0.291 1.130 0.258
## .feltobl 0.024 0.017 1.423 0.155
##
## R-Square:
## Estimate
## prvoice 0.569
## feltobl 0.131The chi-square goodness of fit test was not statistically significant, \(\chi^2\)(4)=2.502, p=.644 (lead to a failure to reject the exact fit hypothesis). The CFI=1 and TLI = 1.032 were above .95, indicating a very good fitting model. Similarly, the RMSEA of approximately 0 indicates a very good fitting model. The SRMR within (.026) was below .05 suggesting good fit at Level 1. The SRMR between (.075) was rather high indicating a fit that was not particularly strong.
At Level 1, team identification (b=.285, p<.001) and psychological safety (b=.365, p<.001) were both significant positive predictors of felt obligation. Felt obligation, in turn, was a positive and significant (b=.430, p<.001) predictor of promotive voice.
At Level 2, neither the compositional effects of team identification
or psychological safety predicted felt obligation (at the team level).
Moreover, felt obligation did not significantly predict variation
between teams in promotive voice.
If you would like to report
on the intraclass correlations for the endogenous level 1 variables use
the lavInspect() function.
lavInspect(fit1, "icc")## prvoice feltobl psysafe teamID
## 0.523 0.051 0.000 0.000We can test for mediation at Level 1 and Level 2 using lavaan. To do this, you will need to (a) assign labels to the paths used in the construction of the indirect effects and then (b) create estimands for the indirect effects. The figure below shows you the paths involved in the construction of the indirect effect estimands.
Specification of our model…
model1a <- '
level: 1 # specification of level 1 model
prvoice ~ l1b*feltobl
feltobl ~ l1a1*psysafe + l1a2*teamID
# Adding estimands to capture indirect effects
IE1 := l1a1*l1b #psysafe->feltobl->prvoice
IE2 := l1a2*l1b #teamID->feltobl->prvoice
level: 2 # specification of level 2 model
prvoice ~ l2b*feltobl
feltobl ~ l2a1*psysafe_M + l2a2*teamID_M
# Adding estimands to capture indirect effects
IE3 := l2a1*l2b #psysafe_M ->feltobl->prvoice
IE4 := l2a2*l2b #teamID_M ->feltobl->prvoice
'
fit1 <- sem(model1a, data = prvoiceV5, cluster = "teamcode")
summary(fit1, fit.measures=T, rsquare=T)## lavaan 0.6.13 ended normally after 51 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 12
##
## Number of observations 162
## Number of clusters [teamcode] 51
##
## Model Test User Model:
##
## Test statistic 2.502
## Degrees of freedom 4
## P-value (Chi-square) 0.644
##
## Model Test Baseline Model:
##
## Test statistic 127.043
## Degrees of freedom 10
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.032
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -343.478
## Loglikelihood unrestricted model (H1) -342.227
##
## Akaike (AIC) 710.956
## Bayesian (BIC) 748.008
## Sample-size adjusted Bayesian (SABIC) 710.018
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent confidence interval - lower 0.000
## 90 Percent confidence interval - upper 0.095
## P-value H_0: RMSEA <= 0.050 0.788
## P-value H_0: RMSEA >= 0.080 0.091
##
## Standardized Root Mean Square Residual (corr metric):
##
## SRMR (within covariance matrix) 0.026
## SRMR (between covariance matrix) 0.075
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Observed
## Observed information based on Hessian
##
##
## Level 1 [within]:
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## prvoice ~
## feltobl (l1b) 0.430 0.109 3.926 0.000
## feltobl ~
## psysafe (l1a1) 0.365 0.051 7.179 0.000
## teamID (l1a2) 0.285 0.051 5.559 0.000
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .prvoice 0.000
## .feltobl 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .prvoice 0.625 0.084 7.434 0.000
## .feltobl 0.229 0.029 7.933 0.000
##
## R-Square:
## Estimate
## prvoice 0.123
## feltobl 0.518
##
##
## Level 2 [teamcode]:
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## prvoice ~
## feltobl (l2b) 3.945 2.372 1.663 0.096
## feltobl ~
## psysf_M (l2a1) -0.041 0.071 -0.575 0.565
## temID_M (l2a2) -0.077 0.092 -0.838 0.402
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .prvoice 0.067 2.699 0.025 0.980
## .feltobl 1.457 0.452 3.223 0.001
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .prvoice 0.328 0.291 1.130 0.258
## .feltobl 0.024 0.017 1.423 0.155
##
## R-Square:
## Estimate
## prvoice 0.569
## feltobl 0.131
##
## Defined Parameters:
## Estimate Std.Err z-value P(>|z|)
## IE1 0.157 0.047 3.345 0.001
## IE2 0.122 0.039 3.104 0.002
## IE3 -0.161 0.290 -0.555 0.579
## IE4 -0.304 0.268 -1.135 0.256We see that the both indirect effects at Level 1 were positive and statistically significant. Both indirect effects at Level 2 were not statistically significant.
References
Xie, X.-Y., Ling, C.-D., Mo, S.-J., Luan, K. (2015). Linking colleague support to employees’ promotive voice: A moderated mediation model. PLoS ONE 10(7): e0132123. doi:10.1371/journal.pone.0132123. Downloaded April 21, 2024 from: https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0132123
Comments
Post a Comment