lavaanif (!require(lavaan)) install.packages("lavaan")
library(lavaan)Example data: 1336 college students self-reporting on 49 items (measuring five factors) assessing childhood maltreatment: Items are answered on a 1–5 scale: 1=Strongly Disagree, 2=Disagree, 3=Neutral, 4=Agree, 5=Strongly Agree. The items are NOT normally distributed, so we’ll use both CFA with MLR and IFA with WLSMV as two options to examine the fit of these models (as an example of how to do each, but NOT to compare between estimators).
1. Spurning: Verbal and nonverbal caregiver acts that reject and degrade a child
2. Terrorizing: Caregiver behaviors that threaten or are likely to physically hurt, kill, abandon, or place the child or the child’s loved ones or objects in recognizably dangerous situations.
3. Isolating: Caregiver acts that consistently deny the child opportunities to meet needs for interacting or communicating with peers or adults inside or outside the home.
4. Corrupting: Caregiver acts that encourage the child to develop inappropriate behaviors (self-destructive, antisocial, criminal, deviant, or other maladaptive behaviors).
5. Ignoring: Emotional unresponsiveness includes caregiver acts that ignore the child’s attempts and needs to interact (failing to express affection, caring, and love for the child) and show no emotion in interactions with the child
abuseData = read.csv(file = "abuse.csv", col.names = c("ID", paste0("p0",1:9), paste0("p",10:57)))First, we separately build each one-factor model:
spurningSyntax = "
spurn =~ p06 + p10 + p14 + p25 + p27 + p29 + p33 + p35 + p48 + p49 + p53 + p54
"
spurningEstimatesMLR = cfa(model = spurningSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "MLR")
fitResultsMLR = data.frame(Model = "Spurning", rbind(inspect(object = spurningEstimatesMLR, what = "fit")), stringsAsFactors = FALSE)
spurningEstimatesWLSMV = cfa(model = spurningSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "WLSMV",
ordered = c("p06", "p10", "p14", "p25", "p27", "p29", "p33", "p35", "p48", "p49", "p53", "p54"),
parameterization = "theta")lavaan WARNING: 22 bivariate tables have empty cells; to see them, use:
lavInspect(fit, "zero.cell.tables")fitResultsWLSMV = data.frame(Model = "Spurning", rbind(inspect(object = spurningEstimatesWLSMV, what = "fit")), stringsAsFactors = FALSE)
spurningParams = cbind(inspect(object = spurningEstimatesMLR, what = "std")$lambda, inspect(object = spurningEstimatesWLSMV, what = "std")$lambda)
colnames(spurningParams) = c("spurningMLR", "spurningWLSMV")
terrorizingSyntax = "
terror =~ p07 + p11 + p13 + p17 + p24 + p26 + p36 + p55 + p56
"
terrorizingEstimatesMLR = cfa(model = terrorizingSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "MLR")
fitResultsMLR = rbind(fitResultsMLR, c("Terrorizing", inspect(object = terrorizingEstimatesMLR, what = "fit")))
terrorizingEstimatesWLSMV = cfa(model = terrorizingSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "WLSMV",
ordered = c("p07", "p11", "p13", "p17", "p24", "p26", "p36", "p55", "p56"), parameterization = "theta")lavaan WARNING: 20 bivariate tables have empty cells; to see them, use:
lavInspect(fit, "zero.cell.tables")fitResultsWLSMV = rbind(fitResultsWLSMV, c("Terrorizing", inspect(object = terrorizingEstimatesWLSMV, what = "fit")))
terrorizingParams = cbind(inspect(object = terrorizingEstimatesMLR, what = "std")$lambda, inspect(object = terrorizingEstimatesWLSMV, what = "std")$lambda)
colnames(terrorizingParams) = c("terrorizingMLR", "terrorizingWLSMV")
isolatingSyntax = "
isolate =~ p01 + p18 + p19 + p23 + p39 + p43
"
isolatingEstimatesMLR = cfa(model = isolatingSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "MLR")
fitResultsMLR = rbind(fitResultsMLR, c("Isolating", inspect(object = isolatingEstimatesMLR, what = "fit")))
isolatingEstimatesWLSMV = cfa(model = isolatingSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "WLSMV",
ordered = c("p01", "p18", "p19", "p23", "p39", "p43"), parameterization = "theta")lavaan WARNING: 11 bivariate tables have empty cells; to see them, use:
lavInspect(fit, "zero.cell.tables")fitResultsWLSMV = rbind(fitResultsWLSMV, c("Isolating", inspect(object = isolatingEstimatesWLSMV, what = "fit")))
isolatingParams = cbind(inspect(object = isolatingEstimatesMLR, what = "std")$lambda, inspect(object = isolatingEstimatesWLSMV, what = "std")$lambda)
colnames(isolatingParams) = c("isolatingMLR", "isolatingWLSMV")
corruptingSyntax = "
corrupt =~ p09 + p12 + p16 + p20 + p28 + p47 + p50
"
corruptingEstimatesMLR = cfa(model = corruptingSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "MLR")
fitResultsMLR = rbind(fitResultsMLR, c("Corrupting", inspect(object = corruptingEstimatesMLR, what = "fit")))
corruptingEstimatesWLSMV = cfa(model = corruptingSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "WLSMV",
ordered = c("p09", "p12", "p16", "p20", "p28", "p47", "p50"), parameterization = "theta")lavaan WARNING: 20 bivariate tables have empty cells; to see them, use:
lavInspect(fit, "zero.cell.tables")fitResultsWLSMV = rbind(fitResultsWLSMV, c("Corrupting", inspect(object = corruptingEstimatesWLSMV, what = "fit")))
corruptingParams = cbind(inspect(object = corruptingEstimatesMLR, what = "std")$lambda, inspect(object = corruptingEstimatesWLSMV, what = "std")$lambda)
colnames(corruptingParams) = c("corruptingMLR", "corruptingWLSMV")
ignoringSyntax = "
ignore =~ p02 + p03 + p04 + p21 + p22 + p30 + p31 + p37 + p40 + p44 + p45 + p46 + p51 + p52 + p57
"
ignoringEstimatesMLR = cfa(model = ignoringSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "MLR")
fitResultsMLR = rbind(fitResultsMLR, c("Ignoring", inspect(object = ignoringEstimatesMLR, what = "fit")))
ignoringEstimatesWLSMV = cfa(model = ignoringSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "WLSMV",
ordered = c("p02", "p03", "p04", "p21", "p22", "p30", "p31", "p37", "p40", "p44", "p45", "p46", "p51", "p52", "p57"),
parameterization = "theta")lavaan WARNING: 66 bivariate tables have empty cells; to see them, use:
lavInspect(fit, "zero.cell.tables")fitResultsWLSMV = rbind(fitResultsWLSMV, c("Ignoring", inspect(object = ignoringEstimatesWLSMV, what = "fit")))
ignoringParams = cbind(inspect(object = ignoringEstimatesMLR, what = "std")$lambda, inspect(object = ignoringEstimatesWLSMV, what = "std")$lambda)
colnames(ignoringParams) = c("ignoringMLR", "ignoringWLSMV")fitResultsMLR[,c("Model", "chisq.scaled", "chisq.scaling.factor", "df.scaled", "pvalue.scaled", "cfi.scaled", "tli.scaled","rmsea.scaled")]fitResultsWLSMV[,c("Model", "chisq.scaled", "chisq.scaling.factor", "df.scaled", "pvalue.scaled", "cfi.scaled", "tli.scaled","rmsea.scaled")]spurningParams
terrorizingParams
isolatingParams
corruptingParams
ignoringParamscfaHigherSyntax = "
spurn =~ p06 + p10 + p14 + p25 + p27 + p29 + p33 + p35 + p48 + p49 + p53 + p54
terror =~ p07 + p11 + p13 + p17 + p24 + p26 + p36 + p55 + p56
isolate =~ p01 + p18 + p19 + p23 + p39 + p43
corrupt =~ p09 + p12 + p16 + p20 + p28 + p47 + p50
ignore =~ p02 + p03 + p04 + p21 + p22 + p30 + p31 + p37 + p40 + p44 + p45 + p46 + p51 + p52 + p57
abuse =~ spurn + terror + isolate + corrupt + ignore
"
cfaHigherEstimates = cfa(model = cfaHigherSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "MLR")
summary(cfaHigherEstimates, fit.measures = TRUE, rsquare = TRUE, standardized = TRUE)lavaan (0.5-23.1097) converged normally after 78 iterations
Number of observations 1335
Number of missing patterns 1
Estimator ML Robust
Minimum Function Test Statistic 6597.050 4489.494
Degrees of freedom 1122 1122
P-value (Chi-square) 0.000 0.000
Scaling correction factor 1.469
for the Yuan-Bentler correction (Mplus variant)
Model test baseline model:
Minimum Function Test Statistic 35067.550 22808.622
Degrees of freedom 1176 1176
P-value 0.000 0.000
User model versus baseline model:
Comparative Fit Index (CFI) 0.838 0.844
Tucker-Lewis Index (TLI) 0.831 0.837
Robust Comparative Fit Index (CFI) 0.851
Robust Tucker-Lewis Index (TLI) 0.844
Loglikelihood and Information Criteria:
Loglikelihood user model (H0) -69010.792 -69010.792
Scaling correction factor 2.505
for the MLR correction
Loglikelihood unrestricted model (H1) -65712.267 -65712.267
Scaling correction factor 1.593
for the MLR correction
Number of free parameters 152 152
Akaike (AIC) 138325.584 138325.584
Bayesian (BIC) 139115.480 139115.480
Sample-size adjusted Bayesian (BIC) 138632.643 138632.643
Root Mean Square Error of Approximation:
RMSEA 0.060 0.047
90 Percent Confidence Interval 0.059 0.062 0.046 0.049
P-value RMSEA <= 0.05 0.000 1.000
Robust RMSEA 0.057
90 Percent Confidence Interval 0.056 0.059
Standardized Root Mean Square Residual:
SRMR 0.058 0.058
Parameter Estimates:
Information Observed
Standard Errors Robust.huber.white
Latent Variables:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
spurn =~
p06 1.000 0.697 0.579
p10 0.792 0.062 12.710 0.000 0.552 0.445
p14 1.065 0.063 17.019 0.000 0.742 0.763
p25 0.897 0.061 14.749 0.000 0.624 0.524
p27 1.015 0.059 17.274 0.000 0.707 0.594
p29 1.319 0.069 19.075 0.000 0.919 0.796
p33 1.142 0.064 17.880 0.000 0.795 0.824
p35 0.747 0.063 11.848 0.000 0.520 0.512
p48 0.545 0.060 9.101 0.000 0.380 0.568
p49 0.927 0.061 15.296 0.000 0.646 0.664
p53 1.041 0.063 16.503 0.000 0.725 0.681
p54 1.098 0.069 15.908 0.000 0.765 0.628
terror =~
p07 1.000 0.483 0.534
p11 1.341 0.097 13.872 0.000 0.648 0.673
p13 0.622 0.065 9.628 0.000 0.301 0.451
p17 1.070 0.088 12.210 0.000 0.517 0.600
p24 0.610 0.058 10.590 0.000 0.295 0.576
p26 1.247 0.111 11.191 0.000 0.602 0.602
p36 1.228 0.098 12.498 0.000 0.594 0.673
p55 1.589 0.130 12.227 0.000 0.768 0.633
p56 1.793 0.134 13.406 0.000 0.866 0.706
isolate =~
p01 1.000 0.358 0.491
p18 2.139 0.219 9.778 0.000 0.766 0.611
p19 1.209 0.117 10.344 0.000 0.433 0.606
p23 1.685 0.168 10.004 0.000 0.603 0.591
p39 0.903 0.088 10.281 0.000 0.323 0.488
p43 1.557 0.134 11.634 0.000 0.558 0.672
corrupt =~
p09 1.000 0.360 0.602
p12 0.961 0.103 9.367 0.000 0.346 0.541
p16 1.014 0.116 8.772 0.000 0.365 0.368
p20 0.645 0.080 8.086 0.000 0.232 0.497
p28 1.177 0.097 12.150 0.000 0.424 0.624
p47 1.347 0.112 12.030 0.000 0.485 0.614
p50 1.041 0.074 14.039 0.000 0.375 0.649
ignore =~
p02 1.000 0.461 0.681
p03 1.318 0.082 16.102 0.000 0.607 0.653
p04 1.139 0.072 15.741 0.000 0.525 0.651
p21 1.317 0.093 14.221 0.000 0.607 0.717
p22 1.046 0.081 12.922 0.000 0.482 0.474
p30 1.504 0.090 16.642 0.000 0.693 0.743
p31 1.437 0.082 17.502 0.000 0.662 0.841
p37 1.161 0.078 14.957 0.000 0.535 0.708
p40 1.431 0.081 17.589 0.000 0.659 0.807
p44 1.302 0.079 16.483 0.000 0.600 0.764
p45 0.915 0.053 17.137 0.000 0.422 0.670
p46 1.439 0.083 17.405 0.000 0.663 0.822
p51 1.484 0.099 14.975 0.000 0.684 0.700
p52 1.673 0.106 15.727 0.000 0.771 0.753
p57 1.302 0.072 18.105 0.000 0.600 0.823
abuse =~
spurn 1.000 0.971 0.971
terror 0.680 0.064 10.572 0.000 0.952 0.952
isolate 0.494 0.056 8.762 0.000 0.934 0.934
corrupt 0.397 0.049 8.189 0.000 0.745 0.745
ignore 0.577 0.054 10.585 0.000 0.846 0.846
Intercepts:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
.p06 2.520 0.033 76.549 0.000 2.520 2.095
.p10 2.208 0.034 65.045 0.000 2.208 1.780
.p14 1.600 0.027 60.165 0.000 1.600 1.647
.p25 2.029 0.033 62.184 0.000 2.029 1.702
.p27 2.229 0.033 68.385 0.000 2.229 1.872
.p29 1.898 0.032 60.059 0.000 1.898 1.644
.p33 1.601 0.026 60.633 0.000 1.601 1.659
.p35 1.776 0.028 63.917 0.000 1.776 1.749
.p48 1.236 0.018 67.548 0.000 1.236 1.849
.p49 1.649 0.027 61.987 0.000 1.649 1.697
.p53 1.844 0.029 63.324 0.000 1.844 1.733
.p54 1.934 0.033 58.053 0.000 1.934 1.589
.p07 1.622 0.025 65.517 0.000 1.622 1.793
.p11 1.586 0.026 60.218 0.000 1.586 1.648
.p13 1.213 0.018 66.573 0.000 1.213 1.822
.p17 1.493 0.024 63.352 0.000 1.493 1.734
.p24 1.196 0.014 85.313 0.000 1.196 2.335
.p26 2.026 0.027 73.948 0.000 2.026 2.024
.p36 1.459 0.024 60.477 0.000 1.459 1.655
.p55 1.837 0.033 55.295 0.000 1.837 1.513
.p56 1.923 0.034 57.270 0.000 1.923 1.567
.p01 1.303 0.020 65.295 0.000 1.303 1.787
.p18 2.318 0.034 67.527 0.000 2.318 1.848
.p19 1.288 0.020 65.846 0.000 1.288 1.802
.p23 2.022 0.028 72.385 0.000 2.022 1.981
.p39 1.311 0.018 72.292 0.000 1.311 1.979
.p43 1.656 0.023 72.927 0.000 1.656 1.996
.p09 1.246 0.016 76.116 0.000 1.246 2.083
.p12 1.338 0.018 76.389 0.000 1.338 2.091
.p16 1.692 0.027 62.260 0.000 1.692 1.704
.p20 1.109 0.013 86.722 0.000 1.109 2.373
.p28 1.205 0.019 64.736 0.000 1.205 1.772
.p47 1.370 0.022 63.288 0.000 1.370 1.732
.p50 1.184 0.016 74.812 0.000 1.184 2.048
.p02 1.298 0.019 70.090 0.000 1.298 1.918
.p03 1.630 0.025 64.022 0.000 1.630 1.752
.p04 1.573 0.022 71.253 0.000 1.573 1.950
.p21 1.562 0.023 67.411 0.000 1.562 1.845
.p22 1.831 0.028 65.796 0.000 1.831 1.801
.p30 1.706 0.026 66.859 0.000 1.706 1.830
.p31 1.514 0.022 70.256 0.000 1.514 1.923
.p37 1.479 0.021 71.457 0.000 1.479 1.956
.p40 1.467 0.022 65.622 0.000 1.467 1.796
.p44 1.599 0.022 74.349 0.000 1.599 2.035
.p45 1.282 0.017 74.467 0.000 1.282 2.038
.p46 1.502 0.022 68.064 0.000 1.502 1.863
.p51 1.619 0.027 60.522 0.000 1.619 1.656
.p52 1.804 0.028 64.384 0.000 1.804 1.762
.p57 1.378 0.020 69.055 0.000 1.378 1.890
spurn 0.000 0.000 0.000
terror 0.000 0.000 0.000
isolate 0.000 0.000 0.000
corrupt 0.000 0.000 0.000
ignore 0.000 0.000 0.000
abuse 0.000 0.000 0.000
Variances:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
.p06 0.961 0.049 19.570 0.000 0.961 0.665
.p10 1.234 0.048 25.501 0.000 1.234 0.802
.p14 0.394 0.025 15.580 0.000 0.394 0.417
.p25 1.032 0.042 24.741 0.000 1.032 0.726
.p27 0.919 0.040 22.915 0.000 0.919 0.648
.p29 0.489 0.027 17.888 0.000 0.489 0.367
.p33 0.298 0.021 14.357 0.000 0.298 0.321
.p35 0.760 0.045 17.037 0.000 0.760 0.738
.p48 0.303 0.031 9.754 0.000 0.303 0.677
.p49 0.528 0.032 16.668 0.000 0.528 0.559
.p53 0.607 0.039 15.389 0.000 0.607 0.536
.p54 0.897 0.043 21.030 0.000 0.897 0.605
.p07 0.584 0.037 15.883 0.000 0.584 0.715
.p11 0.506 0.034 15.099 0.000 0.506 0.547
.p13 0.353 0.041 8.653 0.000 0.353 0.796
.p17 0.474 0.032 14.885 0.000 0.474 0.640
.p24 0.175 0.018 9.605 0.000 0.175 0.669
.p26 0.639 0.043 14.816 0.000 0.639 0.638
.p36 0.425 0.031 13.811 0.000 0.425 0.547
.p55 0.883 0.049 18.120 0.000 0.883 0.600
.p56 0.754 0.044 17.223 0.000 0.754 0.501
.p01 0.404 0.043 9.317 0.000 0.404 0.759
.p18 0.986 0.045 21.986 0.000 0.986 0.627
.p19 0.323 0.029 11.130 0.000 0.323 0.633
.p23 0.678 0.033 20.665 0.000 0.678 0.651
.p39 0.334 0.035 9.588 0.000 0.334 0.762
.p43 0.378 0.026 14.324 0.000 0.378 0.548
.p09 0.228 0.027 8.392 0.000 0.228 0.637
.p12 0.289 0.030 9.810 0.000 0.289 0.707
.p16 0.853 0.047 18.246 0.000 0.853 0.865
.p20 0.164 0.030 5.404 0.000 0.164 0.753
.p28 0.283 0.036 7.880 0.000 0.283 0.611
.p47 0.390 0.036 10.858 0.000 0.390 0.623
.p50 0.193 0.030 6.448 0.000 0.193 0.579
.p02 0.246 0.028 8.724 0.000 0.246 0.536
.p03 0.497 0.036 13.665 0.000 0.497 0.574
.p04 0.375 0.032 11.878 0.000 0.375 0.576
.p21 0.348 0.025 14.074 0.000 0.348 0.486
.p22 0.801 0.039 20.556 0.000 0.801 0.775
.p30 0.389 0.036 10.728 0.000 0.389 0.448
.p31 0.181 0.019 9.457 0.000 0.181 0.292
.p37 0.285 0.027 10.566 0.000 0.285 0.499
.p40 0.232 0.030 7.729 0.000 0.232 0.348
.p44 0.258 0.021 12.304 0.000 0.258 0.417
.p45 0.218 0.020 11.032 0.000 0.218 0.551
.p46 0.210 0.026 8.198 0.000 0.210 0.324
.p51 0.487 0.036 13.473 0.000 0.487 0.510
.p52 0.454 0.028 16.305 0.000 0.454 0.433
.p57 0.172 0.021 8.308 0.000 0.172 0.323
spurn 0.028 0.009 2.984 0.003 0.058 0.058
terror 0.022 0.005 4.189 0.000 0.093 0.093
isolate 0.016 0.005 3.448 0.001 0.129 0.129
corrupt 0.058 0.010 5.778 0.000 0.445 0.445
ignore 0.060 0.008 7.512 0.000 0.284 0.284
abuse 0.457 0.047 9.730 0.000 1.000 1.000
R-Square:
Estimate
p06 0.335
p10 0.198
p14 0.583
p25 0.274
p27 0.352
p29 0.633
p33 0.679
p35 0.262
p48 0.323
p49 0.441
p53 0.464
p54 0.395
p07 0.285
p11 0.453
p13 0.204
p17 0.360
p24 0.331
p26 0.362
p36 0.453
p55 0.400
p56 0.499
p01 0.241
p18 0.373
p19 0.367
p23 0.349
p39 0.238
p43 0.452
p09 0.363
p12 0.293
p16 0.135
p20 0.247
p28 0.389
p47 0.377
p50 0.421
p02 0.464
p03 0.426
p04 0.424
p21 0.514
p22 0.225
p30 0.552
p31 0.708
p37 0.501
p40 0.652
p44 0.583
p45 0.449
p46 0.676
p51 0.490
p52 0.567
p57 0.677
spurn 0.942
terror 0.907
isolate 0.871
corrupt 0.555
ignore 0.716NOTE: With respect to fit of the structural model, we are now fitting a single higher-order factor INSTEAD OF covariances among the 5 factors.
To test the fit against the saturated (all possible factor correlations model), we can do a −2ΔLL scaled difference test.
anova(cfaNoHighEstimates, cfaHigherEstimates)Scaled Chi Square Difference Test (method = "satorra.bentler.2001")
Df AIC BIC Chisq Chisq diff Df diff Pr(>Chisq)
cfaNoHighEstimates 1117 138229 139045 6490.3
cfaHigherEstimates 1122 138326 139115 6597.1 47.083 5 5.465e-09 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1This higher-order factor model uses 5 fewer parameters (5 higher-order loadings to replace the 10 covariances among the factors).
According to the −2ΔLL scaled difference relative to the previous model,
−2ΔLL (5) = 47.083, p < .0001
trying to reproduce the 5 factor covariances with a single higher-order factor results in a significant decrease in fit. Based on the factor correlations we examined earlier and the standardized higher-order loadings, I’d guess the issue lies with the “corrupting”" factor not being as related to the others.
For the sake of illustration, we can try one more alternative – what if the items were measuring a single factor (i.e., a single score)? Syntax for CFA model with MLR including a single factor instead of a higher-order factor (“smallest model” for comparison):
cfaSingleSyntax = "
abuse =~ p06 + p10 + p14 + p25 + p27 + p29 + p33 + p35 + p48 + p49 + p53 + p54 +
p07 + p11 + p13 + p17 + p24 + p26 + p36 + p55 + p56 + p01 + p18 + p19 +
p23 + p39 + p43 + p09 + p12 + p16 + p20 + p28 + p47 + p50 + p02 + p03 +
p04 + p21 + p22 + p30 + p31 + p37 + p40 + p44 + p45 + p46 + p51 + p52 + p57
"
cfaSingleEstimates = cfa(model = cfaSingleSyntax, data = abuseData, std.lv = FALSE, mimic = "mplus", estimator = "MLR")
summary(cfaSingleEstimates, fit.measures = TRUE, rsquare = TRUE, standardized = TRUE)lavaan (0.5-23.1097) converged normally after 47 iterations
Number of observations 1335
Number of missing patterns 1
Estimator ML Robust
Minimum Function Test Statistic 9209.963 6186.390
Degrees of freedom 1127 1127
P-value (Chi-square) 0.000 0.000
Scaling correction factor 1.489
for the Yuan-Bentler correction (Mplus variant)
Model test baseline model:
Minimum Function Test Statistic 35067.550 22808.622
Degrees of freedom 1176 1176
P-value 0.000 0.000
User model versus baseline model:
Comparative Fit Index (CFI) 0.762 0.766
Tucker-Lewis Index (TLI) 0.751 0.756
Robust Comparative Fit Index (CFI) 0.774
Robust Tucker-Lewis Index (TLI) 0.764
Loglikelihood and Information Criteria:
Loglikelihood user model (H0) -70317.248 -70317.248
Scaling correction factor 2.392
for the MLR correction
Loglikelihood unrestricted model (H1) -65712.267 -65712.267
Scaling correction factor 1.593
for the MLR correction
Number of free parameters 147 147
Akaike (AIC) 140928.496 140928.496
Bayesian (BIC) 141692.409 141692.409
Sample-size adjusted Bayesian (BIC) 141225.455 141225.455
Root Mean Square Error of Approximation:
RMSEA 0.073 0.058
90 Percent Confidence Interval 0.072 0.075 0.057 0.059
P-value RMSEA <= 0.05 0.000 0.000
Robust RMSEA 0.071
90 Percent Confidence Interval 0.069 0.072
Standardized Root Mean Square Residual:
SRMR 0.062 0.062
Parameter Estimates:
Information Observed
Standard Errors Robust.huber.white
Latent Variables:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
abuse =~
p06 1.000 0.640 0.532
p10 0.784 0.066 11.916 0.000 0.502 0.405
p14 1.084 0.068 15.992 0.000 0.694 0.714
p25 0.898 0.064 14.002 0.000 0.575 0.482
p27 1.026 0.062 16.490 0.000 0.657 0.551
p29 1.333 0.075 17.690 0.000 0.854 0.739
p33 1.183 0.071 16.751 0.000 0.757 0.785
p35 0.939 0.077 12.203 0.000 0.601 0.592
p48 0.598 0.066 8.990 0.000 0.383 0.573
p49 0.949 0.064 14.723 0.000 0.608 0.625
p53 1.089 0.068 16.023 0.000 0.697 0.655
p54 1.099 0.074 14.765 0.000 0.704 0.578
p07 0.748 0.067 11.082 0.000 0.479 0.529
p11 0.934 0.075 12.466 0.000 0.598 0.622
p13 0.431 0.061 7.022 0.000 0.276 0.414
p17 0.719 0.067 10.725 0.000 0.460 0.535
p24 0.427 0.053 8.106 0.000 0.273 0.534
p26 0.915 0.060 15.249 0.000 0.586 0.585
p36 0.828 0.067 12.339 0.000 0.530 0.602
p55 1.067 0.069 15.406 0.000 0.683 0.563
p56 1.183 0.075 15.714 0.000 0.758 0.618
p01 0.526 0.060 8.737 0.000 0.337 0.462
p18 1.066 0.064 16.776 0.000 0.683 0.544
p19 0.639 0.066 9.710 0.000 0.409 0.573
p23 0.844 0.061 13.757 0.000 0.540 0.529
p39 0.473 0.055 8.581 0.000 0.303 0.457
p43 0.812 0.061 13.257 0.000 0.520 0.627
p09 0.430 0.051 8.357 0.000 0.275 0.460
p12 0.421 0.052 8.132 0.000 0.270 0.422
p16 0.389 0.057 6.819 0.000 0.249 0.251
p20 0.216 0.041 5.257 0.000 0.138 0.295
p28 0.473 0.068 6.991 0.000 0.303 0.445
p47 0.624 0.070 8.949 0.000 0.399 0.505
p50 0.438 0.062 7.082 0.000 0.280 0.485
p02 0.721 0.069 10.525 0.000 0.462 0.682
p03 0.899 0.069 13.001 0.000 0.576 0.619
p04 0.754 0.062 12.191 0.000 0.483 0.598
p21 0.874 0.074 11.734 0.000 0.559 0.661
p22 0.882 0.057 15.442 0.000 0.565 0.556
p30 1.017 0.070 14.548 0.000 0.651 0.698
p31 0.960 0.069 13.984 0.000 0.615 0.781
p37 0.775 0.066 11.728 0.000 0.496 0.657
p40 0.975 0.071 13.805 0.000 0.624 0.765
p44 0.916 0.064 14.284 0.000 0.587 0.746
p45 0.679 0.063 10.774 0.000 0.435 0.691
p46 0.953 0.072 13.227 0.000 0.610 0.757
p51 0.955 0.073 13.161 0.000 0.612 0.626
p52 1.212 0.068 17.943 0.000 0.776 0.758
p57 0.878 0.068 12.944 0.000 0.562 0.771
Intercepts:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
.p06 2.520 0.033 76.549 0.000 2.520 2.095
.p10 2.208 0.034 65.045 0.000 2.208 1.780
.p14 1.600 0.027 60.165 0.000 1.600 1.647
.p25 2.029 0.033 62.184 0.000 2.029 1.702
.p27 2.229 0.033 68.385 0.000 2.229 1.872
.p29 1.898 0.032 60.059 0.000 1.898 1.644
.p33 1.601 0.026 60.633 0.000 1.601 1.659
.p35 1.776 0.028 63.917 0.000 1.776 1.749
.p48 1.236 0.018 67.548 0.000 1.236 1.849
.p49 1.649 0.027 61.987 0.000 1.649 1.697
.p53 1.844 0.029 63.324 0.000 1.844 1.733
.p54 1.934 0.033 58.053 0.000 1.934 1.589
.p07 1.622 0.025 65.517 0.000 1.622 1.793
.p11 1.586 0.026 60.218 0.000 1.586 1.648
.p13 1.213 0.018 66.573 0.000 1.213 1.822
.p17 1.493 0.024 63.352 0.000 1.493 1.734
.p24 1.196 0.014 85.313 0.000 1.196 2.335
.p26 2.026 0.027 73.948 0.000 2.026 2.024
.p36 1.459 0.024 60.477 0.000 1.459 1.655
.p55 1.837 0.033 55.295 0.000 1.837 1.513
.p56 1.923 0.034 57.270 0.000 1.923 1.567
.p01 1.303 0.020 65.295 0.000 1.303 1.787
.p18 2.318 0.034 67.527 0.000 2.318 1.848
.p19 1.288 0.020 65.846 0.000 1.288 1.802
.p23 2.022 0.028 72.385 0.000 2.022 1.981
.p39 1.311 0.018 72.292 0.000 1.311 1.979
.p43 1.656 0.023 72.927 0.000 1.656 1.996
.p09 1.246 0.016 76.116 0.000 1.246 2.083
.p12 1.338 0.018 76.389 0.000 1.338 2.091
.p16 1.692 0.027 62.260 0.000 1.692 1.704
.p20 1.109 0.013 86.722 0.000 1.109 2.373
.p28 1.205 0.019 64.736 0.000 1.205 1.772
.p47 1.370 0.022 63.288 0.000 1.370 1.732
.p50 1.184 0.016 74.812 0.000 1.184 2.048
.p02 1.298 0.019 70.090 0.000 1.298 1.918
.p03 1.630 0.025 64.022 0.000 1.630 1.752
.p04 1.573 0.022 71.253 0.000 1.573 1.950
.p21 1.562 0.023 67.411 0.000 1.562 1.845
.p22 1.831 0.028 65.796 0.000 1.831 1.801
.p30 1.706 0.026 66.859 0.000 1.706 1.830
.p31 1.514 0.022 70.256 0.000 1.514 1.923
.p37 1.479 0.021 71.457 0.000 1.479 1.956
.p40 1.467 0.022 65.622 0.000 1.467 1.796
.p44 1.599 0.022 74.349 0.000 1.599 2.035
.p45 1.282 0.017 74.467 0.000 1.282 2.038
.p46 1.502 0.022 68.064 0.000 1.502 1.863
.p51 1.619 0.027 60.522 0.000 1.619 1.656
.p52 1.804 0.028 64.384 0.000 1.804 1.762
.p57 1.378 0.020 69.055 0.000 1.378 1.890
abuse 0.000 0.000 0.000
Variances:
Estimate Std.Err z-value P(>|z|) Std.lv Std.all
.p06 1.037 0.048 21.813 0.000 1.037 0.717
.p10 1.287 0.048 26.819 0.000 1.287 0.836
.p14 0.462 0.029 16.035 0.000 0.462 0.490
.p25 1.091 0.042 25.785 0.000 1.091 0.767
.p27 0.987 0.041 24.269 0.000 0.987 0.696
.p29 0.605 0.033 18.229 0.000 0.605 0.454
.p33 0.357 0.023 15.700 0.000 0.357 0.384
.p35 0.669 0.041 16.194 0.000 0.669 0.649
.p48 0.300 0.032 9.440 0.000 0.300 0.672
.p49 0.576 0.033 17.198 0.000 0.576 0.609
.p53 0.646 0.040 16.033 0.000 0.646 0.571
.p54 0.986 0.045 21.967 0.000 0.986 0.666
.p07 0.589 0.036 16.149 0.000 0.589 0.720
.p11 0.568 0.036 15.739 0.000 0.568 0.614
.p13 0.368 0.042 8.712 0.000 0.368 0.829
.p17 0.529 0.034 15.407 0.000 0.529 0.714
.p24 0.187 0.020 9.378 0.000 0.187 0.715
.p26 0.659 0.041 16.117 0.000 0.659 0.658
.p36 0.496 0.035 14.348 0.000 0.496 0.638
.p55 1.007 0.051 19.864 0.000 1.007 0.683
.p56 0.931 0.046 20.420 0.000 0.931 0.619
.p01 0.419 0.044 9.581 0.000 0.419 0.787
.p18 1.106 0.043 25.620 0.000 1.106 0.704
.p19 0.343 0.030 11.613 0.000 0.343 0.672
.p23 0.750 0.033 22.901 0.000 0.750 0.720
.p39 0.347 0.035 9.948 0.000 0.347 0.791
.p43 0.418 0.023 17.799 0.000 0.418 0.607
.p09 0.282 0.028 9.966 0.000 0.282 0.788
.p12 0.337 0.031 11.006 0.000 0.337 0.822
.p16 0.924 0.046 19.893 0.000 0.924 0.937
.p20 0.199 0.033 6.011 0.000 0.199 0.913
.p28 0.371 0.039 9.444 0.000 0.371 0.802
.p47 0.466 0.037 12.715 0.000 0.466 0.745
.p50 0.255 0.033 7.825 0.000 0.255 0.765
.p02 0.245 0.027 9.217 0.000 0.245 0.534
.p03 0.534 0.038 14.053 0.000 0.534 0.617
.p04 0.418 0.033 12.790 0.000 0.418 0.642
.p21 0.404 0.028 14.351 0.000 0.404 0.563
.p22 0.714 0.035 20.130 0.000 0.714 0.691
.p30 0.446 0.039 11.520 0.000 0.446 0.512
.p31 0.242 0.023 10.688 0.000 0.242 0.390
.p37 0.325 0.030 10.910 0.000 0.325 0.569
.p40 0.277 0.032 8.742 0.000 0.277 0.415
.p44 0.274 0.021 12.945 0.000 0.274 0.443
.p45 0.207 0.017 12.315 0.000 0.207 0.523
.p46 0.277 0.029 9.723 0.000 0.277 0.427
.p51 0.581 0.040 14.598 0.000 0.581 0.608
.p52 0.447 0.027 16.499 0.000 0.447 0.426
.p57 0.216 0.022 9.615 0.000 0.216 0.406
abuse 0.410 0.045 9.048 0.000 1.000 1.000
R-Square:
Estimate
p06 0.283
p10 0.164
p14 0.510
p25 0.233
p27 0.304
p29 0.546
p33 0.616
p35 0.351
p48 0.328
p49 0.391
p53 0.429
p54 0.334
p07 0.280
p11 0.386
p13 0.171
p17 0.286
p24 0.285
p26 0.342
p36 0.362
p55 0.317
p56 0.381
p01 0.213
p18 0.296
p19 0.328
p23 0.280
p39 0.209
p43 0.393
p09 0.212
p12 0.178
p16 0.063
p20 0.087
p28 0.198
p47 0.255
p50 0.235
p02 0.466
p03 0.383
p04 0.358
p21 0.437
p22 0.309
p30 0.488
p31 0.610
p37 0.431
p40 0.585
p44 0.557
p45 0.477
p46 0.573
p51 0.392
p52 0.574
p57 0.594NOTE: With respect to fit of the structural model, we are now fitting a single factor INSTEAD OF 5 factors and a higher-order factor. This will tell us the extent to which a “total score” is appropriate.
anova(cfaSingleEstimates, cfaNoHighEstimates, cfaHigherEstimates)Scaled Chi Square Difference Test (method = "satorra.bentler.2001")
Df AIC BIC Chisq Chisq diff Df diff Pr(>Chisq)
cfaNoHighEstimates 1117 138229 139045 6490.3
cfaHigherEstimates 1122 138326 139115 6597.1 47.08 5 5.465e-09 ***
cfaSingleEstimates 1127 140928 141692 9210.0 448.91 5 < 2.2e-16 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1According to the −2ΔLL scaled difference relative to the previous model, −2ΔLL (5) = 448.91, p < .0001
Therefore, a single factor fits significantly worse than 5 factors + a higher-order factor, and so one factor does not capture the covariances for these 49 items.
We can try one more alternative – what if the items were measuring a single factor (i.e., a single score)?
NOTE: With respect to fit of the structural model, we are now fitting a single factor INSTEAD OF 5 factors and a higher-order factor. This will tell us the extent to which a single score is appropriate.
To test the fit against the higher-order factor model, we direct DIFFTEST on the ANALYSIS command to use the results from the previous model.
anova(ifaSingleEstimates, ifaNoHighEstimates)Error in lav_test_diff_af_h1(m1 = m1, m0 = m0) :
lavaan ERROR: unconstrained parameter set is not the same in m0 and m1Again, lavaan throws an error. We’ll use the Mplus result in our write up below.
After examining the fit of each of the five factors individually, as described previously, a combined model was estimated in which all five factors were fit simultaneously with covariances estimated freely among them. A total of 49 items were thus included. Each factor was identified by fixing the first item loading on each factor to 1, estimating the factor variance, and then fixing the factor mean to 0, while estimating all possible item intercepts, item residual variances, and remaining item loadings. Robust maximum likelihood (MLR) estimation was used to estimate all higher-order models using the lavaan package (Rosseel, 2012) in R (R Core Team, 2017), and differences in fit between nested models were evaluated using −2* rescaled difference in the model log-likelihood values.
As shown in Table 1, the fit of the model with five correlated factors was acceptable by the RMSEA (.047), but not by the CFI (.844). Standardized model parameters (loadings, intercepts, and residual variances) are shown in Table 2. Correlations of .6 or higher were found among the five factors, suggesting evidence that the five factors may indicate a single higher-order factor. This idea was testing by eliminating the covariances among the factors and instead estimating loadings for the five factors from a single higher-order factor (whose variance was fixed to 1). Although the fit of the higher-order factor model remained marginal (see Table 1), a nested model comparison revealed a significant decrease in fit, −2ΔLL(5) = 47.083, p < .0001, indicating that a single factor did not appear adequate to describe the pattern of correlation amongst the five factors. A further nested model comparison was conducted to examine the extent to which a single factor could describe the covariances among the items rather than five lower-order factors and a single higher-order factor. Fit of the single factor only model was poor, as shown in Table 1, and was significantly worse than the higher-order factor model, −2ΔLL(5) = 448.91, p < .0001, indicating that a single “total score” would not be recommended.
After examining the fit of each of the five factors individually, as described previously, a combined model was estimated in which all five factors were fit simultaneously with covariances estimated freely among them. A total of 49 items were thus included. Each factor was identified by fixing the first item loading on each factor to 1, estimating the factor variance, and then fixing the factor mean to 0, while estimating all possible item thresholds (four for each item given five response options) and remaining item loadings. WLSMV estimation in the lavaan package (Rosseel, 2012) in R (R Core Team, 2017) including a probit link and the THETA parameterization (such that all item residual variances were constrained to 1) was used to estimate all higher-order models. Thus, model fit statistics describe the fit of the item factor model to the polychoric correlation matrix among the items. Nested model comparisons were conducted using the Mplus DIFFTEST procedure.
As shown in Table 1, the fit of the model with five correlated factors was acceptable. Item factor analysis parameters (loadings and thresholds) and their corresponding item response model parameters (discriminations and difficulties) are shown in Table 2. Correlations of .7 or higher were found amongst the five factors, suggesting evidence that the five factors may indicate a single higher-order factor. This idea was testing by eliminating the covariances among the factors and instead estimating loadings for the five factors from a single higher-order factor (whose variance was fixed to 1). Although the fit of the higher-order factor model remained acceptable (see Table 1), a nested model comparison via the DIFFTEST procedure revealed a significant decrease in fit, DIFFTEST(5) = 92.05, p < .0001, indicating that a single factor did not appear adequate to describe the pattern of correlation amongst the five factors. A further nested model comparison was conducted to examine the extent to which a single factor could describe the polychoric correlations among the items rather than five lower-order factors and a single higher-order factor. Fit of the single factor only model was poor, as shown in Table 1, and was significantly worse than the higher-order factor model, DIFFTEST(5) = 611.95, p < .0001, indicating that a single score would not be recommended.
Table 1 = table with fit info per model Table 2 would have actual model parameters…. (unstandardized and standardized estimates and their SEs, so 4 columns)