Show code
library(lavaan)
library(semPlot)
library(semTools)Lab 05 — SEM capstone with the EAT dataset
Goals
This lab is a full SEM workflow exercise (measurement + structure):
- Apply the Two-Step rule (identify/fit measurement first, then structure)
- Evaluate the model using fit indices + diagnostics (global + local)
- Use MI + EPC/SEPC to propose a theory-justified modification
- Compare a latent SEM vs a sum-score path model
You will produce: - a short “model checking log” (what you checked, what you changed, why) - a figure of your final model (R / PowerPoint / hand-drawn)
Setup
Data: the EAT dataset
The dataset includes 13 items measuring:
- peer pressure (PP1–PP4)
- social media use (SM1–SM4)
- social comparison (SC1–SC3)
- eating disorder (ED1–ED2)
Show code
# The object in the .Rdata file is dE4_1 (as in the original slides)
load("../data/Exercise4_1.Rdata")
# Inspect
str(dE4_1)'data.frame': 1423 obs. of 13 variables:
$ PP1: num 0.228 -0.647 0.259 -0.685 -0.064 ...
$ PP2: num 1.909 -3.185 -1.463 -0.289 0.889 ...
$ PP3: num 0.5869 -0.8423 -0.4295 -0.076 0.0449 ...
$ PP4: num -0.7008 -0.4283 -0.0468 -2.2361 -0.4145 ...
$ SM1: num 2.312 0.686 -0.21 -0.642 0.171 ...
$ SM2: num 1.062 0.105 -0.459 0.468 1.024 ...
$ SM3: num 1.4272 1.2365 -0.0203 -0.7571 0.1755 ...
$ SM4: num -0.771 0.509 1.258 0.584 0.609 ...
$ SC1: num 1.18 -1.87 1.39 1.22 3.74 ...
$ SC2: num 2.2 -1.174 1.03 0.967 1.382 ...
$ SC3: num 0.177 -1.233 0.901 2.734 1.852 ...
$ ED1: num -0.385 -1.739 1.074 0.707 1.253 ...
$ ED2: num 1.606 -1.734 2.041 0.32 0.664 ...Show code
round(head(dE4_1), 2) PP1 PP2 PP3 PP4 SM1 SM2 SM3 SM4 SC1 SC2 SC3 ED1 ED2
1 0.23 1.91 0.59 -0.70 2.31 1.06 1.43 -0.77 1.18 2.20 0.18 -0.38 1.61
2 -0.65 -3.18 -0.84 -0.43 0.69 0.10 1.24 0.51 -1.87 -1.17 -1.23 -1.74 -1.73
3 0.26 -1.46 -0.43 -0.05 -0.21 -0.46 -0.02 1.26 1.39 1.03 0.90 1.07 2.04
4 -0.68 -0.29 -0.08 -2.24 -0.64 0.47 -0.76 0.58 1.22 0.97 2.73 0.71 0.32
5 -0.06 0.89 0.04 -0.41 0.17 1.02 0.18 0.61 3.74 1.38 1.85 1.25 0.66
6 -0.29 -2.74 -0.52 2.58 0.15 -1.08 1.08 0.99 1.32 -0.24 0.93 -0.97 -0.68
Theoretical model (concept)
- peer pressure + social media → social comparison → eating disorder
You will test this model as a latent SEM (CFA + structure), then as a sum-score path model.
Exercise 1 — STEP 1 (Two-Step rule): fit the CFA model
1a) Specify the measurement model
Show code
m_cfa <- "
peerPressure =~ PP1 + PP2 + PP3 + PP4
socialMedia =~ SM1 + SM2 + SM3 + SM4
socialComparison =~ SC1 + SC2 + SC3
eatingDisorder =~ ED1 + ED2
"1b) Fit and check convergence
Show code
fit_cfa <- sem(m_cfa, data = dE4_1, std.lv = TRUE)
lavInspect(fit_cfa, "converged")[1] TRUE1c) Inspect measurement results
Show code
summary(fit_cfa, standardized = TRUE, fit.measures = TRUE)lavaan 0.6-19 ended normally after 30 iterations
Estimator ML
Optimization method NLMINB
Number of model parameters 32
Number of observations 1423
Model Test User Model:
Test statistic 421.825
Degrees of freedom 59
P-value (Chi-square) 0.000
Model Test Baseline Model:
Test statistic 3150.947
Degrees of freedom 78
P-value 0.000
User Model versus Baseline Model:
Comparative Fit Index (CFI) 0.882
Tucker-Lewis Index (TLI) 0.844
Loglikelihood and Information Criteria:
Loglikelihood user model (H0) -28616.602
Loglikelihood unrestricted model (H1) -28405.690
Akaike (AIC) 57297.204
Bayesian (BIC) 57465.541
Sample-size adjusted Bayesian (SABIC) 57363.888
Root Mean Square Error of Approximation:
RMSEA 0.066
90 Percent confidence interval - lower 0.060
90 Percent confidence interval - upper 0.072
P-value H_0: RMSEA <= 0.050 0.000
P-value H_0: RMSEA >= 0.080 0.000
Standardized Root Mean Square Residual:
SRMR 0.042
Parameter Estimates:
Standard errors Standard
Information Expected
Information saturated (h1) model Structured
Latent Variables:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
peerPressure =~
PP1 1.016 0.036 27.977 0.000 1.016 0.821
PP2 0.998 0.037 26.980 0.000 0.998 0.786
PP3 0.424 0.033 12.944 0.000 0.424 0.370
PP4 0.468 0.034 13.945 0.000 0.468 0.397
socialMedia =~
SM1 0.592 0.039 15.294 0.000 0.592 0.517
SM2 0.575 0.040 14.449 0.000 0.575 0.486
SM3 0.554 0.039 14.234 0.000 0.554 0.479
SM4 0.645 0.041 15.928 0.000 0.645 0.541
socialComparison =~
SC1 0.881 0.038 23.042 0.000 0.881 0.662
SC2 0.940 0.038 24.865 0.000 0.940 0.718
SC3 0.856 0.039 22.167 0.000 0.856 0.637
eatingDisorder =~
ED1 0.720 0.067 10.709 0.000 0.720 0.596
ED2 0.671 0.064 10.491 0.000 0.671 0.547
Covariances:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
peerPressure ~~
socialMedia -0.013 0.039 -0.345 0.730 -0.013 -0.013
socialComparsn 0.126 0.035 3.636 0.000 0.126 0.126
eatingDisorder -0.035 0.042 -0.822 0.411 -0.035 -0.035
socialMedia ~~
socialComparsn 0.422 0.037 11.347 0.000 0.422 0.422
eatingDisorder 0.052 0.049 1.056 0.291 0.052 0.052
socialComparison ~~
eatingDisorder 0.368 0.042 8.679 0.000 0.368 0.368
Variances:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
.PP1 0.500 0.053 9.365 0.000 0.500 0.326
.PP2 0.615 0.054 11.489 0.000 0.615 0.382
.PP3 1.133 0.044 25.604 0.000 1.133 0.863
.PP4 1.175 0.046 25.410 0.000 1.175 0.843
.SM1 0.958 0.048 20.116 0.000 0.958 0.733
.SM2 1.068 0.051 21.126 0.000 1.068 0.763
.SM3 1.033 0.048 21.355 0.000 1.033 0.771
.SM4 1.005 0.052 19.225 0.000 1.005 0.707
.SC1 0.993 0.054 18.400 0.000 0.993 0.561
.SC2 0.831 0.053 15.544 0.000 0.831 0.485
.SC3 1.075 0.055 19.543 0.000 1.075 0.595
.ED1 0.942 0.094 9.995 0.000 0.942 0.645
.ED2 1.052 0.085 12.310 0.000 1.052 0.700
peerPressure 1.000 1.000 1.000
socialMedia 1.000 1.000 1.000
socialComparsn 1.000 1.000 1.000
eatingDisorder 1.000 1.000 1.000Tasks
- Report the standardized loadings and identify any weak items.
- Check residual variances (Heywood cases? negative variances?).
- Extract the core fit indices: CFI, TLI, RMSEA, SRMR.
Show code
fitMeasures(fit_cfa, c("cfi", "tli", "rmsea", "srmr")) cfi tli rmsea srmr
0.882 0.844 0.066 0.042 Exercise 2 — STEP 2 (Two-Step rule): fit the structural model
2a) Add structural regressions among the latent variables
Show code
m_sem <- "
# CFA model
peerPressure =~ PP1 + PP2 + PP3 + PP4
socialMedia =~ SM1 + SM2 + SM3 + SM4
socialComparison =~ SC1 + SC2 + SC3
eatingDisorder =~ ED1 + ED2
# Structural model
eatingDisorder ~ socialComparison
socialComparison ~ peerPressure + socialMedia
"2b) Fit and inspect
Show code
fit_sem <- sem(m_sem, data = dE4_1, std.lv = TRUE)
summary(fit_sem, standardized = TRUE, fit.measures = TRUE)lavaan 0.6-19 ended normally after 30 iterations
Estimator ML
Optimization method NLMINB
Number of model parameters 30
Number of observations 1423
Model Test User Model:
Test statistic 430.575
Degrees of freedom 61
P-value (Chi-square) 0.000
Model Test Baseline Model:
Test statistic 3150.947
Degrees of freedom 78
P-value 0.000
User Model versus Baseline Model:
Comparative Fit Index (CFI) 0.880
Tucker-Lewis Index (TLI) 0.846
Loglikelihood and Information Criteria:
Loglikelihood user model (H0) -28620.977
Loglikelihood unrestricted model (H1) -28405.690
Akaike (AIC) 57301.954
Bayesian (BIC) 57459.770
Sample-size adjusted Bayesian (SABIC) 57364.470
Root Mean Square Error of Approximation:
RMSEA 0.065
90 Percent confidence interval - lower 0.060
90 Percent confidence interval - upper 0.071
P-value H_0: RMSEA <= 0.050 0.000
P-value H_0: RMSEA >= 0.080 0.000
Standardized Root Mean Square Residual:
SRMR 0.044
Parameter Estimates:
Standard errors Standard
Information Expected
Information saturated (h1) model Structured
Latent Variables:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
peerPressure =~
PP1 1.019 0.036 27.954 0.000 1.019 0.823
PP2 0.996 0.037 26.873 0.000 0.996 0.785
PP3 0.423 0.033 12.909 0.000 0.423 0.369
PP4 0.468 0.034 13.932 0.000 0.468 0.396
socialMedia =~
SM1 0.582 0.039 15.049 0.000 0.582 0.509
SM2 0.576 0.040 14.452 0.000 0.576 0.487
SM3 0.557 0.039 14.279 0.000 0.557 0.481
SM4 0.653 0.041 16.062 0.000 0.653 0.548
socialComparison =~
SC1 0.799 0.037 21.609 0.000 0.885 0.665
SC2 0.847 0.037 22.797 0.000 0.937 0.716
SC3 0.777 0.037 20.925 0.000 0.860 0.640
eatingDisorder =~
ED1 0.654 0.066 9.920 0.000 0.700 0.579
ED2 0.645 0.065 9.945 0.000 0.690 0.563
Regressions:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
eatingDisorder ~
socialComparsn 0.344 0.047 7.322 0.000 0.356 0.356
socialComparison ~
peerPressure 0.139 0.038 3.636 0.000 0.125 0.125
socialMedia 0.455 0.049 9.221 0.000 0.411 0.411
Covariances:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
peerPressure ~~
socialMedia -0.014 0.039 -0.354 0.723 -0.014 -0.014
Variances:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
.PP1 0.495 0.054 9.205 0.000 0.495 0.323
.PP2 0.619 0.054 11.521 0.000 0.619 0.384
.PP3 1.134 0.044 25.612 0.000 1.134 0.864
.PP4 1.175 0.046 25.414 0.000 1.175 0.843
.SM1 0.969 0.048 20.363 0.000 0.969 0.741
.SM2 1.067 0.051 21.069 0.000 1.067 0.762
.SM3 1.030 0.048 21.257 0.000 1.030 0.768
.SM4 0.995 0.053 18.929 0.000 0.995 0.700
.SC1 0.986 0.054 18.198 0.000 0.986 0.558
.SC2 0.836 0.054 15.566 0.000 0.836 0.488
.SC3 1.068 0.055 19.361 0.000 1.068 0.591
.ED1 0.971 0.094 10.279 0.000 0.971 0.665
.ED2 1.026 0.093 11.017 0.000 1.026 0.683
peerPressure 1.000 1.000 1.000
socialMedia 1.000 1.000 1.000
.socialComparsn 1.000 0.816 0.816
.eatingDisorder 1.000 0.873 0.873Show code
fitMeasures(fit_sem, c("cfi", "tli", "rmsea", "srmr")) cfi tli rmsea srmr
0.880 0.846 0.065 0.044 Questions
- Are the key hypotheses supported (sign + magnitude + uncertainty)?
- How does fit compare to the CFA-only model? (Interpret carefully.)
(Optional)
Show code
anova(fit_cfa, fit_sem)
Chi-Squared Difference Test
Df AIC BIC Chisq Chisq diff RMSEA Df diff Pr(>Chisq)
fit_cfa 59 57297 57466 422
fit_sem 61 57302 57460 431 8.75 0.0487 2 0.013 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Exercise 3 — Diagnostics: residuals, MI, EPC/SEPC
3a) Standardized residual covariances (local misfit)
Show code
res_std <- residuals(fit_sem, type = "standardized")$cov
res_std PP1 PP2 PP3 PP4 SM1 SM2 SM3 SM4 SC1 SC2
PP1 0.000
PP2 4.114 0.000
PP3 -0.956 -1.461 0.000
PP4 -1.782 -0.714 3.843 0.000
SM1 -0.025 0.866 0.527 -1.447 0.000
SM2 -0.201 0.135 1.297 -0.472 -0.979 0.000
SM3 0.527 0.490 0.127 -0.211 -3.014 1.973 0.000
SM4 -1.102 -0.389 1.354 -0.710 -2.151 1.660 2.005 0.000
SC1 -0.195 -0.601 0.677 0.508 -0.370 -2.178 -0.679 -0.570 0.000
SC2 0.298 -0.359 1.676 0.657 14.915 -1.812 -0.951 -1.042 -3.833 0.000
SC3 0.244 0.045 1.346 1.211 -3.254 -3.230 -1.506 -2.873 4.203 -1.867
ED1 -0.693 -0.951 -2.649 -0.337 -2.138 -2.846 -1.096 -1.179 1.416 -1.514
ED2 -1.299 -2.440 -0.893 -1.020 -0.471 -0.566 -1.018 0.947 0.108 -1.475
SC3 ED1 ED2
PP1
PP2
PP3
PP4
SM1
SM2
SM3
SM4
SC1
SC2
SC3 0.000
ED1 1.659 0.000
ED2 1.978 0.000 0.000Find the largest absolute residual covariances:
Show code
R <- res_std
diag(R) <- NA
top_res <- as.data.frame(as.table(R))
top_res <- top_res[order(abs(top_res$Freq), decreasing = TRUE), ]
head(top_res, 10) Var1 Var2 Freq
122 SM1 SC2 14.92
62 SC2 SM1 14.92
115 SC3 SC1 4.20
139 SC1 SC3 4.20
2 PP2 PP1 4.11
14 PP1 PP2 4.11
30 PP4 PP3 3.84
42 PP3 PP4 3.84
114 SC2 SC1 -3.83
126 SC1 SC2 -3.833b) Modification indices (MI) + effect size (EPC/SEPC)
Show code
mi <- modificationIndices(fit_sem, sort. = TRUE)
head(mi[, c("lhs","op","rhs","mi","epc","sepc.all")], 15) lhs op rhs mi epc sepc.all
120 SM1 ~~ SC2 328.97 0.565 0.628
49 socialMedia =~ SC2 56.91 0.378 0.289
57 socialComparison =~ SM1 54.67 0.292 0.283
121 SM1 ~~ SC3 51.71 -0.235 -0.231
50 socialMedia =~ SC3 30.54 -0.277 -0.206
119 SM1 ~~ SC1 22.07 -0.151 -0.154
143 SC1 ~~ SC3 19.17 0.293 0.285
74 PP1 ~~ PP2 18.34 0.665 1.201
97 PP3 ~~ PP4 15.33 0.126 0.109
58 socialComparison =~ SM2 13.08 -0.146 -0.137
142 SC1 ~~ SC2 12.91 -0.277 -0.304
138 SM4 ~~ SC2 9.69 -0.101 -0.110
117 SM1 ~~ SM3 8.09 -0.112 -0.112
51 socialMedia =~ ED1 8.02 -0.130 -0.108
105 PP3 ~~ ED1 7.86 -0.091 -0.087Tasks
- Do the top residual pairs match the top MI suggestions?
- Pick one candidate modification that you can justify substantively.
- Write a 1–2 sentence justification (content overlap? method effect? plausible local dependence?).
Exercise 4 — Model modification (one change), refit, and re-check
In the original materials, a key misspecification is a residual correlation:
SM1 ~~ SC2
4a) Add one residual correlation and refit
Show code
m_sem_mod <- "
# CFA model
peerPressure =~ PP1 + PP2 + PP3 + PP4
socialMedia =~ SM1 + SM2 + SM3 + SM4
socialComparison =~ SC1 + SC2 + SC3
eatingDisorder =~ ED1 + ED2
# Structural model
eatingDisorder ~ socialComparison
socialComparison ~ peerPressure + socialMedia
# Residual correlation (one modification)
SM1 ~~ SC2
"
fit_sem_mod <- sem(m_sem_mod, data = dE4_1, std.lv = TRUE)4b) Compare fit + nested test
Show code
rbind(
original = fitMeasures(fit_sem, c("cfi","tli","rmsea","srmr")),
modified = fitMeasures(fit_sem_mod, c("cfi","tli","rmsea","srmr"))
) cfi tli rmsea srmr
original 0.880 0.846 0.0653 0.0440
modified 0.995 0.994 0.0130 0.0233Show code
anova(fit_sem, fit_sem_mod)
Chi-Squared Difference Test
Df AIC BIC Chisq Chisq diff RMSEA Df diff Pr(>Chisq)
fit_sem_mod 60 56948 57111 74.4
fit_sem 61 57302 57460 430.6 356 0.5 1 <2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 14c) Re-check diagnostics
Show code
mi_mod <- modificationIndices(fit_sem_mod, sort. = TRUE)
head(mi_mod[, c("lhs","op","rhs","mi","epc","sepc.all")], 10) lhs op rhs mi epc sepc.all
75 PP1 ~~ PP2 19.27 0.679 1.226
98 PP3 ~~ PP4 15.29 0.126 0.109
106 PP3 ~~ ED1 8.04 -0.092 -0.089
51 socialMedia =~ SC3 7.89 -0.127 -0.094
120 SM1 ~~ SC1 5.37 0.075 0.074
52 socialMedia =~ ED1 4.89 -0.092 -0.076
129 SM2 ~~ ED1 4.21 -0.068 -0.068
160 peerPressure ~ eatingDisorder 3.99 -0.090 -0.097
153 peerPressure ~~ eatingDisorder 3.99 -0.090 -0.090
69 eatingDisorder =~ SM2 3.93 -0.077 -0.070Documentation (write these down in your log)
- What did you change and why?
- What evidence supported it (residuals? MI+EPC?)?
- Does the modification change interpretation of the constructs or paths?
Exercise 5 — Draw the final model
Use one of the following:
- R (semPlot)
- PowerPoint / Keynote
- Hand-drawn
Show code
semPaths(
fit_sem_mod,
style = "lisrel",
layout = "tree2",
rotation = 2,
residuals = FALSE,
nCharNodes = 0
)Exercise 6 — Sum-score path model (manifest approximation)
Now reproduce the structural model using sum scores (no latent variables).
Show code
d_sum <- data.frame(
peerPressure = dE4_1$PP1 + dE4_1$PP2 + dE4_1$PP3 + dE4_1$PP4,
socialMedia = dE4_1$SM1 + dE4_1$SM2 + dE4_1$SM3 + dE4_1$SM4,
socialComparison = dE4_1$SC1 + dE4_1$SC2 + dE4_1$SC3,
eatingDisorder = dE4_1$ED1 + dE4_1$ED2
)
m_path <- "
eatingDisorder ~ socialComparison
socialComparison ~ peerPressure + socialMedia
"
fit_path <- sem(m_path, data = d_sum)
summary(fit_path, standardized = TRUE, fit.measures = TRUE)lavaan 0.6-19 ended normally after 1 iteration
Estimator ML
Optimization method NLMINB
Number of model parameters 5
Number of observations 1423
Model Test User Model:
Test statistic 5.993
Degrees of freedom 2
P-value (Chi-square) 0.050
Model Test Baseline Model:
Test statistic 189.127
Degrees of freedom 5
P-value 0.000
User Model versus Baseline Model:
Comparative Fit Index (CFI) 0.978
Tucker-Lewis Index (TLI) 0.946
Loglikelihood and Information Criteria:
Loglikelihood user model (H0) -6565.599
Loglikelihood unrestricted model (H1) -6562.602
Akaike (AIC) 13141.197
Bayesian (BIC) 13167.500
Sample-size adjusted Bayesian (SABIC) 13151.616
Root Mean Square Error of Approximation:
RMSEA 0.037
90 Percent confidence interval - lower 0.000
90 Percent confidence interval - upper 0.074
P-value H_0: RMSEA <= 0.050 0.663
P-value H_0: RMSEA >= 0.080 0.024
Standardized Root Mean Square Residual:
SRMR 0.019
Parameter Estimates:
Standard errors Standard
Information Expected
Information saturated (h1) model Structured
Regressions:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
eatingDisorder ~
socialComparsn 0.139 0.016 8.623 0.000 0.139 0.223
socialComparison ~
peerPressure 0.097 0.023 4.151 0.000 0.097 0.106
socialMedia 0.258 0.026 9.917 0.000 0.258 0.253
Variances:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
.eatingDisorder 3.736 0.140 26.674 0.000 3.736 0.950
.socialComparsn 9.339 0.350 26.674 0.000 9.339 0.925Show code
fitMeasures(fit_path, c("cfi","tli","rmsea","srmr")) cfi tli rmsea srmr
0.978 0.946 0.037 0.019 Exercise 7 — Compare latent SEM vs sum-score path results
7a) Compare key standardized structural coefficients
Extract standardized regressions from both models:
Show code
pe_lat <- parameterEstimates(fit_sem_mod, standardized = TRUE)
lat_paths <- subset(pe_lat, op == "~" &
lhs %in% c("eatingDisorder","socialComparison"))[,
c("lhs","rhs","est","se","pvalue","std.all")]
lat_paths lhs rhs est se pvalue std.all
14 eatingDisorder socialComparison 0.377 0.049 0 0.371
15 socialComparison peerPressure 0.133 0.036 0 0.126
16 socialComparison socialMedia 0.324 0.045 0 0.306Show code
pe_sum <- parameterEstimates(fit_path, standardized = TRUE)
sum_paths <- subset(pe_sum, op == "~")[, c("lhs","rhs","est","se","pvalue","std.all")]
sum_paths lhs rhs est se pvalue std.all
1 eatingDisorder socialComparison 0.139 0.016 0 0.223
2 socialComparison peerPressure 0.097 0.023 0 0.106
3 socialComparison socialMedia 0.258 0.026 0 0.2537b) Discussion (write a short paragraph)
- Which paths change the most (magnitude / uncertainty)?
- Does the pattern of support for hypotheses change?
- Offer at least two measurement-driven explanations for differences:
- unequal loadings / reliability differences
- correlated residuals / local dependence
- attenuation / error-in-variables logic
- factor correlation vs sum-score correlation differences
Wrap-up
Deliverables
- A 1-page model-checking log:
- CFA fit + key issues
- SEM fit + key issues
- diagnostics used (residuals, MI/EPC)
- what you modified and why
- A diagram of the final SEM
- A short comparison of latent vs sum-score conclusions