Measurement invariance
Dimitrios Zacharatos
library(rwf)CFA measurement invariance
graphics.off()
options(width=500)
options(digits=3)Data simulation
The data are simulated in lavaan. From this simulation it is possible to extract 3 factors and 2 arbitrary groups.
# specify population model
population_model<-'f1=~x1+.6*x2+.7*x3
f2=~x4+.6*x5+.7*x6
f3=~x7+.6*x8+.7*x9'
# generate data
set.seed(1234)
model_data<-lavaan::simulateData(population_model,sample.nobs=20000)
model_data$group<-NA
model_data[1:10000,]$group<-"Group1"
model_data[10001:20000,]$group<-"Group2"Data summary
rwf::cdf(model_data,name_length=5)## $summary ## COLLUMNS ROWS TOTAL EMPTY null NAN na INF FIN FACTOR ## 1 10 20000 200000 0 0 0 0 0 180000 0 ## ## $check ## NAMES EMPTY null na NOT_NA NAN INF FIN RANGE MEAN MEDIAN SD MIN MAX MODE TYPE CLASS FACTOR ## 1 x1 0 0 0 20000 0 0 20000 20000 0 0.01 1.42 numeric double numeric FALSE ## 2 x2 0 0 0 20000 0 0 20000 20000 0 -0.01 1.17 numeric double numeric FALSE ## 3 x3 0 0 0 20000 0 0 20000 20000 0.01 0.01 1.22 numeric double numeric FALSE ## 4 x4 0 0 0 20000 0 0 20000 20000 -0.01 -0.01 1.42 numeric double numeric FALSE ## 5 x5 0 0 0 20000 0 0 20000 20000 -0.02 -0.02 1.16 numeric double numeric FALSE ## 6 x6 0 0 0 20000 0 0 20000 20000 0 0.01 1.22 numeric double numeric FALSE ## 7 x7 0 0 0 20000 0 0 20000 20000 0 -0.01 1.4 numeric double numeric FALSE ## 8 x8 0 0 0 20000 0 0 20000 20000 0 -0.01 1.17 numeric double numeric FALSE ## 9 x9 0 0 0 20000 0 0 20000 20000 0 -0.01 1.22 numeric double numeric FALSE ## 10 group 0 0 0 20000 0 0 0 2 NA NA NA character character character FALSE
Model spesification
The spesified model is shown in the figure bellow.
bf_model<-'f1=~x1+x2+x3
f2=~x4+x5+x6
f3=~x7+x8+x9'
cat(bf_model)## f1=~x1+x2+x3 ## f2=~x4+x5+x6 ## f3=~x7+x8+x9
bfa<-lavaan::cfa(bf_model,data=model_data)
semPlot::semPaths(bfa,"model",layout="tree2",what="est")
Model estimates and model fit
Model estimates are provided by three main matrices:
1. The
lambda matrix representing factor loadings
2. The theta matrix
representing standard errors for items in the diagonal, and error
covariances between items in the off diagonal
3. The psi matrix
representing covariances between factors
Model estimates are derived from a converged CFA model which in turn
can produce a covariance matrix by performing these computations:
covariance=lambda*psi+theta, which are the estimated data. In addition,
a covariance matrix can be computed directly from item responses, which
are the observed data. Because of that, we can compare the covariance
matrix estimated from the model and the covariance matrix computed from
item responses. In a good model, these two matrices should have similar
values. Fit indices, shown bellow quantify model fit. There is
conflicting information and several opinions on critical values of fit
indices. In our case we consider good model fit to be indicated by:
1. RMSEA < .07
2. CFI > .9
3. TLI > .9
4. SRMR < .8
5. GFI > .95
lavaan::fitMeasures(bfa,c("rmsea","cfi","tli","rni","rfi","ifi","srmr","gfi"))## rmsea cfi tli rni rfi ifi srmr gfi ## 0.000 1.000 1.000 1.000 0.999 1.000 0.004 1.000
lavaan::inspect(bfa,"est",add.labels=TRUE,add.class=TRUE,list.by.group=TRUE,drop.list.single.group=TRUE)## Warning: lavaan->inspect.lavaan():
## Unknown argument 'add.labels' for 'inspect'## Warning: lavaan->inspect.lavaan():
## Unknown argument 'add.class' for 'inspect'## Warning: lavaan->inspect.lavaan():
## Unknown argument 'list.by.group' for 'inspect'## Warning: lavaan->inspect.lavaan():
## Unknown argument 'drop.list.single.group' for 'inspect'## $lambda ## f1 f2 f3 ## x1 1.000 0.000 0.000 ## x2 0.598 0.000 0.000 ## x3 0.715 0.000 0.000 ## x4 0.000 1.000 0.000 ## x5 0.000 0.609 0.000 ## x6 0.000 0.693 0.000 ## x7 0.000 0.000 1.000 ## x8 0.000 0.000 0.621 ## x9 0.000 0.000 0.687 ## ## $theta ## x1 x2 x3 x4 x5 x6 x7 x8 x9 ## x1 1.032 ## x2 0.000 1.009 ## x3 0.000 0.000 0.979 ## x4 0.000 0.000 0.000 1.001 ## x5 0.000 0.000 0.000 0.000 0.983 ## x6 0.000 0.000 0.000 0.000 0.000 0.997 ## x7 0.000 0.000 0.000 0.000 0.000 0.000 0.971 ## x8 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.978 ## x9 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 1.012 ## ## $psi ## f1 f2 f3 ## f1 0.982 ## f2 -0.009 1.007 ## f3 -0.008 -0.004 0.996
Testing for invariance
We perform 2 seperate CFA parameter estimations for each group. If these models display significantly different lambda, theta, and psi estimates, that suggests non invariant models.
Test for STRUCTURAL INVARIANCE
For structural invariance, we create 2 models one using Group1 data, the other using the Group2 data, and their respective fit indices are inspected.
# Model with the first group responses
ft<-lavaan::cfa(bf_model,data=model_data[model_data$group %in% "Group1",])
# Model with the second group responses
st<-lavaan::cfa(bf_model,data=model_data[model_data$group %in% "Group2",])
data.frame(Group1=lavaan::fitMeasures(ft,c("rmsea","cfi","tli","rni","rfi","ifi","srmr","gfi")),
Group2=lavaan::fitMeasures(st,c("rmsea","cfi","tli","rni","rfi","ifi","srmr","gfi")))## Group1 Group2 ## rmsea 0.0000 0.00512 ## cfi 1.0000 0.99941 ## tli 1.0012 0.99912 ## rni 1.0008 0.99941 ## rfi 0.9978 0.99575 ## ifi 1.0008 0.99941 ## srmr 0.0049 0.00749 ## gfi 0.9997 0.99933
Test for CONFIGURAL INVARIANCE
For configural invariance we spesify a model with 2 groups. This procedure estimates 2 lambda theta psi and mu (intercepts) matrices, for each group. By comparing the two lambda and theta matrices we notice that their values are not identical. The question is if the differences between these values are significantly different. In any case, the configural invariance model bellow displays similar fit indices compared with the structural invariance models above. But more importantly, according to fit indices the model fits the data and this suggests configural invariance.
# CONFIGURAL INVARIANCE
bf<-lavaan::cfa(bf_model,data=model_data,group="group")
lavaan::fitMeasures(bf,c("rmsea","cfi","tli","rni","rfi","ifi","srmr","gfi"))## rmsea cfi tli rni rfi ifi srmr gfi ## 0.000 1.000 1.000 1.000 0.997 1.000 0.006 0.999
lavaan::inspect(bf,"est",add.labels=TRUE,add.class=TRUE,list.by.group=TRUE,drop.list.single.group=TRUE)## Warning: lavaan->inspect.lavaan():
## Unknown argument 'add.labels' for 'inspect'## Warning: lavaan->inspect.lavaan():
## Unknown argument 'add.class' for 'inspect'## Warning: lavaan->inspect.lavaan():
## Unknown argument 'list.by.group' for 'inspect'## Warning: lavaan->inspect.lavaan():
## Unknown argument 'drop.list.single.group' for 'inspect'## $Group1 ## $Group1$lambda ## f1 f2 f3 ## x1 1.000 0.000 0.000 ## x2 0.607 0.000 0.000 ## x3 0.740 0.000 0.000 ## x4 0.000 1.000 0.000 ## x5 0.000 0.609 0.000 ## x6 0.000 0.704 0.000 ## x7 0.000 0.000 1.000 ## x8 0.000 0.000 0.618 ## x9 0.000 0.000 0.680 ## ## $Group1$theta ## x1 x2 x3 x4 x5 x6 x7 x8 x9 ## x1 1.027 ## x2 0.000 1.002 ## x3 0.000 0.000 0.965 ## x4 0.000 0.000 0.000 0.989 ## x5 0.000 0.000 0.000 0.000 0.983 ## x6 0.000 0.000 0.000 0.000 0.000 0.987 ## x7 0.000 0.000 0.000 0.000 0.000 0.000 0.985 ## x8 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.976 ## x9 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 1.024 ## ## $Group1$psi ## f1 f2 f3 ## f1 0.977 ## f2 -0.019 0.995 ## f3 -0.007 -0.022 0.968 ## ## $Group1$nu ## intrcp ## x1 0.004 ## x2 -0.012 ## x3 0.013 ## x4 -0.005 ## x5 -0.011 ## x6 0.003 ## x7 -0.011 ## x8 -0.002 ## x9 -0.003 ## ## $Group1$alpha ## intrcp ## f1 0 ## f2 0 ## f3 0 ## ## ## $Group2 ## $Group2$lambda ## f1 f2 f3 ## x1 1.000 0.000 0.000 ## x2 0.588 0.000 0.000 ## x3 0.690 0.000 0.000 ## x4 0.000 1.000 0.000 ## x5 0.000 0.609 0.000 ## x6 0.000 0.683 0.000 ## x7 0.000 0.000 1.000 ## x8 0.000 0.000 0.624 ## x9 0.000 0.000 0.692 ## ## $Group2$theta ## x1 x2 x3 x4 x5 x6 x7 x8 x9 ## x1 1.035 ## x2 0.000 1.016 ## x3 0.000 0.000 0.992 ## x4 0.000 0.000 0.000 1.014 ## x5 0.000 0.000 0.000 0.000 0.983 ## x6 0.000 0.000 0.000 0.000 0.000 1.007 ## x7 0.000 0.000 0.000 0.000 0.000 0.000 0.957 ## x8 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.980 ## x9 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 1.001 ## ## $Group2$psi ## f1 f2 f3 ## f1 0.988 ## f2 0.002 1.018 ## f3 -0.009 0.015 1.024 ## ## $Group2$nu ## intrcp ## x1 0.004 ## x2 0.019 ## x3 0.004 ## x4 -0.014 ## x5 -0.022 ## x6 0.003 ## x7 0.012 ## x8 0.003 ## x9 0.006 ## ## $Group2$alpha ## intrcp ## f1 0 ## f2 0 ## f3 0
par(mfrow=c(1,2))
semPlot::semPaths(bf,"model",layout="tree2",what="est")
In order to test if differences between estimated matrices among models are significant, we have to create a series of models for which impose equality constraints for parameters between groups. By imposing constraints we impose the estimator procedure to keep a matrix between groups equal. For instance by imposing an equality constraint on loadings, we force the parameter estimator to produce models that have equal lambda matrices between groups. Finally we compare the unconstrained and constrained models to see if they have any differences. If no significant differences exist, this suggest that the models are equivalent between groups and thus the models display measurement invariance.
Test for WEAK METRIC INVARIANCE
For weak metric invariance we impose a restriction that groups have equal loadings.
# WEAK METRIC INVARIANCE
bf_metric<-lavaan::cfa(bf_model,data=model_data,group="group",group.equal=c("loadings"))
lavaan::anova(bf,bf_metric)## ## Chi-Squared Difference Test ## ## Df AIC BIC Chisq Chisq diff RMSEA Df diff Pr(>Chisq) ## bf 48 573054 573528 46.0 ## bf_metric 54 573045 573472 48.8 2.89 0 6 0.82
Test for STRONG SCALAR INVARIANCE
For strong scalar invariance we impose the restriction that groups have equal loadings and intercepts and then compare the strong scalar invariance model with the weak metric invariance model.
# STRONG SCALAR INVARIANCE
bf_scalar<-lavaan::cfa(bf_model,data=model_data,group="group",group.equal=c("loadings","intercepts"))
lavaan::anova(bf_metric,bf_scalar)## ## Chi-Squared Difference Test ## ## Df AIC BIC Chisq Chisq diff RMSEA Df diff Pr(>Chisq) ## bf_metric 54 573045 573472 48.8 ## bf_scalar 60 573038 573418 54.4 5.53 0 6 0.48
Test for STRICT ERROR VARIANCE INVARIANCE
For strict error variance invariance, we impose the restriction that groups have equal loadings intercepts and residuals.
# STRICT ERROR VARIANCE INVARIANCE
bf_full_strict<-lavaan::cfa(bf_model,data=model_data,group="group",group.equal=c("loadings","intercepts","residuals"))
lavaan::anova(bf_scalar,bf_full_strict)## ## Chi-Squared Difference Test ## ## Df AIC BIC Chisq Chisq diff RMSEA Df diff Pr(>Chisq) ## bf_scalar 60 573038 573418 54.4 ## bf_full_strict 69 573028 573336 61.7 7.36 0 9 0.6
Test for FACTOR VARIANCE COVARIANCE
For factor variance covariance invariance we impose the restriction that groups have equal loadings intercepts residuals latent variable variances and latent variable covariances.
# FACTOR VARIANCE COVARIANCE
bf_factor_vc<-lavaan::cfa(bf_model,data=model_data,group="group",group.equal=c("loadings","intercepts","residuals","lv.variances","lv.covariances"))
lavaan::anova(bf_full_strict,bf_factor_vc)## ## Chi-Squared Difference Test ## ## Df AIC BIC Chisq Chisq diff RMSEA Df diff Pr(>Chisq) ## bf_full_strict 69 573028 573336 61.7 ## bf_factor_vc 75 573025 573285 70.7 8.97 0.00703 6 0.18
Test for FACTOR MEAN INVARIANCE
For factor mean invariance we impose the restriction that groups have equal loadings intercepts and residuals and we impose the restriction that factors have equal variances covariances and means between groups.
# FACTOR MEAN INVARIANCE
bf_factor_mean<-lavaan::cfa(bf_model,data=model_data,group="group",group.equal=c("loadings","intercepts","residuals","lv.variances","lv.covariances","means"))
lavaan::anova(bf_factor_vc,bf_factor_mean)## ## Chi-Squared Difference Test ## ## Df AIC BIC Chisq Chisq diff RMSEA Df diff Pr(>Chisq) ## bf_factor_vc 75 573025 573285 70.7 ## bf_factor_mean 78 573020 573257 72.0 1.33 0 3 0.72
Comparing models
Finally we compare models using an ANOVA procedure where values above 0.05 indicate no significant differences and thus measurement invariance.
lavaan::anova(bf,bf_metric,bf_scalar,bf_factor_vc,bf_factor_mean,SB.classic=TRUE)## ## Chi-Squared Difference Test ## ## Df AIC BIC Chisq Chisq diff RMSEA Df diff Pr(>Chisq) ## bf 48 573054 573528 46.0 ## bf_metric 54 573045 573472 48.8 2.89 0.00000 6 0.82 ## bf_scalar 60 573038 573418 54.4 5.53 0.00000 6 0.48 ## bf_factor_vc 75 573025 573285 70.7 16.33 0.00297 15 0.36 ## bf_factor_mean 78 573020 573257 72.0 1.33 0.00000 3 0.72
In addition, the fit indices for each model above are summarized bellow.
findex<-c("rmsea","cfi","tli","rni","rfi","ifi","srmr","gfi")
fit_indices<-data.frame(rbind(lavaan::fitMeasures(bf,findex),
lavaan::fitMeasures(bf_metric,findex),
lavaan::fitMeasures(bf_scalar,findex),
lavaan::fitMeasures(bf_factor_vc,findex),
lavaan::fitMeasures(bf_factor_mean,findex)))
row.names(fit_indices)<-c("CONFIGURAL INVARIANCE",
"WEAK METRIC INVARIANCE",
"STRONG SCALAR INVARIANCE",
"FACTOR VARIANCE COVARIANCE",
"FACTOR MEAN INVARIANCE")
fit_indices## rmsea cfi tli rni rfi ifi srmr gfi ## CONFIGURAL INVARIANCE 0 1 1 1 0.997 1 0.00566 0.999 ## WEAK METRIC INVARIANCE 0 1 1 1 0.997 1 0.00588 0.999 ## STRONG SCALAR INVARIANCE 0 1 1 1 0.997 1 0.00619 0.999 ## FACTOR VARIANCE COVARIANCE 0 1 1 1 NA 1 0.00853 0.999 ## FACTOR MEAN INVARIANCE 0 1 1 1 NA 1 0.00865 0.999