diff --git a/15-Factor-Analysis-PCA.Rmd b/15-Factor-Analysis-PCA.Rmd index c05e61ec..41c9db9e 100644 --- a/15-Factor-Analysis-PCA.Rmd +++ b/15-Factor-Analysis-PCA.Rmd @@ -51,1499 +51,3 @@ HolzingerSwineford1939 <- as.data.frame(matrix( names(HolzingerSwineford1939) <- varNames vars <- c("x1","x2","x3","x4","x5","x6","x7","x8","x9") ``` - -## Descriptive Statistics and Correlations {#descriptiveStatsCorrelations-factorAnalysis} - -Before conducting a [factor analysis](#factorAnalysis), it is important examine descriptive statistics and correlations among variables.\index{correlation}\index{descriptive statistics} - -### Descriptive Statistics {#descriptiveStats-factorAnalysis} - -Descriptive statistics are presented in Table \@ref(tab:descriptiveStats).\index{descriptive statistics} - -```{r, eval = knitr::is_html_output(excludes = "epub"), echo = knitr::is_html_output(excludes = "epub")} -paged_table(psych::describe(HolzingerSwineford1939[,vars])) -``` - -```{r, eval = knitr::is_latex_output(), echo = knitr::is_latex_output()} -kable(psych::describe(HolzingerSwineford1939[,vars]), - format = "latex", - booktabs = TRUE) %>% - kable_styling(latex_options = "scale_down") -``` - -```{r, eval = knitr::is_html_output(excludes = c("markdown","html","html4","html5","revealjs","s5","slideous","slidy","gfm")), echo = knitr::is_html_output(excludes = c("markdown","html","html4","html5","revealjs","s5","slideous","slidy","gfm"))} -kable(psych::describe(HolzingerSwineford1939[,vars]), - booktabs = TRUE) -``` - -```{r descriptiveStats} -summaryTable <- HolzingerSwineford1939 %>% - dplyr::select(all_of(vars)) %>% - summarise(across( - everything(), - .fns = list( - n = ~ length(na.omit(.)), - missingness = ~ mean(is.na(.)) * 100, - M = ~ mean(., na.rm = TRUE), - SD = ~ sd(., na.rm = TRUE), - min = ~ min(., na.rm = TRUE), - max = ~ max(., na.rm = TRUE), - skewness = ~ psych::skew(., na.rm = TRUE), - kurtosis = ~ kurtosi(., na.rm = TRUE)), - .names = "{.col}.{.fn}")) %>% - pivot_longer( - cols = everything(), - names_to = c("variable","index"), - names_sep = "\\.") %>% - pivot_wider( - names_from = index, - values_from = value) - -summaryTableTransposed <- summaryTable[-1] %>% - t() %>% - as.data.frame() %>% - setNames(summaryTable$variable) %>% - round(., digits = 2) - -summaryTableTransposed -``` - -### Correlations {#correlations-factorAnalysis} - -```{r} -cor( - HolzingerSwineford1939[,vars], - use = "pairwise.complete.obs") -``` - -Correlation matrices of various types using the `cor.table()` function from the [`petersenlab`](https://github.com/DevPsyLab/petersenlab) package [@R-petersenlab] are in Tables \@ref(tab:corTable1), \@ref(tab:corTable2), and \@ref(tab:corTable3).\index{petersenlab package}\index{correlation} - -```{r, eval = FALSE} -cor.table( - HolzingerSwineford1939[,vars], - dig = 2) - -cor.table( - HolzingerSwineford1939[,vars], - type = "manuscript", - dig = 2) - -cor.table( - HolzingerSwineford1939[,vars], - type = "manuscriptBig", - dig = 2) -``` - -```{r, include = FALSE} -corTable1 <- cor.table( - HolzingerSwineford1939[,vars], - dig = 2) - -corTable2 <- cor.table( - HolzingerSwineford1939[,vars], - type = "manuscript", - dig = 2) - -corTable3 <- cor.table( - HolzingerSwineford1939[,vars], - type = "manuscriptBig", - dig = 2) -``` - -```{r corTable1, echo = FALSE} -corTable1 %>% - kable(., - caption = "Correlation Matrix with *r*, *n*, and *p*-values.", - booktabs = TRUE, - linesep = c("", "", "\\addlinespace"), - escape = FALSE) -``` - -```{r corTable2, echo = FALSE} -corTable2 %>% - kable(., - caption = "Correlation Matrix with Asterisks for Significant Associations.", - booktabs = TRUE, - linesep = "", - escape = FALSE) -``` - -```{r corTable3, echo = FALSE} -corTable3 %>% - kable(., - caption = "Correlation Matrix.", - booktabs = TRUE, - linesep = "") -``` - -Pairs panel plots were generated using the `psych` package [@R-psych].\index{correlation} -Correlation plots were generated using the `corrplot` package [@R-corrplot].\index{correlation} - -A pairs panel plot is in Figure \@ref(fig:pairsPanel).\index{correlation} - -```{r pairsPanel, out.width = "100%", fig.align = "center", fig.cap = "Pairs Panel Plot."} -pairs.panels(HolzingerSwineford1939[,vars]) -``` - -A correlation plot is in Figure \@ref(fig:correlationPlot).\index{correlation} - -```{r correlationPlot, out.width = "100%", fig.align = "center", fig.cap = "Correlation Plot."} -corrplot(cor( - HolzingerSwineford1939[,vars], - use = "pairwise.complete.obs")) -``` - -## Factor Analysis {#factorAnalysisExamples} - -### Exploratory Factor Analysis (EFA) {#efaExamples} - -I introduced [exploratory factor analysis](#efa) (EFA) models in Section \@ref(efa).\index{factor analysis!exploratory} - -#### Determine number of factors {#determineNumberOfFactors-efaExamples} - -Determine the number of factors to retain using the [Scree plot](#screePlot) and [Very Simple Structure plot](#vssPlot).\index{scree plot}\index{very simple structure plot} - -##### Scree Plot {#screePlot} - -Scree plots were generated using the `psych` [@R-psych] and `nFactors` [@R-nFactors] packages.\index{scree plot} -The optimal coordinates and the acceleration factor attempt to operationalize the Cattell scree test: i.e., the "elbow" of the scree plot [@Ruscio2012].\index{scree plot} -The optimal coordinators factor is quantified using a series of linear equations to determine whether observed eigenvalues exceed the predicted values.\index{scree plot} -The acceleration factor is quantified using the acceleration of the curve, that is, the second derivative.\index{scree plot} -The Kaiser-Guttman rule states to keep principal components whose eigenvalues are greater than 1.\index{scree plot}\index{eigenvalue} -However, for [exploratory factor analysis](#efa) (as opposed to [PCA](#pca)), the criterion is to keep the factors whose eigenvalues are greater than zero (i.e., not the factors whose eigenvalues are greater than 1) [@Dinno2014].\index{scree plot}\index{factor analysis!exploratory}\index{principal component analysis}\index{eigenvalue} - -The number of factors to keep would depend on which criteria one uses.\index{scree plot} -Based on the rule to keep factors whose eigenvalues are greater than zero and based on the parallel test, we would keep three factors.\index{scree plot}\index{eigenvalue} -However, based on the Cattell scree test (as operationalized by the optimal coordinates and acceleration factor), we would keep one factor.\index{scree plot} -Therefore, interpretability of the factors would be important for deciding between whether to keep one, two, or three factors.\index{scree plot} - -A scree plot from a parallel analysis is in Figure \@ref(fig:screePlotEFAparallel).\index{scree plot}\index{factor analysis!exploratory} - -```{r screePlotEFAparallel, out.width = "100%", fig.align = "center", fig.cap = "Scree Plot from Parallel Analysis in Exploratory Factor Analysis."} -fa.parallel( - x = HolzingerSwineford1939[,vars], - fm = "ml", - fa = "fa") -``` - -A scree plot from EFA is in Figure \@ref(fig:screePlotEFA).\index{scree plot}\index{factor analysis!exploratory} - -```{r screePlotEFA, out.width = "100%", fig.align = "center", fig.cap = "Scree Plot in Exploratory Factor Analysis."} -plot( - nScree(x = cor( - HolzingerSwineford1939[,vars], - use = "pairwise.complete.obs"), - model = "factors")) -``` - -##### Very Simple Structure (VSS) Plot {#vssPlot} - -The very simple structure (VSS) is another criterion that can be used to determine the optimal number of factors or components to retain.\index{very simple structure plot} -Using the VSS criterion, the optimal number of factors to retain is the number of factors that maximizes the VSS criterion [@Revelle1979].\index{very simple structure plot} -The VSS criterion is evaluated with models in which factor loadings for a given item that are less than the maximum factor loading for that item are suppressed to zero, thus forcing simple structure (i.e., no cross-loadings).\index{cross-loading}\index{very simple structure plot} -The goal is finding a factor structure with interpretability so that factors are clearly distinguishable.\index{very simple structure plot} -Thus, we want to identify the number of factors with the highest VSS criterion (i.e., the highest line).\index{very simple structure plot} -Very simple structure (VSS) plots were generated using the `psych` package [@R-psych].\index{very simple structure plot} - -The output also provides additional criteria by which to determine the optimal number of factors, each for which lower values are better, including the Velicer minimum average partial (MAP) test, the Bayesian information criterion (BIC), the sample size-adjusted BIC (SABIC), and the root mean square error of approximation (RMSEA).\index{very simple structure plot} - -###### Orthogonal (Varimax) rotation {#orthogonal-vss} - -In the example with orthogonal rotation below, the VSS criterion is highest with 3 or 4 factors.\index{very simple structure plot}\index{orthogonal rotation} -A three-factor solution is supported by the lowest BIC, whereas a four-factor solution is supported by the lowest SABIC.\index{very simple structure plot}\index{orthogonal rotation} - -A VSS plot is in Figure \@ref(fig:vssOrthogonal).\index{very simple structure plot}\index{orthogonal rotation} - -```{r vssOrthogonal, out.width = "100%", fig.align = "center", fig.cap = "Very Simple Structure Plot With Orthogonal Rotation in Exploratory Factor Analysis."} -vss( - HolzingerSwineford1939[,vars], - rotate = "varimax", - fm = "ml") -``` - -Multiple VSS-related fit indices are in Figure \@ref(fig:vssComplexOrthogonal).\index{very simple structure plot}\index{orthogonal rotation} - -```{r vssComplexOrthogonal, out.width = "100%", fig.align = "center", fig.cap = "Very Simple Structure-Related Indices With Orthogonal Rotation in Exploratory Factor Analysis."} -nfactors( - HolzingerSwineford1939[,vars], - rotate = "varimax", - fm = "ml") -``` - -###### Oblique (Oblimin) rotation {#oblique-vss} - -In the example with oblique rotation below, the VSS criterion is highest with 3 factors.\index{very simple structure plot}\index{oblique rotation} -A three-factor solution is supported by the lowest BIC.\index{very simple structure plot}\index{oblique rotation} - -A VSS plot is in Figure \@ref(fig:vssOblique).\index{very simple structure plot}\index{oblique rotation} - -```{r vssOblique, out.width = "100%", fig.align = "center", fig.cap = "Very Simple Structure Plot With Oblique Rotation in Exploratory Factor Analysis."} -vss( - HolzingerSwineford1939[,vars], - rotate = "oblimin", - fm = "ml") -``` - -Multiple VSS-related fit indices are in Figure \@ref(fig:vssComplexOblique).\index{very simple structure plot}\index{oblique rotation} - -```{r vssComplexOblique, out.width = "100%", fig.align = "center", fig.cap = "Very Simple Structure-Related Indices With Oblique Rotation in Exploratory Factor Analysis."} -nfactors( - HolzingerSwineford1939[,vars], - rotate = "oblimin", - fm = "ml") -``` - -###### No rotation {#noRotation-vss} - -In the example with no rotation below, the VSS criterion is highest with 3 or 4 factors.\index{very simple structure plot} -A three-factor solution is supported by the lowest BIC, whereas a four-factor solution is supported by the lowest SABIC.\index{very simple structure plot} - -A VSS plot is in Figure \@ref(fig:vssUnrotated).\index{very simple structure plot} - -```{r vssUnrotated, out.width = "100%", fig.align = "center", fig.cap = "Very Simple Structure Plot With no Rotation in Exploratory Factor Analysis."} -nfactors( - HolzingerSwineford1939[,vars], - rotate = "none", - fm = "ml") -``` - -Multiple VSS-related fit indices are in Figure \@ref(fig:vssComplexUnrotated).\index{very simple structure plot} - -```{r vssComplexUnrotated, out.width = "100%", fig.align = "center", fig.cap = "Very Simple Structure-Related Indices With no Rotation in Exploratory Factor Analysis."} -nfactors( - HolzingerSwineford1939[,vars], - rotate = "none", - fm = "ml") -``` - -#### Run factor analysis {#runFactorAnalysis-efa} - -[Exploratory factor analysis](#efa) (EFA) models were fit using the `fa()` function of the `psych` package [@R-psych] and the `sem()` and `efa()` functions of the `lavaan` package [@R-lavaan].\index{factor analysis!exploratory} - -##### Orthogonal (Varimax) rotation {#orthogonal-efa} - -###### `psych` {#orthogonalPsych-efa} - -Fit a different model with each number of possible factors:\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r} -efa1factorOrthogonal <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 1, - rotate = "varimax", - fm = "ml") - -efa2factorOrthogonal <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 2, - rotate = "varimax", - fm = "ml") - -efa3factorOrthogonal <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 3, - rotate = "varimax", - fm = "ml") - -efa4factorOrthogonal <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 4, - rotate = "varimax", - fm = "ml") - -efa5factorOrthogonal <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 5, - rotate = "varimax", - fm = "ml") - -efa6factorOrthogonal <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 6, - rotate = "varimax", - fm = "ml") - -efa7factorOrthogonal <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 7, - rotate = "varimax", - fm = "ml") - -efa8factorOrthogonal <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 8, - rotate = "varimax", - fm = "ml") - -efa9factorOrthogonal <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 9, - rotate = "varimax", - fm = "ml") -``` - -###### `lavaan` {#orthogonalLavaan-efa} - -Model syntax is specified below:\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r} -efa1factorLavaan_syntax <- ' - # EFA Factor Loadings - efa("efa1")*f1 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 -' - -efa2factorLavaan_syntax <- ' - # EFA Factor Loadings - efa("efa1")*f1 + - efa("efa1")*f2 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 -' - -efa3factorLavaan_syntax <- ' - # EFA Factor Loadings - efa("efa1")*f1 + - efa("efa1")*f2 + - efa("efa1")*f3 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 -' - -efa4factorLavaan_syntax <- ' - # EFA Factor Loadings - efa("efa1")*f1 + - efa("efa1")*f2 + - efa("efa1")*f3 + - efa("efa1")*f4 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 -' - -efa5factorLavaan_syntax <- ' - # EFA Factor Loadings - efa("efa1")*f1 + - efa("efa1")*f2 + - efa("efa1")*f3 + - efa("efa1")*f4 + - efa("efa1")*f5 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 -' - -efa6factorLavaan_syntax <- ' - # EFA Factor Loadings - efa("efa1")*f1 + - efa("efa1")*f2 + - efa("efa1")*f3 + - efa("efa1")*f4 + - efa("efa1")*f5 + - efa("efa1")*f6 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 -' - -efa7factorLavaan_syntax <- ' - # EFA Factor Loadings - efa("efa1")*f1 + - efa("efa1")*f2 + - efa("efa1")*f3 + - efa("efa1")*f4 + - efa("efa1")*f5 + - efa("efa1")*f6 + - efa("efa1")*f7 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 -' - -efa8factorLavaan_syntax <- ' - # EFA Factor Loadings - efa("efa1")*f1 + - efa("efa1")*f2 + - efa("efa1")*f3 + - efa("efa1")*f4 + - efa("efa1")*f5 + - efa("efa1")*f6 + - efa("efa1")*f7 + - efa("efa1")*f8 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 -' - -efa9factorLavaan_syntax <- ' - # EFA Factor Loadings - efa("efa1")*f1 + - efa("efa1")*f2 + - efa("efa1")*f3 + - efa("efa1")*f4 + - efa("efa1")*f5 + - efa("efa1")*f6 + - efa("efa1")*f7 + - efa("efa1")*f8 + - efa("efa1")*f9 =~ x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 -' -``` - -The models are fitted below:\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r} -efa1factorOrthogonalLavaan_fit <- sem( - efa1factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "varimax", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE)) - -efa2factorOrthogonalLavaan_fit <- sem( - efa2factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "varimax", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE)) - -efa3factorOrthogonalLavaan_fit <- sem( - efa3factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "varimax", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE)) - -efa4factorOrthogonalLavaan_fit <- sem( - efa4factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "varimax", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE)) - -efa5factorOrthogonalLavaan_fit <- sem( - efa5factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "varimax", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE)) - -efa6factorOrthogonalLavaan_fit <- sem( - efa6factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "varimax", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE)) - -efa7factorOrthogonalLavaan_fit <- sem( - efa7factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "varimax", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE)) - -efa8factorOrthogonalLavaan_fit <- sem( - efa8factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "varimax", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE)) - -efa9factorOrthogonalLavaan_fit <- sem( - efa9factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "varimax", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE)) -``` - -The `efa()` wrapper can fit multiple orthogonal [EFA](#efa) models in one function call:\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r} -efaOrthogonalLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 1:5, - rotation = "varimax", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE) -) -``` - -The `efa()` wrapper can also fit individual orthogonal EFA models (with `output = "lavaan"`):\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r} -efaFactor1OrthogonalLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 1, - rotation = "varimax", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE), - output = "lavaan" -) - -efaFactor2OrthogonalLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 2, - rotation = "varimax", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE), - output = "lavaan" -) - -efaFactor3OrthogonalLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 3, - rotation = "varimax", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE), - output = "lavaan" -) - -efaFactor4OrthogonalLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 4, - rotation = "varimax", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE), - output = "lavaan" -) - -efaFactor5OrthogonalLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 5, - rotation = "varimax", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = list(orthogonal = TRUE), - output = "lavaan" -) -``` - -##### Oblique (Oblimin) rotation {#oblique-efa} - -###### `psych` {#obliquePsych-efa} - -Fit a different model with each number of possible factors:\index{factor analysis!exploratory}\index{oblique rotation} - -```{r} -efa1factorOblique <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 1, - rotate = "oblimin", - fm = "ml") - -efa2factorOblique <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 2, - rotate = "oblimin", - fm = "ml") - -efa3factorOblique <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 3, - rotate = "oblimin", - fm = "ml") - -efa4factorOblique <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 4, - rotate = "oblimin", - fm = "ml") - -efa5factorOblique <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 5, - rotate = "oblimin", - fm = "ml") - -efa6factorOblique <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 6, - rotate = "oblimin", - fm = "ml") - -efa7factorOblique <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 7, - rotate = "oblimin", - fm = "ml") - -efa8factorOblique <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 8, - rotate = "oblimin", - fm = "ml") #no convergence - -efa9factorOblique <- fa( - r = HolzingerSwineford1939[,vars], - nfactors = 9, - rotate = "oblimin", - fm = "ml") -``` - -###### `lavaan` {#obliqueLavaan-efa} - -The models are fitted below:\index{factor analysis!exploratory}\index{oblique rotation} - -```{r} -# settings to mimic Mplus -mplusRotationArgs <- - list( - rstarts = 30, - row.weights = "none", - algorithm = "gpa", - orthogonal = FALSE, - jac.init.rot = TRUE, - std.ov = TRUE, # row standard = correlation - geomin.epsilon = 0.0001) - -efa1factorObliqueLavaan_fit <- sem( - efa1factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "geomin", - meanstructure = TRUE, - rotation.args = mplusRotationArgs) - -efa2factorObliqueLavaan_fit <- sem( - efa2factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "geomin", - meanstructure = TRUE, - rotation.args = mplusRotationArgs) - -efa3factorObliqueLavaan_fit <- sem( - efa3factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "geomin", - meanstructure = TRUE, - rotation.args = mplusRotationArgs) - -efa4factorObliqueLavaan_fit <- sem( - efa4factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "geomin", - meanstructure = TRUE, - rotation.args = mplusRotationArgs) - -efa5factorObliqueLavaan_fit <- sem( - efa5factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "geomin", - meanstructure = TRUE, - rotation.args = mplusRotationArgs) - -efa6factorObliqueLavaan_fit <- sem( - efa6factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "geomin", - meanstructure = TRUE, - rotation.args = mplusRotationArgs) - -efa7factorObliqueLavaan_fit <- sem( - efa7factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "geomin", - meanstructure = TRUE, - rotation.args = mplusRotationArgs) - -efa8factorObliqueLavaan_fit <- sem( - efa8factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "geomin", - meanstructure = TRUE, - rotation.args = mplusRotationArgs) - -efa9factorObliqueLavaan_fit <- sem( - efa9factorLavaan_syntax, - data = HolzingerSwineford1939, - information = "observed", - missing = "ML", - estimator = "MLR", - rotation = "geomin", - meanstructure = TRUE, - rotation.args = mplusRotationArgs) -``` - -The `efa()` wrapper can fit multiple oblique [EFA](#efa) models in one function call:\index{factor analysis!exploratory}\index{oblique rotation} - -```{r} -efaObliqueLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 1:5, - rotation = "geomin", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = mplusRotationArgs -) -``` - -The `efa()` wrapper can also fit individual oblique EFA models (with `output = "lavaan"`):\index{factor analysis!exploratory}\index{oblique rotation} - -```{r} -efaFactor1ObliqueLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 1, - rotation = "geomin", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = mplusRotationArgs, - output = "lavaan" -) - -efaFactor2ObliqueLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 2, - rotation = "geomin", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = mplusRotationArgs, - output = "lavaan" -) - -efaFactor3ObliqueLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 3, - rotation = "geomin", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = mplusRotationArgs, - output = "lavaan" -) - -efaFactor4ObliqueLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 4, - rotation = "geomin", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = mplusRotationArgs, - output = "lavaan" -) - -efaFactor5ObliqueLavaan_fit <- efa( - data = HolzingerSwineford1939, - ov.names = vars, - nfactors = 5, - rotation = "geomin", - missing = "ML", - estimator = "MLR", - meanstructure = TRUE, - rotation.args = mplusRotationArgs, - output = "lavaan" -) -``` - -#### Factor Loadings {#factorLoadings-efa} - -##### Orthogonal (Varimax) rotation {#orthogonalFactorLoadings-efa} - -###### `psych` {#orthogonalFactorLoadingsPsych-efa} - -The factor loadings and summaries of the model results are below:\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r} -efa1factorOrthogonal -efa2factorOrthogonal -efa3factorOrthogonal -efa4factorOrthogonal -efa5factorOrthogonal -efa6factorOrthogonal -efa7factorOrthogonal -efa8factorOrthogonal -efa9factorOrthogonal -``` - -###### `lavaan` {#orthogonalFactorLoadingsLavaan-efa} - -The factor loadings and summaries of the model results are below:\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r} -summary( - efaOrthogonalLavaan_fit) - -summary( - efa1factorOrthogonalLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa2factorOrthogonalLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa3factorOrthogonalLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa4factorOrthogonalLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa5factorOrthogonalLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa6factorOrthogonalLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa7factorOrthogonalLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa8factorOrthogonalLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa9factorOrthogonalLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) -``` - -##### Oblique (Oblimin) rotation {#obliqueFactorLoadings-efa} - -###### `psych` {#obliqueFactorLoadingsPsych-efa} - -The factor loadings and summaries of the model results are below:\index{factor analysis!exploratory}\index{oblique rotation} - -```{r} -efa1factorOblique -efa2factorOblique -efa3factorOblique -efa4factorOblique -efa5factorOblique -efa6factorOblique -efa7factorOblique -efa8factorOblique -efa9factorOblique -``` - -###### `lavaan` {#obliqueFactorLoadingsLavaan-efa} - -The factor loadings and summaries of the model results are below:\index{factor analysis!exploratory}\index{oblique rotation} - -```{r} -summary( - efaObliqueLavaan_fit) - -summary( - efa1factorObliqueLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa2factorObliqueLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa3factorObliqueLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa4factorObliqueLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa5factorObliqueLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa6factorObliqueLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa7factorObliqueLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa8factorObliqueLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - efa9factorObliqueLavaan_fit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) -``` - -#### Estimates of Model Fit {#fitEstimates-efa} - -##### Orthogonal (Varimax) rotation {#fitEstimatesOrthogonal-efa} - -Fit indices are generated using the syntax below.\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r efaOrthogonal_fitIndices} -fitMeasures( - efaOrthogonalLavaan_fit, - fit.measures = c( - "chisq", "df", "pvalue", - "chisq.scaled", "df.scaled", "pvalue.scaled", - "chisq.scaling.factor", - "baseline.chisq","baseline.df","baseline.pvalue", - "rmsea", "cfi", "tli", "srmr", - "rmsea.robust", "cfi.robust", "tli.robust")) -``` - -##### Oblique (Oblimin) rotation {#fitEstimatesOblique-efa} - -\index{factor analysis!exploratory}\index{oblique rotation} - -```{r} -fitMeasures( - efaObliqueLavaan_fit, - fit.measures = c( - "chisq", "df", "pvalue", - "chisq.scaled", "df.scaled", "pvalue.scaled", - "chisq.scaling.factor", - "baseline.chisq","baseline.df","baseline.pvalue", - "rmsea", "cfi", "tli", "srmr", - "rmsea.robust", "cfi.robust", "tli.robust")) -``` - -#### Factor Scores {#factorScores-efa} - -##### Orthogonal (Varimax) rotation {#factorScoresOrthogonal-efa} - -###### `psych` {#factorScoresOrthogonalPsych-efa} - -\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r} -fa3Orthogonal <- efa3factorOrthogonal$scores -``` - -###### `lavaan` {#factorScoresOrthogonalLavaan-efa} - -\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r} -fa3Orthogonal_lavaan <- lavPredict(efa3factorOrthogonalLavaan_fit) -``` - -##### Oblique (Oblimin) rotation {#factorScoresOblique-efa} - -###### `psych` {#factorScoresObliquePsych-efa} - -\index{factor analysis!exploratory}\index{oblique rotation} - -```{r} -fa3Oblique <- efa3factorOblique$scores -``` - -###### `lavaan` {#factorScoresObliqueLavaan-efa} - -\index{factor analysis!exploratory}\index{oblique rotation} - -```{r} -fa3Oblique_lavaan <- lavPredict(efa3factorObliqueLavaan_fit) -``` - -#### Plots {#factorScoresPlots-efa} - -Biplots were generated using the `psych` package [@R-psych].\index{correlation} -Pairs panel plots were generated using the `psych` package [@R-psych].\index{correlation} -Correlation plots were generated using the `corrplot` package [@R-corrplot].\index{correlation} - -##### Orthogonal (Varimax) rotation {#factorScoresPlotsOrthogonal-efa} - -A scree plot from a model with orthogonal rotation is in Figure \@ref(fig:screePlotPsychOrthogonal).\index{factor analysis!exploratory}\index{orthogonal rotation}\index{scree plot} - -```{r screePlotPsychOrthogonal, out.width = "100%", fig.align = "center", fig.cap = "Scree Plot With Orthogonal Rotation in Exploratory Factor Analysis: psych."} -plot( - efa9factorOrthogonal$e.values, - xlab = "Factor", - ylab = "Eigenvalue") -``` - -A scree plot based on factor loadings from `lavaan` is in Figure \@ref(fig:screePlotLavaanOrthogonal).\index{factor analysis!exploratory}\index{orthogonal rotation}\index{scree plot} -When the factors are uncorrelated (orthogonal rotation), the eigenvalue for a factor is calculated as the sum of squared standardized factor loadings across all items, as described in Section \@ref(pca).\index{factor analysis!exploratory}\index{orthogonal rotation}\index{scree plot}\index{eigenvalue} - -```{r screePlotLavaanOrthogonal, out.width = "100%", fig.align = "center", class.source = "fold-hide", fig.cap = "Scree Plot With Orthogonal Rotation in Exploratory Factor Analysis: lavaan."} -param1Factor <- parameterEstimates( - efa1factorOrthogonalLavaan_fit, - standardized = TRUE) - -param2Factor <- parameterEstimates( - efa2factorOrthogonalLavaan_fit, - standardized = TRUE) - -param3Factor <- parameterEstimates( - efa3factorOrthogonalLavaan_fit, - standardized = TRUE) - -param4Factor <- parameterEstimates( - efa4factorOrthogonalLavaan_fit, - standardized = TRUE) - -param5Factor <- parameterEstimates( - efa5factorOrthogonalLavaan_fit, - standardized = TRUE) - -param6Factor <- parameterEstimates( - efa6factorOrthogonalLavaan_fit, - standardized = TRUE) - -param7Factor <- parameterEstimates( - efa7factorOrthogonalLavaan_fit, - standardized = TRUE) - -param8Factor <- parameterEstimates( - efa8factorOrthogonalLavaan_fit, - standardized = TRUE) - -param9Factor <- parameterEstimates( - efa9factorOrthogonalLavaan_fit, - standardized = TRUE) - -factorNames <- c("f1","f2","f3","f4","f5","f6","f7","f8","f9") - -loadings1Factor <- param1Factor[which( - param1Factor$lhs %in% factorNames & - param1Factor$rhs %in% vars), - c("lhs","std.all")] - -loadings2Factor <- param2Factor[which( - param2Factor$lhs %in% factorNames & - param2Factor$rhs %in% vars), - c("lhs","std.all")] - -loadings3Factor <- param3Factor[which( - param3Factor$lhs %in% factorNames & - param3Factor$rhs %in% vars), - c("lhs","std.all")] - -loadings4Factor <- param4Factor[which( - param4Factor$lhs %in% factorNames & - param4Factor$rhs %in% vars), - c("lhs","std.all")] - -loadings5Factor <- param5Factor[which( - param5Factor$lhs %in% factorNames & - param5Factor$rhs %in% vars), - c("lhs","std.all")] - -loadings6Factor <- param6Factor[which( - param6Factor$lhs %in% factorNames & - param6Factor$rhs %in% vars), - c("lhs","std.all")] - -loadings7Factor <- param7Factor[which( - param7Factor$lhs %in% factorNames & - param7Factor$rhs %in% vars), - c("lhs","std.all")] - -loadings8Factor <- param8Factor[which( - param8Factor$lhs %in% factorNames & - param8Factor$rhs %in% vars), - c("lhs","std.all")] - -loadings9Factor <- param9Factor[which( - param9Factor$lhs %in% factorNames & - param9Factor$rhs %in% vars), - c("lhs","std.all")] - -eigenData <- data.frame( - Factor = 1:9, - Eigenvalue = NA) - -eigenData$Eigenvalue[which( - eigenData$Factor == 1)] <- sum(loadings1Factor$std.all^2) - -eigenData$Eigenvalue[which(eigenData$Factor == 2)] <- - loadings2Factor %>% - group_by(lhs) %>% - summarise( - Eigenvalue = sum(std.all^2), - .groups = "drop") %>% - summarise(min(Eigenvalue)) - -eigenData$Eigenvalue[which(eigenData$Factor == 3)] <- - loadings3Factor %>% - group_by(lhs) %>% - summarise( - Eigenvalue = sum(std.all^2), - .groups = "drop") %>% - summarise(min(Eigenvalue)) - -eigenData$Eigenvalue[which(eigenData$Factor == 4)] <- - loadings4Factor %>% - group_by(lhs) %>% - summarise( - Eigenvalue = sum(std.all^2), - .groups = "drop") %>% - summarise(min(Eigenvalue)) - -eigenData$Eigenvalue[which(eigenData$Factor == 5)] <- - loadings5Factor %>% - group_by(lhs) %>% - summarise( - Eigenvalue = sum(std.all^2), - .groups = "drop") %>% - summarise(min(Eigenvalue)) - -eigenData$Eigenvalue[which(eigenData$Factor == 6)] <- - loadings6Factor %>% - group_by(lhs) %>% - summarise( - Eigenvalue = sum(std.all^2), - .groups = "drop") %>% - summarise(min(Eigenvalue)) - -eigenData$Eigenvalue[which(eigenData$Factor == 7)] <- - loadings7Factor %>% - group_by(lhs) %>% - summarise( - Eigenvalue = sum(std.all^2), - .groups = "drop") %>% - summarise(min(Eigenvalue)) - -eigenData$Eigenvalue[which(eigenData$Factor == 8)] <- - loadings8Factor %>% - group_by(lhs) %>% - summarise( - Eigenvalue = sum(std.all^2), - .groups = "drop") %>% - summarise(min(Eigenvalue)) - -eigenData$Eigenvalue[which(eigenData$Factor == 9)] <- - loadings9Factor %>% - group_by(lhs) %>% - summarise( - Eigenvalue = sum(std.all^2), - .groups = "drop") %>% - summarise(min(Eigenvalue)) - -plot( - eigenData$Factor, - eigenData$Eigenvalue, - xlab = "Factor", - ylab = "Eigevalue") -``` - -A biplot is in Figure \@ref(fig:biplotOrthogonal).\index{factor analysis!exploratory}\index{orthogonal rotation}\index{correlation} - -```{r biplotOrthogonal, out.width = "100%", fig.align = "center", fig.cap = "Biplot Using Orthogonal Rotation in Exploratory Factor Analysis."} -biplot(efa2factorOrthogonal) -abline(h = 0, v = 0, lty = 2) -``` - -A factor plot is in Figure \@ref(fig:factorPlotOrthogonal).\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r factorPlotOrthogonal, out.width = "100%", fig.align = "center", fig.cap = "Factor Plot With Orthogonal Rotation in Exploratory Factor Analysis."} -factor.plot(efa3factorOrthogonal, cut = 0.5) -``` - -A factor diagram is in Figure \@ref(fig:factorDiagramOrthogonal).\index{factor analysis!exploratory}\index{orthogonal rotation} - -```{r factorDiagramOrthogonal, out.width = "100%", fig.align = "center", fig.cap = "Factor Diagram With Orthogonal Rotation in Exploratory Factor Analysis."} -fa.diagram(efa3factorOrthogonal, digits = 2) -``` - -A pairs panel plot is in Figure \@ref(fig:pairsPlotOrthogonal).\index{factor analysis!exploratory}\index{orthogonal rotation}\index{correlation} - -```{r pairsPlotOrthogonal, out.width = "100%", fig.align = "center", fig.cap = "Pairs Panel Plot With Orthogonal Rotation in Exploratory Factor Analysis."} -pairs.panels(fa3Orthogonal) -``` - -A correlation plot is in Figure \@ref(fig:correlationPlotOrthogonal).\index{factor analysis!exploratory}\index{orthogonal rotation}\index{correlation} - -```{r correlationPlotOrthogonal, out.width = "100%", fig.align = "center", fig.cap = "Correlation Plot With Orthogonal Rotation in Exploratory Factor Analysis."} -corrplot(cor( - fa3Orthogonal, - use = "pairwise.complete.obs")) -``` - -##### Oblique (Oblimin) rotation {#factorScoresPlotsOblique-efa} - -A biplot is in Figure \@ref(fig:biplotOblique).\index{factor analysis!exploratory}\index{oblique rotation}\index{correlation} - -```{r biplotOblique, out.width = "100%", fig.align = "center", fig.cap = "Biplot Using Oblique Rotation in Exploratory Factor Analysis."} -biplot(efa2factorOblique) -``` - -A factor plot is in Figure \@ref(fig:factorPlotOblique).\index{factor analysis!exploratory}\index{oblique rotation} - -```{r factorPlotOblique, out.width = "100%", fig.align = "center", fig.cap = "Factor Plot With Oblique Rotation in Exploratory Factor Analysis."} -factor.plot(efa3factorOblique, cut = 0.5) -``` - -A factor diagram is in Figure \@ref(fig:factorDiagramOblique).\index{factor analysis!exploratory}\index{oblique rotation} - -```{r factorDiagramOblique, out.width = "100%", fig.align = "center", fig.cap = "Factor Diagram With Oblique Rotation in Exploratory Factor Analysis."} -fa.diagram(efa3factorOblique, digits = 2) -``` - -A pairs panel plot is in Figure \@ref(fig:pairsPlotOblique).\index{factor analysis!exploratory}\index{oblique rotation}\index{correlation} - -```{r pairsPlotOblique, out.width = "100%", fig.align = "center", fig.cap = "Pairs Panel Plot With Oblique Rotation in Exploratory Factor Analysis."} -pairs.panels(fa3Oblique) -``` - -A correlation plot is in Figure \@ref(fig:correlationPlotOblique).\index{factor analysis!exploratory}\index{oblique rotation}\index{correlation} - -```{r correlationPlotOblique, out.width = "100%", fig.align = "center", fig.cap = "Correlation Plot With Oblique Rotation in Exploratory Factor Analysis."} -corrplot(cor( - fa3Oblique, - use = "pairwise.complete.obs")) -``` - -### Confirmatory Factor Analysis (CFA) {#cfaExamples} - -I introduced [confirmatory factor analysis](#cfa) (CFA) models in Section \@ref(cfa-sem) in the chapter on [structural equation models](#sem).\index{factor analysis!confirmatory} - -The [confirmatory factor analysis](#cfa) (CFA) models were fit in the `lavaan` package [@R-lavaan].\index{factor analysis!confirmatory} -The examples were adapted from the `lavaan` documentation:\index{factor analysis!confirmatory} https://lavaan.ugent.be/tutorial/cfa.html (archived at https://perma.cc/GKY3-9YE4) - -#### Specify the model {#syntax-cfa} - -```{r} -cfaModel_syntax <- ' - #Factor loadings - visual =~ x1 + x2 + x3 - textual =~ x4 + x5 + x6 - speed =~ x7 + x8 + x9 -' - -cfaModel_fullSyntax <- ' - #Factor loadings (free the factor loading of the first indicator) - visual =~ NA*x1 + x2 + x3 - textual =~ NA*x4 + x5 + x6 - speed =~ NA*x7 + x8 + x9 - - #Fix latent means to zero - visual ~ 0 - textual ~ 0 - speed ~ 0 - - #Fix latent variances to one - visual ~~ 1*visual - textual ~~ 1*textual - speed ~~ 1*speed - - #Estimate covariances among latent variables - visual ~~ textual - visual ~~ speed - textual ~~ speed - - #Estimate residual variances of manifest variables - x1 ~~ x1 - x2 ~~ x2 - x3 ~~ x3 - x4 ~~ x4 - x5 ~~ x5 - x6 ~~ x6 - x7 ~~ x7 - x8 ~~ x8 - x9 ~~ x9 - - #Free intercepts of manifest variables - x1 ~ int1*1 - x2 ~ int2*1 - x3 ~ int3*1 - x4 ~ int4*1 - x5 ~ int5*1 - x6 ~ int6*1 - x7 ~ int7*1 - x8 ~ int8*1 - x9 ~ int9*1 -' -``` - -##### Model syntax in table form: {#tabular-cfa} - -```{r} -lavaanify(cfaModel_syntax) -lavaanify(cfaModel_fullSyntax) -``` - -#### Fit the model {#fit-cfa} - -```{r} -cfaModelFit <- cfa( - cfaModel_syntax, - data = HolzingerSwineford1939, - missing = "ML", - estimator = "MLR", - std.lv = TRUE) - -cfaModelFit_full <- lavaan( - cfaModel_fullSyntax, - data = HolzingerSwineford1939, - missing = "ML", - estimator = "MLR") -``` - -#### Display summary output {#output-cfa} - -```{r} -summary( - cfaModelFit, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) - -summary( - cfaModelFit_full, - fit.measures = TRUE, - standardized = TRUE, - rsquare = TRUE) -``` - -#### Estimates of model fit {#fitCriteria-cfa} - -```{r} -fitMeasures( - cfaModelFit, - fit.measures = c( - "chisq", "df", "pvalue", - "chisq.scaled", "df.scaled", "pvalue.scaled", - "chisq.scaling.factor", - "baseline.chisq","baseline.df","baseline.pvalue", - "rmsea", "cfi", "tli", "srmr", - "rmsea.robust", "cfi.robust", "tli.robust")) -``` - -#### Residuals {#residuals-cfa} - -```{r} -residuals(cfaModelFit, type = "cor") -``` - -#### Modification indices {#modIndices-cfa} - -Modification indices are generated using the syntax below.\index{structural equation modeling!modification indices} - -```{r modIndicesCFA} -modificationindices(cfaModelFit, sort. = TRUE) -``` - -#### Factor scores {#factorScores-cfa} - -```{r} -cfaFactorScores <- lavPredict(cfaModelFit) -``` - -##### Compare CFA factor scores to EFA factor scores {#syntax-cfaEFA} - -As would be expected, the factor scores from the [CFA](#cfa) model are highly correlated with the factor scores from the [EFA](#efa) model.\index{factor analysis!confirmatory}\index{factor analysis!exploratory} - -```{r} -cor.test( - x = cfaFactorScores[,"visual"], - y = fa3Orthogonal[,"ML3"]) -cor.test( - x = cfaFactorScores[,"textual"], - y = fa3Orthogonal[,"ML1"]) -cor.test( - x = cfaFactorScores[,"speed"], - y = fa3Orthogonal[,"ML2"]) -``` - -#### Internal Consistency Reliability {#reliability-cfa} - -[Internal consistency reliability](#internalConsistency-reliability) of items composing the latent factors, as quantified by [omega ($\omega$)](#coefficientOmega) and [average variance extracted](#averageVarianceExtracted) (AVE), was estimated using the `semTools` package [@R-semTools].\index{factor analysis!confirmatory}\index{reliability!internal consistency}\index{reliability!internal consistency!omega}\index{reliability!internal consistency!average variance extracted} - -```{r} -compRelSEM(cfaModelFit) -AVE(cfaModelFit) -``` - -#### Path diagram {#pathDiagram-cfa} - -A path diagram of the model generated using the `semPlot` package [@R-semPlot] is in Figure \@ref(fig:cfaPathDiagram).\index{factor analysis!confirmatory} - -```{r cfaPathDiagram, out.width = "100%", fig.align = "center", fig.cap = "Confirmatory Factor Analysis Model Diagram."} -semPaths( - cfaModelFit, - what = "std", - layout = "tree2", - edge.label.cex = 1.3) -```