# principles and practice of sem, 5th ed. # rex b. kline, guilford press, 2023 # chapter 22, table 22.1, analysis 2 # two-factor cfa model of divergent thinking # analyzed simultaneously over chinese and american samples # to avoid the problem that some R packages share the same name # for different functions, all functions are specified next as # package::function # which prevents masking, or the default hiding of # a function with a redundant name from a package used next date() v <- R.Version() print(paste0(v$language, " version ", v$major, ".", v$minor, " (", v$year, "-", v$month, "-", v$day, ")")) library(lavaan) # get citation information citation("lavaan", auto = TRUE) # name subtests subtest.name <- c("LM1", "LM2", "LM3", "RW1", "RW2") # chinese, n = 316, originality chineseLower.cor <- ' 1.000 .715 1.000 .631 .738 1.000 .377 .408 .319 1.000 .452 .495 .391 .445 1.000 ' # name the variables and convert to full correlation matrix chinese.cor <- lavaan::getCov(chineseLower.cor, names = subtest.name) # add the standard deviations and convert to covariances chinese.sd <- c(1.532, 1.878, 2.012, 1.002, 1.055) chinese.cov <- lavaan::cor2cov(chinese.cor, sds = chinese.sd) # create mean vector chinese.mean = c(1.24, 1.47, 1.86, .77, .89) # american, n = 302, originality americanLower.cor <- ' 1.000 .519 1.000 .543 .639 1.000 .280 .321 .317 1.000 .242 .380 .300 .665 1.000 ' # name the variables and convert to full correlation matrix american.cor <- lavaan::getCov(americanLower.cor, names = subtest.name) # add the standard deviations and convert to covariances american.sd <- c(1.094, 1.350, 1.494, 2.257, 1.541) american.cov <- lavaan::cor2cov(american.cor, sds = american.sd) # create mean vector american.mean = c(.88, 1.35, 1.27, 1.71, 1.37) # combine covariances, group sizes, and means # into single list objects combined.cov <- list(chinese = chinese.cov, american = american.cov) combined.n <- list(chinese = 316, american = 302) combined.mean <- list(chinese = chinese.mean, american = american.mean) # specify configural invariance model configural.model <- ' # turn off (NA) default reference variable scaling for # first indicator of each factor # label loadings and intercepts for both factors # over both groups LineMeaning =~ NA*LM1 + c(a1,a2)*LM1 + c(b1,b2)*LM2 + c(c1,c2)*LM3 LM1 ~ c(d1,d2)*1 LM2 ~ c(e1,e2)*1 LM3 ~ c(f1,f2)*1 RealWorld =~ NA*RW1 + c(g1,g2)*RW1 + c(h1,h2)*RW2 RW1 ~ c(i1,i2)*1 RW2 ~ c(j1,j2)*1 # free factor means in both groups LineMeaning ~ c(m1,m2)*1 RealWorld ~ c(q1,q2)*1 ' # specify configural invariance model constraints # all constraints are for effects coding only # average loading is 1, average intercept is zero # there are no cross-group equality constraints configural.constraints <- ' # LineMeaning, chinese, loadings and intercepts 3 - a1 - b1 - c1 == 0 0 - d1 - e1 - f1 == 0 # LineMeaning, american, loadings and intercepts 3 - a2 - b2 - c2 == 0 0 - d2 - e2 - f2 == 0 # RealWorld, chinese, loadings and intercepts 2 - g1 - h1 == 0 0 - i1 - j1 == 0 # Real World, american, loadings and intercepts 2 - g2 - h2 == 0 0 - i2 - j2 == 0 ' # specify weak invariance model constraints # when loadings are constrained to equality over groups, # release the effects coding scaling constraint that the # average loading is 1.0 in one group (here, american) weak.constraints <- ' # LineMeaning, chinese, loadings and intercepts 3 - a1 - b1 - c1 == 0 0 - d1 - e1 - f1 == 0 # LineMeaning, american, loadings and intercepts # effects coding scaling constraint for # loadings is released in this group by "commenting # out" the original specification # 3 - a2 - b2 - c2 == 0 0 - d2 - e2 - f2 == 0 # RealWorld, chinese, loadings and intercepts 2 - g1 - h1 == 0 0 - i1 - j1 == 0 # RealWorld, american, loadings and intercepts # release constraint for 2nd factor in this group # 2 - g2 - h2 == 0 0 - i2 - j2 == 0 # add equality constraints for loadings a1 == a2 b1 == b2 c1 == c2 g1 == g2 h1 == h2 ' # specify constraints for strong invariance model # when intercepts are also constrained to equality over groups, # release the effects coding scaling constraint that the # average loading is 0 in one group (here, american) strong.constraints <- ' # LineMeaning, chinese, loadings and intercepts 3 - a1 - b1 - c1 == 0 0 - d1 - e1 - f1 == 0 # LineMeaning, american, loadings and intercepts # both effects coding constraints are released in this group # 3 - a2 - b2 - c2 == 0 # 0 - d2 - e2 - f2 == 0 # RealWorld, chinese, loadings and intercepts 2 - g1 - h1 == 0 0 - i1 - j1 == 0 # RealWorld, american, loadings and intercepts # both effects coding constraints are released in this group # 2 - g2 - h2 == 0 # 0 - i2 - j2 == 0 # equality constraints for loadings a1 == a2 b1 == b2 c1 == c2 g1 == g2 h1 == h2 # add equality constraints for intecepts d1 == d2 e1 == e2 f1 == f2 i1 == i2 j1 == j2 ' # specify constraints for partial strong invariance model # release equality constraints for two indicators # with the highest modification indices # LM2 and RW1 partialStrong.constraints <- ' # LineMeaning, chinese, loadings and intercepts 3 - a1 - b1 - c1 == 0 0 - d1 - e1 - f1 == 0 # LineMeaning, american, loadings and intercepts # both effects coding constraints are released in this group # 3 - a2 - b2 - c2 == 0 # 0 - d2 - e2 - f2 == 0 # RealWorld, chinese, loadings and intercepts 2 - g1 - h1 == 0 0 - i1 - j1 == 0 # # RealWorld, american, loadings and intercepts # both effects coding constraints are released in this group # 2 - g2 - h2 == 0 # 0 - i2 - j2 == 0 # equality constraints for loadings a1 == a2 b1 == b2 c1 == c2 g1 == g2 h1 == h2 # add equality constraints for intecepts d1 == d2 # release equality constraint on LM2 intercept # e1 == e2 f1 == f2 # release equality constraint on RW1 intercept # i1 == i2 j1 == j2 ' # fit models to data configural <- lavaan::cfa(configural.model, sample.cov = combined.cov, sample.mean = combined.mean, sample.nobs = combined.n, meanstructure = TRUE, constraints = configural.constraints) weak <- lavaan::cfa(configural.model, sample.cov = combined.cov, sample.mean = combined.mean, sample.nobs = combined.n, meanstructure = TRUE, constraints = weak.constraints) strong <- lavaan::cfa(configural.model, sample.cov = combined.cov, sample.mean = combined.mean, sample.nobs = combined.n, meanstructure = TRUE, constraints = strong.constraints) partialStrong <- lavaan::cfa(configural.model, sample.cov = combined.cov, sample.mean = combined.mean, sample.nobs = combined.n, meanstructure = TRUE, constraints = partialStrong.constraints) # generate chi-square difference tests anova(configural, weak) anova(weak, strong) anova(weak, partialStrong) # configural model estimates, predicted covariances # and means, and residuals lavaan::summary(configural, standardized = TRUE, fit.measures = TRUE) lavaan::fitted(configural) lavaan::residuals(configural, type = "raw") lavaan::residuals(configural, type = "standardized") lavaan::residuals(configural, type = "cor.bollen") # weak model estimates, predicted covariances # and means, and residuals lavaan::summary(weak, standardized = TRUE, fit.measures = TRUE) lavaan::fitted(weak) lavaan::residuals(weak, type = "raw") lavaan::residuals(weak, type = "standardized") lavaan::residuals(weak, type = "cor.bollen") # strong model estimates, predicted covariances # and means, residuals, and modification indices # for equality-constrained parameters lavaan::summary(strong, standardized = TRUE, fit.measures = TRUE) lavaan::fitted (strong) # observed means and standard deviations chinese.mean chinese.sd american.mean american.sd lavaan::residuals(strong, type = "raw") lavaan::residuals(strong, type = "standardized") lavaan::residuals(strong, type = "cor.bollen") # equality constraints for intercepts # are labeled d, e, f, i, j lavaan::lavTestScore(weak, univariate = TRUE) # partial strong model estimates, predicted covariances # and means, residuals, and modification indices # for equality-constrained parameters lavaan::summary(partialStrong, standardized = TRUE, fit.measures = TRUE) # observed means, chinese, american chinese.mean american.mean # predicted covariances and means lavaan::fitted (partialStrong) lavaan::residuals(partialStrong, type = "raw") lavaan::residuals(partialStrong, type = "standardized") lavaan::residuals(partialStrong, type = "cor.bollen") lavaan::lavTestScore(partialStrong, univariate = TRUE)