load1 <- (5:8)/10
load2 <- load1# equal loadings
set.seed(1)
n <- 3000 # xi2 = .5*xi +1
xi <- MASS::mvrnorm(n, c(0,1), Sigma=matrix(c(1,.5, .5, 1),2))
ind1 <- xi[,1 ]%*% matrix(load1, nrow=1)+rnorm(4*n)# unique variance 1
ind2 <- xi[,2]%*% matrix(load2, nrow=1)+rnorm(4*n)# unique variance 1
dat <- cbind(ind1, ind2)
dat <- apply(dat, 2, cut, breaks=c(-Inf, -1,0, 1, Inf), labels=F)
colnames(dat) <- c(paste0("pc", 1:4,"_1"),paste0("pc", 1:4,"_2"))m0 <- "ly1 =~ pc1_1+pc2_1+pc3_1+pc4_1; ly2 =~ pc1_2+pc2_2+pc3_2+pc4_2";
f0 <- sem(m0, dat, ordered=T, std.lv=T, se="robust.sem", test="satorra.bentler", parameterization="theta")
#using semtools
syntax.config <- measEq.syntax(configural.model = f0, data=dat, parameterization="theta",
ID.fac = "std.lv", ID.cat = "Wu.Estabrook.2016",
ordered=T,
longFacNames = list(LY=c("ly1", "ly2")))
cat(as.character(syntax.config))## ## LOADINGS:
##
## ly1 =~ NA*pc1_1 + lambda.1_1*pc1_1
## ly1 =~ NA*pc2_1 + lambda.2_1*pc2_1
## ly1 =~ NA*pc3_1 + lambda.3_1*pc3_1
## ly1 =~ NA*pc4_1 + lambda.4_1*pc4_1
## ly2 =~ NA*pc1_2 + lambda.5_2*pc1_2
## ly2 =~ NA*pc2_2 + lambda.6_2*pc2_2
## ly2 =~ NA*pc3_2 + lambda.7_2*pc3_2
## ly2 =~ NA*pc4_2 + lambda.8_2*pc4_2
##
## ## THRESHOLDS:
##
## pc1_1 | NA*t1 + pc1_1.thr1*t1
## pc1_1 | NA*t2 + pc1_1.thr2*t2
## pc1_1 | NA*t3 + pc1_1.thr3*t3
## pc2_1 | NA*t1 + pc2_1.thr1*t1
## pc2_1 | NA*t2 + pc2_1.thr2*t2
## pc2_1 | NA*t3 + pc2_1.thr3*t3
## pc3_1 | NA*t1 + pc3_1.thr1*t1
## pc3_1 | NA*t2 + pc3_1.thr2*t2
## pc3_1 | NA*t3 + pc3_1.thr3*t3
## pc4_1 | NA*t1 + pc4_1.thr1*t1
## pc4_1 | NA*t2 + pc4_1.thr2*t2
## pc4_1 | NA*t3 + pc4_1.thr3*t3
## pc1_2 | NA*t1 + pc1_2.thr1*t1
## pc1_2 | NA*t2 + pc1_2.thr2*t2
## pc1_2 | NA*t3 + pc1_2.thr3*t3
## pc2_2 | NA*t1 + pc2_2.thr1*t1
## pc2_2 | NA*t2 + pc2_2.thr2*t2
## pc2_2 | NA*t3 + pc2_2.thr3*t3
## pc3_2 | NA*t1 + pc3_2.thr1*t1
## pc3_2 | NA*t2 + pc3_2.thr2*t2
## pc3_2 | NA*t3 + pc3_2.thr3*t3
## pc4_2 | NA*t1 + pc4_2.thr1*t1
## pc4_2 | NA*t2 + pc4_2.thr2*t2
## pc4_2 | NA*t3 + pc4_2.thr3*t3
##
## ## INTERCEPTS:
##
## pc1_1 ~ 0*1 + nu.1*1
## pc2_1 ~ 0*1 + nu.2*1
## pc3_1 ~ 0*1 + nu.3*1
## pc4_1 ~ 0*1 + nu.4*1
## pc1_2 ~ 0*1 + nu.5*1
## pc2_2 ~ 0*1 + nu.6*1
## pc3_2 ~ 0*1 + nu.7*1
## pc4_2 ~ 0*1 + nu.8*1
##
## ## UNIQUE-FACTOR VARIANCES:
##
## pc1_1 ~~ 1*pc1_1 + theta.1_1*pc1_1
## pc2_1 ~~ 1*pc2_1 + theta.2_2*pc2_1
## pc3_1 ~~ 1*pc3_1 + theta.3_3*pc3_1
## pc4_1 ~~ 1*pc4_1 + theta.4_4*pc4_1
## pc1_2 ~~ 1*pc1_2 + theta.5_5*pc1_2
## pc2_2 ~~ 1*pc2_2 + theta.6_6*pc2_2
## pc3_2 ~~ 1*pc3_2 + theta.7_7*pc3_2
## pc4_2 ~~ 1*pc4_2 + theta.8_8*pc4_2
##
## ## UNIQUE-FACTOR COVARIANCES:
##
## pc1_1 ~~ NA*pc1_2 + theta.5_1*pc1_2
## pc2_1 ~~ NA*pc2_2 + theta.6_2*pc2_2
## pc3_1 ~~ NA*pc3_2 + theta.7_3*pc3_2
## pc4_1 ~~ NA*pc4_2 + theta.8_4*pc4_2
##
## ## LATENT MEANS/INTERCEPTS:
##
## ly1 ~ 0*1 + alpha.1*1
## ly2 ~ 0*1 + alpha.2*1
##
## ## COMMON-FACTOR VARIANCES:
##
## ly1 ~~ 1*ly1 + psi.1_1*ly1
## ly2 ~~ 1*ly2 + psi.2_2*ly2
##
## ## COMMON-FACTOR COVARIANCES:
##
## ly1 ~~ NA*ly2 + psi.2_1*ly2f0new <- sem(as.character(syntax.config), dat, ordered=T, std.lv=T, se="robust.sem", test="satorra.bentler", parameterization="theta")
summary(f0new)## lavaan 0.6-7 ended normally after 33 iterations
##
## Estimator DWLS
## Optimization method NLMINB
## Number of free parameters 37
##
## Number of observations 3000
##
## Model Test User Model:
## Standard Robust
## Test Statistic 7.877 10.898
## Degrees of freedom 15 15
## P-value (Chi-square) 0.929 0.760
## Scaling correction factor 0.723
## Satorra-Bentler correction
##
## Parameter Estimates:
##
## Standard errors Robust.sem
## Information Expected
## Information saturated (h1) model Unstructured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## ly1 =~
## pc1_1 (l.1_) 0.487 0.031 15.633 0.000
## pc2_1 (l.2_) 0.608 0.035 17.568 0.000
## pc3_1 (l.3_) 0.766 0.043 17.615 0.000
## pc4_1 (l.4_) 0.826 0.047 17.653 0.000
## ly2 =~
## pc1_2 (l.5_) 0.574 0.035 16.177 0.000
## pc2_2 (l.6_) 0.622 0.038 16.481 0.000
## pc3_2 (l.7_) 0.704 0.041 17.227 0.000
## pc4_2 (l.8_) 0.774 0.046 16.927 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## .pc1_1 ~~
## .pc1_2 (t.5_) -0.018 0.025 -0.734 0.463
## .pc2_1 ~~
## .pc2_2 (t.6_) 0.047 0.026 1.840 0.066
## .pc3_1 ~~
## .pc3_2 (t.7_) -0.005 0.030 -0.152 0.879
## .pc4_1 ~~
## .pc4_2 (t.8_) 0.019 0.032 0.601 0.548
## ly1 ~~
## ly2 (p.2_) 0.479 0.028 16.908 0.000
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .pc1_1 (nu.1) 0.000
## .pc2_1 (nu.2) 0.000
## .pc3_1 (nu.3) 0.000
## .pc4_1 (nu.4) 0.000
## .pc1_2 (nu.5) 0.000
## .pc2_2 (nu.6) 0.000
## .pc3_2 (nu.7) 0.000
## .pc4_2 (nu.8) 0.000
## ly1 (al.1) 0.000
## ly2 (al.2) 0.000
##
## Thresholds:
## Estimate Std.Err z-value P(>|z|)
## p1_1| (p1_1.1) -0.991 0.031 -32.213 0.000
## p1_1| (p1_1.2) 0.011 0.025 0.438 0.661
## p1_1| (p1_1.3) 1.015 0.031 32.710 0.000
## p2_1| (p2_1.1) -0.968 0.032 -29.834 0.000
## p2_1| (p2_1.2) 0.018 0.027 0.657 0.511
## p2_1| (p2_1.3) 1.030 0.033 31.117 0.000
## p3_1| (p3_1.1) -1.007 0.037 -27.261 0.000
## p3_1| (p3_1.2) 0.032 0.029 1.095 0.274
## p3_1| (p3_1.3) 1.048 0.037 28.114 0.000
## p4_1| (p4_1.1) -1.024 0.039 -26.536 0.000
## p4_1| (p4_1.2) -0.036 0.030 -1.205 0.228
## p4_1| (p4_1.3) 1.003 0.038 26.263 0.000
## p1_2| (p1_2.1) -1.527 0.041 -37.064 0.000
## p1_2| (p1_2.2) -0.522 0.028 -18.548 0.000
## p1_2| (p1_2.3) 0.496 0.028 17.812 0.000
## p2_2| (p2_2.1) -1.606 0.044 -36.435 0.000
## p2_2| (p2_2.2) -0.636 0.030 -21.389 0.000
## p2_2| (p2_2.3) 0.392 0.028 14.039 0.000
## p3_2| (p3_2.1) -1.696 0.047 -35.875 0.000
## p3_2| (p3_2.2) -0.679 0.032 -21.537 0.000
## p3_2| (p3_2.3) 0.315 0.029 11.017 0.000
## p4_2| (p4_2.1) -1.703 0.051 -33.229 0.000
## p4_2| (p4_2.2) -0.726 0.034 -21.508 0.000
## p4_2| (p4_2.3) 0.209 0.029 7.162 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .pc1_1 (t.1_) 1.000
## .pc2_1 (t.2_) 1.000
## .pc3_1 (t.3_) 1.000
## .pc4_1 (t.4_) 1.000
## .pc1_2 (t.5_) 1.000
## .pc2_2 (t.6_) 1.000
## .pc3_2 (t.7_) 1.000
## .pc4_2 (t.8_) 1.000
## ly1 (p.1_) 1.000
## ly2 (p.2_) 1.000
##
## Scales y*:
## Estimate Std.Err z-value P(>|z|)
## pc1_1 0.899
## pc2_1 0.855
## pc3_1 0.794
## pc4_1 0.771
## pc1_2 0.867
## pc2_2 0.849
## pc3_2 0.818
## pc4_2 0.791We get estimates close to the population parameter values.
syntax.thr <- measEq.syntax(configural.model = f0, data=dat, parameterization="theta",
ID.fac = "std.lv", ID.cat = "Wu.Estabrook.2016",
longFacNames = list(LY=c("ly1", "ly2")),
ordered=T,
long.equal="thresholds")
cat(as.character(syntax.thr))## ## LOADINGS:
##
## ly1 =~ NA*pc1_1 + lambda.1_1*pc1_1
## ly1 =~ NA*pc2_1 + lambda.2_1*pc2_1
## ly1 =~ NA*pc3_1 + lambda.3_1*pc3_1
## ly1 =~ NA*pc4_1 + lambda.4_1*pc4_1
## ly2 =~ NA*pc1_2 + lambda.5_2*pc1_2
## ly2 =~ NA*pc2_2 + lambda.6_2*pc2_2
## ly2 =~ NA*pc3_2 + lambda.7_2*pc3_2
## ly2 =~ NA*pc4_2 + lambda.8_2*pc4_2
##
## ## THRESHOLDS:
##
## pc1_1 | NA*t1 + pc1_1.thr1*t1
## pc1_1 | NA*t2 + pc1_1.thr2*t2
## pc1_1 | NA*t3 + pc1_1.thr3*t3
## pc2_1 | NA*t1 + pc2_1.thr1*t1
## pc2_1 | NA*t2 + pc2_1.thr2*t2
## pc2_1 | NA*t3 + pc2_1.thr3*t3
## pc3_1 | NA*t1 + pc3_1.thr1*t1
## pc3_1 | NA*t2 + pc3_1.thr2*t2
## pc3_1 | NA*t3 + pc3_1.thr3*t3
## pc4_1 | NA*t1 + pc4_1.thr1*t1
## pc4_1 | NA*t2 + pc4_1.thr2*t2
## pc4_1 | NA*t3 + pc4_1.thr3*t3
## pc1_2 | NA*t1 + pc1_1.thr1*t1
## pc1_2 | NA*t2 + pc1_1.thr2*t2
## pc1_2 | NA*t3 + pc1_1.thr3*t3
## pc2_2 | NA*t1 + pc2_1.thr1*t1
## pc2_2 | NA*t2 + pc2_1.thr2*t2
## pc2_2 | NA*t3 + pc2_1.thr3*t3
## pc3_2 | NA*t1 + pc3_1.thr1*t1
## pc3_2 | NA*t2 + pc3_1.thr2*t2
## pc3_2 | NA*t3 + pc3_1.thr3*t3
## pc4_2 | NA*t1 + pc4_1.thr1*t1
## pc4_2 | NA*t2 + pc4_1.thr2*t2
## pc4_2 | NA*t3 + pc4_1.thr3*t3
##
## ## INTERCEPTS:
##
## pc1_1 ~ 0*1 + nu.1*1
## pc2_1 ~ 0*1 + nu.2*1
## pc3_1 ~ 0*1 + nu.3*1
## pc4_1 ~ 0*1 + nu.4*1
## pc1_2 ~ NA*1 + nu.5*1
## pc2_2 ~ NA*1 + nu.6*1
## pc3_2 ~ NA*1 + nu.7*1
## pc4_2 ~ NA*1 + nu.8*1
##
## ## UNIQUE-FACTOR VARIANCES:
##
## pc1_1 ~~ 1*pc1_1 + theta.1_1*pc1_1
## pc2_1 ~~ 1*pc2_1 + theta.2_2*pc2_1
## pc3_1 ~~ 1*pc3_1 + theta.3_3*pc3_1
## pc4_1 ~~ 1*pc4_1 + theta.4_4*pc4_1
## pc1_2 ~~ NA*pc1_2 + theta.5_5*pc1_2
## pc2_2 ~~ NA*pc2_2 + theta.6_6*pc2_2
## pc3_2 ~~ NA*pc3_2 + theta.7_7*pc3_2
## pc4_2 ~~ NA*pc4_2 + theta.8_8*pc4_2
##
## ## UNIQUE-FACTOR COVARIANCES:
##
## pc1_1 ~~ NA*pc1_2 + theta.5_1*pc1_2
## pc2_1 ~~ NA*pc2_2 + theta.6_2*pc2_2
## pc3_1 ~~ NA*pc3_2 + theta.7_3*pc3_2
## pc4_1 ~~ NA*pc4_2 + theta.8_4*pc4_2
##
## ## LATENT MEANS/INTERCEPTS:
##
## ly1 ~ 0*1 + alpha.1*1
## ly2 ~ 0*1 + alpha.2*1
##
## ## COMMON-FACTOR VARIANCES:
##
## ly1 ~~ 1*ly1 + psi.1_1*ly1
## ly2 ~~ 1*ly2 + psi.2_2*ly2
##
## ## COMMON-FACTOR COVARIANCES:
##
## ly1 ~~ NA*ly2 + psi.2_1*ly2f1 <- sem(as.character(syntax.thr), dat, ordered=T, std.lv=T, se="robust.sem", test="satorra.bentler", parameterization="theta")
summary(f1)## lavaan 0.6-7 ended normally after 51 iterations
##
## Estimator DWLS
## Optimization method NLMINB
## Number of free parameters 45
## Number of equality constraints 12
##
## Number of observations 3000
##
## Model Test User Model:
## Standard Robust
## Test Statistic 8.702 13.236
## Degrees of freedom 19 19
## P-value (Chi-square) 0.978 0.826
## Scaling correction factor 0.657
## Satorra-Bentler correction
##
## Parameter Estimates:
##
## Standard errors Robust.sem
## Information Expected
## Information saturated (h1) model Unstructured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## ly1 =~
## pc1_1 (l.1_) 0.487 0.031 15.633 0.000
## pc2_1 (l.2_) 0.608 0.035 17.568 0.000
## pc3_1 (l.3_) 0.766 0.043 17.615 0.000
## pc4_1 (l.4_) 0.826 0.047 17.653 0.000
## ly2 =~
## pc1_2 (l.5_) 0.569 0.032 17.662 0.000
## pc2_2 (l.6_) 0.621 0.035 17.841 0.000
## pc3_2 (l.7_) 0.719 0.039 18.245 0.000
## pc4_2 (l.8_) 0.823 0.045 18.453 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## .pc1_1 ~~
## .pc1_2 (t.5_) -0.018 0.024 -0.734 0.463
## .pc2_1 ~~
## .pc2_2 (t.6_) 0.047 0.026 1.838 0.066
## .pc3_1 ~~
## .pc3_2 (t.7_) -0.005 0.030 -0.152 0.879
## .pc4_1 ~~
## .pc4_2 (t.8_) 0.020 0.034 0.601 0.548
## ly1 ~~
## ly2 (p.2_) 0.479 0.028 16.908 0.000
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .pc1_1 (nu.1) 0.000
## .pc2_1 (nu.2) 0.000
## .pc3_1 (nu.3) 0.000
## .pc4_1 (nu.4) 0.000
## .pc1_2 (nu.5) 0.526 0.032 16.300 0.000
## .pc2_2 (nu.6) 0.644 0.034 19.118 0.000
## .pc3_2 (nu.7) 0.726 0.038 18.983 0.000
## .pc4_2 (nu.8) 0.765 0.041 18.455 0.000
## ly1 (al.1) 0.000
## ly2 (al.2) 0.000
##
## Thresholds:
## Estimate Std.Err z-value P(>|z|)
## p1_1| (p1_1.1) -0.990 0.030 -32.671 0.000
## p1_1| (p1_1.2) 0.010 0.024 0.404 0.686
## p1_1| (p1_1.3) 1.016 0.030 33.368 0.000
## p2_1| (p2_1.1) -0.965 0.032 -30.170 0.000
## p2_1| (p2_1.2) 0.013 0.025 0.529 0.597
## p2_1| (p2_1.3) 1.033 0.033 31.649 0.000
## p3_1| (p3_1.1) -1.007 0.036 -27.618 0.000
## p3_1| (p3_1.2) 0.032 0.027 1.177 0.239
## p3_1| (p3_1.3) 1.048 0.037 28.429 0.000
## p4_1| (p4_1.1) -1.032 0.038 -26.948 0.000
## p4_1| (p4_1.2) -0.023 0.028 -0.806 0.420
## p4_1| (p4_1.3) 0.995 0.038 26.445 0.000
## p1_2| (p1_1.1) -0.990 0.030 -32.671 0.000
## p1_2| (p1_1.2) 0.010 0.024 0.404 0.686
## p1_2| (p1_1.3) 1.016 0.030 33.368 0.000
## p2_2| (p2_1.1) -0.965 0.032 -30.170 0.000
## p2_2| (p2_1.2) 0.013 0.025 0.529 0.597
## p2_2| (p2_1.3) 1.033 0.033 31.649 0.000
## p3_2| (p3_1.1) -1.007 0.036 -27.618 0.000
## p3_2| (p3_1.2) 0.032 0.027 1.177 0.239
## p3_2| (p3_1.3) 1.048 0.037 28.429 0.000
## p4_2| (p4_1.1) -1.032 0.038 -26.948 0.000
## p4_2| (p4_1.2) -0.023 0.028 -0.806 0.420
## p4_2| (p4_1.3) 0.995 0.038 26.445 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .pc1_1 (t.1_) 1.000
## .pc2_1 (t.2_) 1.000
## .pc3_1 (t.3_) 1.000
## .pc4_1 (t.4_) 1.000
## .pc1_2 (t.5_) 0.983 0.060 16.284 0.000
## .pc2_2 (t.6_) 0.998 0.066 15.219 0.000
## .pc3_2 (t.7_) 1.044 0.075 13.855 0.000
## .pc4_2 (t.8_) 1.131 0.089 12.709 0.000
## ly1 (p.1_) 1.000
## ly2 (p.2_) 1.000
##
## Scales y*:
## Estimate Std.Err z-value P(>|z|)
## pc1_1 0.899
## pc2_1 0.855
## pc3_1 0.794
## pc4_1 0.771
## pc1_2 0.875
## pc2_2 0.850
## pc3_2 0.801
## pc4_2 0.743anova(f0new, f1)## Scaled Chi-Squared Difference Test (method = "satorra.bentler.2001")
##
## lavaan NOTE:
## The "Chisq" column contains standard test statistics, not the
## robust test that should be reported per model. A robust difference
## test is a function of two standard (not robust) statistics.
##
## Df AIC BIC Chisq Chisq diff Df diff Pr(>Chisq)
## f0new 15 7.8768
## f1 19 8.7024 2.0004 4 0.7357compareFit(f0new,f1)## ################### Nested Model Comparison #########################
## Scaled Chi-Squared Difference Test (method = "satorra.bentler.2001")
##
## lavaan NOTE:
## The "Chisq" column contains standard test statistics, not the
## robust test that should be reported per model. A robust difference
## test is a function of two standard (not robust) statistics.
##
## Df AIC BIC Chisq Chisq diff Df diff Pr(>Chisq)
## f0new 15 7.8768
## f1 19 8.7024 2.0004 4 0.7357
##
## ####################### Model Fit Indices ###########################
## chisq.scaled df.scaled pvalue.scaled cfi.robust tli.robust rmsea.robust
## f0new 10.898† 15 .760 1.000† 1.001 .000†
## f1 13.236 19 .826 1.000† 1.001† .000†
## srmr
## f0new .010
## f1 .010†
##
## ################## Differences in Fit Indices #######################
## df.scaled cfi.robust tli.robust rmsea.robust srmr
## f1 - f0new 4 0 0 0 0syntax.load <- measEq.syntax(configural.model = f0, data=dat, parameterization="theta",
ID.fac = "std.lv", ID.cat = "Wu.Estabrook.2016",
longFacNames = list(LY=c("ly1", "ly2")),
ordered=T,
long.equal=c("thresholds", "loadings"))
cat(as.character(syntax.load))## ## LOADINGS:
##
## ly1 =~ NA*pc1_1 + lambda.1_1*pc1_1
## ly1 =~ NA*pc2_1 + lambda.2_1*pc2_1
## ly1 =~ NA*pc3_1 + lambda.3_1*pc3_1
## ly1 =~ NA*pc4_1 + lambda.4_1*pc4_1
## ly2 =~ NA*pc1_2 + lambda.1_1*pc1_2
## ly2 =~ NA*pc2_2 + lambda.2_1*pc2_2
## ly2 =~ NA*pc3_2 + lambda.3_1*pc3_2
## ly2 =~ NA*pc4_2 + lambda.4_1*pc4_2
##
## ## THRESHOLDS:
##
## pc1_1 | NA*t1 + pc1_1.thr1*t1
## pc1_1 | NA*t2 + pc1_1.thr2*t2
## pc1_1 | NA*t3 + pc1_1.thr3*t3
## pc2_1 | NA*t1 + pc2_1.thr1*t1
## pc2_1 | NA*t2 + pc2_1.thr2*t2
## pc2_1 | NA*t3 + pc2_1.thr3*t3
## pc3_1 | NA*t1 + pc3_1.thr1*t1
## pc3_1 | NA*t2 + pc3_1.thr2*t2
## pc3_1 | NA*t3 + pc3_1.thr3*t3
## pc4_1 | NA*t1 + pc4_1.thr1*t1
## pc4_1 | NA*t2 + pc4_1.thr2*t2
## pc4_1 | NA*t3 + pc4_1.thr3*t3
## pc1_2 | NA*t1 + pc1_1.thr1*t1
## pc1_2 | NA*t2 + pc1_1.thr2*t2
## pc1_2 | NA*t3 + pc1_1.thr3*t3
## pc2_2 | NA*t1 + pc2_1.thr1*t1
## pc2_2 | NA*t2 + pc2_1.thr2*t2
## pc2_2 | NA*t3 + pc2_1.thr3*t3
## pc3_2 | NA*t1 + pc3_1.thr1*t1
## pc3_2 | NA*t2 + pc3_1.thr2*t2
## pc3_2 | NA*t3 + pc3_1.thr3*t3
## pc4_2 | NA*t1 + pc4_1.thr1*t1
## pc4_2 | NA*t2 + pc4_1.thr2*t2
## pc4_2 | NA*t3 + pc4_1.thr3*t3
##
## ## INTERCEPTS:
##
## pc1_1 ~ 0*1 + nu.1*1
## pc2_1 ~ 0*1 + nu.2*1
## pc3_1 ~ 0*1 + nu.3*1
## pc4_1 ~ 0*1 + nu.4*1
## pc1_2 ~ NA*1 + nu.5*1
## pc2_2 ~ NA*1 + nu.6*1
## pc3_2 ~ NA*1 + nu.7*1
## pc4_2 ~ NA*1 + nu.8*1
##
## ## UNIQUE-FACTOR VARIANCES:
##
## pc1_1 ~~ 1*pc1_1 + theta.1_1*pc1_1
## pc2_1 ~~ 1*pc2_1 + theta.2_2*pc2_1
## pc3_1 ~~ 1*pc3_1 + theta.3_3*pc3_1
## pc4_1 ~~ 1*pc4_1 + theta.4_4*pc4_1
## pc1_2 ~~ NA*pc1_2 + theta.5_5*pc1_2
## pc2_2 ~~ NA*pc2_2 + theta.6_6*pc2_2
## pc3_2 ~~ NA*pc3_2 + theta.7_7*pc3_2
## pc4_2 ~~ NA*pc4_2 + theta.8_8*pc4_2
##
## ## UNIQUE-FACTOR COVARIANCES:
##
## pc1_1 ~~ NA*pc1_2 + theta.5_1*pc1_2
## pc2_1 ~~ NA*pc2_2 + theta.6_2*pc2_2
## pc3_1 ~~ NA*pc3_2 + theta.7_3*pc3_2
## pc4_1 ~~ NA*pc4_2 + theta.8_4*pc4_2
##
## ## LATENT MEANS/INTERCEPTS:
##
## ly1 ~ 0*1 + alpha.1*1
## ly2 ~ 0*1 + alpha.2*1
##
## ## COMMON-FACTOR VARIANCES:
##
## ly1 ~~ 1*ly1 + psi.1_1*ly1
## ly2 ~~ NA*ly2 + psi.2_2*ly2
##
## ## COMMON-FACTOR COVARIANCES:
##
## ly1 ~~ NA*ly2 + psi.2_1*ly2f2 <- sem(as.character(syntax.load), dat, ordered=T, std.lv=T, se="robust.sem", test="satorra.bentler", parameterization="theta")
summary(f2)## lavaan 0.6-7 ended normally after 52 iterations
##
## Estimator DWLS
## Optimization method NLMINB
## Number of free parameters 46
## Number of equality constraints 16
##
## Number of observations 3000
##
## Model Test User Model:
## Standard Robust
## Test Statistic 15.043 20.668
## Degrees of freedom 22 22
## P-value (Chi-square) 0.860 0.541
## Scaling correction factor 0.728
## Satorra-Bentler correction
##
## Parameter Estimates:
##
## Standard errors Robust.sem
## Information Expected
## Information saturated (h1) model Unstructured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## ly1 =~
## pc1_1 (l.1_) 0.521 0.026 20.110 0.000
## pc2_1 (l.2_) 0.609 0.029 20.953 0.000
## pc3_1 (l.3_) 0.733 0.035 21.179 0.000
## pc4_1 (l.4_) 0.816 0.039 20.864 0.000
## ly2 =~
## pc1_2 (l.1_) 0.521 0.026 20.110 0.000
## pc2_2 (l.2_) 0.609 0.029 20.953 0.000
## pc3_2 (l.3_) 0.733 0.035 21.179 0.000
## pc4_2 (l.4_) 0.816 0.039 20.864 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## .pc1_1 ~~
## .pc1_2 (t.5_) -0.018 0.024 -0.751 0.453
## .pc2_1 ~~
## .pc2_2 (t.6_) 0.047 0.026 1.818 0.069
## .pc3_1 ~~
## .pc3_2 (t.7_) -0.006 0.030 -0.191 0.848
## .pc4_1 ~~
## .pc4_2 (t.8_) 0.020 0.034 0.580 0.562
## ly1 ~~
## ly2 (p.2_) 0.490 0.033 14.979 0.000
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .pc1_1 (nu.1) 0.000
## .pc2_1 (nu.2) 0.000
## .pc3_1 (nu.3) 0.000
## .pc4_1 (nu.4) 0.000
## .pc1_2 (nu.5) 0.521 0.032 16.127 0.000
## .pc2_2 (nu.6) 0.644 0.034 19.064 0.000
## .pc3_2 (nu.7) 0.729 0.038 19.204 0.000
## .pc4_2 (nu.8) 0.767 0.042 18.470 0.000
## ly1 (al.1) 0.000
## ly2 (al.2) 0.000
##
## Thresholds:
## Estimate Std.Err z-value P(>|z|)
## p1_1| (p1_1.1) -0.991 0.031 -31.991 0.000
## p1_1| (p1_1.2) 0.011 0.024 0.456 0.648
## p1_1| (p1_1.3) 1.015 0.031 32.536 0.000
## p2_1| (p2_1.1) -0.966 0.032 -30.252 0.000
## p2_1| (p2_1.2) 0.013 0.025 0.528 0.598
## p2_1| (p2_1.3) 1.034 0.033 31.678 0.000
## p3_1| (p3_1.1) -1.002 0.035 -28.624 0.000
## p3_1| (p3_1.2) 0.030 0.027 1.128 0.259
## p3_1| (p3_1.3) 1.044 0.035 29.509 0.000
## p4_1| (p4_1.1) -1.030 0.037 -27.593 0.000
## p4_1| (p4_1.2) -0.023 0.028 -0.825 0.410
## p4_1| (p4_1.3) 0.994 0.037 26.979 0.000
## p1_2| (p1_1.1) -0.991 0.031 -31.991 0.000
## p1_2| (p1_1.2) 0.011 0.024 0.456 0.648
## p1_2| (p1_1.3) 1.015 0.031 32.536 0.000
## p2_2| (p2_1.1) -0.966 0.032 -30.252 0.000
## p2_2| (p2_1.2) 0.013 0.025 0.528 0.598
## p2_2| (p2_1.3) 1.034 0.033 31.678 0.000
## p3_2| (p3_1.1) -1.002 0.035 -28.624 0.000
## p3_2| (p3_1.2) 0.030 0.027 1.128 0.259
## p3_2| (p3_1.3) 1.044 0.035 29.509 0.000
## p4_2| (p4_1.1) -1.030 0.037 -27.593 0.000
## p4_2| (p4_1.2) -0.023 0.028 -0.825 0.410
## p4_2| (p4_1.3) 0.994 0.037 26.979 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .pc1_1 (t.1_) 1.000
## .pc2_1 (t.2_) 1.000
## .pc3_1 (t.3_) 1.000
## .pc4_1 (t.4_) 1.000
## .pc1_2 (t.5_) 0.985 0.060 16.448 0.000
## .pc2_2 (t.6_) 1.001 0.065 15.378 0.000
## .pc3_2 (t.7_) 1.032 0.072 14.268 0.000
## .pc4_2 (t.8_) 1.128 0.086 13.148 0.000
## ly1 (p.1_) 1.000
## ly2 (p.2_) 1.042 0.067 15.543 0.000
##
## Scales y*:
## Estimate Std.Err z-value P(>|z|)
## pc1_1 0.887
## pc2_1 0.854
## pc3_1 0.807
## pc4_1 0.775
## pc1_2 0.888
## pc2_2 0.849
## pc3_2 0.793
## pc4_2 0.741compareFit(f0new,f1, f2)## ################### Nested Model Comparison #########################
## Scaled Chi-Squared Difference Test (method = "satorra.bentler.2001")
##
## lavaan NOTE:
## The "Chisq" column contains standard test statistics, not the
## robust test that should be reported per model. A robust difference
## test is a function of two standard (not robust) statistics.
##
## Df AIC BIC Chisq Chisq diff Df diff Pr(>Chisq)
## f0new 15 7.8768
## f1 19 8.7024 2.0004 4 0.7357
## f2 22 15.0428 5.4042 3 0.1445
##
## ####################### Model Fit Indices ###########################
## chisq.scaled df.scaled pvalue.scaled cfi.robust tli.robust rmsea.robust
## f0new 10.898† 15 .760 1.000† 1.001 .000†
## f1 13.236 19 .826 1.000† 1.001† .000†
## f2 20.668 22 .541 1.000† 1.000 .000†
## srmr
## f0new .010
## f1 .010†
## f2 .012
##
## ################## Differences in Fit Indices #######################
## df.scaled cfi.robust tli.robust rmsea.robust srmr
## f1 - f0new 4 0 0.000 0 0.000
## f2 - f1 3 0 -0.001 0 0.002