Mediator Analysis in R
Li Zhang
Wednesday, September 03, 2014
Introduction
This article is using two mediator analysis methods in R to analyze two datasets. The two methods are 3-step mediator analysis and structural Equation Model (SEM) respectively. The first dataset to analyze is called ‘AERA Final Dataset.xlsx’, and the second is called ‘Last version wdaysstep-0418.xlsx’. In the following section, we will give a brief introduciton to the two mediator analysis methods, apply the methods to the two datasets, and display the results without further interpretation.
Mediator Anlysis
1 AERA Final Dataset
We perform single-mediator analysis on the AERA Final Dataset. the independent varible (x) is called PA1, the mediator (m) is called ZFITNESSFINAL, and the output variable (y) is HRQOL2_A. We removed missing values from the original dataset and as a result there is a total of 194 observations in the final dataset. Let’s first take a look at these variabes:
## PA1 UpdatedStepmin NEWTOTALFITNESS4 ZFITNESSFINAL
## Min. :0.00 Min. : 53.2 Min. : 48.3 Min. :-14.060
## 1st Qu.:2.93 1st Qu.: 83.3 1st Qu.: 78.8 1st Qu.: -1.638
## Median :3.39 Median : 98.6 Median : 93.5 Median : 0.271
## Mean :3.33 Mean : 96.9 Mean : 94.2 Mean : -0.005
## 3rd Qu.:3.83 3rd Qu.:109.0 3rd Qu.:106.2 3rd Qu.: 1.974
## Max. :4.71 Max. :143.4 Max. :165.4 Max. : 8.727
## HRQOL2_A
## Min. : 43.5
## 1st Qu.: 81.5
## Median : 88.0
## Mean : 85.5
## 3rd Qu.: 93.5
## Max. :100.01.a 3-step mediator analysis
In the 3-step mediator analysis, we do seperate regressions on the following equations to get the regression coefficients.
(1) \(y = cx+d_1+e_1\)
(2) \(m = ax+d_2+e_2\)
(3) \(y = c'x+bm+d_3+e_3\)
We use lm function in R to do the regressions. the details are as below.
# step1: lm(Y~X) regression coefficient for X should be significant
# Y = cX+d1+e1
fit1 <- lm(aera$HRQOL2_A~aera$PA1)
summary(fit1) # significant##
## Call:
## lm(formula = aera$HRQOL2_A ~ aera$PA1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -43.49 -4.85 1.72 6.92 18.35
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 71.973 3.213 22.40 < 2e-16 ***
## aera$PA1 4.068 0.938 4.34 2.3e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 10.6 on 192 degrees of freedom
## Multiple R-squared: 0.0893, Adjusted R-squared: 0.0845
## F-statistic: 18.8 on 1 and 192 DF, p-value: 2.32e-05# step2: lm(M~X) regression coefficient for X should be significant
# M = aX + d2 + e2
fit2 <- lm(aera$ZFITNESSFINAL ~ aera$PA1)
summary(fit2)##
## Call:
## lm(formula = aera$ZFITNESSFINAL ~ aera$PA1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.776 -1.528 0.311 2.016 8.052
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.964 0.924 -2.13 0.035 *
## aera$PA1 0.588 0.270 2.18 0.030 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.04 on 192 degrees of freedom
## Multiple R-squared: 0.0242, Adjusted R-squared: 0.0191
## F-statistic: 4.76 on 1 and 192 DF, p-value: 0.0304# step3: lm(Y ~X+M) regression coefficient for M should be significant
# Y = c'X + bM + d3 + e3
fit3 <- lm(aera$HRQOL2_A~aera$ZFITNESSFINAL+ aera$PA1)
summary(fit3) # both M and X significant##
## Call:
## lm(formula = aera$HRQOL2_A ~ aera$ZFITNESSFINAL + aera$PA1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -40.40 -4.92 1.41 7.72 20.88
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 73.941 3.121 23.69 < 2e-16 ***
## aera$ZFITNESSFINAL 1.002 0.241 4.16 4.8e-05 ***
## aera$PA1 3.479 0.911 3.82 0.00018 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 10.2 on 191 degrees of freedom
## Multiple R-squared: 0.165, Adjusted R-squared: 0.156
## F-statistic: 18.9 on 2 and 191 DF, p-value: 3.34e-08A direct effect relating x to y is quantitifd by c’, that is 3.479 with standard error = 0.911 as shown in summary(fit3).
To evaluate the mediated (indirect) effect, we use the medci function in the RMediation package to calculate \(ab\). PRODCLIN method is chosen in the medci function to evaluate the 95% confidence interval of \(ab\).
# CI for a*b:
a = summary(fit2)$coefficients[2,1]
se.a = summary(fit2)$coefficients[2,2]
b = summary(fit3)$coefficients[2,1]
se.b = summary(fit3)$coefficients[2,2]
library(RMediation)
medci(mu.x=a, mu.y=b, se.x=se.a, se.y=se.b, rho=0, alpha=.05,
type="prodclin", plot=TRUE, plotCI=TRUE)
## 2.5 % 97.5 %
## 0.05409 1.27333As shown in the resulted figure, the mean of ab is 0.59, and standard error 0.312. The corresponding 95% CI is (0.054,1.273).
1.b SEM analysis
In this section, We use the sem function in the lavaan package to perform the mediator analysis. In the definiton of the model variable, we specify the linear regression equations relating y to x and m, and relating m to x. Then we ask the model to calculate the indirect effect as well as the total effect.
library(lavaan)
model <-'
# outcome model
HRQOL2_A ~ c*PA1+b*ZFITNESSFINAL
# mediator models
ZFITNESSFINAL ~ a*PA1
# indirect effects
IDE := a*b
# total effect
total := c+a*b
'Then we call the sem function applying the model to the aera dataset. A summary of the sem fitting is displayed as follows.
fitsem <- sem(model,data=aera,bootstrap = 10000)
summary(fitsem, fit.measures=TRUE, standardize=TRUE, rsquare=TRUE)## lavaan (0.5-16) converged normally after 29 iterations
##
## Number of observations 194
##
## Estimator ML
## Minimum Function Test Statistic 0.000
## Degrees of freedom 0
## P-value (Chi-square) 0.000
##
## Model test baseline model:
##
## Minimum Function Test Statistic 39.715
## Degrees of freedom 3
## P-value 0.000
##
## User model versus baseline model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.000
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -1447.873
## Loglikelihood unrestricted model (H1) -1447.873
##
## Number of free parameters 5
## Akaike (AIC) 2905.747
## Bayesian (BIC) 2922.086
## Sample-size adjusted Bayesian (BIC) 2906.247
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent Confidence Interval 0.000 0.000
## P-value RMSEA <= 0.05 1.000
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.000
##
## Parameter estimates:
##
## Information Expected
## Standard Errors Standard
##
## Estimate Std.err Z-value P(>|z|) Std.lv Std.all
## Regressions:
## HRQOL2_A ~
## PA1 (c) 3.479 0.904 3.847 0.000 3.479 0.255
## ZFITNESSF (b) 1.002 0.239 4.193 0.000 1.002 0.278
## ZFITNESSFINAL ~
## PA1 (a) 0.588 0.268 2.192 0.028 0.588 0.155
##
## Variances:
## HRQOL2_A 101.512 10.307 101.512 0.835
## ZFITNESSFINAL 9.160 0.930 9.160 0.976
##
## Defined parameters:
## IDE 0.590 0.303 1.943 0.052 0.590 0.043
## total 4.068 0.933 4.361 0.000 4.068 0.299
##
## R-Square:
##
## HRQOL2_A 0.165
## ZFITNESSFINAL 0.024Finally, we call the parameterEstimates function in the lavaan package to estimate the indirect effect (ab) and the total effect (c+ab).
boot.fit <- parameterEstimates(fitsem, boot.ci.type="bca.simple",level=0.95, ci=TRUE,standardized = FALSE)
boot.fit## lhs op rhs label est se z pvalue
## 1 HRQOL2_A ~ PA1 c 3.479 0.904 3.847 0.000
## 2 HRQOL2_A ~ ZFITNESSFINAL b 1.002 0.239 4.193 0.000
## 3 ZFITNESSFINAL ~ PA1 a 0.588 0.268 2.192 0.028
## 4 HRQOL2_A ~~ HRQOL2_A 101.512 10.307 9.849 0.000
## 5 ZFITNESSFINAL ~~ ZFITNESSFINAL 9.160 0.930 9.849 0.000
## 6 PA1 ~~ PA1 0.656 0.000 NA NA
## 7 IDE := a*b IDE 0.590 0.303 1.943 0.052
## 8 total := c+a*b total 4.068 0.933 4.361 0.000
## ci.lower ci.upper
## 1 1.706 5.251
## 2 0.534 1.471
## 3 0.062 1.114
## 4 81.310 121.713
## 5 7.337 10.983
## 6 0.656 0.656
## 7 -0.005 1.184
## 8 2.240 5.897The result shows that the 95% CI for ab is (-0.005, 1.184).
1.c Independent variable = UpdatedStepmin
Now we use Updatedstepmin as the independent variable. Firstly, we do the 3-step regressions to get the coefficients a, b, and c.
# step1: lm(Y~X) regression coefficient for X should be significant
# Y = cX+d1+e1
fit1 <- lm(aera$HRQOL2_A~aera$UpdatedStepmin)
summary(fit1) # significant##
## Call:
## lm(formula = aera$HRQOL2_A ~ aera$UpdatedStepmin)
##
## Residuals:
## Min 1Q Median 3Q Max
## -38.45 -4.40 2.04 7.14 17.47
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 73.3579 4.3259 16.96 <2e-16 ***
## aera$UpdatedStepmin 0.1255 0.0439 2.86 0.0047 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 10.9 on 192 degrees of freedom
## Multiple R-squared: 0.0408, Adjusted R-squared: 0.0358
## F-statistic: 8.16 on 1 and 192 DF, p-value: 0.00474# step2: lm(M~X) regression coefficient for X should be significant
# M = aX + d2 + e2
fit2 <- lm(aera$ZFITNESSFINAL ~ aera$UpdatedStepmin)
summary(fit2)##
## Call:
## lm(formula = aera$ZFITNESSFINAL ~ aera$UpdatedStepmin)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.472 -1.736 0.219 2.164 7.831
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.6353 1.1981 -3.03 0.0027 **
## aera$UpdatedStepmin 0.0375 0.0122 3.08 0.0024 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.01 on 192 degrees of freedom
## Multiple R-squared: 0.0471, Adjusted R-squared: 0.0421
## F-statistic: 9.49 on 1 and 192 DF, p-value: 0.00237# step3: lm(Y ~X+M) regression coefficient for M should be significant
# Y = c'X + bM + d3 + e3
fit3 <- lm(aera$HRQOL2_A~aera$ZFITNESSFINAL+ aera$UpdatedStepmin)
summary(fit3) # both M and X significant##
## Call:
## lm(formula = aera$HRQOL2_A ~ aera$ZFITNESSFINAL + aera$UpdatedStepmin)
##
## Residuals:
## Min 1Q Median 3Q Max
## -35.45 -4.33 1.68 7.68 16.92
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 77.1246 4.2532 18.13 < 2e-16 ***
## aera$ZFITNESSFINAL 1.0361 0.2503 4.14 5.2e-05 ***
## aera$UpdatedStepmin 0.0867 0.0432 2.01 0.046 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 10.4 on 191 degrees of freedom
## Multiple R-squared: 0.12, Adjusted R-squared: 0.111
## F-statistic: 13 on 2 and 191 DF, p-value: 5.11e-06Next, we use the medci function to evaluate the 95% confidence interval of \(ab\).
# CI for a*b:
a = summary(fit2)$coefficients[2,1]
se.a = summary(fit2)$coefficients[2,2]
b = summary(fit3)$coefficients[2,1]
se.b = summary(fit3)$coefficients[2,2]
library(RMediation)
medci(mu.x=a, mu.y=b, se.x=se.a, se.y=se.b, rho=0, alpha=.05,
type="prodclin", plot=TRUE, plotCI=TRUE)
## 2.5 % 97.5 %
## 0.01173 0.07400Finally, we apply the sem function.
library(lavaan)
model <-'
# outcome model
HRQOL2_A ~ c*UpdatedStepmin+b*ZFITNESSFINAL
# mediator models
ZFITNESSFINAL ~ a*UpdatedStepmin
# indirect effects
IDE := a*b
# total effect
total := c+a*b
'
fitsem <- sem(model,data=aera,bootstrap = 10000)
summary(fitsem, fit.measures=TRUE, standardize=TRUE, rsquare=TRUE)## lavaan (0.5-16) converged normally after 23 iterations
##
## Number of observations 194
##
## Estimator ML
## Minimum Function Test Statistic 0.000
## Degrees of freedom 0
## P-value (Chi-square) 1.000
##
## Model test baseline model:
##
## Minimum Function Test Statistic 34.108
## Degrees of freedom 3
## P-value 0.000
##
## User model versus baseline model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.000
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -2049.577
## Loglikelihood unrestricted model (H1) -2049.577
##
## Number of free parameters 5
## Akaike (AIC) 4109.153
## Bayesian (BIC) 4125.492
## Sample-size adjusted Bayesian (BIC) 4109.653
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent Confidence Interval 0.000 0.000
## P-value RMSEA <= 0.05 1.000
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.000
##
## Parameter estimates:
##
## Information Expected
## Standard Errors Standard
##
## Estimate Std.err Z-value P(>|z|) Std.lv Std.all
## Regressions:
## HRQOL2_A ~
## UpdtdStpm (c) 0.087 0.043 2.022 0.043 0.087 0.139
## ZFITNESSF (b) 1.036 0.248 4.173 0.000 1.036 0.288
## ZFITNESSFINAL ~
## UpdtdStpm (a) 0.037 0.012 3.096 0.002 0.037 0.217
##
## Variances:
## HRQOL2_A 107.000 10.864 107.000 0.880
## ZFITNESSFINAL 8.945 0.908 8.945 0.953
##
## Defined parameters:
## IDE 0.039 0.016 2.486 0.013 0.039 0.062
## total 0.125 0.044 2.872 0.004 0.125 0.202
##
## R-Square:
##
## HRQOL2_A 0.120
## ZFITNESSFINAL 0.047boot.fit <- parameterEstimates(fitsem, boot.ci.type="bca.simple",level=0.95, ci=TRUE,standardized = FALSE)
boot.fit## lhs op rhs label est se z pvalue
## 1 HRQOL2_A ~ UpdatedStepmin c 0.087 0.043 2.022 0.043
## 2 HRQOL2_A ~ ZFITNESSFINAL b 1.036 0.248 4.173 0.000
## 3 ZFITNESSFINAL ~ UpdatedStepmin a 0.037 0.012 3.096 0.002
## 4 HRQOL2_A ~~ HRQOL2_A 107.000 10.864 9.849 0.000
## 5 ZFITNESSFINAL ~~ ZFITNESSFINAL 8.945 0.908 9.849 0.000
## 6 UpdatedStepmin ~~ UpdatedStepmin 314.861 0.000 NA NA
## 7 IDE := a*b IDE 0.039 0.016 2.486 0.013
## 8 total := c+a*b total 0.125 0.044 2.872 0.004
## ci.lower ci.upper
## 1 0.003 0.171
## 2 0.549 1.523
## 3 0.014 0.061
## 4 85.707 128.293
## 5 7.165 10.725
## 6 314.861 314.861
## 7 0.008 0.069
## 8 0.040 0.2112 Last version wdaysstep-0418
The mediator analysis for the second dataset is complicated relative to the first datasets. There are 2 independent variables - Mastery and Performance, 2 mediators - expectancybelief and taskvalue, and 3 output variables - HRQOL, PA, and Av.Stepsmin.
First, let’s take a look at the data. One missing value is removed from the original dataset and the following is a summary of a total of 336 observations.
## Mastery Performance expectancybelief taskvalue
## Min. :1.00 Min. :1.00 Min. :2.00 Min. :1.83
## 1st Qu.:4.00 1st Qu.:2.00 1st Qu.:4.00 1st Qu.:3.83
## Median :4.40 Median :2.75 Median :4.40 Median :4.33
## Mean :4.27 Mean :2.76 Mean :4.27 Mean :4.21
## 3rd Qu.:4.80 3rd Qu.:3.25 3rd Qu.:4.60 3rd Qu.:4.67
## Max. :5.00 Max. :5.00 Max. :5.00 Max. :5.00
## HRQOL PA Av.Stepsmin Grade Teacher
## Min. : 13.8 Min. :1.44 Min. : 22.8 4:172 1: 51
## 1st Qu.: 73.0 1st Qu.:2.96 1st Qu.: 51.7 5:164 2: 49
## Median : 83.1 Median :3.43 Median : 78.9 3:100
## Mean : 80.2 Mean :3.42 Mean : 77.3 4: 30
## 3rd Qu.: 90.4 3rd Qu.:3.93 3rd Qu.:100.9 5: 34
## Max. :100.0 Max. :4.86 Max. :143.4 6: 72In this multiple-mediator multivariate example, we define each of the mediators in relation to both the independent variables. We define each of the output variable as a linear relationship with both the independent variables and both the mediator variables. The detailed equations are shown in the followin subsections. It is noted that the model is so complicated that we avoid the 3-step mediation analysis and directly use the sem function to perform the analysis.
2.a output variable = HRQOL
The following is the SEM analysis for the output variable HRQOL. As shown in the definiton of model1 below, HRQOL is linearly related to the 2 independent variables and the 2 mediators.
library(lavaan)
# 2 x variables + 2 mediators
model1 <-'
# outcome model
HRQOL ~ c1*Mastery+c2*Performance+b1*expectancybelief+b2*taskvalue
# mediator models
expectancybelief ~ a11*Mastery+a12*Performance
taskvalue ~ a21*Mastery+a22*Performance
# indirect effects
MasteryExpectancy := b1*a11
Masterytaskvalue := b2*a21
MasteryIDE := b1*a11 + b2*a21
PerformancyExpectancy := b1*a12
PerformancyTaskvalue := b2*a22
PerfIDE := b1*a12+b2*a22
expectancybelief ~~ taskvalue # model correlation between mediators
'
# SEM analysis
fitsem1 <- sem(model1,data=xym,group=NULL)
summary(fitsem1, fit.measures=TRUE, standardize=TRUE, rsquare=TRUE)## lavaan (0.5-16) converged normally after 53 iterations
##
## Number of observations 336
##
## Estimator ML
## Minimum Function Test Statistic 0.000
## Degrees of freedom 0
## P-value (Chi-square) 1.000
##
## Model test baseline model:
##
## Minimum Function Test Statistic 231.330
## Degrees of freedom 9
## P-value 0.000
##
## User model versus baseline model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.000
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -2552.508
## Loglikelihood unrestricted model (H1) -2552.508
##
## Number of free parameters 12
## Akaike (AIC) 5129.016
## Bayesian (BIC) 5174.821
## Sample-size adjusted Bayesian (BIC) 5136.756
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent Confidence Interval 0.000 0.000
## P-value RMSEA <= 0.05 1.000
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.000
##
## Parameter estimates:
##
## Information Expected
## Standard Errors Standard
##
## Estimate Std.err Z-value P(>|z|) Std.lv Std.all
## Regressions:
## HRQOL ~
## Mastery (c1) 4.101 1.186 3.457 0.001 4.101 0.191
## Prfrmnc (c2) -1.898 0.788 -2.408 0.016 -1.898 -0.128
## expctnc (b1) 5.318 1.576 3.374 0.001 5.318 0.218
## taskval (b2) -0.436 1.595 -0.273 0.785 -0.436 -0.018
## expectancybelief ~
## Mastery (a11) 0.262 0.045 5.779 0.000 0.262 0.297
## Prfrmnc (a12) 0.102 0.031 3.262 0.001 0.102 0.168
## taskvalue ~
## Mastery (a21) 0.245 0.045 5.471 0.000 0.245 0.282
## Prfrmnc (a22) 0.102 0.031 3.305 0.001 0.102 0.170
##
## Covariances:
## expectancybelief ~~
## taskvalue 0.137 0.016 8.424 0.000 0.137 0.517
##
## Variances:
## HRQOL 163.833 12.640 163.833 0.896
## expectancyblf 0.268 0.021 0.268 0.871
## taskvalue 0.262 0.020 0.262 0.879
##
## Defined parameters:
## MstryExpctncy 1.391 0.477 2.914 0.004 1.391 0.065
## Masterytaskvl -0.107 0.391 -0.273 0.785 -0.107 -0.005
## MasteryIDE 1.284 0.457 2.808 0.005 1.284 0.060
## PrfrmncyExpct 0.542 0.231 2.345 0.019 0.542 0.037
## PrfrmncyTskvl -0.044 0.163 -0.273 0.785 -0.044 -0.003
## PerfIDE 0.497 0.225 2.210 0.027 0.497 0.034
##
## R-Square:
##
## HRQOL 0.104
## expectancybelief 0.129
## taskvalue 0.121# mediator effect
boot.fit1 <- parameterEstimates(fitsem1, boot.ci.type="bca.simple")
boot.fit1## lhs op rhs label est
## 1 HRQOL ~ Mastery c1 4.101
## 2 HRQOL ~ Performance c2 -1.898
## 3 HRQOL ~ expectancybelief b1 5.318
## 4 HRQOL ~ taskvalue b2 -0.436
## 5 expectancybelief ~ Mastery a11 0.262
## 6 expectancybelief ~ Performance a12 0.102
## 7 taskvalue ~ Mastery a21 0.245
## 8 taskvalue ~ Performance a22 0.102
## 9 expectancybelief ~~ taskvalue 0.137
## 10 HRQOL ~~ HRQOL 163.833
## 11 expectancybelief ~~ expectancybelief 0.268
## 12 taskvalue ~~ taskvalue 0.262
## 13 Mastery ~~ Mastery 0.396
## 14 Mastery ~~ Performance 0.074
## 15 Performance ~~ Performance 0.832
## 16 MasteryExpectancy := b1*a11 MasteryExpectancy 1.391
## 17 Masterytaskvalue := b2*a21 Masterytaskvalue -0.107
## 18 MasteryIDE := b1*a11+b2*a21 MasteryIDE 1.284
## 19 PerformancyExpectancy := b1*a12 PerformancyExpectancy 0.542
## 20 PerformancyTaskvalue := b2*a22 PerformancyTaskvalue -0.044
## 21 PerfIDE := b1*a12+b2*a22 PerfIDE 0.497
## se z pvalue ci.lower ci.upper
## 1 1.186 3.457 0.001 1.776 6.426
## 2 0.788 -2.408 0.016 -3.443 -0.353
## 3 1.576 3.374 0.001 2.229 8.407
## 4 1.595 -0.273 0.785 -3.561 2.689
## 5 0.045 5.779 0.000 0.173 0.350
## 6 0.031 3.262 0.001 0.041 0.163
## 7 0.045 5.471 0.000 0.157 0.332
## 8 0.031 3.305 0.001 0.042 0.163
## 9 0.016 8.424 0.000 0.105 0.169
## 10 12.640 12.961 0.000 139.059 188.606
## 11 0.021 12.961 0.000 0.228 0.309
## 12 0.020 12.961 0.000 0.222 0.302
## 13 0.000 NA NA 0.396 0.396
## 14 0.000 NA NA 0.074 0.074
## 15 0.000 NA NA 0.832 0.832
## 16 0.477 2.914 0.004 0.455 2.326
## 17 0.391 -0.273 0.785 -0.873 0.659
## 18 0.457 2.808 0.005 0.388 2.180
## 19 0.231 2.345 0.019 0.089 0.995
## 20 0.163 -0.273 0.785 -0.364 0.275
## 21 0.225 2.210 0.027 0.056 0.9392.b output variable = PA
The floowing is the SEM analysis for the output variable PA.
library(lavaan)
# 2 x variables + 2 mediators
model2 <-'
# outcome model
PA ~ c1*Mastery+c2*Performance+b1*expectancybelief+b2*taskvalue
# mediator models
expectancybelief ~ a11*Mastery+a12*Performance
taskvalue ~ a21*Mastery+a22*Performance
# indirect effects
MasteryExpectancy := b1*a11
Masterytaskvalue := b2*a21
MasteryIDE := b1*a11 + b2*a21
PerformancyExpectancy := b1*a12
PerformancyTaskvalue := b2*a22
PerfIDE := b1*a12+b2*a22
expectancybelief ~~ taskvalue # model correlation between mediators
'
# SEM analysis
fitsem2 <- sem(model2,data=xym,group=NULL)
summary(fitsem2, fit.measures=TRUE, standardize=TRUE, rsquare=TRUE)## lavaan (0.5-16) converged normally after 21 iterations
##
## Number of observations 336
##
## Estimator ML
## Minimum Function Test Statistic 0.000
## Degrees of freedom 0
## P-value (Chi-square) 1.000
##
## Model test baseline model:
##
## Minimum Function Test Statistic 241.504
## Degrees of freedom 9
## P-value 0.000
##
## User model versus baseline model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.000
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -1553.011
## Loglikelihood unrestricted model (H1) -1553.011
##
## Number of free parameters 12
## Akaike (AIC) 3130.022
## Bayesian (BIC) 3175.827
## Sample-size adjusted Bayesian (BIC) 3137.761
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent Confidence Interval 0.000 0.000
## P-value RMSEA <= 0.05 1.000
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.000
##
## Parameter estimates:
##
## Information Expected
## Standard Errors Standard
##
## Estimate Std.err Z-value P(>|z|) Std.lv Std.all
## Regressions:
## PA ~
## Mastery (c1) -0.039 0.061 -0.648 0.517 -0.039 -0.035
## Prfrmnc (c2) 0.020 0.040 0.503 0.615 0.020 0.026
## expctnc (b1) 0.356 0.080 4.421 0.000 0.356 0.282
## taskval (b2) 0.156 0.081 1.914 0.056 0.156 0.121
## expectancybelief ~
## Mastery (a11) 0.262 0.045 5.779 0.000 0.262 0.297
## Prfrmnc (a12) 0.102 0.031 3.262 0.001 0.102 0.168
## taskvalue ~
## Mastery (a21) 0.245 0.045 5.471 0.000 0.245 0.282
## Prfrmnc (a22) 0.102 0.031 3.305 0.001 0.102 0.170
##
## Covariances:
## expectancybelief ~~
## taskvalue 0.137 0.016 8.424 0.000 0.137 0.517
##
## Variances:
## PA 0.427 0.033 0.427 0.869
## expectancyblf 0.268 0.021 0.268 0.871
## taskvalue 0.262 0.020 0.262 0.879
##
## Defined parameters:
## MstryExpctncy 0.093 0.026 3.511 0.000 0.093 0.084
## Masterytaskvl 0.038 0.021 1.807 0.071 0.038 0.034
## MasteryIDE 0.131 0.029 4.553 0.000 0.131 0.118
## PrfrmncyExpct 0.036 0.014 2.625 0.009 0.036 0.047
## PrfrmncyTskvl 0.016 0.010 1.657 0.098 0.016 0.021
## PerfIDE 0.052 0.016 3.187 0.001 0.052 0.068
##
## R-Square:
##
## PA 0.131
## expectancybelief 0.129
## taskvalue 0.121# mediator effect
boot.fit2 <- parameterEstimates(fitsem2, boot.ci.type="bca.simple")
boot.fit2## lhs op rhs label est
## 1 PA ~ Mastery c1 -0.039
## 2 PA ~ Performance c2 0.020
## 3 PA ~ expectancybelief b1 0.356
## 4 PA ~ taskvalue b2 0.156
## 5 expectancybelief ~ Mastery a11 0.262
## 6 expectancybelief ~ Performance a12 0.102
## 7 taskvalue ~ Mastery a21 0.245
## 8 taskvalue ~ Performance a22 0.102
## 9 expectancybelief ~~ taskvalue 0.137
## 10 PA ~~ PA 0.427
## 11 expectancybelief ~~ expectancybelief 0.268
## 12 taskvalue ~~ taskvalue 0.262
## 13 Mastery ~~ Mastery 0.396
## 14 Mastery ~~ Performance 0.074
## 15 Performance ~~ Performance 0.832
## 16 MasteryExpectancy := b1*a11 MasteryExpectancy 0.093
## 17 Masterytaskvalue := b2*a21 Masterytaskvalue 0.038
## 18 MasteryIDE := b1*a11+b2*a21 MasteryIDE 0.131
## 19 PerformancyExpectancy := b1*a12 PerformancyExpectancy 0.036
## 20 PerformancyTaskvalue := b2*a22 PerformancyTaskvalue 0.016
## 21 PerfIDE := b1*a12+b2*a22 PerfIDE 0.052
## se z pvalue ci.lower ci.upper
## 1 0.061 -0.648 0.517 -0.158 0.079
## 2 0.040 0.503 0.615 -0.059 0.099
## 3 0.080 4.421 0.000 0.198 0.513
## 4 0.081 1.914 0.056 -0.004 0.315
## 5 0.045 5.779 0.000 0.173 0.350
## 6 0.031 3.262 0.001 0.041 0.163
## 7 0.045 5.471 0.000 0.157 0.332
## 8 0.031 3.305 0.001 0.042 0.163
## 9 0.016 8.424 0.000 0.105 0.169
## 10 0.033 12.961 0.000 0.363 0.492
## 11 0.021 12.961 0.000 0.228 0.309
## 12 0.020 12.961 0.000 0.222 0.302
## 13 0.000 NA NA 0.396 0.396
## 14 0.000 NA NA 0.074 0.074
## 15 0.000 NA NA 0.832 0.832
## 16 0.026 3.511 0.000 0.041 0.145
## 17 0.021 1.807 0.071 -0.003 0.080
## 18 0.029 4.553 0.000 0.075 0.188
## 19 0.014 2.625 0.009 0.009 0.063
## 20 0.010 1.657 0.098 -0.003 0.035
## 21 0.016 3.187 0.001 0.020 0.0842.c output variable = Av.Stepsmin
The floowing is the SEM analysis for the output variable Av.Stepsmin.
library(lavaan)
# 2 x variables + 2 mediators
model3 <-'
# outcome model
Av.Stepsmin ~ c1*Mastery+c2*Performance+b1*expectancybelief+b2*taskvalue
# mediator models
expectancybelief ~ a11*Mastery+a12*Performance
taskvalue ~ a21*Mastery+a22*Performance
# indirect effects
MasteryExpectancy := b1*a11
Masterytaskvalue := b2*a21
MasteryIDE := b1*a11 + b2*a21
PerformancyExpectancy := b1*a12
PerformancyTaskvalue := b2*a22
PerfIDE := b1*a12+b2*a22
expectancybelief ~~ taskvalue # model correlation between mediators
'
# SEM analysis
fitsem3 <- sem(model3,data=xym,group=NULL)
summary(fitsem3, fit.measures=TRUE, standardize=TRUE, rsquare=TRUE)## lavaan (0.5-16) converged normally after 60 iterations
##
## Number of observations 336
##
## Estimator ML
## Minimum Function Test Statistic 0.000
## Degrees of freedom 0
## P-value (Chi-square) 0.000
##
## Model test baseline model:
##
## Minimum Function Test Statistic 233.489
## Degrees of freedom 9
## P-value 0.000
##
## User model versus baseline model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.000
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -2799.616
## Loglikelihood unrestricted model (H1) -2799.616
##
## Number of free parameters 12
## Akaike (AIC) 5623.232
## Bayesian (BIC) 5669.038
## Sample-size adjusted Bayesian (BIC) 5630.972
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent Confidence Interval 0.000 0.000
## P-value RMSEA <= 0.05 1.000
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.000
##
## Parameter estimates:
##
## Information Expected
## Standard Errors Standard
##
## Estimate Std.err Z-value P(>|z|) Std.lv Std.all
## Regressions:
## Av.Stepsmin ~
## Mastery (c1) 8.063 2.475 3.258 0.001 8.063 0.179
## Prfrmnc (c2) -7.598 1.645 -4.620 0.000 -7.598 -0.245
## expctnc (b1) -1.818 3.288 -0.553 0.580 -1.818 -0.036
## taskval (b2) -7.866 3.327 -2.364 0.018 -7.866 -0.152
## expectancybelief ~
## Mastery (a11) 0.262 0.045 5.779 0.000 0.262 0.297
## Prfrmnc (a12) 0.102 0.031 3.262 0.001 0.102 0.168
## taskvalue ~
## Mastery (a21) 0.245 0.045 5.471 0.000 0.245 0.282
## Prfrmnc (a22) 0.102 0.031 3.305 0.001 0.102 0.170
##
## Covariances:
## expectancybelief ~~
## taskvalue 0.137 0.016 8.424 0.000 0.137 0.517
##
## Variances:
## Av.Stepsmin 713.175 55.023 713.175 0.890
## expectancyblf 0.268 0.021 0.268 0.871
## taskvalue 0.262 0.020 0.262 0.879
##
## Defined parameters:
## MstryExpctncy -0.475 0.864 -0.550 0.582 -0.475 -0.011
## Masterytaskvl -1.925 0.887 -2.170 0.030 -1.925 -0.043
## MasteryIDE -2.400 0.916 -2.622 0.009 -2.400 -0.053
## PrfrmncyExpct -0.185 0.340 -0.545 0.586 -0.185 -0.006
## PrfrmncyTskvl -0.803 0.417 -1.923 0.054 -0.803 -0.026
## PerfIDE -0.988 0.432 -2.289 0.022 -0.988 -0.032
##
## R-Square:
##
## Av.Stepsmin 0.110
## expectancybelief 0.129
## taskvalue 0.121# mediator effect
boot.fit3 <- parameterEstimates(fitsem3, boot.ci.type="bca.simple")
boot.fit3## lhs op rhs label est
## 1 Av.Stepsmin ~ Mastery c1 8.063
## 2 Av.Stepsmin ~ Performance c2 -7.598
## 3 Av.Stepsmin ~ expectancybelief b1 -1.818
## 4 Av.Stepsmin ~ taskvalue b2 -7.866
## 5 expectancybelief ~ Mastery a11 0.262
## 6 expectancybelief ~ Performance a12 0.102
## 7 taskvalue ~ Mastery a21 0.245
## 8 taskvalue ~ Performance a22 0.102
## 9 expectancybelief ~~ taskvalue 0.137
## 10 Av.Stepsmin ~~ Av.Stepsmin 713.175
## 11 expectancybelief ~~ expectancybelief 0.268
## 12 taskvalue ~~ taskvalue 0.262
## 13 Mastery ~~ Mastery 0.396
## 14 Mastery ~~ Performance 0.074
## 15 Performance ~~ Performance 0.832
## 16 MasteryExpectancy := b1*a11 MasteryExpectancy -0.475
## 17 Masterytaskvalue := b2*a21 Masterytaskvalue -1.925
## 18 MasteryIDE := b1*a11+b2*a21 MasteryIDE -2.400
## 19 PerformancyExpectancy := b1*a12 PerformancyExpectancy -0.185
## 20 PerformancyTaskvalue := b2*a22 PerformancyTaskvalue -0.803
## 21 PerfIDE := b1*a12+b2*a22 PerfIDE -0.988
## se z pvalue ci.lower ci.upper
## 1 2.475 3.258 0.001 3.212 12.914
## 2 1.645 -4.620 0.000 -10.822 -4.375
## 3 3.288 -0.553 0.580 -8.263 4.627
## 4 3.327 -2.364 0.018 -14.386 -1.345
## 5 0.045 5.779 0.000 0.173 0.350
## 6 0.031 3.262 0.001 0.041 0.163
## 7 0.045 5.471 0.000 0.157 0.332
## 8 0.031 3.305 0.001 0.042 0.163
## 9 0.016 8.424 0.000 0.105 0.169
## 10 55.023 12.961 0.000 605.333 821.017
## 11 0.021 12.961 0.000 0.228 0.309
## 12 0.020 12.961 0.000 0.222 0.302
## 13 0.000 NA NA 0.396 0.396
## 14 0.000 NA NA 0.074 0.074
## 15 0.000 NA NA 0.832 0.832
## 16 0.864 -0.550 0.582 -2.169 1.218
## 17 0.887 -2.170 0.030 -3.663 -0.187
## 18 0.916 -2.622 0.009 -4.195 -0.606
## 19 0.340 -0.545 0.586 -0.851 0.481
## 20 0.417 -1.923 0.054 -1.620 0.015
## 21 0.432 -2.289 0.022 -1.834 -0.142