Package {EFAtools}

Title:	Fast and Flexible Implementations of Exploratory Factor Analysis Tools
Version:	0.8.0
Description:	Provides functions to perform exploratory factor analysis (EFA) procedures and compare their solutions. The goal is to provide state-of-the-art factor retention methods and a high degree of flexibility in the EFA procedures. This way, for example, implementations from R 'psych' and 'SPSS' can be compared. Moreover, functions for Schmid-Leiman transformation and the computation of omegas are provided. To speed up the analyses, some of the iterative procedures, like principal axis factoring (PAF), are implemented in C++.
Depends:	R (≥ 4.1.0)
License:	GPL-3
Encoding:	UTF-8
LazyData:	true
URL:	https://github.com/mdsteiner/EFAtools
BugReports:	https://github.com/mdsteiner/EFAtools/issues
Imports:	psych, stats, ggplot2 (≥ 3.4.0), cli, Rcpp, future.apply, future, checkmate, progressr, rlang, clue
LinkingTo:	Rcpp, RcppArmadillo, roptim
Suggests:	testthat (≥ 3.0.0), GPArotation (≥ 2022.4-1), lavaan, lavaan.mi, MASS, knitr, rmarkdown, microbenchmark, vdiffr, polycor, mnormt, fungible, withr
VignetteBuilder:	knitr
Config/testthat/edition:	3
Config/roxygen2/version:	8.0.0
NeedsCompilation:	yes
Packaged:	2026-07-07 08:31:44 UTC; ma1159938
Author:	Markus Steiner [aut, cre], Silvia Steiner [aut], William Revelle [ctb], Max Auerswald [ctb], Morten Moshagen [ctb], John Ruscio [ctb], Brendan Roche [ctb], Urbano Lorenzo-Seva [ctb], David Navarro-Gonzalez [ctb], Johan Braeken [ctb], Andreas Soteriades [ctb]
Maintainer:	Markus Steiner <markus.d.steiner@gmail.com>
Repository:	CRAN
Date/Publication:	2026-07-07 09:20:02 UTC

EFAtools: Fast and Flexible Implementations of Exploratory Factor Analysis Tools

Description

Provides functions to perform exploratory factor analysis (EFA) procedures and compare their solutions. The goal is to provide state-of-the-art factor retention methods and a high degree of flexibility in the EFA procedures. This way, for example, implementations from R 'psych' and 'SPSS' can be compared. Moreover, functions for Schmid-Leiman transformation and the computation of omegas are provided. To speed up the analyses, some of the iterative procedures, like principal axis factoring (PAF), are implemented in C++.

Author(s)

Maintainer: Markus Steiner markus.d.steiner@gmail.com (ORCID)

Authors:

• Markus Steiner markus.d.steiner@gmail.com (ORCID)

• Silvia Steiner silvia.steiner.grieder@gmail.com (ORCID)

Other contributors:

• William Revelle [contributor]

• Max Auerswald [contributor]

• Morten Moshagen [contributor]

• John Ruscio [contributor]

• Brendan Roche [contributor]

• Urbano Lorenzo-Seva [contributor]

• David Navarro-Gonzalez [contributor]

• Johan Braeken [contributor]

• Andreas Soteriades [contributor]

See Also

Useful links:

• https://github.com/mdsteiner/EFAtools

• Report bugs at https://github.com/mdsteiner/EFAtools/issues

Average a list of matrices elementwise

Description

Average a list of matrices elementwise

Usage

.average_matrices(x)

Arguments

Value

A matrix with the same dimensions as the inputs.

Covert a "LOADINGS" table to matrix or a matrix to "LOADINGS"

Description

Covert a "LOADINGS" table to matrix or a matrix to "LOADINGS"

Usage

.change_class(x, cl = "matrix")

Arguments

Value

A table with the loadings, of class either "LOADINGS" or "matrix".

Author(s)

Andreas Soteriades

The loadings tables returned by EFA() are of class "LOADINGS", which prevents applying functions on them. This function allows to change their class to "matrix", and to change back to "LOADINGS" when done.

Compute explained variances from loadings

Description

From unrotated loadings compute the communalities and uniquenesses for total variance. Compute explained variances per factor from rotated loadings (and factor intercorrelations Phi if oblique rotation was used).

Usage

.compute_vars(L_unrot, L_rot, Phi = NULL)

Arguments

Value

A matrix with sum of squared loadings, proportion explained variance from total variance per factor, same as previous but cumulative, Proportion of explained variance from total explained variance, and same as previous but cumulative.

Mean squared discrepancy to a consensus target

Description

Mean squared discrepancy to a consensus target

Usage

.consensus_loss(aligned_loadings, target)

Arguments

Value

Mean sum of squared deviations from the target across matrices.

Internal single-start GPA-consensus engine

Description

Performs a single GPA-consensus run from one starting target. The multi-start wrapper .gpa_consensus_target() dispatches here.

Usage

.consensus_target_procrustes_single(
  unrotated_list,
  init_targets = NULL,
  rotation = c("orthogonal", "oblique"),
  start = 1,
  tol = 0.001,
  loss_tol = 1e-07,
  loss_patience = 5,
  convergence = c("either", "target", "loss", "both"),
  min_iter = 2,
  max_iter = 200,
  alpha = 1,
  match_target = TRUE,
  hyper_cutoff = 0.15,
  verbose = FALSE
)

Arguments

Extract a list object by its name

Description

Extract a list object by its name

Usage

.extract_list_object(alist, object)

Arguments

Value

A list of length m, with each element containing the extracted object for the kth element (k = 1,..., m).

Author(s)

Andreas Soteriades

Consider a list of named sub-lists. This function extracts, for each sub-list, the sub-list element that is specified by the user. This function is useful for extracting results from EFA() for each permutation run in EFA_POOLED().

Compute number of non-matching indicator-to-factor correspondences

Description

Compute number of non-matching indicator-to-factor correspondences

Usage

.factor_corres(x, y, thresh = 0.3)

Arguments

Generalized Procrustes Analysis consensus target across loading matrices

Description

Internal helper that constructs a Generalized Procrustes Analysis (GPA) consensus target across a list of loading matrices and returns the aligned loadings, the centroid target, and convergence diagnostics. Used by EFA_POOLED() under target_method = "consensus" to build a common rotation target across imputations. Oblique rotations are not supported here: the iteration is degenerate for oblique transforms with more than one factor (cf. Lorenzo-Seva & Van Ginkel 2016, who use a Promin step on top of the centroid rather than iterated oblique Procrustes); callers should pass the unrotated solutions of an orthogonal rotation, or use target_method = "first_target".

Usage

.gpa_consensus_target(
  unrotated_list,
  init_targets = NULL,
  rotation = c("orthogonal", "oblique"),
  start = 1,
  multi_start = FALSE,
  starts = NULL,
  tol = 0.001,
  loss_tol = 1e-06,
  loss_patience = 5,
  convergence = c("either", "target", "loss", "both"),
  min_iter = 2,
  max_iter = 200,
  alpha = 1,
  match_target = TRUE,
  hyper_cutoff = 0.15,
  verbose = FALSE
)

Arguments

Details

The iteration alternates two steps:

1. each loading matrix is aligned to the current target with PROCRUSTES();

2. the target is updated to the elementwise centroid of the aligned matrices.

The outer loop stops when the target stabilises, when the consensus loss stabilises, or when both criteria are satisfied.

If multi_start = FALSE, one consensus run is performed. If multi_start = TRUE, the same engine is repeated for the selected starting targets and the run with the smallest final mean loss is returned as the main result; all runs and a between-run congruence summary are retained in the multi_start component.

Value

A list with the converged target, aligned matrices, pooled loadings, pooled Phi, convergence history, inner-alignment diagnostics, and hyperplane summaries. If multi_start = TRUE, the multi_start element also contains the per-start losses, convergence indicators, run summaries, all run objects, and between-run Tucker congruence matrices.

References

Gower, J. C. (1975). Generalized Procrustes analysis. Psychometrika, 40, 33-51.

Van Ginkel, J. R., & Kroonenberg, P. M. (2014). Using Generalized Procrustes Analysis for Multiple Imputation in Principal Component Analysis. Journal of Classification, 31, 242-269.

Lorenzo-Seva, U., & Van Ginkel, J. R. (2016). Multiple Imputation of missing values in exploratory factor analysis of multidimensional scales: estimating latent trait scores. Anales de Psicologia, 32, 596-608.

Count near-zero loadings

Description

Hyperplane count is the number of loadings with absolute value smaller than a user-specified cutoff.

Usage

.hyperplane_count(L, cutoff = 0.15)

Arguments

Value

A list with the total hyperplane count and counts by factor and item.

Oblique Procrustes target rotation using a k x k inner objective

Description

Compute an oblique target rotation for a loading matrix using a targetQ-compatible parameterization and a ⁠k x k⁠ objective.

Usage

.oblique_procrustes(
  A,
  B,
  S_r = NULL,
  T_init_r = NULL,
  eps = 1e-05,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  normalize = FALSE,
  random_starts = 0L,
  screen_keep = 2L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

The rotated loading matrix is defined as L = A %*% solve(t(T)), and the corresponding factor correlation matrix is Phi = t(T) %*% T. The optimization is carried out over the transformation matrix T under the oblique normalization constraint diag(t(T) %*% T) = 1.

Non-invertible candidate transformations are rejected rather than evaluated through a pseudo-inverse.

Additional random starts may be requested. To reduce runtime, the solver uses a two-stage strategy for extra starts: cheap objective screening, followed by short triage optimization, followed by full optimization only for starts that improve on the current incumbent by at least triage_improve_tol.

The routine is intended for repeated oblique target rotations in workflows such as bootstrap alignment or consensus alignment of exploratory factor solutions across multiply imputed datasets. It follows the same oblique transformation convention as GPArotation::targetQ().

Value

A named list containing the rotated loadings, transformation matrix, factor correlation matrix, target criterion value, convergence diagnostics, line-search diagnostics, and multi-start summaries.

References

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Browne, M. W. (2001). An overview of analytic rotation in exploratory factor analysis. Multivariate Behavioral Research, 36, 111-150.

Batched oblique Procrustes target rotation over a cube of loading matrices

Description

Align each slice of a loading-matrix cube to a single shared target using the same oblique target rotation as .oblique_procrustes(), in one call. This removes the per-replicate marshalling overhead of looping PROCRUSTES() in R over bootstrap or multiple-imputation arrays.

Usage

.oblique_procrustes_batch(
  A,
  B,
  eps = 1e-05,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  normalize = FALSE,
  random_starts = 0L,
  screen_keep = 2L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

Each slice A[, , i] is aligned to B. For a single-factor cube the alignment reduces to the closed-form sign match T = sign(crossprod(A_i, B)) with factor correlation 1, matching the one-factor short-circuit in PROCRUSTES(). For two or more factors the slice is warm-started from the closed-form orthogonal Procrustes solution (mirroring PROCRUSTES()) and optimized with the same multi-start oblique solver as .oblique_procrustes(). Random starts are drawn serially with R::rnorm in the calling process.

Slices are aligned independently. A slice that cannot be aligned (a non-finite loading matrix, a failed warm-start decomposition, an invalid fit, or any linear-algebra exception) is reported with valid = FALSE and NA for the loadings, factor correlations, and all other per-slice diagnostics, rather than aborting the whole call, so one degenerate replicate does not discard the rest.

Value

A named list with the aligned-loadings array loadings (⁠n x m x b⁠), the factor-correlation array Phi (⁠m x m x b⁠), and the per-slice diagnostics valid, convergence, value, iterations, and line_search_failed.

Closed-form orthogonal Procrustes rotation

Description

Rotate A to the orthogonal target B by minimizing ⁠||A %*% T - B||_F^2⁠ subject to t(T) %*% T = I.

Usage

.orthogonal_procrustes(A, B)

Arguments

Value

A list with the rotated loadings, orthogonal transformation matrix, target criterion value, and basic diagnostics.

References

Schoenemann, P. H. (1966). A generalized solution of the orthogonal Procrustes problem. Psychometrika, 31, 1-10.

Perform the iterative PAF procedure

Description

Function called from within PAF so usually no call to this is needed by the user. Provides a C++ implementation of the PAF procedure

Usage

.paf_iter(h2, criterion, R, n_fac, abs_eig, crit_type, max_iter)

Arguments

Parallel analysis on simulated data.

Description

Function called from within PARALLEL so usually no call to this is needed by the user. Provides a C++ implementation of the PARALLEL simulation procedure

Usage

.parallel_sim(n_datasets, n_vars, N, eigen_type, maxit = 10000L)

Arguments

Oblique Bentler factor rotation

Description

Rotate a loading matrix obliquely under Bentler's invariant pattern simplicity criterion using a gradient-projection optimizer along the oblique (column-normalized) manifold.

Usage

.rotate_bentler_oblq(
  L,
  eps = 1e-05,
  normalize = TRUE,
  random_starts = 0L,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  screen_keep = 5L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

The criterion value f and its gradient dQ/dL at the rotated loadings L = A %*% solve(t(T)) define the search; the engine maps the gradient to the transformation T on the manifold diag(t(T) %*% T) = 1, projects it onto the tangent space, performs a non-monotone line search, and retracts back onto the manifold by column normalization. The Bentler criterion measures the departure of the cross-products of squared loadings from a diagonal pattern; it is prone to local minima, so additional random starts are recommended.

Additional random starts may be requested. To bound runtime the solver screens each random start by its objective, runs a short triage optimization on the best-screened starts, and fully optimizes only those that improve on the current incumbent by at least triage_improve_tol.

Value

A named list with the rotated loadings, the transformation matrix Th (with L %*% t(solve(Th)) reproducing the rotated loadings), the factor correlation matrix Phi (t(Th) %*% Th), the attained criterion value, and the convergence and validity flags. The list additionally reports the criterion value reached at each optimized start in all_values, with a per-start convergence flag in all_converged.

References

Bentler, P. M. (1977). Factor simplicity index and transformations. Psychometrika, 42, 277-295.

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Orthogonal Bentler factor rotation

Description

Rotate a loading matrix orthogonally under Bentler's invariant pattern simplicity criterion using a gradient-projection optimizer along the orthogonal (Stiefel) manifold.

Usage

.rotate_bentler_orth(
  L,
  eps = 1e-05,
  normalize = TRUE,
  random_starts = 0L,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  screen_keep = 5L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

The criterion value f and its gradient dQ/dL at the rotated loadings L = A %*% T define the search; the engine maps the gradient to the orthogonal transformation T, projects it onto the tangent space, performs a non-monotone line search, and retracts back onto the orthogonal group via a polar (singular value) projection. The Bentler criterion measures the departure of the cross-products of squared loadings from a diagonal pattern; it is prone to local minima, so additional random starts are recommended.

Additional random orthogonal starts may be requested. To bound runtime the solver screens each random start by its objective, runs a short triage optimization on the best-screened starts, and fully optimizes only those that improve on the current incumbent by at least triage_improve_tol.

Value

A named list with the rotated loadings, the orthogonal rotation matrix Th (with L %*% Th reproducing the rotated loadings), the attained criterion value, and the convergence and validity flags. The list additionally reports the criterion value reached at each optimized start in all_values, with a per-start convergence flag in all_converged.

References

Bentler, P. M. (1977). Factor simplicity index and transformations. Psychometrika, 42, 277-295.

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Oblique bifactor factor rotation

Description

Rotate a loading matrix obliquely under the Jennrich-Bentler bifactor criterion using a gradient-projection optimizer along the oblique (column-normalized) manifold.

Usage

.rotate_bifactor_oblq(
  L,
  eps = 1e-05,
  normalize = TRUE,
  random_starts = 0L,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  screen_keep = 5L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

The criterion value f and its gradient dQ/dL at the rotated loadings L = A %*% solve(t(T)) define the search; the engine maps the gradient to the transformation T on the manifold diag(t(T) %*% T) = 1, projects it onto the tangent space, performs a non-monotone line search, and retracts back onto the manifold by column normalization. The first factor is treated as a general factor and is exempt from the penalty; the criterion measures the between-group-factor cross-products of the squared loadings, so it is minimized when each variable loads on the general factor plus at most one group factor. The criterion is prone to local minima, so additional random starts are recommended.

Additional random starts may be requested. To bound runtime the solver screens each random start by its objective, runs a short triage optimization on the best-screened starts, and fully optimizes only those that improve on the current incumbent by at least triage_improve_tol.

Value

A named list with the rotated loadings, the transformation matrix Th (with L %*% t(solve(Th)) reproducing the rotated loadings), the factor correlation matrix Phi (t(Th) %*% Th), the attained criterion value, and the convergence and validity flags. The list additionally reports the criterion value reached at each optimized start in all_values, with a per-start convergence flag in all_converged.

References

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Jennrich, R. I., & Bentler, P. M. (2011). Exploratory bi-factor analysis. Psychometrika, 76, 537-549.

Orthogonal bifactor factor rotation

Description

Rotate a loading matrix orthogonally under the Jennrich-Bentler bifactor criterion using a gradient-projection optimizer along the orthogonal (Stiefel) manifold.

Usage

.rotate_bifactor_orth(
  L,
  eps = 1e-05,
  normalize = TRUE,
  random_starts = 0L,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  screen_keep = 5L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

The criterion value f and its gradient dQ/dL at the rotated loadings L = A %*% T define the search; the engine maps the gradient to the orthogonal transformation T, projects it onto the tangent space, performs a non-monotone line search, and retracts back onto the orthogonal group via a polar (singular value) projection. The first factor is treated as a general factor and is exempt from the penalty; the criterion measures the between-group-factor cross-products of the squared loadings, so it is minimized when each variable loads on the general factor plus at most one group factor. The criterion is prone to local minima, so additional random starts are recommended.

Additional random orthogonal starts may be requested. To bound runtime the solver screens each random start by its objective, runs a short triage optimization on the best-screened starts, and fully optimizes only those that improve on the current incumbent by at least triage_improve_tol.

Value

A named list with the rotated loadings, the orthogonal rotation matrix Th (with L %*% Th reproducing the rotated loadings), the attained criterion value, and the convergence and validity flags. The list additionally reports the criterion value reached at each optimized start in all_values, with a per-start convergence flag in all_converged.

References

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Jennrich, R. I., & Bentler, P. M. (2011). Exploratory bi-factor analysis. Psychometrika, 76, 537-549.

Orthogonal Crawford-Ferguson factor rotation

Description

Rotate a loading matrix orthogonally under the Crawford-Ferguson criterion using a gradient-projection optimizer along the orthogonal (Stiefel) manifold.

Usage

.rotate_cf_orth(
  L,
  kappa,
  eps = 1e-05,
  normalize = TRUE,
  random_starts = 0L,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  screen_keep = 5L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

The criterion value f and its gradient dQ/dL at the rotated loadings L = A %*% T define the search; the engine maps the gradient to the orthogonal transformation T, projects it onto the tangent space, performs a non-monotone line search, and retracts back onto the orthogonal group via a polar (singular value) projection. kappa = 0 is the quartimax criterion and kappa = ncol(A) / (2 * nrow(A)) is the equamax criterion.

Additional random orthogonal starts may be requested. To bound runtime the solver screens each random start by its objective, runs a short triage optimization on the best-screened starts, and fully optimizes only those that improve on the current incumbent by at least triage_improve_tol.

Value

A named list with the rotated loadings, the orthogonal rotation matrix Th (with L %*% Th reproducing the rotated loadings), the attained criterion value, and the convergence and validity flags. The list additionally reports the criterion value reached at each optimized start in all_values, with a per-start convergence flag in all_converged.

References

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Crawford, C. B., & Ferguson, G. A. (1970). A general rotation criterion and its use in orthogonal rotation. Psychometrika, 35, 321-332.

Oblique geomin factor rotation

Description

Rotate a loading matrix obliquely under the geomin criterion using a gradient-projection optimizer along the oblique (column-normalized) manifold.

Usage

.rotate_geomin_oblq(
  L,
  delta = 0.01,
  eps = 1e-05,
  normalize = TRUE,
  random_starts = 0L,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  screen_keep = 5L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

The criterion value f and its gradient dQ/dL at the rotated loadings L = A %*% solve(t(T)) define the search; the engine maps the gradient to the transformation T on the manifold diag(t(T) %*% T) = 1, projects it onto the tangent space, performs a non-monotone line search, and retracts back onto the manifold by column normalization. The geomin criterion sums the per-variable geometric mean of the squared loadings offset by delta; it is prone to local minima, so additional random starts are recommended.

Additional random starts may be requested. To bound runtime the solver screens each random start by its objective, runs a short triage optimization on the best-screened starts, and fully optimizes only those that improve on the current incumbent by at least triage_improve_tol.

Value

A named list with the rotated loadings, the transformation matrix Th (with L %*% t(solve(Th)) reproducing the rotated loadings), the factor correlation matrix Phi (t(Th) %*% Th), the attained criterion value, and the convergence and validity flags. The list additionally reports the criterion value reached at each optimized start in all_values, with a per-start convergence flag in all_converged.

References

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Browne, M. W. (2001). An overview of analytic rotation in exploratory factor analysis. Multivariate Behavioral Research, 36, 111-150.

Orthogonal geomin factor rotation

Description

Rotate a loading matrix orthogonally under the geomin criterion using a gradient-projection optimizer along the orthogonal (Stiefel) manifold.

Usage

.rotate_geomin_orth(
  L,
  delta = 0.01,
  eps = 1e-05,
  normalize = TRUE,
  random_starts = 0L,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  screen_keep = 5L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

The criterion value f and its gradient dQ/dL at the rotated loadings L = A %*% T define the search; the engine maps the gradient to the orthogonal transformation T, projects it onto the tangent space, performs a non-monotone line search, and retracts back onto the orthogonal group via a polar (singular value) projection. The geomin criterion sums the per-variable geometric mean of the squared loadings offset by delta; it is prone to local minima, so additional random starts are recommended.

Additional random orthogonal starts may be requested. To bound runtime the solver screens each random start by its objective, runs a short triage optimization on the best-screened starts, and fully optimizes only those that improve on the current incumbent by at least triage_improve_tol.

Value

A named list with the rotated loadings, the orthogonal rotation matrix Th (with L %*% Th reproducing the rotated loadings), the attained criterion value, and the convergence and validity flags. The list additionally reports the criterion value reached at each optimized start in all_values, with a per-start convergence flag in all_converged.

References

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Browne, M. W. (2001). An overview of analytic rotation in exploratory factor analysis. Multivariate Behavioral Research, 36, 111-150.

Oblique oblimin factor rotation

Description

Rotate a loading matrix obliquely under the oblimin criterion using a gradient-projection optimizer along the oblique (column-normalized) manifold.

Usage

.rotate_oblimin(
  L,
  gam = 0,
  eps = 1e-05,
  normalize = TRUE,
  random_starts = 0L,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1,
  screen_keep = 5L,
  triage_maxit = 25L,
  triage_improve_tol = 0
)

Arguments

Details

The criterion value f and its gradient dQ/dL at the rotated loadings L = A %*% solve(t(T)) define the search; the engine maps the gradient to the transformation T on the manifold diag(t(T) %*% T) = 1, projects it onto the tangent space, performs a non-monotone line search, and retracts back onto the manifold by column normalization. gam = 0 is the quartimin criterion.

Additional random starts may be requested. To bound runtime the solver screens each random start by its objective, runs a short triage optimization on the best-screened starts, and fully optimizes only those that improve on the current incumbent by at least triage_improve_tol.

Value

A named list with the rotated loadings, the transformation matrix Th (with L %*% t(solve(Th)) reproducing the rotated loadings), the factor correlation matrix Phi (t(Th) %*% Th), the attained criterion value, and the convergence and validity flags. The list additionally reports the criterion value reached at each optimized start in all_values, with a per-start convergence flag in all_converged.

References

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Jennrich, R. I., & Sampson, P. F. (1966). Rotation for simple loadings. Psychometrika, 31, 313-323.

Oblique simplimax factor rotation

Description

Rotate a loading matrix obliquely under the simplimax criterion using a gradient-projection optimizer along the oblique (column-normalized) manifold.

Usage

.rotate_simplimax_oblq(
  L,
  k,
  eps = 1e-05,
  normalize = TRUE,
  random_starts = 0L,
  maxit = 1000L,
  max_line_search = 10L,
  step0 = 1
)

Arguments

Details

The criterion value f and its gradient dQ/dL at the rotated loadings L = A %*% solve(t(T)) define the search; the engine maps the gradient to the transformation T on the manifold diag(t(T) %*% T) = 1, projects it onto the tangent space, performs a non-monotone line search, and retracts back onto the manifold by column normalization. The simplimax criterion sums the k smallest squared loadings, so it is minimized when the k "close-to-zero" loadings are driven toward zero; the count k is a tuning parameter. Because the set of k smallest loadings is reselected at every evaluation, the criterion is only piecewise smooth: its gradient jumps as loadings cross the kth-smallest threshold, so the line search accepts a step whenever it decreases the largest objective over a short window of recent iterations (a non-monotone test; Grippo, Lampariello, & Lucidi, 1986), letting the optimizer step across the kinks where a strictly monotone descent would stall.

The criterion is strongly prone to local minima, so the solver fully optimizes the identity start together with random_starts random orthogonal starts and keeps the solution with the lowest criterion value. Fully optimizing every start – rather than the screen-and-triage strategy used for the smooth criteria, which assumes the rational start lies in the global basin – is the standard remedy for the local minima of complexity-based rotation criteria (Kiers, 1994; Browne, 2001).

Value

A named list with the rotated loadings, the transformation matrix Th (with L %*% t(solve(Th)) reproducing the rotated loadings), the factor correlation matrix Phi (t(Th) %*% Th), the attained criterion value, and the convergence and validity flags. The list additionally reports the criterion value reached at each optimized start in all_values, with a per-start convergence flag in all_converged.

References

Bernaards, C. A., & Jennrich, R. I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676-696.

Browne, M. W. (2001). An overview of analytic rotation in exploratory factor analysis. Multivariate Behavioral Research, 36, 111-150.

Grippo, L., Lampariello, F., & Lucidi, S. (1986). A nonmonotone line search technique for Newton's method. SIAM Journal on Numerical Analysis, 23, 707-716.

Kiers, H. A. L. (1994). Simplimax: Oblique rotation to an optimal target with simple structure. Psychometrika, 59, 567-579.

Rotation Jacobians for analytic rotation standard errors

Description

Forward-difference the warm-started re-rotation map ⁠A -> (rotated loadings, Phi)⁠ over the unrotated loadings A to obtain the rotation Jacobians used by the analytic standard errors for rotated loadings (se = "information" in EFA()). The full nrow(A) * ncol(A) finite- difference loop runs in compiled code, re-solving the rotation from the converged transformation T_init at each perturbation; the caller forms ⁠J V J'⁠ in R.

Usage

.rotation_se_jacobian(
  A,
  T_init,
  method,
  param,
  normalize,
  oblique,
  eps,
  general_col = 0L
)

Arguments

Value

A named list with the Jacobian J_L (⁠pk x pk⁠), the re-rotated base_loadings, a validity flag, and – when oblique – the Jacobian J_Phi (⁠k^2 x pk⁠) and base_Phi.

References

Jennrich, R. I. (1973). Standard errors for obliquely rotated factor loadings. Psychometrika, 38, 593-604.

Reference eigenvalues for the NEST simulation via the shared kernel.

Description

Internal helper called from NEST(). Simulates nreps datasets from an (nf - 1)-factor reference model, given that model's loadings Lambda and uniquenesses Psi, and returns the nf-th largest eigenvalue of each simulated correlation matrix. The data are drawn with the shared Z * M rule (see .simulate_cfm_mvn()) using the factor-score square root ⁠M = t([Lambda | diag(sqrt(Psi))])⁠, so a row randn(1, nf - 1 + p) * M is N(0, Lambda Lambda' + diag(Psi)). Drawing nf - 1 + p standard normals and post-multiplying by the factor-score matrix is faster than forming the model- implied correlation matrix and drawing from it, and matches the position at which NEST() reads the reference eigenvalue.

Usage

.simulate_cfm_eigen(nf, N, Lambda, Psi, nreps = 1000L)

Arguments

Draw multivariate-normal data from a population correlation matrix.

Description

Internal helper called from efa_simulate(). Draws N cases from a p-variate normal with correlation (or covariance) R by post-multiplying a matrix of standard normal deviates by a matrix square root M of R (with M' M = R, so the rows of Z * M are N(0, R)). This is the same Z * M rule used by the NEST reference simulation (.simulate_cfm_eigen): there M is the transposed factor-score matrix, here it is a Cholesky or eigen square root. A positive-definite R is factored by Cholesky; a positive-semidefinite but singular R (which makes the Cholesky fail although it is still a valid covariance, e.g. a no-factor block or a smoothed factor intercorrelation matrix) falls back to a symmetric eigen square root.

Usage

.simulate_cfm_mvn(R, N, tol = 1e-08)

Arguments

Tucker congruence between factors

Description

Compute the Tucker congruence matrix between the columns of two loading matrices.

Usage

.tucker_congruence(L1, L2)

Arguments

Value

A square matrix whose ⁠(i, j)⁠ entry is the Tucker congruence between column i of L1 and column j of L2.

References

Lorenzo-Seva, U., and ten Berge, J. M. F. (2006). Tucker's congruence coefficient as a meaningful index of factor similarity. Methodology, 2, 57-64.

Bartlett's test of sphericity

Description

This function tests whether a correlation matrix is significantly different from an identity matrix (Bartlett, 1951). If the Bartlett's test is not significant, the correlation matrix is not suitable for factor analysis because the variables show too little covariance.

Usage

BARTLETT(
  x,
  N = NA,
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra")
)

Arguments

Details

Bartlett (1951) proposed this statistic to determine a correlation matrix' suitability for factor analysis. The statistic is approximately chi square distributed with df = \frac{p(p - 1)}{2} and is given by

chi^2 = -log(det(R)) (N - 1 - (2 * p + 5)/6)

where det(R) is the determinant of the correlation matrix, N is the sample size, and p is the number of variables.

This tests requires multivariate normality. If this condition is not met, the Kaiser-Meyer-Olkin criterion (KMO()) can still be used.

This function was heavily influenced by the psych::cortest.bartlett() function from the psych package.

The BARTLETT function can also be called together with the (KMO()) function and with factor retention criteria in the N_FACTORS() function.

Value

A list containing

Source

Bartlett, M. S. (1951). The effect of standardization on a Chi-square approximation in factor analysis. Biometrika, 38, 337-344.

See Also

KMO() for another measure to determine suitability for factor analysis.

N_FACTORS() as a wrapper function for this function, KMO() and several factor retention criteria.

Examples

BARTLETT(test_models$baseline$cormat, N = 500)

Comparison Data

Description

Factor retention method introduced by Ruscio and Roche (2012). The code was adapted from the CD code by Auerswald and Moshagen (2019) available at https://osf.io/x5cz2/?view_only=d03efba1fd0f4c849a87db82e6705668

Usage

CD(
  x,
  n_factors_max = NA,
  N_pop = 10000,
  N_samples = 500,
  alpha = 0.3,
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra"),
  max_iter = 50
)

Arguments

Details

Comparison data (CD) extends parallel analysis by reproducing the observed correlation matrix rather than generating random data: datasets with a known factor structure are generated with an increasing number of factors, and the smallest number for which adding a further factor no longer significantly improves the reproduction of the observed eigenvalues is retained (Ruscio & Roche, 2012).

Note that if the data contains missing values, these will be removed for the comparison data procedure using stats::na.omit(). If missing data should be treated differently, e.g., by imputation, do this outside CD and then pass the complete data.

Value

An object of class efa_retention (see print.efa_retention() and plot.efa_retention() for the print and plot methods). Its main fields are:

Source

Auerswald, M., & Moshagen, M. (2019). How to determine the number of factors to retain in exploratory factor analysis: A comparison of extraction methods under realistic conditions. Psychological Methods, 24(4), 468–491. https://doi.org/10.1037/met0000200

Ruscio, J., & Roche, B. (2012). Determining the number of factors to retain in an exploratory factor analysis using comparison data of known factorial structure. Psychological Assessment, 24, 282–292. doi: 10.1037/a0025697

See Also

N_FACTORS() as a wrapper function for this and the other factor retention criteria.

Other factor retention criteria: EKC(), HULL(), KGC(), MAP(), NEST(), PARALLEL(), SCREE(), SMT()

Examples

# determine n factors of the GRiPS
CD(GRiPS_raw)

# determine n factors of the DOSPERT risk subscale
CD(DOSPERT_raw)

Compare two vectors or matrices (communalities or loadings)

Description

The function takes two objects of the same dimensions containing numeric information (loadings or communalities) and returns a list of class COMPARE containing summary information of the differences of the objects.

Usage

COMPARE(
  x,
  y,
  reorder = c("congruence", "names", "none"),
  corres = TRUE,
  thresh = 0.3,
  digits = 4,
  m_red = 0.001,
  range_red = 0.001,
  round_red = 3,
  print_diff = TRUE,
  na.rm = FALSE,
  x_labels = c("x", "y"),
  plot = TRUE,
  plot_red = 0.01
)

Arguments

Value

A list of class COMPARE containing summary statistics on the differences of x and y.

Examples

# A type SPSS EFA to mimick the SPSS implementation
EFA_SPSS_6 <- EFA(test_models$case_11b$cormat, n_factors = 6, type = "SPSS")

# A type psych EFA to mimick the psych::fa() implementation
EFA_psych_6 <- EFA(test_models$case_11b$cormat, n_factors = 6, type = "psych")

# compare the two
COMPARE(EFA_SPSS_6$unrot_loadings, EFA_psych_6$unrot_loadings,
        x_labels = c("SPSS", "psych"))

DOSPERT

Description

A list containing the the bivariate correlations (cormat) of the 40 items of the Domain Specific Risk Taking Scale (DOSPERT; Weber, Blais, & Betz, 2002) and the sample size (N) based on the publicly available dataset at (https://osf.io/rce7g) of the Basel-Berlin Risk Study (Frey et al., 2017). The items measure risk-taking propensity on six different domains: social, recreational, gambling, health/ safety, investment, and ethical.

Usage

DOSPERT

Format

An object of class list of length 2.

Source

Weber, E. U., Blais, A.-R., & Betz, N. E. (2002). A domain specific risk-attitude scale: Measuring risk perceptions and risk behaviors. Journal of Behavioral Decision Making, 15(4), 263–290. doi: 10.1002/bdm.414

Frey, R., Pedroni, A., Mata, R., Rieskamp, J., & Hertwig, R. (2017). Risk preference shares the psychometric structure of major psychological traits. Science Advances, 3, e1701381.

https://osf.io/rce7g

DOSPERT_raw

Description

A data.frame containing responses to the risk subscale of the Domain Specific Risk Taking Scale (DOSPERT; Weber, Blais, & Betz, 2002) based on the publicly available dataset (at https://osf.io/pjt57/) by Frey, Duncan, and Weber (2020). The items measure risk-taking propensity on six different domains: social, recreational, gambling, health/ safety, investment, and ethical.

Usage

DOSPERT_raw

Format

An object of class data.frame with 3123 rows and 30 columns.

Source

Blais, A.-R., & Weber, E. U. (2002). A domain-specific risk-taking (DOSPERT) scale for adult populations. Judgment and Decision Making, 15(4), 263–290. doi: 10.1002/bdm.414

Frey, R., Duncan, S. M., & Weber, E. U. (2020). Towards a typology of risk preference: Four risk profiles describe two thirds of individuals in a large sample of the U.S. population. PsyArXiv Preprint. doi:10.31234/osf.io/yjwr9

Exploratory factor analysis (EFA)

Description

This function does an EFA with either PAF, ML, ULS/MINRES, or DWLS with or without subsequent rotation. All arguments with default value NA can be left to default if type is set to one of "EFAtools", "SPSS", or "psych". The respective specifications are then handled according to the specified type (see details).

Usage

EFA(
  x,
  n_factors,
  N = NA,
  method = c("PAF", "ML", "ULS", "MINRES", "DWLS"),
  rotation = c("none", "varimax", "equamax", "quartimax", "geominT", "bentlerT",
    "bifactorT", "promax", "oblimin", "quartimin", "simplimax", "bentlerQ", "geominQ",
    "bifactorQ"),
  se = c("none", "information", "sandwich", "np-boot"),
  type = c("EFAtools", "psych", "SPSS", "none"),
  max_iter = NA,
  init_comm = NA,
  criterion = NA,
  criterion_type = NA,
  abs_eigen = NA,
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  varimax_type = NA,
  k = NA,
  normalize = TRUE,
  P_type = NA,
  precision = 1e-05,
  order_type = NA,
  start_method = "psych",
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra", "fiml"),
  b_boot = 1000,
  ci = 0.95,
  randomStarts = 100,
  seed = NULL,
  ...
)

Arguments

Details

There are two main ways to use this function. The easiest way is to use it with a specified type (see above), which sets most of the other arguments accordingly. Another way is to use it more flexibly by explicitly specifying all arguments used and set type to "none" (see examples). A mix of the two can also be done by specifying a type as well as additional arguments. However, this will throw warnings to avoid unintentional deviations from the implementations according to the specified type.

Estimators

The estimator is chosen with method.

• PAF (principal axis factoring) iteratively estimates the communalities and makes no distributional assumptions, which makes it robust and a good general-purpose default. Because it minimises no likelihood or weighted discrepancy it provides no model chi-square, and hence no chi-square-based fit indices (see Fit indices). The PAF iteration is governed by init_comm, criterion, criterion_type, max_iter, and abs_eigen (set by type; see Using the type presets).

• ML (maximum likelihood) maximises the normal-theory likelihood. It yields the full set of fit indices and is the only estimator with analytic expected-information standard errors (se = "information"), but it assumes multivariate normality and is the most prone to Heywood (improper) cases. Its starting values are set by start_method.

• ULS (unweighted least squares) minimises the sum of squared correlation residuals. "MINRES" (minimum residual) is the same estimator under a different name and returns identical results. It makes no normality assumption, is robust to mild non-normality, and yields the full set of fit indices.

• DWLS (diagonally weighted least squares) is the recommended estimator for ordinal data. It weights each off-diagonal correlation residual by the inverse asymptotic variance of the corresponding polychoric correlation (Muthén, 1984), reproducing the loadings of a diagonally weighted least squares fit (e.g. lavaan::efa(..., estimator = "DWLS")). It therefore requires raw ordinal data with cor_method = "poly" or "tetra" and has no fallback for a supplied correlation matrix or a continuous cor_method. Because the weighting follows the polychoric asymptotic covariance, the matrix and the weights are estimated on the listwise-complete cases. Its fit-index behaviour is described under Fit indices.

Correlation methods

When raw data are supplied, cor_method selects how the correlation matrix is computed (it is ignored when a correlation matrix is entered directly).

• "pearson" (default), "spearman", and "kendall" are passed to stats::cor() for continuous or rank data.

• "poly" / "tetra" compute polychoric / tetrachoric correlations for ordinal / binary data, assuming an underlying bivariate-normal latent variable. They use a two-step estimator with no empty-cell continuity correction, matching polycor::polychor() and lavaan. The polychoric asymptotic covariance that underlies both the DWLS weights and the scaled (sandwich) statistic relies on large-sample theory that degrades for empty or near-empty response-category combinations; with very sparse cells the resulting weights and standard errors can be unreliable (a warning is issued when empty cells are present), so interpret them with caution and consider collapsing rare categories.

• "fiml" estimates a two-stage full-information maximum-likelihood correlation. The saturated multivariate-normal mean and covariance are estimated from raw data with missing values by an EM algorithm assuming the data are missing at random (Yuan, Marshall, & Bentler, 2002; Little & Rubin, 2002), and the standardized covariance is then analysed. This reproduces psych::corFiml() followed by psych::fa() and lavaan(missing = "two.stage"), not lavaan::efa(missing = "ml"), so the point estimates are not expected to match the latter. The model fit indices are corrected two-stage statistics (see Fit indices). "fiml" uses every case and handles the missingness itself, so use is ignored; it supplies a continuous (Pearson-type) correlation only and is therefore not compatible with method = "DWLS". Standard errors are available analytically for method = "ML" or "ULS" and, for any method, by the non-parametric bootstrap (see Standard errors). For multiply imputed data, EFA_POOLED() is the alternative route to handling missingness.

Rotations

A rotation transforms the unrotated loadings toward a simpler, more interpretable pattern; all rotations are performed by rotation engines built into the package. Orthogonal rotations keep the factors uncorrelated, whereas oblique rotations let them correlate (returning a pattern matrix, a structure matrix, and the factor intercorrelations Phi) and are usually more realistic for psychological constructs.

Orthogonal rotations:

• varimax maximises the variance of the squared loadings within each factor (column simplicity). It is the most widely used orthogonal rotation and spreads variance across factors rather than concentrating it in a general factor.

• quartimax simplifies the variables (rows) so that each loads mainly on one factor; it tends to produce a strong general factor.

• equamax is a Crawford-Ferguson compromise between varimax (column) and quartimax (row) simplicity.

• geominT uses a geometric-mean criterion that rewards a sparse pattern and tolerates variables with cross-loadings; a smaller offset delta gives a sparser solution but sharper local minima.

• bentlerT uses Bentler's invariant pattern simplicity criterion.

• bifactorT is the Jennrich-Bentler orthogonal bifactor criterion: a general factor plus group factors (bifactor simple structure).

Oblique rotations:

• promax is a fast two-step rotation: a varimax solution is raised to a power (controlled by k and P_type) to form a target that is then fitted obliquely. It is the common, inexpensive oblique default.

• oblimin is a flexible oblique family controlled by gam (default 0); a good general-purpose criterion.

• quartimin is oblimin pinned at gam = 0; a robust default oblique criterion.

• simplimax drives the k smallest loadings toward zero. Its criterion is only piecewise smooth, so it is the most prone to local minima and relies on several randomStarts.

• bentlerQ is the oblique Bentler invariant pattern simplicity criterion.

• geominQ is the oblique geomin criterion; it handles complex (cross-loading) structure well but is multimodal, so it benefits from more randomStarts (and uses a more thorough multi-start search internally).

• bifactorQ is the oblique (correlated) Jennrich-Bentler bifactor criterion.

The criterion-based rotations (all except varimax and promax) are fitted by gradient projection with randomStarts random starts to guard against local minima; the complexity criteria (simplimax and geominQ in particular) are the most multimodal. The type argument changes the varimax and promax settings (see Using the type presets) and, for every rotation, the factor order_type. A single factor cannot be rotated.

Standard errors

se selects whether and how standard errors (and matching confidence intervals) are computed. They cover the unrotated loadings and uniquenesses and, when a rotation is applied, the rotated loadings, the communalities, and – for oblique rotations – the factor correlations and the structure coefficients (see the SE and CI entries in Value).

• "none" (default) computes no standard errors.

• "information" returns analytic standard errors from the expected (Fisher) information matrix of the maximum-likelihood solution, and therefore requires method = "ML". The rotated standard errors are obtained by propagating the unrotated-loading covariance through the rotation by the delta method (Jennrich, 1973); because rotated quantities are identification-invariant they are directly comparable across programs. Unlike the bootstrap it also works from a correlation matrix as long as N is supplied. The covariance is the inverse expected information under the identification constraint that \Lambda' \Psi^{-1} \Lambda is diagonal, scaled by 1 / (N - 1); the confidence intervals are Wald intervals (estimate \pm z * SE). These standard errors assume multivariate normality and a correctly specified model; under heavy-tailed data or model misfit they can understate the sampling variability, where a bootstrap is more robust. The rotated structure-coefficient intervals are somewhat conservative for high-communality variables, where "sandwich" or "np-boot" give sharper intervals.

• "sandwich" returns robust (Godambe sandwich) standard errors from raw data, combining the estimator weight with an asymptotic-distribution-free covariance of the correlations, so it stays valid under non-normality and weight misspecification (Browne, 1984; Satorra & Bentler, 1994). It is available either for ordinal data with cor_method = "poly" or "tetra" and method one of "ML", "ULS", or "DWLS" (the meat is the polychoric / tetrachoric asymptotic covariance), or for continuous data with cor_method = "pearson" and method = "ML" or "ULS" (the meat is the fourth-moment ADF covariance of the sample correlations, the basis of the MLM / MLR robust statistics). It reports the same standard errors as "information", propagated by the same delta method, and additionally fills the model fit's chi-square block with a scaled (Satorra-Bentler / scaled-and-shifted) chi-square (see Fit indices). Because the asymptotic covariance must describe the same cases as the correlation matrix, the sandwich (like method = "DWLS") is computed on the listwise-complete cases; on data with missing values the reported N, the correlation matrix, and the point estimate therefore reflect the complete cases regardless of use.

• "np-boot" draws a non-parametric (case-resampling) bootstrap and needs raw data. It is the most general method – available for any method, rotation, and cor_method – and the most robust to non-normality and misfit, at the cost of speed; its intervals are bootstrap percentile intervals. The replicate fits are run across replicates with the future framework. By default they run sequentially; to run them in parallel, register a plan with future::plan() (e.g. future::plan(future::multisession, workers = 2); see examples). With a fixed seed the bootstrap is reproducible and yields the same result regardless of the number of workers. Under cor_method = "fiml" each resample also re-runs the EM moment estimation and is therefore slow, so a smaller b_boot may be advisable.

The analytic methods ("information" and "sandwich") are not available with the "promax" or "simplimax" rotations, which have no usable analytic rotation Jacobian; use "np-boot" there. Under cor_method = "fiml", "information" and "sandwich" instead return, for method = "ML" or "ULS", the corrected two-stage (Yuan & Bentler, 2000; Savalei & Bentler, 2009) sandwich standard errors, built on the saturated FIML asymptotic covariance with the estimator's own Stage-2 weight: the model is fitted to the EM-estimated correlation, so the naive Stage-2 standard errors (treating that correlation as complete data) are inconsistent under missingness and are not reported (method = "PAF" carries no Stage-2 weight, so use se = "np-boot" there).

Fit indices

For ML and ULS, EFA() returns the model chi-square (with its p-value and degrees of freedom), the Comparative Fit Index (CFI; Bentler, 1990), the Tucker-Lewis Index (TLI, also called the non-normed fit index; Tucker & Lewis, 1973), the Root Mean Square Error of Approximation (RMSEA) with its 90% confidence interval (Browne & Cudeck, 1992), the Akaike and Bayesian Information Criteria (AIC, BIC), the Expected Cross-Validation Index (ECVI; Browne & Cudeck, 1989), the Root Mean Squared Residual (RMSR), the Standardized Root Mean Squared Residual (SRMR; Bentler, 1995), and the common-part-accounted-for (CAF) index (Lorenzo-Seva, Timmerman, & Kiers, 2011). The print and summary methods show SRMR, not RMSR, because the two residual summaries differ only by the fixed scaling \sqrt{(p - 1) / (p + 1)} for a fixed number of variables; RMSR remains in the returned object. The model chi-square is the Bartlett-corrected discrepancy (matching stats::factanal() for ML); the AIC, BIC, and ECVI are the minimum-fit-function (chi-square-based) forms (\chi^2 - 2\,df and \chi^2 - \log(N)\,df for AIC and BIC, as in psych::fa()) and can therefore be negative. The RMSEA, CFI, and TLI place the model and baseline noncentrality on the uncorrected N - 1 discrepancy scale on which these approximate-fit indices are defined, so the Bartlett small-sample correction enters only the chi-square test, not the approximate-fit indices.

Which indices are reported depends on the estimator:

• ML and ULS compute the full set above.

• PAF returns only the descriptive residual indices (RMSR, SRMR, CAF) and df; the printed model-fit block shows CAF and SRMR. The chi-square-based indices are NA, because PAF minimises no discrepancy.

• DWLS by default returns only RMSR, SRMR, CAF, and df, because the ordinary maximum-likelihood discrepancy is not its fit function. When se = "sandwich", a scaled (Satorra & Bentler, 1994; Asparouhov & Muthén, 2010) chi-square and the CFI, TLI, and RMSEA derived from it are reported (AIC and BIC remain NA). That scaled statistic is a two-stage correction applied to the polychoric-correlation residuals (Browne, 1984), so it is not identical to the full WLSMV test of lavaan or Mplus, which also projects the thresholds.

• cor_method = "fiml" (with ML or ULS) reports Satorra-Bentler-corrected two-stage statistics (Yuan, Marshall, & Bentler, 2002): the normal-theory discrepancy on the EM-estimated correlation, rescaled by the saturated FIML asymptotic covariance, because the plain two-stage likelihood-ratio statistic is not asymptotically \chi^2(df). The CFI, TLI, and RMSEA follow from the scaled statistics; AIC, BIC, and ECVI are left NA, as for any scaled (moment-adjusted) chi-square.

Whenever the chi-square is a scaled one (se = "sandwich", or any cor_method = "fiml" fit), the AIC, BIC, and ECVI are NA and the returned fit_indices additionally carry the scaled-statistic components (see the fit_indices entry in Value). Note that Lorenzo-Seva, Timmerman, and Kiers (2011) introduce the CAF as ranging between 0 and 1, with values close to 1 indicating close fit; this does not match the formula they apply, 1 - KMO(residuals), which only works if the diagonal of the residual matrix is set to 1s and then approximates 0.5 with close fit.

Available combinations

Not every estimator, rotation, standard-error, and correlation method can be combined:

• Estimator and correlation method. method = "DWLS" requires ordinal data with cor_method = "poly" or "tetra". cor_method = "fiml" works with PAF, ML, and ULS (not DWLS) and needs raw data with missing values.

• Standard errors. se = "information" requires method = "ML" and can be computed from a correlation matrix when N is supplied. se = "sandwich" requires raw data, with either a polychoric/tetrachoric cor_method (ML, ULS, or DWLS) or a Pearson cor_method (ML or ULS); it is not available for PAF. Under cor_method = "fiml", "information" and "sandwich" are available for ML and ULS only and both return the corrected two-stage sandwich. se = "np-boot" requires raw data and works with any estimator, rotation, and correlation method. Neither "information" nor "sandwich" is available with the "promax" or "simplimax" rotations.

• Fit indices. The chi-square-based indices are available for ML and ULS (and, as scaled statistics, for cor_method = "fiml" and for DWLS with se = "sandwich"); PAF and DWLS otherwise report only the descriptive residual indices.

Using the type presets

The type argument is evaluated for PAF and for all rotations (mainly important for the varimax and promax rotations). The type-specific settings for these functions are detailed below.

For PAF, the values of init_comm, criterion, criterion_type, max_iter, and abs_eigen depend on the type argument.

type = "EFAtools" will use the following argument specification: ⁠init_comm = "smc", criterion = .001, criterion_type = "sum", max_iter = 300, abs_eigen = TRUE⁠.

type = "psych" will use the following argument specification: ⁠init_comm = "smc", criterion = .001, criterion_type = "sum", max_iter = 50, abs_eigen = FALSE⁠.

type = "SPSS" will use the following argument specification: ⁠init_comm = "smc", criterion = .001, criterion_type = "max_individual", max_iter = 25, abs_eigen = TRUE⁠.

If SMCs fail, SPSS takes "mac". However, as SPSS takes absolute eigenvalues, this is hardly ever the case. Psych, on the other hand, takes "unity" if SMCs fail, but uses the Moore-Penrose Psudo Inverse of a matrix, thus, taking "unity" is only necessary if negative eigenvalues occur afterwards in the iterative PAF procedure. The EFAtools type setting combination was the best in terms of accuracy and number of Heywood cases compared to all the other setting combinations tested in simulation studies in Grieder & Steiner (2022), which is why this type is used as a default here.

For varimax, the values of varimax_type and order_type depend on the type argument.

type = "EFAtools" will use the following argument specification: ⁠varimax_type = "kaiser", order_type = "eigen"⁠.

type = "psych" will use the following argument specification: ⁠varimax_type = "svd", order_type = "eigen"⁠.

type = "SPSS" will use the following argument specification: ⁠varimax_type = "kaiser", order_type = "ss_factors"⁠.

For promax, the values of P_type, order_type, and k depend on the type argument.

type = "EFAtools" will use the following argument specification: ⁠P_type = "norm", order_type = "eigen", k = 4⁠.

type = "psych" will use the following argument specification: ⁠P_type = "unnorm", order_type = "eigen", k = 4⁠.

type = "SPSS" will use the following argument specification: ⁠P_type = "norm", order_type = "ss_factors", k = 4⁠.

The P_type argument can take two values, "unnorm" and "norm". It controls which formula is used to compute the target matrix P in the promax rotation. "unnorm" uses the formula from Hendrickson and White (1964), specifically: P = abs(A^(k + 1)) / A, where A is the unnormalized matrix containing varimax rotated loadings. "norm" uses the normalized varimax rotated loadings. Specifically it used the following formula, which can be found in the SPSS 23 and SPSS 27 Algorithms manuals: P = abs(A / sqrt(rowSums(A^2))) ^(k + 1) * (sqrt(rowSums(A^2)) / A). As for PAF, the EFAtools type setting combination for promax was the best compared to the other setting combinations tested in simulation studies in Grieder & Steiner (2022). Note that all type presets keep the EFAtools default Kaiser normalization (normalize = TRUE), whereas psych::fa() does not normalize before its promax target rotation; set normalize = FALSE to reproduce the psych::fa() promax result exactly.

The varimax_type argument can take two values, "svd", and "kaiser". "svd" uses singular value decomposition, by calling stats::varimax(). "kaiser" performs the varimax procedure as described in the SPSS Algorithms manual and by Kaiser (1958). The varimax simplicity criterion monitored for convergence is sum(n*colSums(lambda ^ 4) - colSums(lambda ^ 2) ^ 2) / n ^ 2, where n is the number of indicators, and lambda is the Kaiser-normalized rotated loadings matrix.

For all other rotations except varimax and promax, the type argument only controls the order_type argument with the same values as stated above for the varimax and promax rotations. Additional arguments can also be specified and will be passed to the rotation procedure (e.g., maxit to change the maximum number of iterations).

The type argument has no effect on ULS and ML. For ULS, no additional arguments are needed. For ML, an additional argument start_method is needed to determine the starting values for the optimization procedure. Default for this argument is "psych" which takes the starting values specified in psych::fa().

Value

A list of class EFA containing (a subset of) the following:

Source

Grieder, S., & Steiner, M. D. (2022). Algorithmic jingle jungle: A comparison of implementations of principal axis factoring and promax rotation in R and SPSS. Behavior Research Methods, 54, 54–74. doi: 10.3758/s13428-021-01581-x

Hendrickson, A. E., & White, P. O. (1964). Promax: A quick method for rotation to oblique simple structure. British Journal of Statistical Psychology, 17 , 65–70. doi: 10.1111/j.2044-8317.1964.tb00244.x

Lorenzo-Seva, U., Timmerman, M. E., & Kiers, H. A. L. (2011). The Hull Method for Selecting the Number of Common Factors, Multivariate Behavioral Research, 46, 340-364, doi: 10.1080/00273171.2011.564527

Kaiser, H. F. (1958). The varimax criterion for analytic rotation in factor analysis. Psychometrika, 23, 187–200. doi: 10.1007/BF02289233

Lawley, D. N., & Maxwell, A. E. (1971). Factor analysis as a statistical method (2nd ed.). Butterworths.

Jennrich, R. I., & Thayer, D. T. (1973). A note on Lawley's formulas for standard errors in maximum likelihood factor analysis. Psychometrika, 38, 571–580. doi: 10.1007/BF02291495

Jennrich, R. I. (1973). Standard errors for obliquely rotated factor loadings. Psychometrika, 38, 593–604. doi: 10.1007/BF02291497

Zhang, G., & Preacher, K. J. (2015). Factor rotation and standard errors in exploratory factor analysis. Journal of Educational and Behavioral Statistics, 40, 579–603. doi: 10.3102/1076998615606098

Browne, M. W. (1984). Asymptotically distribution-free methods for the analysis of covariance structures. British Journal of Mathematical and Statistical Psychology, 37, 62–83. doi: 10.1111/j.2044-8317.1984.tb00789.x

Satorra, A., & Bentler, P. M. (1994). Corrections to test statistics and standard errors in covariance structure analysis. In A. von Eye & C. C. Clogg (Eds.), Latent variables analysis (pp. 399–419). Sage.

Asparouhov, T., & Muthén, B. (2010). Simple second order chi-square correction. Mplus Technical Appendix.

Muthén, B. (1984). A general structural equation model with dichotomous, ordered categorical, and continuous latent variable indicators. Psychometrika, 49, 115–132. doi: 10.1007/BF02294210

Yuan, K.-H., & Bentler, P. M. (2000). Three likelihood-based methods for mean and covariance structure analysis with nonnormal missing data. Sociological Methodology, 30, 165–200. doi: 10.1111/0081-1750.00078

Yuan, K.-H., Marshall, L. L., & Bentler, P. M. (2002). A unified approach to exploratory factor analysis with missing data, nonnormal data, and in the presence of outliers. Psychometrika, 67, 95–121. doi: 10.1007/BF02294711

Savalei, V., & Bentler, P. M. (2009). A two-stage approach to missing data: Theory and application to auxiliary variables. Structural Equation Modeling, 16, 477–497. doi: 10.1080/10705510903008238

Little, R. J. A., & Rubin, D. B. (2002). Statistical analysis with missing data (2nd ed.). Wiley.

Bartlett, M. S. (1951). The effect of standardization on approximation in factor analysis. Biometrika, 38, 337–344.

Bentler, P. M. (1990). Comparative fit indexes in structural models. Psychological Bulletin, 107, 238–246. doi: 10.1037/0033-2909.107.2.238

Tucker, L. R., & Lewis, C. (1973). A reliability coefficient for maximum likelihood factor analysis. Psychometrika, 38, 1–10. doi: 10.1007/BF02291170

Browne, M. W., & Cudeck, R. (1989). Single sample cross-validation indices for covariance structures. Multivariate Behavioral Research, 24, 445–455. doi: 10.1207/s15327906mbr2404_4

Browne, M. W., & Cudeck, R. (1992). Alternative ways of assessing model fit. Sociological Methods & Research, 21, 230–258. doi: 10.1177/0049124192021002005

Bentler, P. M. (1995). EQS structural equations program manual. Multivariate Software.

Examples

# Principal axis factoring with oblimin rotation
mod_oblimin <- EFA(test_models$baseline$cormat, n_factors = 3, N = 500,
                   rotation = "oblimin")
mod_oblimin
summary(mod_oblimin)

# ML estimation with oblimin rotation
mod_oblimin <- EFA(test_models$baseline$cormat, n_factors = 3, N = 500,
                   method = "ML", rotation = "oblimin")
mod_oblimin
summary(mod_oblimin)

# Analytic (expected-information) standard errors for the above
ML_info <- EFA(test_models$baseline$cormat, n_factors = 3, N = 500,
               method = "ML", rotation = "oblimin", se = "information")
ML_info
summary(ML_info)


# Robust (sandwich) standard errors and a scaled chi-square for ordinal raw data.
# These need a polychoric/tetrachoric correlation method and method ML, ULS, or DWLS.
DWLS_rob <- EFA(DOSPERT_raw, n_factors = 6, cor_method = "poly",
                method = "DWLS", rotation = "oblimin", se = "sandwich")
DWLS_rob
summary(DWLS_rob)

# The same robust SEs and scaled chi-square for continuous data: a Pearson
# correlation with method ML or ULS (the fourth-moment ADF covariance).
ML_rob <- EFA(GRiPS_raw, n_factors = 1, cor_method = "pearson",
              method = "ML", rotation = "none", se = "sandwich")
ML_rob
summary(ML_rob)



# Two-stage FIML correlations from raw data with missing values: the saturated
# multivariate-normal moments are EM-estimated (assuming the data are missing at
# random) and the standardized covariance is analysed.
x_miss <- GRiPS_raw
x_miss[cbind(1:20, 1)] <- NA
EFA_fiml <- EFA(x_miss, n_factors = 1, method = "ML", cor_method = "fiml")
EFA_fiml


## Not run: 
# Bootstrap standard errors from raw data, reproducible via a fixed seed and run
# in parallel across replicates.
future::plan(future::multisession, workers = 2)
EFA_boot <- EFA(GRiPS_raw, n_factors = 1, method = "PAF", rotation = "none",
                se = "np-boot", b_boot = 1000, seed = 42)
future::plan(future::sequential)

## End(Not run)

Model averaging across different EFA methods and types

Description

Not all EFA procedures always arrive at the same solution. This function allows you perform a number of EFAs from different methods (e.g., Maximum Likelihood and Principal Axis Factoring), with different implementations (e.g., the SPSS and psych implementations of Principal Axis Factoring), and across different rotations of the same type (e.g., multiple oblique rotations, like promax and oblimin). EFA_AVERAGE will then run all these EFAs (using the EFA() function) and provide a summary across the different solutions.

Usage

EFA_AVERAGE(
  x,
  n_factors,
  N = NA,
  method = "PAF",
  rotation = "promax",
  type = "none",
  averaging = c("mean", "median"),
  trim = 0,
  salience_threshold = 0.3,
  max_iter = 10000,
  init_comm = c("smc", "mac", "unity"),
  criterion = c(0.001),
  criterion_type = c("sum", "max_individual"),
  abs_eigen = c(TRUE),
  varimax_type = c("svd", "kaiser"),
  normalize = TRUE,
  k_promax = 2:4,
  k_simplimax = ncol(x),
  P_type = c("norm", "unnorm"),
  precision = 1e-05,
  start_method = c("psych", "factanal"),
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra", "fiml"),
  show_progress = TRUE
)

Arguments

Details

As a first step in this function, a grid is produced containing the setting combinations for the to-be-performed EFAs. These settings are then entered as arguments to the EFA() function and the EFAs are run in a second step. After all EFAs are run, the factor solutions are averaged and their variability determined in a third step.

When raw data are supplied, the correlation matrix is computed once before the grid is run and reused for every EFA in it. Under cor_method = "fiml" this means the saturated multivariate-normal moments are EM-estimated a single time (from the raw data with missing values, assuming the data are missing at random) and the resulting two-stage correlation is analysed by every solution in the grid; the EM is not re-run per solution. Under cor_method = "fiml", use does not select cases (every case contributes to the EM). The averaged loadings and communalities are then the two-stage FIML estimates, but the averaged Chi-Square and the indices derived from it (CFI, TLI, RMSEA, AIC, BIC, ECVI) are the ordinary ML/ULS discrepancy statistics on the EM correlation, not the corrected two-stage (Satorra-Bentler) statistics that a standalone EFA() with cor_method = "fiml" reports; in particular the averaged AIC and BIC are finite here rather than NA.

The grid containing the setting combinations is produced based on the entries to the respective arguments. To this end, all possible combinations resulting in unique EFA models are considered. That is, if, for example, the type argument was set to c("none", "SPSS") and one combination of the specific settings entered was identical to the SPSS combination, this combination would be included in the grid and run only once. We include here a list of arguments that are only evaluated under specific conditions:

The arguments init_comm, criterion, criterion_type, abs_eigen are only evaluated if "PAF" is included in method and "none" is included in type.

The argument varimax_type is only evaluated if "varimax", "promax", "oblique", or "orthogonal" is included in rotation and "none" is included in type.

The argument normalize is only evaluated if rotation is not set to "none" and "none" is included in type.

The argument k_simplimax is only evaluated if "simplimax" or "oblique" is included in rotation.

The arguments k_promax and P_type are only evaluated if "promax" or "oblique" is included in rotation and "none" is included in type.

The argument start_method is only evaluated if "ML" is included in method.

To avoid a bias in the averaged factor solutions from problematic solutions, these are excluded prior to averaging. A solution is deemed problematic if at least one of the following is true: an error occurred, the model did not converge, or there is at least one Heywood (improper) case (a communality at or above 1, or, for ML/ULS, a uniqueness pinned at the estimator's lower bound). Information on errors, convergence, and Heywood cases are returned in the implementations_grid and a summary of these is given when printing the output. In addition to these, information on the admissibility of the factor solutions is also included. A solution was deemed admissible if (1) no error occurred, (2) the model converged, (3) no Heywood cases are present, and (4) there are at least two salient loadings (i.e., loadings exceeding the specified salience_threshold) for each factor. So, solutions failing one of the first three of these criteria of admissibility are also deemed problematic and therefore excluded from averaging. However, solutions failing only the fourth criterion of admissibility are still included for averaging. Finally, if all solutions are problematic (e.g., all solutions contain Heywood cases), no averaging is performed and the respective outputs are NA. In this case, the implementations_grid should be inspected to see if there are any error messages, and the separate EFA solutions that are also included in the output can be inspected as well, for example, to see where Heywood cases occurred.

A core output of this function includes the average, minimum, and maximum loadings derived from all non-problematic (see above) factor solutions. Please note that these are not entire solutions, but the matrices include the average, minimum, or maximum value for each cell (i.e., each loading separately). This means that, for example, the matrix with the minimum loadings will contain the minimum value in any of the factor solutions for each specific loading, and therefore most likely contains loadings from different factor solutions. The matrices containing the minimum and maximum factor solutions can therefore not be interpreted as whole factor solutions.

The output also includes information on the average, minimum, maximum, and variability of the fit indices across the non-problematic factor solutions. It is important to note that not all fit indices are computed for all fit methods: For ML and ULS, all fit indices can be computed, while for PAF the chi-square-based indices (the chi-square statistic and its significance, CFI, TLI, RMSEA, AIC, BIC, and ECVI) are NA. The common part accounted for (CAF) index (Lorenzo-Seva, Timmerman, & Kiers, 2011) and the residual-based SRMR and RMSR are still computed for PAF. As a consequence, if only "PAF" is included in the method argument, averaging is performed for the CAF, SRMR, and RMSR, while the chi-square-based indices are NA. If a combination of "PAF" and "ML" and/or "ULS" are included in the method argument, the CAF, SRMR, and RMSR are averaged across all non-problematic factor solutions, while the chi-square-based indices are only averaged across the ML and ULS solutions. The user should therefore keep in mind that the number of EFAs across which the fit indices are averaged can diverge for the CAF, SRMR, and RMSR compared to the chi-square-based indices.

Each reported fit index is summarised across the (non-problematic) solutions in the same descriptive way: the average, standard deviation, minimum, and maximum of the per-solution values. This includes the chi-square significance level (p_chi), which is therefore the mean (or median) of the per-solution p-values and is purely descriptive; it is not the p-value of any pooled chi-square test.

Value

A list of class EFA_AVERAGE containing

Source

Grieder, S., & Steiner, M. D. (2022). Algorithmic jingle jungle: A comparison of implementations of principal axis factoring and promax rotation in R and SPSS. Behavior Research Methods, 54, 54–74. doi: 10.3758/s13428-021-01581-x

Hendrickson, A. E., & White, P. O. (1964). Promax: A quick method for rotation to oblique simple structure. British Journal of Statistical Psychology, 17 , 65–70. doi: 10.1111/j.2044-8317.1964.tb00244.x

Lorenzo-Seva, U., Timmerman, M. E., & Kiers, H. A. L. (2011). The Hull Method for Selecting the Number of Common Factors, Multivariate Behavioral Research, 46, 340-364, doi: 10.1080/00273171.2011.564527

Kaiser, H. F. (1958). The varimax criterion for analytic rotation in factor analysis. Psychometrika, 23, 187–200. doi: 10.1007/BF02289233

Examples

## Not run: 
# Averaging across different implementations of PAF and promax rotation (72 EFAs)
Aver_PAF <- EFA_AVERAGE(test_models$baseline$cormat, n_factors = 3, N = 500)

# Use median instead of mean for averaging (72 EFAs)
Aver_PAF_md <- EFA_AVERAGE(test_models$baseline$cormat, n_factors = 3, N = 500,
                           averaging = "median")

# Averaging across different implementations of PAF and promax rotation,
# and across ULS and different versions of ML (108 EFAs)
Aver_meth_ext <- EFA_AVERAGE(test_models$baseline$cormat, n_factors = 3, N = 500,
                             method = c("PAF", "ULS", "ML"))

# Averaging across one implementation each of PAF (EFAtools type), ULS, and
# ML with one implementation of promax (EFAtools type) (3 EFAs)
Aver_meth <- EFA_AVERAGE(test_models$baseline$cormat, n_factors = 3, N = 500,
                         method = c("PAF", "ULS", "ML"), type = "EFAtools",
                         start_method = "psych")

# Averaging across different oblique rotation methods, using one implementation
# of ML and one implementation of promax (EFAtools type) (7 EFAs)
Aver_rot <- EFA_AVERAGE(test_models$baseline$cormat, n_factors = 3, N = 500,
                         method = "ML", rotation = "oblique", type = "EFAtools",
                         start_method = "psych")

## End(Not run)

## Not run: 
# Two-stage FIML correlations from raw data with missing values: the EM
# saturated moments are estimated once and the resulting correlation is
# averaged across the grid of EFAs.
x_miss <- GRiPS_raw
x_miss[cbind(1:20, 1)] <- NA
Aver_fiml <- EFA_AVERAGE(x_miss, n_factors = 1, method = c("PAF", "ML"),
                         cor_method = "fiml")

## End(Not run)

Exploratory factor analysis on multiple data imputations

Description

Fits EFA() to each of several imputed datasets, aligns the factor solutions to a common factor space, and pools the resulting estimates and selected fit quantities across imputations.

Usage

EFA_POOLED(
  data_list,
  p = 0.05,
  target_method = c("first_target", "consensus"),
  align_unrotated = c("signed_tucker_congruence", "none", "procrustes"),
  fit_pool_method = c("D2"),
  consensus_args = list(),
  procrustes_args = list(),
  rmsea_ci_level = 0.9,
  rmsr_upper = TRUE,
  ...
)

Arguments

Details

EFA_POOLED() is the multiple-imputation route to handling missing data: several imputed datasets are each fitted with EFA() and the solutions pooled. A single-fit alternative is full-information maximum likelihood, available directly in EFA() as cor_method = "fiml", which EM-estimates a two-stage correlation from one raw dataset with missing values. Both feed the same correlation-scale EFA core and differ only in how the missingness is handled; FIML is intentionally not routed through EFA_POOLED(), which is a multi-fit pooler by construction.

Standard-error pooling routes

The pooling pathway is selected automatically from the se method recorded on the component EFA() fits, which must be identical across imputations:

• se = "none": no standard errors are pooled.

• se = "information": the per-imputation expected-information standard errors are pooled with Rubin's (1987) rules (Wald intervals).

• se = "sandwich": the two-stage pooled-inputs (MI2S) approach fits a single model on the Rubin-pooled correlation matrix and asymptotic covariance.

• se = "np-boot": the non-parametric bootstrap replicates are re-aligned to the multiple-imputation target and Rubin-pooled.

On the information and np-boot routes, if pooled standard errors cannot be produced (for example an unreliable analytic covariance or too few bootstrap replicates) the pool falls back to point-estimate-only pooling and downgrades settings$se to "none". The MI2S route is the exception: its single fit fuses the point estimates and standard errors through the pooled asymptotic covariance, so a structural failure aborts directly rather than falling back.

Aligning solutions across imputations

The same EFA() model is fitted to each imputed dataset and the solutions are put into a common factor space before averaging. For oblique solutions the factor intercorrelations are aligned together with the loadings so the model stays internally consistent.

target_method controls how rotated solutions are aligned. "first_target" (the default) aligns every imputation to the first imputation's rotated solution by one Procrustes rotation each. "consensus" instead refines a centroid target by Generalized Procrustes Analysis (Gower 1975; van Ginkel & Kroonenberg 2014; Lorenzo-Seva & Van Ginkel 2016). The two give the same pooled estimate for orthogonal rotations (consensus is just more expensive), and "consensus" is only supported there. Anchoring on the first imputation can understate the between-imputation variability when the imputations disagree substantially, whereas "consensus" is more robust to an atypical first imputation (van Ginkel & Kroonenberg 2014).

align_unrotated controls how unrotated loadings are aligned before pooling: "signed_tucker_congruence" (the default) matches them up to factor reordering and sign changes, "procrustes" aligns them to the first imputation by orthogonal Procrustes rotation, and "none" averages them as returned by EFA().

Pooling point estimates

Point estimates are pooled by arithmetic averaging after alignment. For oblique rotations the structure matrix is recomputed from the pooled pattern matrix and pooled factor correlations, Structure = \Lambda \Phi, and communalities are the diagonal of the reproduced correlation matrix, diag(\Lambda \Phi \Lambda') for oblique rotations and diag(\Lambda \Lambda') otherwise. Residuals are not averaged across imputations; they are the pooled observed correlation matrix minus the model-implied correlation of the pooled solution, so RMSR/SRMR are based on these pooled residuals. Both are returned, though the print and summary methods show SRMR only.

Pooling the model chi-square and fit indices

The model chi-square and the indices derived from it (RMSEA, ECVI, and the descriptive AIC/BIC) are pooled with the D2 rule (Li, Meng, Raghunathan & Rubin, 1991), not arithmetically averaged. Because D2 shrinks the pooled chi-square in proportion to the between-imputation variability, the pooled RMSEA can fall below the mean of the per-imputation RMSEAs (as it does in lavaan.mi); read it together with the per-imputation fit. The incremental indices CFI (Bentler, 1990) and TLI (Tucker & Lewis, 1973) are instead the average of the per-imputation indices, which keeps them consistent with the component fits and avoids the out-of-range values that separately pooling the model and baseline noncentralities (as lavaan.mi/semTools do) can produce; those separately pooled noncentralities remain available in mi_diagnostics. AIC and BIC, if returned, are chi-square-derived descriptive quantities and are not likelihood-based MI information criteria. On the sandwich/MI2S route the chi-square is the single fit's scaled statistic rather than a D2 pool.

Bootstrap pooling (np-boot)

If each component EFA call was run with se = "np-boot", pooled bootstrap SEs and Wald-type MI confidence intervals are computed for loadings, communalities, residuals, and, when applicable, factor correlations and structure coefficients. The unrotated bootstrap replicates are re-aligned to the final MI target before the within-imputation covariance is estimated, and Rubin pooling is applied with T = Ubar + (1 + 1/m) B. The confidence level of the pooled intervals is set by p, not by the component EFA calls' ci.

Analytic pooling (information)

With se = "information", the analytic unrotated-loading and uniqueness SEs returned by each fit are pooled element-wise with Rubin's rules (T = Ubar + (1 + 1/m) B), with Wald intervals on the plain Rubin (1987) degrees of freedom (the analytic loadings are asymptotically normal, so the Barnard-Rubin (1999) adjustment reduces to this form, matching lavaan.mi). NA propagation is fail-closed: if any imputation is NA at an element, all pooled outputs for that element are NA. When a rotation was requested, the rotated loadings, communalities, and (for oblique rotations) factor correlations and structure coefficients are pooled as well; residual SE pooling is available only on the bootstrap path. Under align_unrotated = "procrustes" the full unrotated covariance vcov_unrot_loadings (populated by se = "information") is propagated through the alignment, so it must be present and reliable on every fit.

A rotated-loading standard error is conditional on the rotation criterion (Jennrich 1973, 1974; Browne 2001; Zhang & Preacher 2015). For both orthogonal and oblique rotations the within-imputation variance is therefore each fit's own criterion-aware delta-method rotated SE (the quantity EFA() returns), reused after a signed-permutation alignment to the MI target, and the between-imputation variance is the sample variance of the aligned rotated loadings. This is a deliberate approximation – each SE is conditional on its own fit's rotation optimum rather than on a common gauge – and is flagged by ⁠MI$<param>$method = "signed_permutation_approx"⁠. Communalities are rotation-invariant and pool element-wise. For a fully gauge-consistent rotated uncertainty, cross-check with se = "np-boot".

Two-stage pooling (sandwich / MI2S)

With se = "sandwich" (robust SEs from a polychoric/tetrachoric or continuous-Pearson asymptotic covariance), pooling follows the two-stage, pooled-inputs approach (Chung & Cai 2019; Sriutaisuk, Liu, Chung, Kim & Gu 2025): the correlation matrix and the asymptotic covariance of its off-diagonal entries are Rubin-pooled across imputations,

\tilde\Gamma = \Gamma_W + \left(1 + \frac{1}{m}\right)\Gamma_B,

and a single EFA model is fitted to the pooled correlation with \tilde\Gamma as the robust meat (its diagonal as the weights for method = "DWLS"). Because there is only one fit and one rotational gauge, this route bypasses the per-imputation alignment: target_method and align_unrotated do not apply. The fitted object carries native scaled-shifted chi-square statistics and sandwich SEs that already reflect the multiple-imputation uncertainty, so the chi-square is not D2-pooled and the likelihood-ratio-based AIC/BIC/ECVI are NA; it is returned in the mi_fit slot, with the per-imputation fits retained for diagnostics. At least 20 imputations are recommended for the scaled-shifted statistic, and more (around 100) at higher rates of missingness (Sriutaisuk et al. 2025). The polychoric/tetrachoric (ordinal) case is the primary, best-evaluated target; the continuous-Pearson case uses the same recipe but is less benchmarked.

Value

A list of class "EFA_POOLED" containing pooled estimates, residuals, fit indices, the individual fits, and MI diagnostics. In addition to the slots inherited from EFA() (including SE, CI, and, on the bootstrap path, replicates), the object carries:

MI

Multiple-imputation diagnostics for each pooled parameter family. On the bootstrap path: unrot_loadings, h2, residuals, optionally rot_loadings, Phi, Structure, and fit_indices_descriptive, plus integer vectors bootstrap_source_failures (replicates the component EFA could not fit), bootstrap_rotation_failures (replicates whose Procrustes alignment to the target was invalid), and bootstrap_rotation_valid (those that entered the pool, B - source - rotation failures). Both paths use the plain Rubin (1987) df. On the analytic path (se = "information"): unrot_loadings and uniquenesses, plus, when a rotation was requested, rot_loadings, h2, and (oblique) Phi and Structure. Each per-family entry is a list with RIV (relative increase in variance), FMI (the fraction of missing information, reported as Rubin's asymptotic \lambda = RIV / (1 + RIV), equal to lavaan.mi's fmi), and df; the rotated families on the analytic path additionally carry a method string recording the gauge alignment used ("gauge_invariant" for communalities and "signed_permutation_approx" for rotated loadings and, for oblique rotations, factor correlations and structure coefficients). fit_indices_descriptive, on the bootstrap path, pools every per-imputation fit index, so the structural constants among them (df, df_null) appear with a standard error of 0.

mi_fit

On the se = "sandwich" (MI2S) path only: the single EFA() fit on the pooled correlation matrix \bar r and pooled asymptotic covariance \tilde\Gamma. Its orig_R is \bar r and its Gamma is \tilde\Gamma; the pooled SE, CI, and fit_indices are taken from it. MI is NULL on this path because the imputation uncertainty is carried by \tilde\Gamma rather than by per-parameter Rubin pooling.

Conditions

All conditions are classed (prefix efa_pooled_, or efa_consensus_ for the consensus target) so they can be caught by class. The ones most likely to be encountered:

• Inputs. efa_pooled_min_fits (at least two fits are required); efa_pooled_mixed_se (every imputation must use the same se).

• Alignment. efa_consensus_oblique_unsupported (target_method = "consensus" is orthogonal-only).

• Standard errors. efa_pooled_se_unavailable (a warning: pooled SEs could not be produced, so only point estimates are returned); efa_pooled_no_vcov and efa_pooled_unreliable_vcov (the analytic "procrustes" path needs a reliable vcov_unrot_loadings on every fit).

• Two-stage (se = "sandwich"). efa_pooled_mi2s_inputs_inconsistent (every imputation must use se = "sandwich" with the same cor_method); efa_pooled_mi2s_n_too_small (a warning below 20 imputations); efa_pooled_mi2s_acov_not_psd (the pooled covariance is indefinite – use more imputations); efa_pooled_mi2s_alignment_ignored (a warning that the alignment arguments do not apply here).

The remaining conditions concern partial or insufficient bootstrap replicates and unequal sample sizes across imputations.

Author(s)

Andreas Soteriades, Markus Steiner

References

Barnard, J., & Rubin, D. B. (1999). Small-sample degrees of freedom with multiple imputation. Biometrika, 86(4), 948-955.

Bentler, P. M. (1990). Comparative fit indexes in structural models. Psychological Bulletin, 107(2), 238-246.

Browne, M. W. (2001). An overview of analytic rotation in exploratory factor analysis. Multivariate Behavioral Research, 36(1), 111-150.

Chan, K. W., & Meng, X.-L. (2022). Multiple improvements of multiple imputation likelihood ratio tests. Statistica Sinica, 32, 1489-1514.

Chung, S., & Cai, L. (2019). Alternative multiple imputation inference for categorical structural equation modeling. Multivariate Behavioral Research, 54(3), 323-337.

Gower, J. C. (1975). Generalized Procrustes analysis. Psychometrika, 40(1), 33-51.

Li, K. H., Meng, X.-L., Raghunathan, T. E., & Rubin, D. B. (1991). Significance levels from repeated p-values with multiply-imputed data. Statistica Sinica, 1(1), 65-92.

Jennrich, R. I. (1973). Standard errors for obliquely rotated factor loadings. Psychometrika, 38(4), 593-604.

Jennrich, R. I. (1974). Simplified formulae for standard errors in maximum-likelihood factor analysis. British Journal of Mathematical and Statistical Psychology, 27(1), 122-131.

Lorenzo-Seva, U., & Van Ginkel, J. R. (2016). Multiple imputation of missing values in exploratory factor analysis of multidimensional scales. Anales de Psicologia, 32(2), 596-608.

Rubin, D. B. (1987). Multiple imputation for nonresponse in surveys. Wiley.

Schoenemann, P. H. (1966). A generalized solution of the orthogonal Procrustes problem. Psychometrika, 31(1), 1-10.

Sriutaisuk, S., Liu, Y., Chung, S., Kim, H., & Gu, F. (2025). Evaluating imputation-based fit statistics in structural equation modeling with ordinal data: The MI2S approach. Educational and Psychological Measurement, 85(1), 82-113.

Tucker, L. R., & Lewis, C. (1973). A reliability coefficient for maximum likelihood factor analysis. Psychometrika, 38(1), 1-10.

van Ginkel, J. R., & Kroonenberg, P. M. (2014). Using generalized Procrustes analysis for multiple imputation in principal component analysis. Journal of Classification, 31(2), 242-269.

Zhang, G., & Preacher, K. J. (2015). Factor rotation and standard errors in exploratory factor analysis. Journal of Educational and Behavioral Statistics, 40(6), 579-603.

Examples

# create a list of three datasets, mimicking a list you would obtain from
# e.g. mice.
dat_list <- lapply(1:3, function(x) GRiPS_raw[sample(1:nrow(GRiPS_raw), replace = TRUE),])
mod <- EFA_POOLED(dat_list, n_factors = 1, method = "ML")
mod


# add computation of standard errors and CIs
mod <- EFA_POOLED(dat_list, n_factors = 1, method = "ML", se = "np-boot")
mod

Empirical Kaiser Criterion

Description

The empirical Kaiser criterion incorporates random sampling variations of the eigenvalues from the Kaiser-Guttman criterion (KGC(); see Auerswald & Moshagen , 2019; Braeken & van Assen, 2017). The code is based on Braeken & van Assen, (2017) and on Auerswald and Moshagen (2019).

Usage

EKC(
  x,
  N = NA,
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra"),
  type = "BvA2017"
)

Arguments

Details

The Kaiser-Guttman criterion was defined with the intend that a factor should only be extracted if it explains at least as much variance as a single factor (see KGC()). However, this only applies to population-level correlation matrices. Due to sampling variation, the KGC strongly overestimates the number of factors to retrieve (e.g., Zwick & Velicer, 1986). To account for this and to introduce a factor retention method that performs well with small number of indicators and correlated factors (cases where the performance of parallel analysis, see PARALLEL(), is known to deteriorate) Braeken and van Assen (2017) introduced the empirical Kaiser criterion in which a series of reference eigenvalues is created as a function of the variables-to-sample-size ratio and the observed eigenvalues.

Braeken and van Assen (2017) showed that "(a) EKC performs about as well as parallel analysis for data arising from the null, 1-factor, or orthogonal factors model; and (b) clearly outperforms parallel analysis for the specific case of oblique factors, particularly whenever factor intercorrelation is moderate to high and the number of variables per factor is small, which is characteristic of many applications these days" (p.463-464).

In EFAtools version <= 0.5.0 only the implementation of Auerswald and Moshagen (2019) was implemented (now available with type = "AM2019"). However, this implementation, that was probably also used in Caron (2025), differs from the original implementation by Braeken and van Assen (2017) in that it corrects by the reference values, i.e., without using the empirical eigenvalues used in the original implementation. Thanks to Luis Eduardo Garrido for pointing this out and to Johan Braeken for sharing sample code, based on which the original version is now implemented and used by default with type = "BvA2017".

While the adapted version performed relatively well in the simulation studies by Auerswald and Moshagen (2019) and Caron (2025), the theoretical derivations of the EKC as introduced by Braeken and van Assen (2017) may no longer hold. Currently we are unaware of studies comparing the two implementations, but based on our own brief comparisons across multiple datasets, the two implementations appear to often differ substantially regarding the number of factors suggested.

As both implementations exist in different packages and studies, we provide both versions here. Be sure to state clearly which version you use when reporting your results to avoid confusion and ensure reproducibility.

Value

An object of class efa_retention (see print.efa_retention() and plot.efa_retention() for the print and plot methods). Its main fields are:

Source

Auerswald, M., & Moshagen, M. (2019). How to determine the number of factors to retain in exploratory factor analysis: A comparison of extraction methods under realistic conditions. Psychological Methods, 24(4), 468–491. https://doi.org/10.1037/met0000200

Braeken, J., & van Assen, M. A. (2017). An empirical Kaiser criterion. Psychological Methods, 22, 450 – 466. http://dx.doi.org/10.1037/met0000074

Caron, P.-O. (2025). A Comparison of the Next Eigenvalue Sufficiency Test to Other Stopping Rules for the Number of Factors in Factor Analysis. Educational and Psychological Measurement, Online-first publication. https://doi.org/10.1177/00131644241308528

Zwick, W. R., & Velicer, W. F. (1986). Comparison of five rules for determining the number of components to retain. Psychological Bulletin, 99, 432–442. http://dx.doi.org/10.1037/0033-2909.99.3.432

See Also

N_FACTORS() as a wrapper function for this and the other factor retention criteria.

Other factor retention criteria: CD(), HULL(), KGC(), MAP(), NEST(), PARALLEL(), SCREE(), SMT()

Examples

# original implementation
EKC(test_models$baseline$cormat, N = 500)

Estimate factor scores for an EFA model

Description

A convenience wrapper around efa_scores() that returns factor scores and weights in a compact list. Factor scores are calculated according to the specified method if raw data are provided, and only factor weights if a correlation matrix is provided.

Usage

FACTOR_SCORES(
  x,
  f,
  Phi = NULL,
  rho = NULL,
  method = c("Thurstone", "tenBerge", "Anderson", "Bartlett", "Harman", "components")
)

Arguments

Value

A list of class FACTOR_SCORES containing the following:

See Also

efa_scores() for the factor-score weights together with the full set of score-quality diagnostics (determinacy, univocality, and Guttman indeterminacy index) and a print/summary method.

Other factor scoring: efa_scores(), print.efa_scores()

Examples

# Example with raw data with method "Bartlett"
EFA_raw <- EFA(DOSPERT_raw, n_factors = 10, type = "EFAtools", method = "PAF",
               rotation = "oblimin", randomStarts = 0)
fac_scores_raw <- FACTOR_SCORES(DOSPERT_raw, f = EFA_raw, method = "Bartlett")

# Same as above, but with raw data AND a correlation matrix
cor_pearson <- cor(DOSPERT_raw)
EFA_cor_pearson <- EFA(cor_pearson, n_factors = 10, N = nrow(DOSPERT_raw),
                       type = "EFAtools", method = "PAF",
                       rotation = "oblimin", randomStarts = 0)
fac_scores_cor_pearson <- FACTOR_SCORES(DOSPERT_raw, f = EFA_cor_pearson,
                                        rho = cor_pearson,
                                        method = "Bartlett")

# Scores between two alternatives above are identical
isTRUE(all.equal(fac_scores_raw$scores, fac_scores_cor_pearson$scores,
                 check.attributes = FALSE))

# Example with a correlation matrix only (does not return factor scores)
EFA_cor <- EFA(test_models$baseline$cormat, n_factors = 3, N = 500,
               type = "EFAtools", method = "PAF", rotation = "oblimin")
fac_scores_cor <- FACTOR_SCORES(test_models$baseline$cormat, f = EFA_cor)

GRiPS_raw

Description

A data.frame containing responses to the General Risk Propensity Scale (GRiPS, Zhang, Highhouse & Nye, 2018) of 810 participants of Study 1 of Steiner and Frey (2020). The original data can be accessed via https://osf.io/kxp8t/.

Usage

GRiPS_raw

Format

An object of class data.frame with 810 rows and 8 columns.

Source

Zhang, D. C., Highhouse, S., & Nye, C. D. (2018). Development and validation of the general risk propensity scale (GRiPS).Journal of Behavioral Decision Making, 32, 152–167. doi: 10.1002/bdm.2102

Steiner, M., & Frey, R. (2020). Representative design in psychological assessment: A case study using the Balloon Analogue Risk Task (BART). PsyArXiv Preprint. doi:10.31234/osf.io/dg4ks

Hull method for determining the number of factors to retain

Description

Implementation of the Hull method suggested by Lorenzo-Seva, Timmerman, and Kiers (2011), with an extension to principal axis factoring. See details for parallelization.

Usage

HULL(
  x,
  N = NA,
  n_fac_theor = NA,
  method = c("PAF", "ULS", "ML"),
  gof = c("CAF", "CFI", "RMSEA"),
  eigen_type = c("SMC", "PCA", "EFA"),
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra"),
  n_datasets = 1000,
  percent = 95,
  decision_rule = c("means", "percentile", "crawford"),
  n_factors = 1,
  ...
)

Arguments

Details

The Hull method aims to find a model with an optimal balance between model fit and number of parameters, retaining only major factors (Lorenzo-Seva, Timmerman, & Kiers, 2011). It fits 0 to J factors – where J is the number of factors suggested by parallel analysis (or n_fac_theor, if that is larger), plus one – keeps the solutions on the upper boundary of the convex hull of goodness-of-fit against degrees of freedom, and selects the one at the sharpest elbow, i.e. with the highest st value.

The PARALLEL function and the principal axis factoring of the different number of factors can be parallelized using the future framework, by calling the future::plan() function. The examples provide example code on how to enable parallel processing.

Note that if gof = "RMSEA" is used, 1 - RMSEA is actually used to compare the different solutions. Thus, the threshold of .05 is then .95. This is necessary due to how the heuristic to locate the elbow of the hull works.

The ML estimation method uses the psych::fa() starting values. See also the EFA documentation.

Value

An object of class efa_retention (see print.efa_retention() and plot.efa_retention() for the print and plot methods). Its main fields are:

Source

Lorenzo-Seva, U., Timmerman, M. E., & Kiers, H. A. (2011). The Hull method for selecting the number of common factors. Multivariate Behavioral Research, 46(2), 340-364.

See Also

N_FACTORS() as a wrapper function for this and the other factor retention criteria.

Other factor retention criteria: CD(), EKC(), KGC(), MAP(), NEST(), PARALLEL(), SCREE(), SMT()

Examples

# using PAF (this will throw a warning if gof is not specified manually
# and CAF will be used automatically)
HULL(test_models$baseline$cormat, N = 500, gof = "CAF")

# using ML with all available fit indices (CAF, CFI, and RMSEA)
HULL(test_models$baseline$cormat, N = 500, method = "ML")

# using ULS with only RMSEA
HULL(test_models$baseline$cormat, N = 500, method = "ULS", gof = "RMSEA")


## Not run: 
# using parallel processing (Note: plans can be adapted, see the future
# package for details)
future::plan(future::multisession)
HULL(test_models$baseline$cormat, N = 500, gof = "CAF")

## End(Not run)

Intelligence subtests from the Intelligence and Development Scales–2

Description

A matrix containing the bivariate correlations of the 14 intelligence subtests from the Intelligence and Development Scales–2 (IDS-2; Grob & Hagmann-von Arx, 2018), an intelligence and development test battery for children and adolescents aged 5 to 20 years, for the standardization and validation sample (N = 1,991). Details can be found in Grieder & Grob (2019).

Usage

IDS2_R

Format

A 14 x 14 matrix of bivariate correlations

GS

(numeric) - Geometric shapes.

PL

(numeric) - Plates.

TC

(numeric) - Two characteristics.

CB

(numeric) - Crossing out boxes.

NL

(numeric) - Numbers / letters.

NLM

(numeric) - Numbers / letter mixed.

GF

(numeric) - Geometric figures.

RGF

(numeric) - Rotated geometric figures.

CM

(numeric) - Completing matrices.

EP

(numeric) - Excluding pictures.

CA

(numeric) - Categories.

OP

(numeric) - Opposites.

RS

(numeric) - Retelling a story.

DP

(numeric) - Describing pictures.

Source

Grieder, S., & Grob, A. (2019). Exploratory factor analyses of the intelligence and development scales–2: Implications for theory and practice. Assessment. Advance online publication. doi:10.1177/10731911198450

Grob, A., & Hagmann-von Arx, P. (2018). Intelligence and Development Scales–2 (IDS-2). Intelligenz- und Entwicklungsskalen für Kinder und Jugendliche. [Intelligence and Development Scales for Children and Adolescents.]. Bern, Switzerland: Hogrefe.

Kaiser-Guttman Criterion

Description

Probably the most popular factor retention criterion. Kaiser and Guttman suggested to retain as many factors as there are sample eigenvalues greater than 1. This is why the criterion is also known as eigenvalues-greater-than-one rule.

Usage

KGC(
  x,
  eigen_type = c("PCA", "SMC", "EFA"),
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra"),
  n_factors = 1,
  ...
)

Arguments

Details

Originally, the Kaiser-Guttman criterion was intended for the use with principal components, hence with eigenvalues derived from the original correlation matrix. This can be done here by setting eigen_type to "PCA". However, it is well-known that this criterion is often inaccurate and that it tends to overestimate the number of factors, especially for unidimensional or orthogonal factor structures (e.g., Zwick & Velicer, 1986).

The criterion's inaccuracy in these cases is somewhat addressed if it is applied on the correlation matrix with communalities in the diagonal, either initial communalities estimated from SMCs (done setting eigen_type to "SMC") or final communality estimates from an EFA (done setting eigen_type to "EFA"; see Auerswald & Moshagen, 2019). However, although this variant of the KGC is more accurate in some cases compared to the traditional KGC, it is at the same time less accurate than the PCA-variant in other cases, and it is still often less accurate than other factor retention methods, for example parallel analysis (PARALLEL()), the Hull method HULL(), or sequential chi^2 model tests (SMT(); see Auerswald & Moshagen, 2019).

The KGC function can also be called together with other factor retention criteria in the N_FACTORS() function.

Value

An object of class efa_retention (see print.efa_retention() and plot.efa_retention() for the print and plot methods). Its main fields are:

Source

Auerswald, M., & Moshagen, M. (2019). How to determine the number of factors to retain in exploratory factor analysis: A comparison of extraction methods under realistic conditions. Psychological Methods, 24(4), 468–491. https://doi.org/10.1037/met0000200

Guttman, L. (1954). Some necessary conditions for common-factor analysis. Psychometrika, 19, 149 –161. http://dx.doi.org/10.1007/BF02289162

Kaiser, H. F. (1960). The application of electronic computers to factor analysis. Educational and Psychological Measurement, 20, 141–151. http://dx.doi.org/10.1177/001316446002000116

Zwick, W. R., & Velicer, W. F. (1986). Comparison of five rules for determining the number of components to retain. Psychological Bulletin, 99, 432–442. http://dx.doi.org/10.1037/0033-2909.99.3.432

See Also

N_FACTORS() as a wrapper function for this and the other factor retention criteria.

Other factor retention criteria: CD(), EKC(), HULL(), MAP(), NEST(), PARALLEL(), SCREE(), SMT()

Examples

KGC(test_models$baseline$cormat, eigen_type = c("PCA", "SMC"))

Kaiser-Meyer-Olkin criterion

Description

This function computes the Kaiser-Meyer-Olkin (KMO) criterion overall and for each variable in a correlation matrix. The KMO represents the degree to which each observed variable is predicted by the other variables in the dataset and with this indicates the suitability for factor analysis.

Usage

KMO(
  x,
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra")
)

Arguments

Details

Kaiser (1970) proposed this index, originally called measure of sampling adequacy (MSA), that indicates how near the inverted correlation matrix R^{-1} is to a diagonal matrix S to determine a given correlation matrix's (R) suitability for factor analysis. The index is

KMO = \frac{\sum_{i \neq j} r_{ij}^2}{\sum_{i \neq j} r_{ij}^2 + \sum_{i \neq j} q_{ij}^2}

with Q = SR^{-1}S and S = (diag R^{-1})^{-1/2} where \sum_{i \neq j} r_{ij}^2 is the sum of squares of the off-diagonal elements of R and \sum_{i \neq j} q_{ij}^2 is the sum of squares of the off-diagonal elements of Q (see also Cureton & D'Agostino, 1983).

So KMO varies between 0 and 1, with larger values indicating higher suitability for factor analysis. Kaiser and Rice (1974) suggest that KMO should at least exceed .50 for a correlation matrix to be suitable for factor analysis.

This function was heavily influenced by the psych::KMO() function.

See also BARTLETT() for another test of suitability for factor analysis.

The KMO function can also be called together with the BARTLETT() function and with factor retention criteria in the N_FACTORS() function.

Value

A list containing

Source

Kaiser, H. F. (1970). A second generation little jiffy. Psychometrika, 35, 401-415.

Kaiser, H. F. & Rice, J. (1974). Little jiffy, mark IV. Educational and Psychological Measurement, 34, 111-117.

Cureton, E. E. & D'Agostino, R. B. (1983). Factor analysis: An applied approach. Hillsdale, N.J.: Lawrence Erlbaum Associates, Inc.

See Also

BARTLETT() for another measure to determine suitability for factor analysis.

N_FACTORS() as a wrapper function for this function, BARTLETT() and several factor retention criteria.

Examples

KMO(test_models$baseline$cormat)

Velicer's Minimum Average Partial (MAP) Criterion

Description

Computes Velicer's Minimum Average Partial (MAP) criterion for determining the number of factors/components to retain. The function implements the original MAP criterion (Velicer, 1976), expressed via the \mathrm{TR2} representation, and the revised \mathrm{TR4} variant proposed by Velicer, Eaton, and Fava (2000).

Usage

MAP(
  x,
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra")
)

Arguments

Details

MAP partials successive principal components out of the correlation matrix and, after removing m components, summarizes the off-diagonal partial correlations r^*_{ij} that remain; the suggested number of factors is the m that minimizes the criterion. Two criteria are returned:

• TR2 (original MAP; Velicer, 1976):

  \mathrm{MAP}_m = \frac{\sum_{i \neq j} (r^*_{ij})^2}{p(p-1)},

  the mean squared off-diagonal partial correlation.

• TR4 (revised MAP; Velicer, Eaton, & Fava, 2000):

  \mathrm{MAP4}_m = \frac{\sum_{i \neq j} (r^*_{ij})^4}{p(p-1)},

  the mean fourth-power off-diagonal partial correlation, which downweights small partial correlations.

A non-positive-definite input correlation matrix (e.g. from sampling error) is smoothed with psych::cor.smooth().

Value

An object of class efa_retention (see print.efa_retention() for the print method) with the following main elements:

• n_factors: A named numeric vector ("TR2", "TR4") with the index m that minimizes the original (TR2) and revised (TR4) MAP criterion.

• results: A list with one record per criterion, each holding the criterion values over m.

• settings: A list containing use and cor_method.

Source

Velicer, W. F. (1976). Determining the number of components from the matrix of partial correlations. Psychometrika, 41, 321–327.

Velicer, W. F., Eaton, C. A., & Fava, J. L. (2000). Construct explication through factor or component analysis: A review and evaluation of alternative procedures for determining the number of factors or components. In Goffin, R. D. & Helmes, E. (Eds.), Problems and Solutions in Human Assessment: Honoring Douglas N. Jackson at Seventy (pp. 41–71). Boston: Kluwer.

See Also

N_FACTORS() as a wrapper function for this and the other factor retention criteria.

Other factor retention criteria: CD(), EKC(), HULL(), KGC(), NEST(), PARALLEL(), SCREE(), SMT()

Examples

## Example with raw data
res <- MAP(GRiPS_raw)
res

## Example with a correlation matrix
res2 <- MAP(test_models$baseline$cormat)
res2

Next eigenvalue sufficiency test (NEST)

Description

NEST uses many synthetic datasets to generate reference eigenvalues against which to compare the empirical eigenvalues. This is similar to parallel analysis, but other than parallel analysis, NEST does not just rely on synthetic eigenvalues based on an identity matrix as null model. It was introduced by Achim (2017), see also Brandenburg and Papenberg (2024) and Caron (2025) for further simulation studies including NEST.

Usage

NEST(
  x,
  N = NA,
  alpha = 0.05,
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra"),
  n_datasets = 1000,
  ...
)

Arguments

Details

NEST compares the first empirical eigenvalue against the first eigenvalues of n_dataset synthetic datasets based on a null model (i.e., with uncorrelated variables; same as in parallel analysis, see PARALLEL()). The following eigenvalues are compared against synthetic datasets based on an EFA-model with one fewer factors than the position of the respective empirical eigenvalue. E.g, the second empirical eigenvalue is compared against synthetic data based on a one-factor model. In each comparison the k-th empirical eigenvalue is tested against the k-th largest eigenvalue of the synthetic datasets. The alpha-level defines against which percentile of the synthetic eigenvalue distribution to compare the empirical eigenvalues against, i.e., an alpha of .05 (the default) uses the 95th percentile as reference value.

The number of factors tested is capped at \lfloor 0.8 \times p \rfloor (with p the number of variables; Achim, 2017) and additionally limited so that the (k - 1)-factor reference model used at each step stays over-identified. If no empirical eigenvalue falls at or below its reference within this range, every tested factor is accepted and this capped number is returned.

For details on the method, including simulation studies, see Achim (2017), Brandenburg and Papenberg (2024), and Caron (2025).

The NEST function can also be called together with other factor retention criteria in the N_FACTORS() function.

Value

An object of class efa_retention (see print.efa_retention() for the print method). Its main fields are:

Source

Achim, A. (2017). Testing the number of required dimensions in exploratory factor analysis. The Quantitative Methods for Psychology, 13(1), 64–74. https://doi.org/10.20982/tqmp.13.1.p064

Brandenburg, N., & Papenberg, M. (2024). Reassessment of innovative methods to determine the number of factors: A simulation-based comparison of exploratory graph analysis and Next Eigenvalue Sufficiency Test. Psychological Methods, 29(1), 21–47. https://doi.org/10.1037/met0000527

Caron, P.-O. (2025). A Comparison of the Next Eigenvalue Sufficiency Test to Other Stopping Rules for the Number of Factors in Factor Analysis. Educational and Psychological Measurement, Online-first publication. https://doi.org/10.1177/00131644241308528

See Also

N_FACTORS() as a wrapper function for this and the other factor retention criteria.

Other factor retention criteria: CD(), EKC(), HULL(), KGC(), MAP(), PARALLEL(), SCREE(), SMT()

Examples

# with correlation matrix
NEST(test_models$baseline$cormat, N = 500)

# with raw data
NEST(GRiPS_raw)

Various Factor Retention Criteria

Description

Among the most important decisions for an exploratory factor analysis (EFA) is the choice of the number of factors to retain. Several factor retention criteria have been developed for this. With this function, various factor retention criteria can be performed simultaneously. Additionally, the data can be checked for their suitability for factor analysis.

Usage

N_FACTORS(
  x,
  criteria = c("CD", "EKC", "HULL", "MAP", "NEST", "PARALLEL"),
  suitability = TRUE,
  N = NA,
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra"),
  n_factors_max = NA,
  N_pop = 10000,
  N_samples = 500,
  alpha = 0.3,
  max_iter_CD = 50,
  n_fac_theor = NA,
  method = c("ML", "PAF", "ULS"),
  gof = c("CAF", "CFI", "RMSEA"),
  eigen_type_HULL = c("SMC", "PCA", "EFA"),
  eigen_type_other = c("SMC"),
  n_factors = 1,
  n_datasets = 1000,
  percent = 95,
  decision_rule = c("means", "percentile", "crawford"),
  ekc_type = c("BvA2017"),
  n_datasets_nest = 1000,
  alpha_nest = 0.05,
  show_progress = FALSE,
  ...
)

Arguments

Details

By default, the entered data are checked for suitability for factor analysis using the following methods (see respective documentations for details):

• Bartlett's test of sphericity (see BARTLETT())

• Kaiser-Meyer-Olkin criterion (see KMO())

The available factor retention criteria are the following (see respective documentations for details):

• Comparison data (see CD())

• Empirical Kaiser criterion (see EKC())

• Hull method (see HULL())

• Kaiser-Guttman criterion (see KGC())

• Parallel analysis (see PARALLEL())

• Next Eigenvalue Sufficiency Test, NEST (see NEST())

• Scree plot (see SCREE())

• Sequential chi-square model tests, RMSEA lower bound, and AIC (see SMT())

Value

A list of class N_FACTORS containing

Examples

# Default criteria, with correlation matrix and fit method "ML" (where needed)
# This will throw a warning for CD, as no raw data were specified
nfac_all <- N_FACTORS(test_models$baseline$cormat, N = 500, method = "ML")

# The same as above, but without "CD"
nfac_wo_CD <- N_FACTORS(test_models$baseline$cormat, criteria = c("EKC",
                        "HULL", "PARALLEL", "NEST"), N = 500,
                        method = "ML")

# Use PAF instead of ML (this will take longer). For this, gof has
# to be set to "CAF" for the Hull method.
nfac_PAF <- N_FACTORS(test_models$baseline$cormat, criteria = c("EKC",
                        "HULL", "PARALLEL", "NEST"), N = 500,
                      method = "PAF", gof = "CAF")

# Do KGC and PARALLEL with only "PCA" type of eigenvalues
nfac_PCA <- N_FACTORS(test_models$baseline$cormat, criteria = c("EKC",
                      "HULL", "PARALLEL", "NEST"), N = 500,
                      method = "ML", eigen_type_other = "PCA")

# Use raw data, such that CD can also be performed
nfac_raw <- N_FACTORS(GRiPS_raw, method = "ML")

McDonald's omega

Description

This function finds omega total, hierarchical, and subscale, as well as additional model-based indices of interpretive relevance (H index, ECV, PUC) from a Schmid-Leiman (SL) solution or lavaan single factor, second-order (see below), or bifactor solution. The SL-based omegas can either be found from a psych::schmid(), SL(), or, in a more flexible way, by leaving model = NULL and specifying additional arguments. By setting the type argument, results from psych::omega() can be reproduced.

Usage

OMEGA(
  model = NULL,
  type = c("EFAtools", "psych"),
  g_name = "g",
  group_names = NULL,
  add_ind = TRUE,
  factor_corres = NULL,
  var_names = NULL,
  fac_names = NULL,
  g_load = NULL,
  s_load = NULL,
  u2 = NULL,
  cormat = NULL,
  pattern = NULL,
  Phi = NULL,
  variance = c("correlation", "sums_load")
)

Arguments

Details

What this function does

This function calculates McDonald's omegas (McDonald, 1978, 1985, 1999), the H index (Hancock & Mueller, 2001), the explained common variance (ECV; Rodriguez et al., 2016a, 2016b), and the percent of uncontaminated correlations (PUC; Bonifay et al., 2015; Reise et al., 2013).

All types of omegas (total, hierarchical, and subscale) are calculated for the general factor as well as for the subscales / group factors (see, e.g., Gignac, 2014; Rodriguez et al., 2016a, 2016b). Omegas refer to the correlation between a factor and a unit-weighted composite score and thus the true score variance in a unit-weighted composite based on the respective indicators. Omega total is the total true score variance in a composite. Omega hierarchical is the true score variance in a composite that is attributable to the general factor, and omega subscale is the true score variance in a composite attributable to all subscales / group factors (for the whole scale) or to the specific subscale / group factor (for subscale composites).

Accordingly, on a subscale row the hier column reports the share of that subscale's composite variance due to the general factor and the sub column the share due to the subscale-specific factor; the latter corresponds to the omega hierarchical subscale of Rodriguez et al. (2016a, 2016b).

The H index (also construct reliability or replicability index) is the correlation between an optimally-weighted composite score and a factor (Hancock & Mueller, 2001; Rodriguez et al., 2016a, 2016b). It, too, can be calculated for the whole scale / general factor as well as for the subscales / group factors. Low values indicate that a latent variable is not well defined by its indicators.

The ECV (Rodriguez et al., 2016a, 2016b) is the ratio of the variance explained by the general factor and the variance explained by the general factor and the group factors.

The PUC (Bonifay et al., 2015; Reise et al., 2013, Rodriguez et al., 2016a, 2016b) refers to the proportion of correlations in the underlying correlation matrix that is not contaminated by variance of both the general factor and the group factors (i.e., correlations between indicators from different group factors, which reflect only general factor variance). The higher the PUC, the more similar a general factor from a multidimensional model will be to the single factor from a unidimensional model.

How to use this function

If model is a lavaan second-order or bifactor solution, only the name of the general factor from the lavaan model needs to be specified additionally with the g_name argument. It is then determined whether this general factor is a second-order factor (second-order model with one second-order factor assumed) or a breadth factor (bifactor model assumed). Please note that this function only works for second-order models if they contain no more than one second-order factor. In case of a second-order solution, a Schmid-Leiman transformation is performed on the first- and second-order loadings and omega coefficients are obtained from the transformed (orthogonalized) solution (see SL() for more information on Schmid-Leiman transformation). There is also the possibility to enter a lavaan single factor solution. In this case, g_name is not needed. Finally, if a solution from a lavaan multiple group analysis is entered, the indices are computed for each group. For lavaan input the composite variances entering the omegas are model-implied (computed from the fitted loadings and residual variances), so the coefficients coincide with the observed-score versions when the model fits perfectly. The type argument is not evaluated if model is of class lavaan.

If model is of class SL() or psych::schmid() only the type and, depending on the type (see below), the factor_corres arguments need to be specified additionally. If model is of class psych::schmid() and variance = "correlation" (default), it is recommended to also provide the original correlation matrix in cormat to get more accurate results. Otherwise, the correlation matrix will be found based on the pattern matrix and Phi from the psych::schmid() output using the psych::factor.model() function.

If model = NULL, the arguments type, factor_corres (depending on the type, see below), var_names, g_load, s_load, and u2 and either cormat (recommended) or Phi and pattern need to be specified. If Phi and pattern are specified instead of cormat, the correlation matrix is found using the psych::factor.model() function.

The only difference between type = "EFAtools" and type = "psych" is the determination of variable-to-factor correspondences. type = "psych" reproduces the psych::omega() results, where variable-to-factor correspondences are found by taking the highest group factor loading for each variable as the relevant group factor loading. To do this, factor_corres must be left NULL.

The variance argument controls how the total variances (for the whole scale and the subscale composites) are computed. "correlation" (default) finds them from the correlation matrix – the observed-variance form of omega, which reproduces psych::omega(). "sums_load" instead uses the model-implied variances from the squared sums of the general and group factor loadings and the sum of the uniquenesses, so the whole-scale omega total is McDonald's model-implied total and partitions exactly into omega hierarchical plus omega subscale. The two settings agree on the whole-scale omega total and omega hierarchical up to model misfit, and differ mainly in the whole-scale omega subscale, which counts all group-factor variance under "sums_load" but only the assigned subscale composites under "correlation". On the subscale rows the two conventions agree when simple structure is well-achieved.

Value

If found for an SL or lavaan second-order or bifactor solution without multiple groups: A matrix with omegas for the whole scale and for the subscales and (only if add_ind = TRUE) with the H index, ECV, and PUC.

If found for a lavaan single factor solution without multiple groups: A (named) vector with omega total and (if add_ind = TRUE) the H index for the single factor.

If found for a lavaan output from a multiple group analysis: A list containing the output described above for each group.

Source

McDonald, R. P. (1978). Generalizability in factorable domains: ‘‘Domain validity and generalizability’’. Educational and Psychological Measurement, 38, 75–79.

McDonald, R. P. (1985). Factor analysis and related methods. Hillsdale, NJ: Erlbaum.

McDonald, R. P. (1999). Test theory: A unified treatment. Mahwah, NJ: Erlbaum.

Rodriguez, A., Reise, S. P., & Haviland, M. G. (2016a). Applying bifactor statistical indices in the evaluation of psychological measures. Journal of Personality Assessment, 98, 223-237.

Rodriguez, A., Reise, S. P., & Haviland, M. G. (2016b). Evaluating bifactor models: Calculating and interpreting statistical indices. Psychological Methods, 21, 137-150.

Hancock, G. R., & Mueller, R. O. (2001). Rethinking construct reliability within latent variable systems. In R. Cudeck, S. du Toit, & D. Sörbom (Eds.), Structural equation modeling: Present and future—A Festschrift in honor of Karl Jöreskog (pp. 195–216). Lincolnwood, IL: Scientific Software International.

Reise, S. P., Scheines, R., Widaman, K. F., & Haviland, M. G. (2013). Multidimensionality and structural coefficient bias in structural equation modeling: A bifactor perspective. Educational and Psychological Measurement, 73, 5–26.

Bonifay, W. E., Reise, S. P., Scheines, R., & Meijer, R. R. (2015). When are multidimensional data unidimensional enough for structural equation modeling?: An evaluation of the DETECT multidimensionality index. Structural Equation Modeling, 22, 504—516.

Gignac, G. E. (2014). On the Inappropriateness of Using Items to Calculate Total Scale Score Reliability via Coefficient Alpha for Multidimensional Scales. European Journal of Psychological Assessment, 30, 130-139.

Examples

## Use with lavaan outputs
if (requireNamespace("lavaan", quietly = TRUE)) {

# Create and fit bifactor model in lavaan (assume all variables have SDs of 1)
mod <- 'F1 =~ V1 + V2 + V3 + V4 + V5 + V6
        F2 =~ V7 + V8 + V9 + V10 + V11 + V12
        F3 =~ V13 + V14 + V15 + V16 + V17 + V18
        g =~ V1 + V2 + V3 + V4 + V5 + V6 + V7 + V8 + V9 + V10 + V11 + V12 +
             V13 + V14 + V15 + V16 + V17 + V18'
fit_bi <- lavaan::cfa(mod, sample.cov = test_models$baseline$cormat,
                      sample.nobs = 500, estimator = "ml", orthogonal = TRUE)

# Compute omegas and additional indices for bifactor solution
OMEGA(fit_bi, g_name = "g")

# Compute only omegas
OMEGA(fit_bi, g_name = "g", add_ind = FALSE)

# Create and fit second-order model in lavaan (assume all variables have SDs of 1)
mod <- 'F1 =~ V1 + V2 + V3 + V4 + V5 + V6
        F2 =~ V7 + V8 + V9 + V10 + V11 + V12
        F3 =~ V13 + V14 + V15 + V16 + V17 + V18
        g =~ F1 + F2 + F3'
fit_ho <- lavaan::cfa(mod, sample.cov = test_models$baseline$cormat,
                      sample.nobs = 500, estimator = "ml")

# Compute omegas and additional indices for second-order solution
OMEGA(fit_ho, g_name = "g")
}


## Use with an output from the SL function, with type EFAtools
efa_mod <- EFA(test_models$baseline$cormat, N = 500, n_factors = 3,
               type = "EFAtools", method = "PAF", rotation = "promax")
sl_mod <- SL(efa_mod, type = "EFAtools", method = "PAF")

# Two examples how to specify the indicator-to-factor correspondences:

# Based on a specific salience threshold for the loadings (here: .20):
factor_corres_1 <- sl_mod$sl[, c("F1", "F2", "F3")] >= .2

# Or more flexibly (could also be TRUE and FALSE instead of 0 and 1):
factor_corres_2 <- matrix(c(rep(0, 12), rep(1, 6), rep(0, 6), rep(1, 6),
                         rep(0, 6), rep(1, 6), rep(0, 12)), ncol = 3,
                         byrow = FALSE)

OMEGA(sl_mod, type = "EFAtools", factor_corres = factor_corres_1)

## Use with an output from the psych::schmid function, with type psych for
## OMEGA
schmid_mod <- psych::schmid(test_models$baseline$cormat, nfactors = 3,
                            n.obs = 500, fm = "pa", rotate = "Promax")
# Find correlation matrix from phi and pattern matrix from psych::schmid output
OMEGA(schmid_mod, type = "psych")
# Use specified correlation matrix
OMEGA(schmid_mod, type = "psych", cormat = test_models$baseline$cormat)

## Manually specify components (useful if omegas should be computed for a SL
## or bifactor solution found with another program)
## As an example, we extract the elements from an SL output here. This gives
## the same results as in the second example above.

efa_mod <- EFA(test_models$baseline$cormat, N = 500, n_factors = 3,
               type = "EFAtools", method = "PAF", rotation = "promax")
sl_mod <- SL(efa_mod, type = "EFAtools", method = "PAF")

factor_corres <- matrix(c(rep(0, 12), rep(1, 6), rep(0, 6), rep(1, 6),
                        rep(0, 6), rep(1, 6), rep(0, 12)), ncol = 3,
                        byrow = FALSE)

OMEGA(model = NULL, type = "EFAtools", var_names = rownames(sl_mod$sl),
      g_load = sl_mod$sl[, "g"], s_load = sl_mod$sl[, c("F1", "F2", "F3")],
      u2 = sl_mod$sl[, "u2"], cormat = test_models$baseline$cormat,
      factor_corres = factor_corres)

Parallel analysis

Description

Various methods for performing parallel analysis. This function uses future_lapply() for which a parallel processing plan can be selected. To do so, call library(future) and, for example, plan(multisession); see examples.

Usage

PARALLEL(
  x = NULL,
  N = NA,
  n_vars = NA,
  n_datasets = 1000,
  percent = 95,
  eigen_type = c("PCA", "SMC", "EFA"),
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra"),
  decision_rule = c("means", "percentile", "crawford"),
  n_factors = 1,
  ...
)

Arguments

Details

Parallel analysis (Horn, 1965) compares the eigenvalues obtained from the sample correlation matrix against those of null model correlation matrices (i.e., with uncorrelated variables) of the same sample size. This way, it accounts for the variation in eigenvalues introduced by sampling error and thus eliminates the main problem inherent in the Kaiser-Guttman criterion (KGC()).

Parallel analysis is often argued to be one of the most accurate factor retention criteria. However, for highly correlated factor structures it has been shown to underestimate the correct number of factors. The reason for this is that a null model (uncorrelated variables) is used as reference. However, when factors are highly correlated, the first eigenvalue will be much larger compared to the following ones, as later eigenvalues are conditional on the earlier ones in the sequence and thus the shared variance is already accounted in the first eigenvalue (e.g., Braeken & van Assen, 2017).

The PARALLEL function can also be called together with other factor retention criteria in the N_FACTORS() function.

Value

An object of class efa_retention (see print.efa_retention() and plot.efa_retention() for the print and plot methods). Its main fields are:

Source

Braeken, J., & van Assen, M. A. (2017). An empirical Kaiser criterion. Psychological Methods, 22, 450 – 466. http://dx.doi.org/10.1037/ met0000074

Crawford, A. V., Green, S. B., Levy, R., Lo, W. J., Scott, L., Svetina, D., & Thompson, M. S. (2010). Evaluation of parallel analysis methods for determining the number of factors. Educational and Psychological Measurement, 70(6), 885-901.

Glorfeld, L. W. (1995). An improvement on Horn's parallel analysis methodology for selecting the correct number of factors to retain. Educational and Psychological Measurement, 55(3), 377-393.

Horn, J. L. (1965). A rationale and test for the number of factors in factor analysis. Psychometrika, 30(2), 179–185. doi: 10.1007/BF02289447

See Also

N_FACTORS() as a wrapper function for this and the other factor retention criteria.

Other factor retention criteria: CD(), EKC(), HULL(), KGC(), MAP(), NEST(), SCREE(), SMT()

Examples

# example without real data
pa_unreal <- PARALLEL(N = 500, n_vars = 10)

# example with correlation matrix with all eigen_types and PAF estimation
pa_paf <- PARALLEL(test_models$case_11b$cormat, N = 500)

# example with correlation matrix with all eigen_types and ML estimation
# this will be faster than the above with PAF)
pa_ml <- PARALLEL(test_models$case_11b$cormat, N = 500, method = "ML")


## Not run: 
# for parallel computation
future::plan(future::multisession)
pa_faster <- PARALLEL(test_models$case_11b$cormat, N = 500)

## End(Not run)

Rotate a loading matrix to a target using Procrustes alignment

Description

PROCRUSTES() aligns one loading matrix to a target loading matrix with the same dimensions. It is used internally by EFA_POOLED(), but can also be used directly when factor columns must be brought into a common orientation before averaging or comparing solutions.

Usage

PROCRUSTES(
  A,
  Target,
  rotation = c("orthogonal", "oblique"),
  S = NULL,
  T_init = NULL,
  oblique_eps = 1e-05,
  oblique_maxit = 1000,
  oblique_max_line_search = 10,
  oblique_step0 = 1,
  oblique_normalize = FALSE,
  oblique_random_starts = 0,
  oblique_screen_keep = 2,
  oblique_triage_maxit = 25,
  oblique_triage_improve_tol = 0
)

Arguments

Details

For rotation = "orthogonal", the function solves the closed-form orthogonal Procrustes problem

\min_T \frac{1}{2}\|A T - B\|_F^2 \quad \textrm{subject to}\quad T'T = I,

where A is the loading matrix and B is Target.

For rotation = "oblique", the function calls the compiled .oblique_procrustes() optimizer. The oblique convention is the same as in GPArotation::targetQ():

L = A T^{-T}, \qquad \Phi = T'T, \qquad diag(\Phi) = 1.

By default the oblique solver is warm-started from the closed-form orthogonal Procrustes solution, which resolves the factor permutation and sign indeterminacy and avoids the poor local minima an identity start can fall into. Supply T_init to override this start. Random starts are only used for oblique alignment. For one-factor models, oblique and orthogonal alignment are equivalent, so the function uses the stable one-factor orthogonal solution instead of calling the oblique optimizer.

Value

A list containing aligned loadings, transformation matrix T, factor intercorrelation matrix Phi, target criterion value, convergence diagnostics, line-search diagnostics, and multi-start summaries. Row and column names are preserved where possible. When oblique_normalize = TRUE the returned loadings are back-transformed to the original scale, but value is the criterion on the Kaiser-normalized loadings, so it is not 0.5 * sum((loadings - Target)^2).

RiskDimensions

Description

A list containing the bivariate correlations (cormat) of the 9 dimensions on which participants in Fischhoff et al. (1978) rated different activities and technologies as well as the sample size (N). This was then analyzed together with ratings of the risks and benefits of these activities and technologies.

Usage

RiskDimensions

Format

An object of class list of length 2.

Source

Fischhoff, B, Slovic, P, Lichtenstein, S, Read, S, and Combs, B. (1978). How safe is safe enough? A psychometric study of attitudes towards technological risks and benefits. Policy Sciences, 9, 127-152. doi: 10.1007/BF00143739

Scree Plot

Description

The scree plot was originally introduced by Cattell (1966) to perform the scree test. In a scree plot, the eigenvalues of the factors / components are plotted against the index of the factors / components, ordered from 1 to N factors components, hence from largest to smallest eigenvalue. According to the scree test, the number of factors / components to retain is the number of factors / components to the left of the "elbow" (where the curve starts to level off) in the scree plot.

Usage

SCREE(
  x,
  eigen_type = c("PCA", "SMC", "EFA"),
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra"),
  n_factors = 1,
  ...
)

Arguments

Details

As the scree test requires visual examination, the test has been especially criticized for its subjectivity and with this low inter-rater reliability. Moreover, a scree plot can be ambiguous if there are either no clear "elbow" or multiple "elbows", making it difficult to judge just where the eigenvalues do level off. Finally, the scree test has also been found to be less accurate than other factor retention criteria. For all these reasons, the scree test has been recommended against, at least for exclusive use as a factor retention criterion (Zwick & Velicer, 1986)

The SCREE function can also be called together with other factor retention criteria in the N_FACTORS() function.

Value

An object of class efa_retention (see print.efa_retention() and plot.efa_retention() for the print and plot methods). The scree plot is a visual criterion, so it returns no numeric suggestion. Its main fields are:

Source

Cattell, R. B. (1966). The scree test for the number of factors. Multivariate Behavioral Research, 1(2), 245–276. https://doi.org/10.1207/s15327906mbr0102_10

Zwick, W. R., & Velicer, W. F. (1986). Comparison of five rules for determining the number of components to retain. Psychological Bulletin, 99, 432–442. http://dx.doi.org/10.1037/0033-2909.99.3.432

See Also

N_FACTORS() as a wrapper function for this and the other factor retention criteria.

Other factor retention criteria: CD(), EKC(), HULL(), KGC(), MAP(), NEST(), PARALLEL(), SMT()

Examples

SCREE(test_models$baseline$cormat, eigen_type = c("PCA", "SMC"))

Schmid-Leiman Transformation

Description

This function implements the Schmid-Leiman (SL) transformation (Schmid & Leiman, 1957). It takes the pattern coefficients and factor intercorrelations from an oblique factor solution as input and can reproduce the results from psych::schmid() and from the SPSS implementation from Wolff & Preising (2005). Other arguments from EFA() can be used to control the procedure to find the second-order loadings more flexibly. The function can also be used on a second-order confirmatory factor analysis (CFA) solution from lavaan.

Usage

SL(
  x,
  Phi = NULL,
  type = c("EFAtools", "psych", "SPSS", "none"),
  method = c("PAF", "ML", "ULS", "MINRES"),
  g_name = "g",
  ...
)

Arguments

Details

The SL transformation (also called SL orthogonalization) is a procedure with which an oblique factor solution is transformed into a hierarchical, orthogonalized solution. As a first step, the factor intercorrelations are factor analyzed to extract a single second-order (general) factor, yielding a two-level hierarchical structure. The first-order factor and the second-order factor are then orthogonalized, resulting in an orthogonalized factor solution with proportionality constraints. The procedure thus makes a suggested hierarchical data structure based on factor intercorrelations explicit. One major advantage of SL transformation is that it enables variance partitioning between higher-order and first-order factors, including the calculation of McDonald's omegas (see OMEGA()).

Value

A list of class SL containing the following

Source

Schmid, J. & Leiman, J. M. (1957). The development of hierarchical factor solutions. Psychometrika, 22(1), 53–61. doi:10.1007/BF02289209

Wolff, H.-G., & Preising, K. (2005). Exploring item and higher order factor structure with the Schmid-Leiman solution: Syntax codes for SPSS and SAS. Behavior Research Methods, 37 , 48–58. doi:10.3758/BF03206397

Examples

## Use with an output from the EFAtools::EFA function, both with type EFAtools
EFA_mod <- EFA(test_models$baseline$cormat, N = 500, n_factors = 3,
               type = "EFAtools", method = "PAF", rotation = "promax")
SL_EFAtools <- SL(EFA_mod, type = "EFAtools", method = "PAF")


## Use with an output from the psych::fa function with type psych in SL
fa_mod <- psych::fa(test_models$baseline$cormat, nfactors = 3, n.obs = 500,
                    fm = "pa", rotate = "Promax")
SL_psych <- SL(fa_mod, type = "psych", method = "PAF")


## Use more flexibly by entering a pattern matrix and phi directly (useful if
## a factor solution found with another program should be subjected to SL
## transformation)

## For demonstration, take pattern matrix and phi from an EFA output
## This gives the same solution as the first example
EFA_mod <- EFA(test_models$baseline$cormat, N = 500, n_factors = 3,
               type = "EFAtools", method = "PAF", rotation = "promax")
SL_flex <- SL(EFA_mod$rot_loadings, Phi = EFA_mod$Phi, type = "EFAtools",
              method = "PAF")


## Use with a lavaan second-order CFA output
if (requireNamespace("lavaan", quietly = TRUE)) {

# Create and fit model in lavaan (assume all variables have SDs of 1)
mod <- 'F1 =~ V1 + V2 + V3 + V4 + V5 + V6
        F2 =~ V7 + V8 + V9 + V10 + V11 + V12
        F3 =~ V13 + V14 + V15 + V16 + V17 + V18
        g =~ F1 + F2 + F3'
fit <- lavaan::cfa(mod, sample.cov = test_models$baseline$cormat,
                   sample.nobs = 500, estimator = "ml")

SL_lav <- SL(fit, g_name = "g")

}

Sequential Chi Square Model Tests, RMSEA lower bound, and AIC

Description

Sequential Chi Square Model Tests (SMT) are a factor retention method where multiple EFAs with increasing numbers of factors are fitted and the number of factors for which the Chi Square value first becomes non-significant is taken as the suggested number of factors. Preacher, Zhang, Kim, & Mels (2013) suggested a similar approach with the lower bound of the 90% confidence interval of the Root Mean Square Error of Approximation (RMSEA; Browne & Cudeck, 1992; Steiger & Lind, 1980), and with the Akaike Information Criterion (AIC). For the RMSEA, the number of factors for which this lower bound first falls below .05 is the suggested number of factors to retain. For the AIC, it is the number of factors where the AIC is lowest.

Usage

SMT(
  x,
  N = NA,
  use = c("pairwise.complete.obs", "all.obs", "complete.obs", "everything",
    "na.or.complete"),
  cor_method = c("pearson", "spearman", "kendall", "poly", "tetra")
)

Arguments

Details

As a first step in the procedure, a maximum number of factors to extract is determined for which the model is still over-identified (df > 0).

Then, EFAs with increasing numbers of factors from 1 to the maximum number are fitted with maximum likelihood estimation.

For the SMT, first the significance of the chi square value for a model with 0 factors is determined. If this value is not significant, 0 factors are suggested to retain. If it is significant, a model with 1 factor is estimated and the significance of its chi square value is determined, and so on, until a non-significant result is obtained. The suggested number of factors is the number of factors for the model where the chi square value first becomes non-significant.

Regarding the RMSEA, the suggested number of factors is the number of factors for the model where the lower bound of the 90% confidence interval of the RMSEA first falls below the .05 threshold.

Regarding the AIC, the suggested number of factors is the number of factors for the model with the lowest AIC.

In comparison with other prominent factor retention criteria, SMT performed well at determining the number of factors to extract in EFA (Auerswald & Moshagen, 2019). The RMSEA lower bound also performed well at determining the true number of factors, while the AIC performed well at determining the most generalizable model (Preacher, Zhang, Kim, & Mels, 2013).

Value

An object of class efa_retention (see print.efa_retention() for the print method). Its main fields are:

Source

Auerswald, M., & Moshagen, M. (2019). How to determine the number of factors to retain in exploratory factor analysis: A comparison of extraction methods under realistic conditions. Psychological Methods, 24(4), 468–491. https://doi.org/10.1037/met0000200

Browne, M.W., & Cudeck, R. (1992). Alternative ways of assessing model fit. Sociological Methods and Research, 21, 230–258.

Preacher, K. J., Zhang G., Kim, C., & Mels, G. (2013). Choosing the Optimal Number of Factors in Exploratory Factor Analysis: A Model Selection Perspective, Multivariate Behavioral Research, 48(1), 28-56, doi:10.1080/00273171.2012.710386

Steiger, J. H., & Lind, J. C. (1980, May). Statistically based tests for the number of common factors. Paper presented at the annual meeting of the Psychometric Society, Iowa City, IA.

See Also

N_FACTORS() as a wrapper function for this and the other factor retention criteria.

Other factor retention criteria: CD(), EKC(), HULL(), KGC(), MAP(), NEST(), PARALLEL(), SCREE()

Examples

SMT_base <- SMT(test_models$baseline$cormat, N = 500)
SMT_base

Various outputs from SPSS (version 23) FACTOR

Description

Various outputs from SPSS (version 23) FACTOR for the IDS-2 (Grob & Hagmann-von Arx, 2018), the WJIV (3 to 5 and 20 to 39 years; McGrew, LaForte, & Schrank, 2014), the DOSPERT (Frey et al., 2017; Weber, Blais, & Betz, 2002), the NEO-PI-R (Costa, & McCrae, 1992), and four simulated datasets (baseline, case_1a, case_6b, and case_11b, see test_models and population_models) used in Grieder and Steiner (2022).

Usage

SPSS_23

Format

A list of 9 containing EFA results for each of the data sets mentioned above. Each of these nine entries is a list of 4 or 8 (see details), of the following structure:

paf_comm

(vector) - The final communalities obtained with the FACTOR algorithm with PAF and no rotation. For details, see Grieder and Grob (2019).

paf_load

(matrix) - F1 to FN = unrotated factor loadings obtained with the FACTOR algorithm with PAF. Rownames are the abbreviated subtest names.

paf_iter

(numeric) - Number of iterations needed for the principal axis factoring to converge.

var_load

(matrix) - F1 to FN = varimax rotated factor loadings obtained with the FACTOR algorithm with PAF. Rownames are the abbreviated subtest names.

pro_load

(matrix) - F1 to FN = promax rotated factor loadings obtained with the FACTOR algorithm with PAF. Rownames are the abbreviated subtest names.

pro_phi

(matrix) - F1 to FN = intercorrelations of the promax rotated loadings.

sl

(matrix) - g = General / second order factor of the Schmid-Leiman solution. F1 to FN = First order factors of the Schmid-Leiman solution. h2 = Communalities of the Schmid-Leiman solution. This Schmid-Leiman solution was found using the SPSS Syntax provided by Wolff and Preising (2005).

L2

(matrix) - Second order loadings used for the Schmid-Leiman transformation. This Schmid-Leiman solution was found using the SPSS Syntax provided by Wolff and Preising (2005).

Details

The IDS-2, the two WJIV, the DOSPERT, and the NEO-PI-R contain all the above entries, while the four simulated datasets contain only paf_load, var_load, pro_load, and pro_phi.

Source

Grieder, S., & Steiner, M. D. (2022). Algorithmic jingle jungle: A comparison of implementations of principal axis factoring and promax rotation in R and SPSS. Behavior Research Methods, 54, 54–74. doi: 10.3758/s13428-021-01581-x

Wolff, H.G., & Preising, K. (2005). Exploring item and higher order factor structure with the Schmid-Leiman solution: Syntax codes for SPSS and SAS. Behavior Research Methods, 37, 48–58. doi: 10.3758/BF03206397

Grieder, S., & Grob, A. (2019). Exploratory factor analyses of the intelligence and development scales–2: Implications for theory and practice. Assessment. Advance online publication. doi:10.1177/10731911198450

Grob, A., & Hagmann-von Arx, P. (2018). Intelligence and Development Scales–2 (IDS-2). Intelligenz- und Entwicklungsskalen für Kinder und Jugendliche. [Intelligence and Development Scales for Children and Adolescents.]. Bern, Switzerland: Hogrefe.

Frey, R., Pedroni, A., Mata, R., Rieskamp, J., & Hertwig, R. (2017). Risk preference shares the psychometric structure of major psychological traits. Science Advances, 3, e1701381.

McGrew, K. S., LaForte, E. M., & Schrank, F. A. (2014). Technical Manual. Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Schrank, F. A., McGrew, K. S., & Mather, N. (2014). Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Costa, P. T., & McCrae, R. R. (1992). NEO PI-R professional manual. Odessa, FL: Psychological Assessment Resources, Inc.

Various outputs from SPSS (version 27) FACTOR

Description

Various outputs from SPSS (version 27) FACTOR for the IDS-2 (Grob & Hagmann-von Arx, 2018), the WJIV (3 to 5 and 20 to 39 years; McGrew, LaForte, & Schrank, 2014), the DOSPERT (Frey et al., 2017; Weber, Blais, & Betz, 2002), the NEO-PI-R (Costa, & McCrae, 1992), and four simulated datasets (baseline, case_1a, case_6b, and case_11b, see test_models and population_models) used in Grieder and Steiner (2022).

Usage

SPSS_27

Format

A list of 9 containing EFA results for each of the data sets mentioned above. Each of these nine entries is a list of 4 or 8 (see details), of the following structure:

paf_comm

(vector) - The final communalities obtained with the FACTOR algorithm with PAF and no rotation. For details, see Grieder and Grob (2019).

paf_load

(matrix) - F1 to FN = unrotated factor loadings obtained with the FACTOR algorithm with PAF. Rownames are the abbreviated subtest names.

paf_iter

(numeric) - Number of iterations needed for the principal axis factoring to converge.

var_load

(matrix) - F1 to FN = varimax rotated factor loadings obtained with the FACTOR algorithm with PAF. Rownames are the abbreviated subtest names.

pro_load

(matrix) - F1 to FN = promax rotated factor loadings obtained with the FACTOR algorithm with PAF. Rownames are the abbreviated subtest names.

pro_phi

(matrix) - F1 to FN = intercorrelations of the promax rotated loadings.

sl

(matrix) - g = General / second order factor of the Schmid-Leiman solution. F1 to FN = First order factors of the Schmid-Leiman solution. h2 = Communalities of the Schmid-Leiman solution. This Schmid-Leiman solution was found using the SPSS Syntax provided by Wolff and Preising (2005).

L2

(matrix) - Second order loadings used for the Schmid-Leiman transformation. This Schmid-Leiman solution was found using the SPSS Syntax provided by Wolff and Preising (2005).

Details

The IDS-2, the two WJIV, the DOSPERT, and the NEO-PI-R contain all the above entries, while the four simulated datasets contain only paf_load, var_load, pro_load, and pro_phi.

Source

Grieder, S., & Steiner, M. D. (2022). Algorithmic jingle jungle: A comparison of implementations of principal axis factoring and promax rotation in R and SPSS. Behavior Research Methods, 54, 54–74. doi: 10.3758/s13428-021-01581-x

Wolff, H.G., & Preising, K. (2005). Exploring item and higher order factor structure with the Schmid-Leiman solution: Syntax codes for SPSS and SAS. Behavior Research Methods, 37, 48–58. doi: 10.3758/BF03206397

Grieder, S., & Grob, A. (2019). Exploratory factor analyses of the intelligence and development scales–2: Implications for theory and practice. Assessment. Advance online publication. doi:10.1177/10731911198450

Grob, A., & Hagmann-von Arx, P. (2018). Intelligence and Development Scales–2 (IDS-2). Intelligenz- und Entwicklungsskalen für Kinder und Jugendliche. [Intelligence and Development Scales for Children and Adolescents.]. Bern, Switzerland: Hogrefe.

Frey, R., Pedroni, A., Mata, R., Rieskamp, J., & Hertwig, R. (2017). Risk preference shares the psychometric structure of major psychological traits. Science Advances, 3, e1701381.

McGrew, K. S., LaForte, E. M., & Schrank, F. A. (2014). Technical Manual. Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Schrank, F. A., McGrew, K. S., & Mather, N. (2014). Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Costa, P. T., & McCrae, R. R. (1992). NEO PI-R professional manual. Odessa, FL: Psychological Assessment Resources, Inc.

UPPS_raw

Description

A dataframe containing responses to the UPPS personality scale (Whiteside & Lynam, 2005) of 645 participants of Study 2 of Steiner and Frey (2020). Each column are the ratings to one of 45 items to assess urgency, premeditation, perseverance, and sensation seeking. The original data can be accessed via https://osf.io/kxp8t/.

Usage

UPPS_raw

Format

An object of class data.frame with 645 rows and 45 columns.

Source

Whiteside, S. P., Lynam, D. R., Miller, J. D., & Reynolds, S. K. (2005). Validation of the UPPS impulsive behaviour scale: A four-factor model of impulsivity. European Journal of Personality, 19 (7), 559–574.

Steiner, M., & Frey, R. (2020). Representative design in psychological assessment: A case study using the Balloon Analogue Risk Task (BART). PsyArXiv Preprint. doi:10.31234/osf.io/dg4ks

Woodcock Johnson IV: ages 14 to 19

Description

A list containing the bivariate correlations (N = 1,685) of the 47 cognitive and achievement subtests from the WJ IV for 14- to 19-year-olds from the standardization sample obtained from the WJ-IV technical manual (McGrew, LaForte, & Schrank, 2014). Tables are reproduced with permission from the publisher.

Usage

WJIV_ages_14_19

Format

A list of 2 with elements "cormat" (47 x 47 matrix of bivariate correlations) and "N" (scalar). The correlation matrix contains the following variables:

ORLVOC

(numeric) - Oral Vocabulary.

NUMSER

(numeric) - Number Series.

VRBATN

(numeric) - Verbal Attention.

LETPAT

(numeric) - Letter-Pattern Matching.

PHNPRO

(numeric) - Phonological Processing.

STYREC

(numeric) - Story Recall.

VISUAL

(numeric) - Visualization.

GENINF

(numeric) - General Information.

CONFRM

(numeric) - Concept Formation.

NUMREV

(numeric) - Numbers Reversed.

NUMPAT

(numeric) - Number-Pattern Matching.

NWDREP

(numeric) - Nonword Repetition.

VAL

(numeric) - Visual-Auditory Learning.

PICREC

(numeric) - Picture Recognition.

ANLSYN

(numeric) - Analysis-Synthesis.

OBJNUM

(numeric) - Object-Number Sequencing.

PAIRCN

(numeric) - Pair Cancellation.

MEMWRD

(numeric) - Memory for Words.

PICVOC

(numeric) - Picture Vocabulary.

ORLCMP

(numeric) - Oral Comprehension.

SEGMNT

(numeric) - Segmentation.

RPCNAM

(numeric) - Rapid Picture Naming.

SENREP

(numeric) - Sentence Repetition.

UNDDIR

(numeric) - Understanding Directions.

SNDBLN

(numeric) - Sound Blending.

RETFLU

(numeric) - Retrieval Fluency.

SNDAWR

(numeric) - Sound Awareness.

LWIDNT

(numeric) - Letter-Word Identification.

APPROB

(numeric) - Applied Problems.

SPELL

(numeric) - Spelling.

PSGCMP

(numeric) - Passage Comprehension.

CALC

(numeric) - Calculation.

WRTSMP

(numeric) - Writing Samples.

WRDATK

(numeric) - Word Attack.

ORLRDG

(numeric) - Oral Reading.

SNRDFL

(numeric) - Sentence Reading Fluency.

MTHFLU

(numeric) - Math Facts Fluency.

SNWRFL

(numeric) - Sentence Writing Fluency.

RDGREC

(numeric) - Reading Recall.

NUMMAT

(numeric) - Number Matrices.

EDIT

(numeric) - Editing.

WRDFLU

(numeric) - Word Reading Fluency.

SPLSND

(numeric) - Spelling of Sounds.

RDGVOC

(numeric) - Reading Vocabulary.

SCI

(numeric) - Science.

SOC

(numeric) - Social Studies.

HUM

(numeric) - Humanities.

Source

McGrew, K. S., LaForte, E. M., & Schrank, F. A. (2014). Technical Manual. Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Schrank, F. A., McGrew, K. S., & Mather, N. (2014). Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Woodcock Johnson IV: ages 20 to 39

Description

A list containing the bivariate correlations (N = 1,251) of the 47 cognitive and achievement subtests from the WJ IV for the 20- to 39-year-olds from the standardization sample obtained from the WJ-IV technical manual (McGrew, LaForte, & Schrank, 2014). Tables are reproduced with permission from the publisher.

Usage

WJIV_ages_20_39

Format

A list of 2 with elements "cormat" (47 x 47 matrix of bivariate correlations) and "N" (scalar). The correlation matrix contains the following variables:

ORLVOC

(numeric) - Oral Vocabulary.

NUMSER

(numeric) - Number Series.

VRBATN

(numeric) - Verbal Attention.

LETPAT

(numeric) - Letter-Pattern Matching.

PHNPRO

(numeric) - Phonological Processing.

STYREC

(numeric) - Story Recall.

VISUAL

(numeric) - Visualization.

GENINF

(numeric) - General Information.

CONFRM

(numeric) - Concept Formation.

NUMREV

(numeric) - Numbers Reversed.

NUMPAT

(numeric) - Number-Pattern Matching.

NWDREP

(numeric) - Nonword Repetition.

VAL

(numeric) - Visual-Auditory Learning.

PICREC

(numeric) - Picture Recognition.

ANLSYN

(numeric) - Analysis-Synthesis.

OBJNUM

(numeric) - Object-Number Sequencing.

PAIRCN

(numeric) - Pair Cancellation.

MEMWRD

(numeric) - Memory for Words.

PICVOC

(numeric) - Picture Vocabulary.

ORLCMP

(numeric) - Oral Comprehension.

SEGMNT

(numeric) - Segmentation.

RPCNAM

(numeric) - Rapid Picture Naming.

SENREP

(numeric) - Sentence Repetition.

UNDDIR

(numeric) - Understanding Directions.

SNDBLN

(numeric) - Sound Blending.

RETFLU

(numeric) - Retrieval Fluency.

SNDAWR

(numeric) - Sound Awareness.

LWIDNT

(numeric) - Letter-Word Identification.

APPROB

(numeric) - Applied Problems.

SPELL

(numeric) - Spelling.

PSGCMP

(numeric) - Passage Comprehension.

CALC

(numeric) - Calculation.

WRTSMP

(numeric) - Writing Samples.

WRDATK

(numeric) - Word Attack.

ORLRDG

(numeric) - Oral Reading.

SNRDFL

(numeric) - Sentence Reading Fluency.

MTHFLU

(numeric) - Math Facts Fluency.

SNWRFL

(numeric) - Sentence Writing Fluency.

RDGREC

(numeric) - Reading Recall.

NUMMAT

(numeric) - Number Matrices.

EDIT

(numeric) - Editing.

WRDFLU

(numeric) - Word Reading Fluency.

SPLSND

(numeric) - Spelling of Sounds.

RDGVOC

(numeric) - Reading Vocabulary.

SCI

(numeric) - Science.

SOC

(numeric) - Social Studies.

HUM

(numeric) - Humanities.

Source

McGrew, K. S., LaForte, E. M., & Schrank, F. A. (2014). Technical Manual. Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Schrank, F. A., McGrew, K. S., & Mather, N. (2014). Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Woodcock Johnson IV: ages 3 to 5

Description

A list containing the bivariate correlations (N = 435) of the 29 cognitive and achievement subtests from the WJ IV for 3- to 5-year-olds from the standardization sample obtained from the WJ IV technical Manual (McGrew, LaForte, & Schrank, 2014). Tables are reproduced with permission from the publisher.

Usage

WJIV_ages_3_5

Format

A list of 2 with elements "cormat" (29 x 29 matrix of bivariate correlations) and "N" (scalar). The correlation matrix contains the following variables:

ORLVOC

(numeric) - Oral Vocabulary.

VRBATN

(numeric) - Verbal Attention.

LETPAT

(numeric) - Phonological Processing.

STYREC

(numeric) - Story Recall.

VISUAL

(numeric) - Visualization.

GENINF

(numeric) - General Information.

CONFRM

(numeric) - Concept Formation.

NUMREV

(numeric) - Numbers Reversed.

NUMPAT

(numeric) - Number-Pattern Matching.

NWDREP

(numeric) - Nonword Repetition.

VAL

(numeric) - Visual-Auditory Learning.

PICREC

(numeric) - Picture Recognition.

MEMWRD

(numeric) - Memory for Words.

PICVOC

(numeric) - Picture Vocabulary.

ORLCMP

(numeric) - Oral Comprehension.

SEGMNT

(numeric) - Segmentation.

RPCNAM

(numeric) - Rapid Picture Naming.

SENREP

(numeric) - Sentence Repetition.

UNDDIR

(numeric) - Understanding Directions.

SNDBLN

(numeric) - Sound Blending.

RETFLU

(numeric) - Retrieval Fluency.

SNDAWR

(numeric) - Sound Awareness.

LWIDNT

(numeric) - Letter-Word Identification.

APPROB

(numeric) - Applied Problems.

SPELL

(numeric) - Spelling.

PSGCMP

(numeric) - Passage Comprehension.

SCI

(numeric) - Science.

SOC

(numeric) - Social Studies.

HUM

(numeric) - Humanities.

Source

McGrew, K. S., LaForte, E. M., & Schrank, F. A. (2014). Technical Manual. Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Schrank, F. A., McGrew, K. S., & Mather, N. (2014). Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Woodcock Johnson IV: ages 40 to 90 plus

Description

A list containing the bivariate correlations (N = 1,146) of the 47 cognitive and achievement subtests from the WJ IV for 40- to 90+-year-olds from the standardization sample obtained from the WJ-IV technical manual (McGrew, LaForte, & Schrank, 2014). Tables are reproduced with permission from the publisher.

Usage

WJIV_ages_40_90

Format

A list of 2 with elements "cormat" (47 x 47 matrix of bivariate correlations) and "N". The correlation matrix contains the following variables:

ORLVOC

(numeric) - Oral Vocabulary.

NUMSER

(numeric) - Number Series.

VRBATN

(numeric) - Verbal Attention.

LETPAT

(numeric) - Letter-Pattern Matching.

PHNPRO

(numeric) - Phonological Processing.

STYREC

(numeric) - Story Recall.

VISUAL

(numeric) - Visualization.

GENINF

(numeric) - General Information.

CONFRM

(numeric) - Concept Formation.

NUMREV

(numeric) - Numbers Reversed.

NUMPAT

(numeric) - Number-Pattern Matching.

NWDREP

(numeric) - Nonword Repetition.

VAL

(numeric) - Visual-Auditory Learning.

PICREC

(numeric) - Picture Recognition.

ANLSYN

(numeric) - Analysis-Synthesis.

OBJNUM

(numeric) - Object-Number Sequencing.

PAIRCN

(numeric) - Pair Cancellation.

MEMWRD

(numeric) - Memory for Words.

PICVOC

(numeric) - Picture Vocabulary.

ORLCMP

(numeric) - Oral Comprehension.

SEGMNT

(numeric) - Segmentation.

RPCNAM

(numeric) - Rapid Picture Naming.

SENREP

(numeric) - Sentence Repetition.

UNDDIR

(numeric) - Understanding Directions.

SNDBLN

(numeric) - Sound Blending.

RETFLU

(numeric) - Retrieval Fluency.

SNDAWR

(numeric) - Sound Awareness.

LWIDNT

(numeric) - Letter-Word Identification.

APPROB

(numeric) - Applied Problems.

SPELL

(numeric) - Spelling.

PSGCMP

(numeric) - Passage Comprehension.

CALC

(numeric) - Calculation.

WRTSMP

(numeric) - Writing Samples.

WRDATK

(numeric) - Word Attack.

ORLRDG

(numeric) - Oral Reading.

SNRDFL

(numeric) - Sentence Reading Fluency.

MTHFLU

(numeric) - Math Facts Fluency.

SNWRFL

(numeric) - Sentence Writing Fluency.

RDGREC

(numeric) - Reading Recall.

NUMMAT

(numeric) - Number Matrices.

EDIT

(numeric) - Editing.

WRDFLU

(numeric) - Word Reading Fluency.

SPLSND

(numeric) - Spelling of Sounds.

RDGVOC

(numeric) - Reading Vocabulary.

SCI

(numeric) - Science.

SOC

(numeric) - Social Studies.

HUM

(numeric) - Humanities.

Source

McGrew, K. S., LaForte, E. M., & Schrank, F. A. (2014). Technical Manual. Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Schrank, F. A., McGrew, K. S., & Mather, N. (2014). Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Woodcock Johnson IV: ages 6 to 8

Description

A list containing the bivariate correlations (N = 825) of the 47 cognitive and achievement subtests from the WJ IV for 6- to 8-year-olds from the standardization sample obtained from the WJ-IV technical manual (McGrew, LaForte, & Schrank, 2014). Tables are reproduced with permission from the publisher.

Usage

WJIV_ages_6_8

Format

A list of 2 with elements "cormat" (47 x 47 matrix of bivariate correlations) and "N". The correlation matrix contains the following variables:

ORLVOC

(numeric) - Oral Vocabulary.

NUMSER

(numeric) - Number Series.

VRBATN

(numeric) - Verbal Attention.

LETPAT

(numeric) - Letter-Pattern Matching.

PHNPRO

(numeric) - Phonological Processing.

STYREC

(numeric) - Story Recall.

VISUAL

(numeric) - Visualization.

GENINF

(numeric) - General Information.

CONFRM

(numeric) - Concept Formation.

NUMREV

(numeric) - Numbers Reversed.

NUMPAT

(numeric) - Number-Pattern Matching.

NWDREP

(numeric) - Nonword Repetition.

VAL

(numeric) - Visual-Auditory Learning.

PICREC

(numeric) - Picture Recognition.

ANLSYN

(numeric) - Analysis-Synthesis.

OBJNUM

(numeric) - Object-Number Sequencing.

PAIRCN

(numeric) - Pair Cancellation.

MEMWRD

(numeric) - Memory for Words.

PICVOC

(numeric) - Picture Vocabulary.

ORLCMP

(numeric) - Oral Comprehension.

SEGMNT

(numeric) - Segmentation.

RPCNAM

(numeric) - Rapid Picture Naming.

SENREP

(numeric) - Sentence Repetition.

UNDDIR

(numeric) - Understanding Directions.

SNDBLN

(numeric) - Sound Blending.

RETFLU

(numeric) - Retrieval Fluency.

SNDAWR

(numeric) - Sound Awareness.

LWIDNT

(numeric) - Letter-Word Identification.

APPROB

(numeric) - Applied Problems.

SPELL

(numeric) - Spelling.

PSGCMP

(numeric) - Passage Comprehension.

CALC

(numeric) - Calculation.

WRTSMP

(numeric) - Writing Samples.

WRDATK

(numeric) - Word Attack.

ORLRDG

(numeric) - Oral Reading.

SNRDFL

(numeric) - Sentence Reading Fluency.

MTHFLU

(numeric) - Math Facts Fluency.

SNWRFL

(numeric) - Sentence Writing Fluency.

RDGREC

(numeric) - Reading Recall.

NUMMAT

(numeric) - Number Matrices.

EDIT

(numeric) - Editing.

WRDFLU

(numeric) - Word Reading Fluency.

SPLSND

(numeric) - Spelling of Sounds.

RDGVOC

(numeric) - Reading Vocabulary.

SCI

(numeric) - Science.

SOC

(numeric) - Social Studies.

HUM

(numeric) - Humanities.

Source

McGrew, K. S., LaForte, E. M., & Schrank, F. A. (2014). Technical Manual. Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Schrank, F. A., McGrew, K. S., & Mather, N. (2014). Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Woodcock Johnson IV: ages 9 to 13

Description

A list containing the bivariate correlations (N = 1,572) of the 47 cognitive and achievement subtests from the WJ IV for 9- to 13-year-olds from the standardization sample obtained from the WJ-IV technical manual (McGrew, LaForte, & Schrank, 2014). Tables are reproduced with permission from the publisher.

Usage

WJIV_ages_9_13

Format

A list of 2 with elements "cormat" (47 x 47 matrix of bivariate correlations) and "N". The correlation matrix contains the following variables:

ORLVOC

(numeric) - Oral Vocabulary.

NUMSER

(numeric) - Number Series.

VRBATN

(numeric) - Verbal Attention.

LETPAT

(numeric) - Letter-Pattern Matching.

PHNPRO

(numeric) - Phonological Processing.

STYREC

(numeric) - Story Recall.

VISUAL

(numeric) - Visualization.

GENINF

(numeric) - General Information.

CONFRM

(numeric) - Concept Formation.

NUMREV

(numeric) - Numbers Reversed.

NUMPAT

(numeric) - Number-Pattern Matching.

NWDREP

(numeric) - Nonword Repetition.

VAL

(numeric) - Visual-Auditory Learning.

PICREC

(numeric) - Picture Recognition.

ANLSYN

(numeric) - Analysis-Synthesis.

OBJNUM

(numeric) - Object-Number Sequencing.

PAIRCN

(numeric) - Pair Cancellation.

MEMWRD

(numeric) - Memory for Words.

PICVOC

(numeric) - Picture Vocabulary.

ORLCMP

(numeric) - Oral Comprehension.

SEGMNT

(numeric) - Segmentation.

RPCNAM

(numeric) - Rapid Picture Naming.

SENREP

(numeric) - Sentence Repetition.

UNDDIR

(numeric) - Understanding Directions.

SNDBLN

(numeric) - Sound Blending.

RETFLU

(numeric) - Retrieval Fluency.

SNDAWR

(numeric) - Sound Awareness.

LWIDNT

(numeric) - Letter-Word Identification.

APPROB

(numeric) - Applied Problems.

SPELL

(numeric) - Spelling.

PSGCMP

(numeric) - Passage Comprehension.

CALC

(numeric) - Calculation.

WRTSMP

(numeric) - Writing Samples.

WRDATK

(numeric) - Word Attack.

ORLRDG

(numeric) - Oral Reading.

SNRDFL

(numeric) - Sentence Reading Fluency.

MTHFLU

(numeric) - Math Facts Fluency.

SNWRFL

(numeric) - Sentence Writing Fluency.

RDGREC

(numeric) - Reading Recall.

NUMMAT

(numeric) - Number Matrices.

EDIT

(numeric) - Editing.

WRDFLU

(numeric) - Word Reading Fluency.

SPLSND

(numeric) - Spelling of Sounds.

RDGVOC

(numeric) - Reading Vocabulary.

SCI

(numeric) - Science.

SOC

(numeric) - Social Studies.

HUM

(numeric) - Humanities.

Source

McGrew, K. S., LaForte, E. M., & Schrank, F. A. (2014). Technical Manual. Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Schrank, F. A., McGrew, K. S., & Mather, N. (2014). Woodcock-Johnson IV. Rolling Meadows, IL: Riverside.

Multigroup exploratory factor analysis

Description

Fit an exploratory factor analysis in each of several groups at a common number of factors and bring the per-group solutions into one shared orientation so their loadings can be compared. Each group is fitted with EFA(); the solutions are then aligned either to a symmetric consensus target or to a chosen reference group (see Alignment).

Usage

efa_group(
  x,
  groups = NULL,
  n_factors,
  N = NA,
  reference_group = NULL,
  b_boot = 0L,
  ci = 0.95,
  seed = NULL,
  delta = 0.1,
  invariance = FALSE,
  ...
)

Arguments

Details

Input

Groups can be supplied in two ways: raw data together with a grouping vector (x a data frame or matrix, groups one value per row), or a named list of per-group data sets in x (with groups left NULL). The list may hold raw data frames or correlation matrices (supply N), but not a mix of the two. All groups must contain the same items in the same order; a different item set or order is an error rather than being silently reordered.

Every group is fitted at the same n_factors (a common-k multigroup model). Extra arguments in ... (for example method, rotation, cor_method, or type) are forwarded unchanged to each EFA() call, so the estimator and rotation are common to all groups.

The k-factor model must be identified for the shared item set: its degrees of freedom ((p - k)^2 - (p + k)) / 2 must be non-negative. Unlike a single EFA() fit – which only warns on an under-identified model – a multigroup fit aborts, because a shared alignment target across an under-identified group is not interpretable.

Alignment

A factor solution is identified only up to a rotation of its factors, so the per-group solutions must be brought into a common orientation before their loadings can be compared. Two strategies are available and are chosen automatically:

• Consensus (the default for orthogonal rotations and for unrotated solutions): a symmetric Generalized Procrustes Analysis target is built across all groups (Gower, 1975), and every group's loadings are rotated to it. No group is privileged.

• Reference: every group's loadings are aligned by Procrustes rotation to one reference group's loadings, which are kept fixed. This path is used when reference_group is given, and is used automatically for oblique rotations because the consensus iteration is not defined for oblique transforms with more than one factor. When an oblique rotation triggers the reference path without an explicit reference_group, the first group is used and a message reports this; the requested rotation is never silently changed.

In both cases the returned per-group loadings share the column order and sign of the returned target.

Comparing the aligned loadings

Because the per-group loadings share one orientation, they can be compared cell by cell. efa_group() reports a per-pair summary of their differences (diffs) and a per-item, per-factor flag table (flags) marking cells whose absolute difference reaches delta. delta is a descriptive salience heuristic, not a significance test; a bootstrap (b_boot > 0) additionally reports whether each difference's confidence interval excludes zero. With invariance = TRUE, a per-factor verdict grades the matched Tucker congruence by the Lorenzo-Seva and ten Berge (2006) similarity bands – "equal" (⁠>= .95⁠) and "fair" (⁠[.85, .95)⁠) – labelling weaker congruences "incongruent"; when a bootstrap is available the verdict uses the congruence CI lower bound, so a factor is judged "equal" only if even the lower bound clears the band.

Value

An object of class efa_group, a list containing:

References

Efron, B., & Tibshirani, R. J. (1993). An Introduction to the Bootstrap. Chapman & Hall.

Gower, J. C. (1975). Generalized Procrustes analysis. Psychometrika, 40, 33-51. doi: 10.1007/BF02291478

Lorenzo-Seva, U., and ten Berge, J. M. F. (2006). Tucker's congruence coefficient as a meaningful index of factor similarity. Methodology, 2, 57-64. doi: 10.1027/1614-2241.2.2.57

See Also

Other factor analysis: plot.efa_group(), print.efa_group()

Examples

# Raw data split by a grouping vector (unrotated, consensus alignment)
g <- rep(c("g1", "g2"), length.out = nrow(GRiPS_raw))
mg <- efa_group(GRiPS_raw, groups = g, n_factors = 1)
mg$loadings

# Per-pair difference summary and the per-item salience-flag table
mg$diffs
mg$flags


# Percentile bootstrap confidence intervals for the between-group congruences, with an
# approximate-invariance verdict read conservatively off the congruence CI lower bound
mg_ci <- efa_group(GRiPS_raw, groups = g, n_factors = 1, b_boot = 100, seed = 42,
                   invariance = TRUE)
mg_ci$congruence$matched_ci
mg_ci$invariance

# A named list of correlation matrices sharing the same items, common
# three-factor model, orthogonal rotation -> symmetric consensus target
bands <- list(age_6_8 = WJIV_ages_6_8$cormat, age_14_19 = WJIV_ages_14_19$cormat)
Ns <- c(WJIV_ages_6_8$N, WJIV_ages_14_19$N)
efa_group(bands, n_factors = 3, N = Ns, rotation = "varimax")

# An oblique rotation aligns to a reference group (reported via a message)
efa_group(bands, n_factors = 3, N = Ns, rotation = "promax")

Power analysis for exploratory factor analysis

Description

Power analysis for exploratory factor analysis, in two modes selected by mode.

mode = "rmsea" (the default) gives the analytic statistical power of the root mean square error of approximation (RMSEA) tests of close and not-close fit (MacCallum, Browne, & Sugawara, 1996): given the model degrees of freedom and a sample size it returns the power to reject the null hypothesis, and given a target power instead of a sample size it returns the sample size required to reach it.

mode = "simulation" runs a Monte-Carlo study of a known population: it draws n_datasets samples of size N (via efa_simulate()) and, over the replicates, reports how often each factor-retention criterion recovers the true number of factors (the hit-rate), how often the fitted loadings recover the population structure (structure recovery, by Tucker congruence), and the convergence and Heywood-case rate of the fit.

Usage

efa_power(
  mode = c("rmsea", "simulation"),
  type = c("close", "notclose"),
  eps0 = NULL,
  eps1 = NULL,
  N = NULL,
  p = NULL,
  k = NULL,
  df = NULL,
  alpha = 0.05,
  power = NULL,
  group = 1,
  Lambda = NULL,
  Phi = NULL,
  Psi = NULL,
  R = NULL,
  n_datasets = 500,
  criteria = c("EKC", "MAP"),
  method = "PAF",
  rotation = NULL,
  recovery_threshold = 0.95,
  model_error = c("TKL", "CB", "WB", "none"),
  target_rmsea = NULL,
  target_cfi = NULL,
  seed = NULL
)

Arguments

Value

An object of class efa_power. For mode = "rmsea", a list containing:

For mode = "simulation", a list containing:

RMSEA mode

The RMSEA fit statistic is referred to a noncentral chi-square distribution. For a model with df degrees of freedom and population RMSEA \varepsilon, the noncentrality parameter is

\lambda = (N - 1)\, df\, \varepsilon^2 / G,

where N is the sample size and G the number of groups (group).

Two tests are supported, selected by type (never by the ordering of eps0 and eps1):

"close"

The test of close fit (MacCallum et al., 1996): the null hypothesis is that the fit is close (RMSEA \le eps0, conventionally 0.05), and power is computed against a worse alternative (eps1, conventionally 0.08, so eps0 < eps1), in the upper tail.

"notclose"

The test of not-close fit: the null hypothesis is that the fit is not close (RMSEA \ge eps0), and power is computed against a better alternative (eps1, conventionally 0.01, so eps0 > eps1), in the lower tail.

When eps0 and eps1 are ordered the wrong way round for the chosen type a message is emitted but the requested test is still computed. Equal eps0 and eps1 leave nothing to detect and are an error.

Power increases monotonically with N, so the required sample size (the smallest N reaching power) is found by bisection. It is a per-group sample size: with group > 1 the noncentrality carries the 1 / group factor above.

Simulation mode

The population is supplied either as a factor model (Lambda, with optional Phi and Psi) or as a ready correlation matrix R, and passed to efa_simulate(), which draws n_datasets samples of size N. By default the population fits the factor model exactly, which overstates how well the criteria and the fit recover the structure; supplying a misfit target (target_rmsea and/or target_cfi) perturbs it with model error (model_error, "TKL" by default) so the model fits only approximately – a more realistic target (MacCallum, 2003). The true number of factors k_true is ncol(Lambda) for a factor-model population, or k for a bare R.

Each replicate is analysed three ways. Hit-rate: every criterion in criteria is run and its suggested number of factors compared with k_true; the hit-rate is the proportion of agreements over the replicates on which the criterion returned a definite factor count (replicates where it errored or was undecided are excluded, not counted as misses). Structure recovery (only for a factor-model population): the k_true-factor model is fitted with EFA(), its rotated loadings are matched to the population loadings, and the matched-factor Tucker congruences (Lorenzo-Seva & ten Berge, 2006) are compared with recovery_threshold; a replicate counts as a success when the smallest (min) or the average (mean) matched congruence reaches it. Convergence: the same fit supplies the proportion of replicates whose fit completed and, among those, the proportion that converged and that produced a Heywood case.

Replicates are analysed in parallel with future.apply; a plan is selected with future::plan(). Each replicate is bound to its own reproducible random-number stream, so with a fixed seed the result is identical regardless of the number of workers, and the caller's random-number state is left unchanged.

References

MacCallum, R. C., Browne, M. W., & Sugawara, H. M. (1996). Power analysis and determination of sample size for covariance structure modeling. Psychological Methods, 1(2), 130-149. doi:10.1037/1082-989X.1.2.130

MacCallum, R. C. (2003). 2001 Presidential Address: Working with imperfect models. Multivariate Behavioral Research, 38(1), 113-139. doi:10.1207/S15327906MBR3801_5

Lorenzo-Seva, U., & ten Berge, J. M. F. (2006). Tucker's congruence coefficient as a meaningful index of factor similarity. Methodology, 2(2), 57-64. doi:10.1027/1614-2241.2.2.57

See Also

Other power analysis: plot.efa_power(), print.efa_power()

Examples

# Power of the test of close fit at N = 200 for a 100-df model
efa_power(df = 100, N = 200)

# Deriving df from the model dimensions instead of giving it directly
efa_power(p = 20, k = 3, N = 200)

# Required (per-group) sample size for 80% power
efa_power(df = 100, power = 0.80)

# Test of not-close fit
efa_power(df = 100, N = 200, type = "notclose")


# Simulation mode: retention hit-rate and structure recovery for a known
# three-factor population at N = 300
efa_power("simulation", Lambda = population_models$loadings$baseline,
          Phi = population_models$phis_3$moderate, N = 300,
          n_datasets = 50, criteria = c("EKC", "MAP"), seed = 42)

Reliability and common-variance coefficients for a factor solution

Description

Computes model-based reliability coefficients (McDonald's omega total, hierarchical, and subscale; standardized Cronbach's alpha; and the H index) together with the bifactor common-variance indices (explained common variance, ECV, and percent of uncontaminated correlations, PUC) for the general factor and each group factor, and returns them as a tidy, long-format table. The coefficients can be obtained from a Schmid-Leiman solution (SL() or psych::schmid()), an oblique EFA() (correlated-factors) solution, a lavaan single-factor, second-order, or bifactor fit, a raw bifactor loading matrix, or manually supplied components.

Usage

efa_reliability(
  model = NULL,
  coefficients = NULL,
  g_name = "g",
  group_names = NULL,
  factor_map = NULL,
  variance = c("correlation", "sums_load"),
  var_names = NULL,
  fac_names = NULL,
  g_load = NULL,
  s_load = NULL,
  u2 = NULL,
  cormat = NULL,
  pattern = NULL,
  Phi = NULL
)

Arguments

Details

Coefficients

The reliability coefficients are McDonald's omegas (McDonald, 1978, 1985, 1999), standardized Cronbach's alpha (Cronbach, 1951), and the H index (construct replicability; Hancock & Mueller, 2001). The common-variance indices are the ECV and PUC (Bonifay et al., 2015; Reise et al., 2013; Rodriguez et al., 2016a, 2016b); they describe the general factor and so are reported for the general factor only. See OMEGA() for the definition and interpretation of each coefficient.

Input

The dispatch on model mirrors OMEGA(), with two additions: an oblique EFA() object is scored as the correlated-factors model it is (having no general factor, it omits the bifactor indices – omega hierarchical, ECV, and PUC), and a bare loading matrix is read as a raw bifactor solution (general factor in the first column). For a correlated-factors EFA() solution variance is always "correlation". The indicator-to-factor correspondences come from factor_map when it is supplied; otherwise each variable is assigned to the group factor on which it loads most strongly. For lavaan input the composite variances are model-implied (variance is not used), and the coefficients are computed per group.

Value

An object of class efa_reliability: a long-format data frame with one row per computed coefficient, with columns

Structurally undefined cells (for example ECV and PUC on a group factor) are omitted. The object carries a settings attribute (the total-variance convention used) and a kind attribute tagging each coefficient as a reliability coefficient or a common-variance index, and has a print() method.

Source

McDonald, R. P. (1978). Generalizability in factorable domains: Domain validity and generalizability. Educational and Psychological Measurement, 38, 75-79.

McDonald, R. P. (1985). Factor analysis and related methods. Hillsdale, NJ: Erlbaum.

McDonald, R. P. (1999). Test theory: A unified treatment. Mahwah, NJ: Erlbaum.

Cronbach, L. J. (1951). Coefficient alpha and the internal structure of tests. Psychometrika, 16, 297-334.

Hancock, G. R., & Mueller, R. O. (2001). Rethinking construct reliability within latent variable systems. In R. Cudeck, S. du Toit, & D. Sörbom (Eds.), Structural equation modeling: Present and future (pp. 195-216). Lincolnwood, IL: Scientific Software International.

Bonifay, W. E., Reise, S. P., Scheines, R., & Meijer, R. R. (2015). When are multidimensional data unidimensional enough for structural equation modeling? An evaluation of the DETECT multidimensionality index. Structural Equation Modeling, 22, 504-516.

Reise, S. P., Scheines, R., Widaman, K. F., & Haviland, M. G. (2013). Multidimensionality and structural coefficient bias in structural equation modeling: A bifactor perspective. Educational and Psychological Measurement, 73, 5-26.

Rodriguez, A., Reise, S. P., & Haviland, M. G. (2016a). Applying bifactor statistical indices in the evaluation of psychological measures. Journal of Personality Assessment, 98, 223-237.

Rodriguez, A., Reise, S. P., & Haviland, M. G. (2016b). Evaluating bifactor models: Calculating and interpreting statistical indices. Psychological Methods, 21, 137-150.

See Also

OMEGA() for the same coefficients in a wide, per-factor layout.

Other reliability coefficients: print.efa_reliability()

Examples

## From an oblique EFA (correlated-factors) solution. With no factor_map, each
## item is auto-assigned to its highest-loading factor.
efa_mod <- EFA(test_models$baseline$cormat, N = 500, n_factors = 3,
               type = "EFAtools", method = "PAF", rotation = "promax")
efa_reliability(efa_mod)

## From a Schmid-Leiman solution, with an explicit indicator-to-factor map.
sl_mod <- SL(efa_mod, type = "EFAtools", method = "PAF")
fc <- sl_mod$sl[, c("F1", "F2", "F3")] >= .2
efa_reliability(sl_mod, factor_map = fc)

## Request a subset of the coefficients only.
efa_reliability(sl_mod, factor_map = fc,
                coefficients = c("omega_total", "alpha"))

## From a lavaan bifactor solution.

if (requireNamespace("lavaan", quietly = TRUE)) {
mod <- 'F1 =~ V1 + V2 + V3 + V4 + V5 + V6
        F2 =~ V7 + V8 + V9 + V10 + V11 + V12
        F3 =~ V13 + V14 + V15 + V16 + V17 + V18
        g =~ V1 + V2 + V3 + V4 + V5 + V6 + V7 + V8 + V9 + V10 + V11 + V12 +
             V13 + V14 + V15 + V16 + V17 + V18'
fit <- lavaan::cfa(mod, sample.cov = test_models$baseline$cormat,
                   sample.nobs = 500, estimator = "ml", orthogonal = TRUE)
efa_reliability(fit, g_name = "g")
}

Estimate factor scores and score-quality diagnostics for an EFA model

Description

Computes factor-score weights (and, from raw data, the factor scores themselves) natively for an EFA() solution or a directly supplied loading matrix, together with the score-quality diagnostics that describe how well the estimated scores represent the factors: the score intercorrelations, the determinacy (validity) and univocality of each score, and Guttman's indeterminacy index. Factor scores are returned only when raw data are supplied; a correlation matrix yields the weights and diagnostics alone.

Usage

efa_scores(
  x,
  f,
  Phi = NULL,
  rho = NULL,
  method = c("regression", "Bartlett", "Anderson", "tenBerge", "Harman", "components")
)

Arguments

Details

The p by m weight matrix W (standardized scores are scale(X) %*% W) is computed from the structure matrix S = Lambda %*% Phi, the model uniquenesses Psi = diag(1 - h2), and the scoring correlation matrix R according to method:

"regression"

Thurstone's (1935) regression scores, ⁠W = R^-1 S⁠.

"Bartlett"

Bartlett's (1937) conditionally unbiased scores.

"Anderson"

Anderson & Rubin's (1956) uncorrelated, unit-variance scores; defined for orthogonal factors only.

"tenBerge"

ten Berge, Krijnen, Wansbeek & Shapiro's (1999) scores, which preserve the factor intercorrelations.

"Harman"

Harman's (1976) idealized-variable scores.

"components"

component scores, W = Lambda.

The determinacy (validity) of a score is its correlation with the factor it estimates, computed from the returned weights; for regression scores it is the multiple correlation between the factor and the observed variables (Guttman, 1955; Grice, 2001). The off-diagonal score-factor correlations give the univocality (the correlation of a score with the other factors), and ⁠2 rho^2 - 1⁠ is Guttman's (1955) indeterminacy index, the minimum correlation between two equally valid sets of scores. For a method other than "regression" both quantities are specific to those scores: the determinacy is the method's own score-factor correlation (never larger than the regression value), and the reported guttman follows from it.

Value

An object of class efa_scores, a list containing: