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).
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,
...
)A list of class EFA containing (a subset of) the following:
Original correlation matrix.
Initial communality estimates from PAF.
Final communality estimates from the unrotated solution.
Eigen values of the original correlation matrix.
Initial eigenvalues, obtained from the correlation matrix with the initial communality estimates as diagonal in PAF.
Eigenvalues obtained from the correlation matrix with the final communality estimates as diagonal.
For PAF, the number of iterations until convergence. For ML, ULS, and DWLS, the number of objective-function evaluations used by the optimiser (not the number of optimiser iterations).
Integer convergence code (0 = converged). For ML, ULS, and
DWLS this is the convergence code from the bounded optimiser (the same codes as
stats::optim()'s "L-BFGS-B"); for PAF it is 1 if the maximum
number of iterations was reached without meeting the convergence criterion and 0
otherwise. A non-zero code is also reported with a warning.
A named integer vector indicating which variables have a Heywood (improper) case in the unrotated solution; empty if there are none.
Loading matrix containing the final unrotated loadings.
Matrix of explained variances and sums of squared loadings. Based on the unrotated loadings.
A named list of fit indices computed from the unrotated
loadings. For ML and ULS it holds the model Chi Square (with its p-value and
df), CFI, TLI, RMSEA with its 90% confidence interval, AIC, BIC, ECVI, RMSR,
SRMR, and CAF; for PAF and DWLS only RMSR, SRMR, CAF, and df are populated and
the Chi-Square-based indices are NA (for DWLS with se = "sandwich" the full
block is filled from a scaled Chi Square instead). Whenever the Chi Square is a
scaled statistic (se = "sandwich", or any cor_method = "fiml" fit) AIC,
BIC, and ECVI are NA and the list additionally carries the scaling
components: chi_scaling (the multiplier a in the scaled-and-shifted statistic
\(aT + b\), i.e. the reciprocal of lavaan's chisq.scaling.factor),
chi_shift (b), chi_unscaled (the unscaled statistic T), and the alternative
chi_mean_adjusted and chi_mean_var statistics with their df_mean_var.
RMSR is retained for programmatic use and backward compatibility, although the
print and summary methods display SRMR. See the Fit indices section in
Details for how each index is defined, scaled, and referenced.
The model implied correlation matrix.
Residual correlations, i.e., orig_R - model_implied_R
Residual correlations standardized by their
bootstrap standard errors. Only returned, if se = "np-boot".
Loading matrix containing the final rotated loadings (pattern matrix).
The factor intercorrelations (only for oblique rotations).
The structure matrix (only for oblique rotations).
The rotation matrix. The rotated loadings are recovered from the
unrotated loadings as unrot_loadings %*% rotmat for orthogonal rotations and
for promax, and as unrot_loadings %*% t(solve(rotmat)) for the other oblique
rotations.
Matrix of explained variances and sums of squared loadings. Based on rotated loadings and, for oblique rotations, the factor intercorrelations.
A list of the settings used.
A named list of standard error matrices. For se = "np-boot": bootstrap standard deviations of the unrotated and (when a rotation is applied) rotated loadings, the residuals, and the fit indices, plus -- for oblique rotations -- the factor correlations (Phi) and the structure coefficients. For se = "information": Wald standard errors from the expected (Fisher) information matrix for the unrotated loadings and the uniquenesses and, when a rotation is applied, the rotated loadings and the communalities (and, for oblique rotations, Phi and the structure coefficients). For se = "sandwich": robust Godambe sandwich standard errors with the same coverage as "information", robust to non-normality and weight misspecification. Only returned if se is not "none".
A named list of confidence intervals of width ci. For se = "np-boot": percentile intervals matching the components of SE. For se = "information" and se = "sandwich": Wald intervals matching the components of SE. Only returned if se is not "none".
A named list of bootstrap replicate cubes for the aligned unrotated and (where applicable) rotated loadings, structure coefficients, factor correlations (Phi), residuals, and fit indices. Each cube's last dimension indexes the replicate. Populated only for se = "np-boot"; NULL for the analytic SE methods.
The full unrotated loading covariance matrix the marginal SE$unrot_loadings were derived from: a p * n_factors by p * n_factors numeric matrix in column-major vec(Lambda) order. Populated for se = "information" (expected-information block) and se = "sandwich" (robust V_AA), even when a rotation is applied (the persisted block is always the unrotated one); NA-filled if the analytic covariance is unreliable (a Heywood case or a singular bordered information matrix); NULL for se = "np-boot" and se = "none".
The asymptotic covariance of the off-diagonal sample correlations -- the meat of the robust sandwich SEs -- on the variance scale (Var(rho-hat); lavaan's correlation NACOV is N * Gamma). A p (p - 1) / 2 by p (p - 1) / 2 numeric matrix; rows and columns ordered by utils::combn() over the column pairs and labelled "<var_i>-<var_j>". Populated only for se = "sandwich"; NULL otherwise.
data.frame or matrix. Dataframe or matrix of raw data or matrix with correlations. If raw data is entered, the correlation matrix is found from the data.
numeric. Number of factors to extract.
numeric. The number of observations. Needs only be specified if a
correlation matrix is used. If input is a correlation matrix and N = NA
(default), not all fit indices can be computed. When raw data with missing
values are entered and use is "complete.obs" or "na.or.complete", rows
are deleted listwise, so N is taken as the number of complete cases.
character. The estimator used to fit the EFA: "PAF" (principal axis factoring), "ML" (maximum likelihood), "ULS" (unweighted least squares; "MINRES" is an accepted alias returning identical results), or "DWLS" (diagonally weighted least squares, for ordinal data). See the Estimators section in Details for their properties and data requirements.
character. Either perform no rotation ("none"; default), an orthogonal rotation ("varimax", "equamax", "quartimax", "geominT", "bentlerT", or "bifactorT"), or an oblique rotation ("promax", "oblimin", "quartimin", "simplimax", "bentlerQ", "geominQ", or "bifactorQ"). See the Rotations section in Details for their properties and known issues.
character. Whether and how to compute standard errors (and matching
confidence intervals): "none" (default, no standard errors), "information" (analytic
standard errors from the expected Fisher information of the ML solution), "sandwich"
(robust Godambe sandwich standard errors from raw data), or "np-boot" (non-parametric
bootstrap). The methods differ in their assumptions, their data requirements, and which
estimator, rotation, and cor_method combinations they support; see the Standard
errors section in Details.
character. If one of "EFAtools" (default), "psych", or "SPSS" is
used, and the following arguments with default NA are left with
NA, these implementations are executed according to the respective program
("psych" and "SPSS") or according to the best solution found in Grieder &
Steiner (2022; "EFAtools"). Individual properties can be adapted using one of
the three types and specifying some of the following arguments. If set to
"none" additional arguments must be specified depending on the method
and rotation used (see details).
numeric. The maximum number of iterations to perform after which
the iterative PAF procedure is halted with a warning. If type is one of
"EFAtools", "SPSS", or "psych", this is automatically specified if max_iter is
left to be NA, but can be overridden by entering a number. Default is
NA.
character. The method to estimate the initial communalities
in PAF. "smc" will use squared multiple correlations, "mac" will use
maximum absolute correlations, "unity" will use 1s (see details).
Default is NA.
numeric. The convergence criterion used for PAF.
If the change in communalities from one iteration to the next is smaller than
this criterion the solution is accepted and the procedure ends.
Default is NA.
character. Type of convergence criterion used for
PAF. "max_individual" selects the maximum change in any of the
communalities from one iteration to the next and tests it against the
specified criterion. This is also used by SPSS. "sum" takes the difference of
the sum of all communalities in one iteration and the sum of all communalities
in the next iteration and tests this against the criterion. This procedure is
used by the psych::fa() function. Default is NA.
logical. Which algorithm to use in the PAF
iterations. If FALSE, the loadings are computed from the eigenvalues. This is
also used by the psych::fa() function. If TRUE the
loadings are computed with the absolute eigenvalues as done by SPSS.
Default is NA.
character. Passed to stats::cor() if raw data
is given as input. Default is "pairwise.complete.obs".
character. The type of the varimax rotation performed.
If "svd", singular value decomposition is used, as stats::varimax() does. If
"kaiser", the varimax procedure performed in SPSS is used, following the original
procedure from Kaiser (1958) (see details). Default is NA.
numeric. Either the power used for computing the target matrix P in
the promax rotation or the number of 'close to zero loadings' for the simplimax
rotation. If left to NA (default), the value for promax depends on the specified type.
For simplimax, nrow(L), where L is the matrix of unrotated loadings,
is used by default.
logical. If TRUE, a kaiser normalization is
performed before the specified rotation. Default is TRUE.
character. This specifies how the target
matrix P is computed in promax rotation. If "unnorm" it will use the
unnormalized target matrix as originally done in Hendrickson and White (1964).
This is also used in the psych and stats packages. If "norm" it will use the
normalized target matrix as used in SPSS. Default is NA.
numeric. The tolerance for stopping in the rotation procedure. Default is 10^-5 for all rotation methods.
character. How to order the factors. "eigen" reorders the
factors by descending explained variance, i.e. by their reported sums of squared
loadings ("SS loadings"): the column sums of squares for orthogonal solutions and
the factor-intercorrelation-weighted sums of squares for oblique solutions, so the
reported variances decrease monotonically (as in the psych package). "ss_factors"
reorders the factors by descending (unweighted) sum of squared factor loadings per
factor; for oblique solutions this can differ from "eigen", whereas for orthogonal
solutions the two coincide. Default is NA.
character. How to specify the starting values for the
optimization procedure for ML. Default is "psych" which takes the
starting values specified in psych::fa(). "factanal" takes the
starting values specified in the stats::factanal() function.
Solutions are very similar.
character. How the correlation is computed from raw data:
"pearson", "spearman", or "kendall" (passed to stats::cor()); "poly" /
"tetra" for polychoric / tetrachoric correlations of ordinal / binary data; or
"fiml" for a two-stage full-information maximum-likelihood correlation from raw data
with missing values. See the Correlation methods section in Details for their
properties and the combinations they support. Default is "pearson".
numeric. The number of bootstrap samples to draw. Default is 1000.
Under cor_method = "fiml" each bootstrap sample re-runs the EM moment
estimation, so a smaller value may be advisable.
numeric. The confidence interval to create from the bootstrap samples. Must be between 0 and 1. Default is .95 for 95% CIs.
numeric. The number of random starts to use in the rotation. Some rotation criteria are prone to produce local minima, and several random starts are usually needed to locate the best solution. The rotation screens the random starts cheaply and fully optimises only the most promising ones, so a large value adds little cost for most criteria. The complexity criteria (simplimax and, to a lesser extent, geomin) are the most multimodal and may need a larger value on difficult data. Default is 100.
numeric. An optional seed for the random-number generator used by the
non-parametric bootstrap (se = "np-boot"), i.e. for the case resampling, the
rotation random starts, and the Procrustes random starts. Setting it makes the
bootstrap reproducible and independent of the number of parallel
workers (see Details); the caller's random-number stream is restored afterwards,
so supplying a seed leaves no lasting effect on it. Default is NULL, which uses
(and advances) the current state of the generator.
Additional arguments passed to the rotation procedure (e.g., maxit for the maximum number of iterations).
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.
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.
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.
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.
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).
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.
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.
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().
# 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)
# \donttest{
# 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)
# }
# \donttest{
# 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
# }
if (FALSE) {
# 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)
}
Run the code above in your browser using DataLab