library(sem)
## Warning: package 'sem' was built under R version 3.5.3
# data from Votaw (1948)
votaw <- readMoments(diag=TRUE,
names=c('orig1', 'hcpy1', 'ccpy1', 'orig2'), text="
25.0704
12.4363 28.2021
11.7257 9.2281 22.7390
20.7510 11.9732 12.0692 21.8707
")
congeneric model
votaw.mod1 <- specifyEquations(text="
orig1 = lam1 * Ability
hcpy1 = lam2 * Ability
ccpy1 = lam3 * Ability
orig2 = lam4 * Ability
V(Ability) = 1
")
## NOTE: adding 4 variances to the model
votaw.sem1 <- sem(votaw.mod1, votaw, N=126)
summary(votaw.sem1)
##
## Model Chisquare = 2.280074 Df = 2 Pr(>Chisq) = 0.3198071
## AIC = 18.28007
## BIC = -7.392489
##
## Normalized Residuals
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.165244 -0.016739 0.002136 0.096890 0.034388 0.906931
##
## R-square for Endogenous Variables
## orig1 hcpy1 ccpy1 orig2
## 0.8341 0.2540 0.3090 0.9405
##
## Parameter Estimates
## Estimate Std Error z value Pr(>|z|)
## lam1 4.572765 0.3622384 12.623637 1.564264e-36 orig1 <--- Ability
## lam2 2.676477 0.4540865 5.894200 3.765007e-09 hcpy1 <--- Ability
## lam3 2.650794 0.4013294 6.605034 3.974274e-11 ccpy1 <--- Ability
## lam4 4.535340 0.3266465 13.884552 7.858988e-44 orig2 <--- Ability
## V[orig1] 4.160215 1.2134631 3.428382 6.071898e-04 orig1 <--> orig1
## V[hcpy1] 21.038578 2.7104389 7.762056 8.356344e-15 hcpy1 <--> hcpy1
## V[ccpy1] 15.712295 2.0381554 7.709076 1.267322e-14 ccpy1 <--> ccpy1
## V[orig2] 1.301390 1.0861703 1.198146 2.308603e-01 orig2 <--> orig2
##
## Iterations = 29
# try with cfa
votaw.mod1a <- cfa(text="
Ability: orig1, hcpy1, ccpy1, orig2
")
## NOTE: adding 4 variances to the model
votaw.sem1a <- sem(votaw.mod1a, votaw, N=126)
summary(votaw.sem1a)
##
## Model Chisquare = 2.280074 Df = 2 Pr(>Chisq) = 0.3198071
## AIC = 18.28007
## BIC = -7.392489
##
## Normalized Residuals
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.165244 -0.016738 0.002136 0.096890 0.034388 0.906931
##
## R-square for Endogenous Variables
## orig1 hcpy1 ccpy1 orig2
## 0.8341 0.2540 0.3090 0.9405
##
## Parameter Estimates
## Estimate Std Error z value Pr(>|z|)
## lam[hcpy1:Ability] 0.5853082 0.09520200 6.148066 7.843352e-10
## lam[ccpy1:Ability] 0.5796916 0.08317646 6.969420 3.182511e-12
## lam[orig2:Ability] 0.9918155 0.06818800 14.545308 6.254826e-48
## V[Ability] 20.9101866 3.31286217 6.311819 2.757745e-10
## V[orig1] 4.1602136 1.21346295 3.428381 6.071921e-04
## V[hcpy1] 21.0385728 2.71043817 7.762056 8.356351e-15
## V[ccpy1] 15.7122919 2.03815509 7.709076 1.267322e-14
## V[orig2] 1.3013913 1.08617001 1.198147 2.308598e-01
##
## lam[hcpy1:Ability] hcpy1 <--- Ability
## lam[ccpy1:Ability] ccpy1 <--- Ability
## lam[orig2:Ability] orig2 <--- Ability
## V[Ability] Ability <--> Ability
## V[orig1] orig1 <--> orig1
## V[hcpy1] hcpy1 <--> hcpy1
## V[ccpy1] ccpy1 <--> ccpy1
## V[orig2] orig2 <--> orig2
##
## Iterations = 55
tau-equivalent model
votaw.mod2 <- specifyEquations(
text="
orig1 = lam * Ability
hcpy1 = lam * Ability
ccpy1 = lam * Ability
orig2 = lam * Ability
V(Ability) = 1
")
## NOTE: adding 4 variances to the model
votaw.sem2 <- sem(votaw.mod2, votaw, N=126)
summary(votaw.sem2)
##
## Model Chisquare = 40.42317 Df = 5 Pr(>Chisq) = 1.226875e-07
## AIC = 50.42317
## BIC = 16.24176
##
## Normalized Residuals
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -2.8941 -2.2406 -2.1109 -1.4566 -1.2594 0.9564
##
## R-square for Endogenous Variables
## orig1 hcpy1 ccpy1 orig2
## 0.8165 0.4486 0.5112 0.8970
##
## Parameter Estimates
## Estimate Std Error z value Pr(>|z|)
## lam 4.282332 0.2898063 14.776531 2.075851e-49 orig1 <--- Ability
## V[orig1] 4.120652 0.8122052 5.073412 3.907453e-07 orig1 <--> orig1
## V[hcpy1] 22.541754 3.0101503 7.488581 6.962217e-14 hcpy1 <--> hcpy1
## V[ccpy1] 17.531813 2.3809960 7.363227 1.795171e-13 ccpy1 <--> ccpy1
## V[orig2] 2.105483 0.6619380 3.180786 1.468759e-03 orig2 <--> orig2
##
## Iterations = 20
anova(votaw.sem1, votaw.sem2)
## LR Test for Difference Between Models
##
## Model Df Model Chisq Df LR Chisq Pr(>Chisq)
## votaw.sem1 2 2.280
## votaw.sem2 5 40.423 3 38.143 2.636e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
try with cfa function
votaw.mod2a <- cfa(reference.indicators=FALSE, text="
Ability: orig1 = hcpy1 = ccpy1 = orig2
")
## NOTE: adding 4 variances to the model
votaw.sem2a <- sem(votaw.mod2a, votaw, N=126)
summary(votaw.sem2a)
##
## Model Chisquare = 40.42317 Df = 5 Pr(>Chisq) = 1.226875e-07
## AIC = 50.42317
## BIC = 16.24176
##
## Normalized Residuals
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -2.8941 -2.2406 -2.1109 -1.4566 -1.2594 0.9564
##
## R-square for Endogenous Variables
## orig1 hcpy1 ccpy1 orig2
## 0.8165 0.4486 0.5112 0.8970
##
## Parameter Estimates
## Estimate Std Error z value
## lam[orig1.hcpy1.ccpy1.orig2:Ability] 4.282332 0.2898063 14.776531
## V[orig1] 4.120652 0.8122052 5.073412
## V[hcpy1] 22.541754 3.0101503 7.488581
## V[ccpy1] 17.531813 2.3809960 7.363227
## V[orig2] 2.105483 0.6619380 3.180786
## Pr(>|z|)
## lam[orig1.hcpy1.ccpy1.orig2:Ability] 2.075851e-49 orig1 <--- Ability
## V[orig1] 3.907453e-07 orig1 <--> orig1
## V[hcpy1] 6.962217e-14 hcpy1 <--> hcpy1
## V[ccpy1] 1.795171e-13 ccpy1 <--> ccpy1
## V[orig2] 1.468759e-03 orig2 <--> orig2
##
## Iterations = 20
votaw.mod2b <- cfa(reference.indicators=TRUE, text="
Ability: orig1 = hcpy1 = ccpy1 = orig2
")
## NOTE: adding 4 variances to the model
votaw.sem2b <- sem(votaw.mod2b, votaw, N=126)
summary(votaw.sem2b)
##
## Model Chisquare = 40.42317 Df = 5 Pr(>Chisq) = 1.226875e-07
## AIC = 50.42317
## BIC = 16.24176
##
## Normalized Residuals
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -2.8941 -2.2406 -2.1109 -1.4566 -1.2594 0.9564
##
## R-square for Endogenous Variables
## orig1 hcpy1 ccpy1 orig2
## 0.8165 0.4486 0.5112 0.8970
##
## Parameter Estimates
## Estimate Std Error z value Pr(>|z|)
## V[Ability] 18.338363 2.4820930 7.388266 1.487562e-13 Ability <--> Ability
## V[orig1] 4.120650 0.8122049 5.073411 3.907473e-07 orig1 <--> orig1
## V[hcpy1] 22.541770 3.0101521 7.488581 6.962195e-14 hcpy1 <--> hcpy1
## V[ccpy1] 17.531814 2.3809960 7.363227 1.795169e-13 ccpy1 <--> ccpy1
## V[orig2] 2.105483 0.6619379 3.180786 1.468761e-03 orig2 <--> orig2
##
## Iterations = 43
parallel model
votaw.mod3 <- specifyEquations(text="
orig1 = lam * Ability
hcpy1 = lam * Ability
ccpy1 = lam * Ability
orig2 = lam * Ability
V(Ability) = 1
V(orig1) = error
V(hcpy1) = error
V(ccpy1) = error
V(orig2) = error
")
votaw.sem3 <- sem(votaw.mod3, votaw, N=126)
summary(votaw.sem3)
##
## Model Chisquare = 109.123 Df = 8 Pr(>Chisq) = 5.765522e-20
## AIC = 113.123
## BIC = 70.43275
##
## Normalized Residuals
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -1.5334554 -0.5345257 -0.4070619 -0.0000006 -0.1312979 3.1134675
##
## R-square for Endogenous Variables
## orig1 hcpy1 ccpy1 orig2
## 0.5325 0.5325 0.5325 0.5325
##
## Parameter Estimates
## Estimate Std Error z value Pr(>|z|)
## lam 3.60979 0.2799109 12.89621 4.727746e-38 orig1 <--- Ability
## error 11.43997 0.8354571 13.69306 1.116981e-42 orig1 <--> orig1
##
## Iterations = 10
Compare models
anova(votaw.sem2, votaw.sem3)
## LR Test for Difference Between Models
##
## Model Df Model Chisq Df LR Chisq Pr(>Chisq)
## votaw.sem2 5 40.423
## votaw.sem3 8 109.123 3 68.7 8.103e-15 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# try with cfa
votaw.mod3a <- cfa(reference.indicators=FALSE, text="
Ability: orig1 = hcpy1 = ccpy1 = orig2
Var: orig1 = hcpy1 = ccpy1 = orig2
")
votaw.mod3a
## Path Parameter StartValue
## 1 Ability -> orig1 lam[orig1.hcpy1.ccpy1.orig2:Ability]
## 2 Ability -> hcpy1 lam[orig1.hcpy1.ccpy1.orig2:Ability]
## 3 Ability -> ccpy1 lam[orig1.hcpy1.ccpy1.orig2:Ability]
## 4 Ability -> orig2 lam[orig1.hcpy1.ccpy1.orig2:Ability]
## 5 Ability <-> Ability <fixed> 1
## 6 orig1 <-> orig1 V[orig1.hcpy1.ccpy1.orig2]
## 7 hcpy1 <-> hcpy1 V[orig1.hcpy1.ccpy1.orig2]
## 8 ccpy1 <-> ccpy1 V[orig1.hcpy1.ccpy1.orig2]
## 9 orig2 <-> orig2 V[orig1.hcpy1.ccpy1.orig2]
# try with cfa
votaw.mod3b <- cfa(reference.indicators=FALSE, subscript="number", text="
Ability: orig1 = hcpy1 = ccpy1 = orig2
Var: orig1 = hcpy1 = ccpy1 = orig2
")
votaw.mod3b
## Path Parameter StartValue
## 1 Ability -> orig1 lam[1:Ability]
## 2 Ability -> hcpy1 lam[1:Ability]
## 3 Ability -> ccpy1 lam[1:Ability]
## 4 Ability -> orig2 lam[1:Ability]
## 5 Ability <-> Ability <fixed> 1
## 6 orig1 <-> orig1 V[1]
## 7 hcpy1 <-> hcpy1 V[1]
## 8 ccpy1 <-> ccpy1 V[1]
## 9 orig2 <-> orig2 V[1]
# back to congeneric model, but test whether orig1 and orig2 are equally reliable
# semi-parallel model: originals vs copies, each equally reliable
semi-parallel: congeneric model, with some equality constraints
votaw.mod4 <- specifyEquations(
text="
orig1 = lam1 * Ability
hcpy1 = lam2 * Ability
ccpy1 = lam2 * Ability
orig2 = lam1 * Ability
V(Ability) = 1
V(orig1) = err1
V(hcpy1) = err2
V(ccpy1) = err2
V(orig2) = err1
")
votaw.sem4 <- sem(votaw.mod4, votaw, N=126)
summary(votaw.sem4)
##
## Model Chisquare = 8.986011 Df = 6 Pr(>Chisq) = 0.1743663
## AIC = 16.98601
## BIC = -20.03168
##
## Normalized Residuals
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.84784 -0.09506 -0.01233 0.09365 0.22467 0.89471
##
## R-square for Endogenous Variables
## orig1 hcpy1 ccpy1 orig2
## 0.8834 0.2792 0.2792 0.8834
##
## Parameter Estimates
## Estimate Std Error z value Pr(>|z|)
## lam1 4.553558 0.3075749 14.804714 1.365670e-49 orig1 <--- Ability
## lam2 2.666804 0.3266719 8.163556 3.253036e-16 hcpy1 <--- Ability
## err1 2.735667 0.3458230 7.910598 2.561549e-15 orig1 <--> orig1
## err2 18.358709 1.6833907 10.905792 1.081476e-27 hcpy1 <--> hcpy1
##
## Iterations = 18
# try with cfa
votaw.mod4a <- cfa(reference.indicators=FALSE, text="
Ability: orig1 = orig2, hcpy1 = ccpy1
Var: orig1 = orig2, hcpy1 = ccpy1
")
# try with cfa
votaw.mod4b <- cfa(reference.indicators=FALSE, subscript="number", text="
Ability: orig1 = orig2, hcpy1 = ccpy1
Var: orig1 = orig2, hcpy1 = ccpy1
")
votaw.mod3b
## Path Parameter StartValue
## 1 Ability -> orig1 lam[1:Ability]
## 2 Ability -> hcpy1 lam[1:Ability]
## 3 Ability -> ccpy1 lam[1:Ability]
## 4 Ability -> orig2 lam[1:Ability]
## 5 Ability <-> Ability <fixed> 1
## 6 orig1 <-> orig1 V[1]
## 7 hcpy1 <-> hcpy1 V[1]
## 8 ccpy1 <-> ccpy1 V[1]
## 9 orig2 <-> orig2 V[1]