Package {GPArotation}

Gradient Projection Algorithms for Factor Rotation

Description

GPA Rotation for Factor Analysis

The GPArotation package contains functions for the rotation of factor loadings matrices. The functions implement Gradient Projection (GP) algorithms for orthogonal and oblique rotation. Additionally, a number of rotation criteria are provided. The GP algorithms minimize the rotation criterion function and provide the corresponding rotation matrix. For oblique rotation, the covariance/correlation matrix of the factors is also provided. The rotation criteria implemented in this package are described in Bernaards and Jennrich (2005) in addition to a number of others. Theory of the GP algorithm is described in Jennrich (2001, 2002).

Vignettes are provided covering general usage, local minima diagnostics, bifactor rotation and reliability, and a visual walkthrough of the gradient projection algorithm. Access them via browseVignettes("GPArotation").

Index of functions:

Rotations using gradient projection algorithms

Core gradient projection algorithms

Random-start wrappers and internal engine

2D orthogonal trajectory plot

S3 methods

Data sets

Other rotations not using gradient projection algorithms

Legacy gradient projection algorithms (code unchanged since 2008)

Utility functions

Utility functions (not exported)

Value, gradient, rotation criterion functions (not exported)

Vignettes

Author(s)

Coen A. Bernaards and Robert I. Jennrich with some R modifications by Paul Gilbert.

References

The software reference is:

Bernaards, C.A. and Jennrich, R.I. (2005). Gradient projection algorithms and software for arbitrary rotation criteria in factor analysis. Educational and Psychological Measurement, 65, 676–696. doi:10.1177/0013164404272507

Theory of gradient projection algorithms:

Jennrich, R.I. (2001). A simple general procedure for orthogonal rotation. Psychometrika, 66, 289–306. doi:10.1007/BF02294840

Jennrich, R.I. (2002). A simple general method for oblique rotation. Psychometrika, 67, 7–19. doi:10.1007/BF02294706

A clear and accessible introduction to gradient projection algorithms for factor rotation is provided in:

Mansolf, M. and Reise, S.P. (2016). Exploratory bifactor analysis: The Schmid-Leiman orthogonalization and Jennrich-Bentler analytic rotations. Multivariate Behavioral Research, 51(5), 698–717. doi:10.1080/00273171.2016.1215898

Barzilai-Borwein step size:

Barzilai, J. and Borwein, J.M. (1988). Two-point step size gradient methods. IMA Journal of Numerical Analysis, 8, 141–148. doi:10.1093/imanum/8.1.141

Cayley transform retraction:

Wen, Z. and Yin, W. (2013). A feasible method for optimization with orthogonality constraints. Mathematical Programming, 142, 397–434. doi:10.1007/s10107-012-0584-1

Non-monotone line search:

Grippo, L., Lampariello, F., and Lucidi, S. (1986). A nonmonotone line search technique for Newton's method. SIAM Journal on Numerical Analysis, 23(4), 707–716. doi:10.1137/0723046

Zhang, H. and Hager, W.W. (2004). A nonmonotone line search technique and its application to unconstrained optimization. SIAM Journal on Optimization, 14(4), 1043–1056. doi:10.1137/S1052623403428208

Simple structure measures:

Liu, X., Wallin, G., Chen, Y., and Moustaki, I. (2023). Rotation to sparse loadings using L^p losses and related inference problems. Psychometrika, 88(2), 527–553. doi:10.1007/s11336-023-09911-y

Lorenzo-Seva, U. (2003). A factor simplicity index. Psychometrika, 68(1), 49–60. doi:10.1007/BF02296652

See Also

GPFRSorth, GPFRSoblq, rotations, vgQ, simple_structure, plot.GPArotation browseVignettes("GPArotation")

The CCAI Climate-Friendly Purchasing Choices Domain

Description

Correlation matrix, factor pattern matrix, and factor intercorrelations for the Climate Change Action Inventory (CCAI) (Barchard et al., 2021) Climate-Friendly Purchasing Choices domain, analyzed in Bi and Barchard (2024). The scale measures the frequency with which individuals make purchasing choices aimed at reducing climate change. Data were collected from 500 United States MTurk workers. After 15 climate change deniers and 24 multivariate outliers were removed, 461 participants remained. Climate change deniers were excluded because the scale measures behaviors intended to reduce climate change — including individuals who do not believe climate change exists would be inconsistent with the measure's intent.

The Climate Change Action Inventory contains eight domains; the Climate-Friendly Purchasing Choices domain analyzed here is one of them.

Usage

data(CCAI, package = "GPArotation")

Format

Three objects are loaded by data(CCAI):

• CCAI_R: a 14 \times 14 numeric correlation matrix. Row and column names are item labels CCAI1 through CCAI14, ordered by factor to match CCAI_pattern.

• CCAI_pattern: a 14 \times 3 numeric matrix of factor pattern coefficients. Row names are item labels CCAI1 through CCAI14, ordered by factor to facilitate interpretation. Column names are factor labels SustainableOptions, CollectiveAction, and AvoidBuyingNew.

• CCAI_Phi: a 3 \times 3 numeric matrix of factor intercorrelations. Row and column names are the three factor labels.

Details

The CCAI Climate-Friendly Purchasing Choices domain consists of 14 items. Items used a nine-point frequency scale ranging from 1 (less than once a year) to 9 (at least 14 times a week). For full details on study design, data collection, analysis, and interpretation see Bi and Barchard (2024).

CCAI_R is the observed 14 \times 14 correlation matrix for the 14 items of the Climate-Friendly Purchasing Choices domain. The correlation matrix was kindly provided by the authors and has not been published separately.

CCAI_pattern is a 14 \times 3 factor pattern matrix from principal components extraction followed by direct oblimin rotation, reported in Table 2 of Bi and Barchard (2024) (the Table is labeled “Factor Structure” in the paper but contains pattern coefficients). The three factors are: Choosing Sustainable Options (F1), Supporting Collective Action (F2), and Avoiding Buying New (F3).

CCAI_Phi is a 3 \times 3 factor intercorrelation matrix. All three intercorrelations exceed 0.50.

The Tandem Paradigm on Correlated Data

The CCAI Purchasing Choices domain illustrates the Tandem criteria on highly correlated factors. Because the underlying dimensions are strongly intercorrelated (\phi_{ij} > 0.50), TandemI pulls shared variance into a single dominant general factor (approximately 59 percent of total variance), alerting the researcher to a unified behavioral core.

TandemII forces an orthogonal simple structure, but the 90-degree constraint creates an artificial compromise when true dimensions are correlated — smearing loadings across factors. Comparing TandemII to oblique Oblimin using the package summary diagnostics (calc_AUC, calc_hyperplane) provides empirical evidence that an oblique structure is conceptually superior for these data.

See vignette("GPA3bifactor", package = "GPArotation") for a bifactor analysis of these data.

Note

The raw data are publicly available on the Open Science Framework. As the data are subject to a license, users should consult the OSF page for terms of use before using the raw data directly: https://osf.io/h38yb/overview.

References

Barchard, K.A., Okagawa, K., Hoffman, C.K., and Odents, O. (2021). Climate Change Action Inventory. Unpublished psychological test. Available from kim.barchard@unlv.edu.

Bi, Y. and Barchard, K.A. (2024). Purchasing choices that reduce climate change: An exploratory factor analysis. Spectra Undergraduate Research Journal, 3(2), 8–14. doi:10.9741/2766-7227.1028.

See Also

rotations, bifactorT, eigen, factanal, Harman, Thurstone, WansbeekMeijer, GriffithMulaik

Examples

  data(CCAI, package = "GPArotation")

  # Observed correlation matrix
  round(CCAI_R, 2)

  # Published pattern matrix and factor intercorrelations
  round(CCAI_pattern, 2)
  round(CCAI_Phi, 2)

  # Reproduce published analysis: PCA extraction via eigendecomposition
  # followed by oblimin rotation. No additional packages required.
  # Gives the same result as psych::principal used by Bi and Barchard (2024).
  ev          <- eigen(CCAI_R)
  k           <- 3
  L_unrotated <- ev$vectors[, 1:k] %*% diag(sqrt(ev$values[1:k]))
  rownames(L_unrotated) <- colnames(CCAI_R)
  res.ob <- oblimin(L_unrotated, randomStarts = 100)

  # Sort and extract loadings for comparison
  res_sorted <- GPArotation:::.sortGPALoadings(res.ob)
  L_repro    <- unclass(res_sorted$loadings)

  # Compare reproduced vs published side by side, alternating by factor
  comparison <- cbind(round(L_repro[, 1], 2), round(CCAI_pattern[, 1], 2),
                      round(L_repro[, 2], 2), round(CCAI_pattern[, 2], 2),
                      round(L_repro[, 3], 2), round(CCAI_pattern[, 3], 2))
  colnames(comparison) <- c("F1.repro", "F1.pub",
                            "F2.repro", "F2.pub",
                            "F3.repro", "F3.pub")
  print(comparison)

  # --- Tandem criteria ---
  # Stage 1: TandemI reveals factor redundancy structure.
  # One dominant general factor (SS = 8.3, 59% variance) suggests
  # highly correlated underlying constructs.
  res.t1 <- tandemI(L_unrotated, randomStarts = 100)
  print(res.t1)

  # Stage 2: TandemII for simple structure (orthogonal constraint).
  res.t2 <- tandemII(L_unrotated, randomStarts = 100)
  summary(res.t2)

  # For heavily correlated constructs, oblique rotation achieves
  # substantially cleaner simple structure than TandemII.
  # Compare AUC and hyperplane counts across methods.
  cat("           AUC    Hyperplane\n")
  cat("TandemII  ",
      round(GPArotation:::calc_AUC(res.t2)$AUC_mean, 3), "  ",
      sum(abs(unclass(res.t2$loadings)) < 0.1), "loadings\n")
  cat("Oblimin   ",
      round(GPArotation:::calc_AUC(res.ob)$AUC_mean, 3), "  ",
      sum(abs(unclass(res.ob$loadings)) < 0.1), "loadings\n")

  # Orthogonal bifactor rotation using MLE extraction
  fa.un <- factanal(factors = 3, covmat = CCAI_R, n.obs = 461,
                    rotation = "none")
  bif <- bifactorT(fa.un)
  print(bif, sortLoadings = FALSE, digits = 3, cutoff = 0.05)

Core Algorithms and Random-Start Wrappers

Description

Gradient projection rotation optimization routines for orthogonal and oblique factor rotation. These functions can be used directly to rotate a loadings matrix, or indirectly through a rotation objective passed to a factor estimation routine such as factanal.

Usage

    GPForth(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=2000, 
       method="varimax", methodArgs=NULL, algorithm = "bb", fwindow = 10)
    GPFoblq(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=2000, 
       method="quartimin", methodArgs=NULL, algorithm = "bb", fwindow = 10)

    GPFRSorth(A, Tmat=NULL, normalize=FALSE, eps=1e-5, maxit=2000, method="varimax", 
       methodArgs=NULL, randomStarts=0, algorithm = "bb", fwindow = 10, ...)
    GPFRSoblq(A, Tmat=NULL, normalize=FALSE, eps=1e-5, maxit=2000, method="quartimin", 
       methodArgs=NULL, randomStarts=0, algorithm = "bb", fwindow = 10, ...)
    

Arguments

Details

The GPFRSorth and GPFRSoblq functions serve as the primary user interfaces for orthogonal and oblique rotations, respectively. They act as wrappers for the core GP algorithms (GPForth for orthogonal rotation and GPFoblq for oblique rotation), extending them with the ability to perform multiple random starts. Any additional arguments provided to these wrappers are passed directly down to the underlying GP algorithms. While the wrappers are generally recommended, the core functions GPForth and GPFoblq can also be invoked directly.

All of these functions require an initial loadings matrix, A, which fixes the equivalence class over which the optimization is performed. This matrix must be the solution to an orthogonal factor analysis problem, such as one obtained from factanal or another factor estimation routine.

Mathematically, a general rotation of a matrix A is defined as A %*% solve(t(Th)). In the case of orthogonal rotation, the initial rotation matrix Tmat is orthonormal, which simplifies the rotation formula to A %*% Th. In all scenarios, the final rotation matrix Th is computed by the GP rotation algorithm.

An accessible introduction to gradient projection algorithms for factor rotation is provided in Mansolf and Reise (2016).

The normalize argument

The normalize argument specifies whether and how the loadings matrix should be normalized prior to rotation, and subsequently denormalized after rotation.

• If FALSE (the default), no normalization is performed.

• If TRUE, Kaiser normalization is applied so that the squared row entries of the normalized matrix A sum to 1.0. This procedure is sometimes referred to as Horst normalization.

• If provided as a vector (which must have a length equal to the number of indicators, i.e., the number of rows in A), the columns of A are divided by this vector before rotation and multiplied by it afterward.

• If provided as a function, it can be used to apply a custom normalization scheme. The function must take A as an argument and return a vector, which is then applied in the same manner as the vector input described above. See NormalizingWeight for an example implementing Cureton-Mulaik normalization.

For a detailed investigation into how normalization affects factor rotations, including its potential impact on the qualitative interpretation of loadings, see Nguyen and Waller (2022).

The method argument

The method argument takes a string specifying the rotation objective function. By default, oblique rotations use "quartimin", while orthogonal rotations default to "varimax". The package supports a comprehensive suite of rotation objectives: "oblimin", "quartimin", "target", "pst", "oblimax", "entropy", "quartimax", "Varimax", "simplimax", "bentler", "tandemI", "tandemII", "geomin", "cf", "infomax", "mccammon", "bifactor", "lp", and "varimin".

Internally, this string is prefixed with "vgQ." to invoke the actual calculation function (see vgQ for underlying mathematical details). It is important to note that several rotation criteria — specifically "oblimin", "target", "pst", "simplimax", "geomin", "cf", and "lp" — require one or more supplementary arguments. These additional arguments can be seamlessly passed via the methodArgs list in the wrapper functions. Default values and direct usage examples for these arguments can be found in the rotations documentation.

The randomStarts argument

Because factor rotation criteria frequently suffer from local minima, trying multiple starting configurations can help identify a superior solution. The randomStarts argument, available exclusively in the GPFRSorth and GPFRSoblq wrappers, facilitates this robust search approach.

• By default, randomStarts = 0, which defaults to using the identity matrix as the initial rotation matrix Tmat. The initial rotation matrix Tmat can also be set by the user.

• Setting randomStarts = 1 initializes Tmat with a single random matrix.

• Setting randomStarts > 1 attempts multiple random starts and returns the rotated loadings matrix that achieved the lowest criterion value f across all attempts. Note that this returned solution is technically still a local minimum, and is not guaranteed to be the global minimum. Users are encouraged to review the random start diagnostics detailed in the package examples.

Under the hood, an internal, unexported engine named .GPA_RS_engine safely manages the random start loop, tracks convergence diagnostics, and handles factor correlation matrix naming.

While the core algorithms GPForth and GPFoblq do not support the randomStarts argument directly, users can manually supply a single random initial rotation matrix to them using Tmat = Random.Start(ncol(A)).

The algorithm argument

The algorithm argument controls how the step size \alpha is initialized at each iteration of the gradient projection loop. Three options are available:

"legacy"

The original Bernaards and Jennrich (2005) heuristic: \alpha is doubled at the start of each iteration and halved until the Armijo sufficient decrease condition is satisfied. Simple and reliable, particularly for smooth orthogonal criteria. This is the default for all orthogonal rotations and for target and PST rotations.

"bb"

Barzilai-Borwein (BB) step size initialization. Uses the ratio of successive changes in the rotation matrix and projected gradient to estimate a good starting step size automatically:

\alpha_{BB} = \frac{\|\Delta T\|^2}{|\langle \Delta T, \Delta G_p \rangle|}

where \Delta T = T_k - T_{k-1} and \Delta G_p = G_{p,k} - G_{p,k-1} are the changes in the rotation matrix and projected gradient between consecutive iterations. The BB estimate is then used as the starting point for the Armijo line search. Benchmarking shows meaningful speedups for oblique rotation criteria on larger matrices (50 variables or more) with many random starts, particularly for non-smooth criteria such as simplimax. This is the default for all oblique rotations except target and PST.

"cayley"

The Cayley transform retraction for orthogonal rotations (Wen and Yin, 2013). Instead of projecting back onto the orthogonal manifold via SVD, the Cayley transform produces an exactly orthogonal update in closed form:

T_{k+1} = \left(I + \frac{\alpha}{2} W\right)^{-1} \left(I - \frac{\alpha}{2} W\right) T_k

where W = G_p T_k^T - T_k G_p^T is a skew-symmetric matrix derived from the projected gradient. The Cayley transform avoids the SVD decomposition at each iteration and can be faster for larger orthogonal problems. Available for orthogonal rotations only.

• algorithm = "bb" with fwindow = 10 is the default for all rotation criteria.

• algorithm = "cayley" is available for orthogonal rotations and can be specified explicitly by the user.

• algorithm = "legacy" is available for exact reproducibility of pre-2027 results.

Users can override the default by passing algorithm explicitly to any rotation function or directly to GPForth and GPFoblq.

The fwindow argument

The fwindow argument controls the width of the non-monotone line search window (Zhang and Hager, 2004). At each iteration, the Armijo sufficient decrease condition is evaluated against the maximum criterion value over the last fwindow iterations rather than just the immediately preceding value:

f(T_{k+1}) < \max_{j \in W_k} f(T_j) - 0.5\, s^2\, \alpha

where W_k is the window of the last fwindow iterates and s is the projected gradient norm.

• fwindow = 1 (default for "legacy") reduces to the standard monotone Armijo condition — the criterion must decrease at every iteration.

• fwindow > 1 allows temporary increases in the criterion value, which can help the algorithm escape flat regions and shallow local minima. This is particularly beneficial in combination with algorithm = "bb".

• fwindow = 10 is the default for algorithm = "bb" and algorithm = "cayley".

Benchmarking across matrix sizes and rotation criteria suggests that fwindow = 10 with algorithm = "bb" provides the best balance of speed and solution quality for oblique rotations on larger models. For orthogonal rotations with algorithm = "legacy", fwindow = 1 is recommended as larger windows can slow convergence.

Legacy functions

The original implementations authored by Bernaards and Jennrich (2005) have been retained as GPForth.legacy and GPFoblq.legacy. These functions are kept purely for historical reference and backward compatibility for reproducibility. They are not exported into the package namespace, meaning they must be explicitly invoked using the triple-colon operator:

GPArotation:::GPForth.legacy(A, method = "varimax")
GPArotation:::GPFoblq.legacy(A, method = "quartimin")

The results generated by these legacy functions should be numerically identical to those produced by the current implementations when algorithm = "legacy" and fwindow = 1 are used.

Value

A GPArotation object which is a list with elements:

Author(s)

Coen A. Bernaards and Robert I. Jennrich with some R modifications by Paul Gilbert

References

Barzilai, J. and Borwein, J.M. (1988) Two-point step size gradient methods. IMA Journal of Numerical Analysis, 8, 141–148. doi:10.1093/imanum/8.1.141

Bernaards, C.A. and Jennrich, R.I. (2005) Gradient Projection Algorithms and Software for Arbitrary Rotation Criteria in Factor Analysis. Educational and Psychological Measurement, 65, 676–696. doi:10.1177/0013164404272507

Jennrich, R.I. (2001) A simple general procedure for orthogonal rotation. Psychometrika, 66, 289–306. doi:10.1007/BF02294840

Jennrich, R.I. (2002) A simple general method for oblique rotation. Psychometrika, 67, 7–19. doi:10.1007/BF02294706

Mansolf, M. and Reise, S.P. (2016) Exploratory Bifactor Analysis: The Schmid-Leiman Orthogonalization and Jennrich-Bentler Analytic Rotations. Multivariate Behavioral Research, 51(5), 698–717. doi:10.1080/00273171.2016.1215898

Nguyen, H.V. and Waller, N.G. (2023) Local minima and factor rotations in exploratory factor analysis. Psychological Methods, 28(5), 1122–1141. doi:10.1037/met0000467

Wen, Z. and Yin, W. (2013) A feasible method for optimization with orthogonality constraints. Mathematical Programming, 142, 397–434. doi:10.1007/s10107-012-0584-1

Zhang, H. and Hager, W.W. (2004) A nonmonotone line search technique and its application to unconstrained optimization. SIAM Journal on Optimization, 14(4), 1043–1056. doi:10.1137/S1052623403428208

See Also

rotations, Random.Start, factanal, Harman8, box26, CCAI, NetherlandsTV

Examples

  # --- Basic rotation calls ---
  data(Harman, package = "GPArotation")           # 8 physical variables
  quartimax(Harman8)                              # direct rotation call
  GPFRSorth(Harman8, method = "quartimax")        # equivalent via wrapper
  GPFRSoblq(Harman8, method = "quartimin", normalize = TRUE)
  loadings(quartimin(Harman8, normalize = TRUE))  # extract loadings directly

  # --- Passing criterion arguments via methodArgs ---
  # Crawford-Ferguson family: kappa selects the criterion.
  # For box26: p = 26 variables, m = 3 factors.
  # Equamax:  kappa = m / (2 * p) = 3 / 52
  # Parsimax: kappa = (m - 1) / (p + m - 2) = 2 / 27
  data(Thurstone, package = "GPArotation") # 26 variable box problem 
  GPFRSoblq(box26, method = "cf", methodArgs = list(kappa = 3/52))  # Equamax
  GPFRSoblq(box26, method = "cf", methodArgs = list(kappa = 2/27))  # Parsimax
  
  # --- Two-step vs single-step factanal for oblique rotation ---
  #
  # The recommended approach for oblique rotation is the two-step procedure:
  # (1) obtain unrotated loadings from factanal, then
  # (2) rotate separately using GPArotation.
  # This gives full control over the rotation, including random starts.
  #
  # Prior to R 4.5.1, the single-step approach (rotation inside factanal)
  # had a bug in factor reordering after oblique rotation. 
  # This was fixed by the R core team in R 4.5.1.

  data("WansbeekMeijer", package = "GPArotation")

  # Step 1: unrotated 3-factor solution
  fa.unrotated <- factanal(factors = 3, covmat = NetherlandsTV,
                           normalize = TRUE, rotation = "none")

  # Step 2: oblique Crawford-Ferguson rotation with kappa = 0.3
  # (non-standard kappa, not corresponding to any named special case)
  set.seed(44)
  fa.cf <- cfQ(fa.unrotated, kappa = 0.3, normalize = TRUE, randomStarts = 100)
  fa.cf

  # Single-step via factanal - correct in R >= 4.5.1
  if (getRversion() >= "4.5.1") {
    set.seed(44)
    fa.factanal <- factanal(factors = 3, covmat = NetherlandsTV, rotation = "cfQ",
      control = list(rotate = list(normalize = TRUE, kappa = 0.3, randomStarts = 100)))
    # The two approaches should agree after sorting
    fa.sorted <- print(fa.cf, sortLoadings = TRUE)
    cat("Maximum difference in loadings between two-step and single-step:\n")
    print(max(abs(abs(fa.sorted$loadings) - abs(fa.factanal$loadings))))
  } else {
    cat("Single-step factanal oblique rotation requires R >= 4.5.1.\n")
    cat("Use the two-step procedure above for correct results.\n")
  }

  # --- Displaying rotation output ---
  origdigits <- options("digits")
  data("CCAI", package = "GPArotation")
  fa.unrotated <- factanal(factors = 3, covmat = CCAI_R, n.obs = 461, rotation = "none")
  res <- oblimin(fa.unrotated, gam = -0.5, randomStarts = 20)
                                    # gam = -0.5: more orthogonal than quartimin
  res                               # default print
  print(res)                        # equivalent to above
  print(res, Table = TRUE)          # include iteration table
  print(res, rotateMat = TRUE)      # include rotating matrix
  print(res, digits = 2)            # rounded to 2 decimal places
  summary(res)                      # pattern and structure matrices for oblique rotation
  summary(res, Structure = FALSE)   # pattern matrix only
  options(digits = origdigits$digits)
    
  # --- Random start diagnostics ---
  # When randomStarts > 1, the output includes randStartChar which summarizes
  # the random start results:
  #   randomStarts : number of random starts attempted
  #   Converged    : number of starts that converged
  #   atMinimum    : number of starts at the same lowest minimum
  #   localMins    : number of distinct local minima found
  data(Thurstone, package = "GPArotation")
  res <- GPFRSoblq(box26, method = "geomin", normalize = TRUE, randomStarts = 50)
  res$randStartChar

  # --- Factor ordering ---
  # Raw GPArotation output is unsorted — factors may appear in any order
  # depending on the starting matrix. Use print() to obtain sorted loadings.
  # Once sorted, repeated calls to print() are stable.
  set.seed(334)
  xusl <- quartimin(Harman8, normalize = TRUE, randomStarts = 100)

  loadings(xusl)                                        # unsorted raw output
  max(abs(print(xusl)$loadings - xusl$loadings)) == 0  # FALSE: print() reorders
  xsl <- print(xusl)                                    # capture sorted result
  max(abs(print(xsl)$loadings - xsl$loadings)) == 0    # TRUE: already sorted
  
  # --- Normalization ---
  # Kaiser normalization
  data("CCAI", package = "GPArotation")
  fa.unrotated <- factanal(factors = 3, covmat = CCAI_R, n.obs = 461, rotation = "none")
  oblimin(fa.unrotated, normalize = TRUE, randomStarts = 100)
    data("CCAI", package = "GPArotation")
  res.u <- factanal(covmat = CCAI_R, factors = 4, rotation = "none",
                    n.obs = 461)

  # TandemI --- non-smooth orthogonal, BB and Cayley both faster
  res.t1l <- tandemI(res.u, algorithm = "legacy", fwindow =  1, maxit = 10000)
  res.t1b <- tandemI(res.u, algorithm = "bb",     fwindow = 10)
  res.t1c <- tandemI(res.u, algorithm = "cayley", fwindow = 10)
  cat("TandemI   legacy:", nrow(res.t1l$Table) - 1, "iterations\n")
  cat("TandemI   bb:    ", nrow(res.t1b$Table) - 1, "iterations",
      " max diff:", format(max(abs(res.t1l$loadings - res.t1b$loadings)),
                           scientific = TRUE, digits = 2), "\n")
  cat("TandemI   cayley:", nrow(res.t1c$Table) - 1, "iterations",
      " max diff:", format(max(abs(res.t1l$loadings - res.t1c$loadings)),
                           scientific = TRUE, digits = 2), "\n")

  # Entropy --- orthogonal, BB and Cayley both faster
  res.t2l <- entropy(res.u, algorithm = "legacy", fwindow =  1, maxit = 10000)
  res.t2b <- entropy(res.u, algorithm = "bb",     fwindow = 10)
  res.t2c <- entropy(res.u, algorithm = "cayley", fwindow = 10)
  cat("Entropy   legacy:", nrow(res.t2l$Table) - 1, "iterations\n")
  cat("Entropy   bb:    ", nrow(res.t2b$Table) - 1, "iterations",
      " max diff:", format(max(abs(res.t2l$loadings - res.t2b$loadings)),
                           scientific = TRUE, digits = 2), "\n")
  cat("Entropy   cayley:", nrow(res.t2c$Table) - 1, "iterations",
      " max diff:", format(max(abs(res.t2l$loadings - res.t2c$loadings)),
                           scientific = TRUE, digits = 2), "\n")

  # Oblimin --- smooth oblique, BB faster, same solution
  res.ol <- oblimin(res.u, algorithm = "legacy", fwindow =  1, maxit = 10000)
  res.ob <- oblimin(res.u, algorithm = "bb",     fwindow = 10)
  cat("Oblimin   legacy:", nrow(res.ol$Table) - 1, "iterations\n")
  cat("Oblimin   bb:    ", nrow(res.ob$Table) - 1, "iterations",
      " max diff:", format(max(abs(res.ol$loadings - res.ob$loadings)),
                           scientific = TRUE, digits = 2), "\n")

  # Simplimax --- non-smooth oblique with many local minima.
  # BB finds a substantially better solution in fewer iterations.
  # Difference in loadings reflects different local minima, not noise.
  res.sl <- simplimax(res.u, algorithm = "legacy", fwindow =  1, maxit = 10000)
  res.sb <- simplimax(res.u, algorithm = "bb",     fwindow = 10)
  cat("Simplimax legacy:", nrow(res.sl$Table) - 1, "iterations",
      " f =", round(res.sl$Table[nrow(res.sl$Table), 2], 6), "\n")
  cat("Simplimax bb:    ", nrow(res.sb$Table) - 1, "iterations",
      " f =", round(res.sb$Table[nrow(res.sb$Table), 2], 6), "\n")
  cat("BB found a lower criterion value --- better solution quality.\n")

Griffith and Mulaik Interpersonal Personality Traits

Description

Correlation matrix for 24 bipolar seven-point trait-adjective scales representing six hypothetical factors from the interpersonal domain of personality. Data from an unpublished senior thesis conducted under the supervision of Stanley Mulaik at the School of Psychology, Georgia Institute of Technology, 1998.

Usage

data(GriffithMulaik, package = "GPArotation")

Format

A 24 \times 24 numeric correlation matrix with n.obs = 523. Row and column names are abbreviated variable labels for 24 interpersonal trait scales: ALLOW, OUTGO, NONTALK, PERMISS, FAIR, UNFRIEND, FORCELSS, RESTRAIN, NONDOMIN, UNLOVING, NONCNFRM, COMPLINT, IMPARTL, LEADER, UNKIND, UNCNVENT, DIRECTNG, PROHIBTV, UNGREGAR, JUST, NEUTRAL, SOCIABLE, UNAFFECT, REBELLOS.

Details

Participants were 523 Georgia Tech introductory psychology students who each rated three stereotypical persons selected at random from a large list of stereotypes on 52 bipolar seven-point trait-adjective scales. The 24 variables analyzed here represent six hypothetical factors from the interpersonal domain of personality: four synonym trait items per factor.

The study was intended to test Mulaik's theory of personality trait factors. Six simple structure factors were obtained as expected, providing tentative support for the theory.

The correlation matrix is reproduced as Table 8.5 in Mulaik (2018). The data are used in that chapter to illustrate maximum likelihood factor analysis with a reduced correlation matrix scree plot.

The six-factor oblimin solution obtained by GPArotation closely replicates the Direct Oblimin solution reported in Mulaik (2018) Table 8.6, with minor differences attributable to factor ordering and sign conventions.

References

Griffith, D. and Mulaik, S.A. (1998). Personality trait factors from the interpersonal domain. Poster presented at the Society for Multivariate Experimental Psychology (SMEP).

Mulaik, S.A. (2018). Fundamentals of common factor analysis. In P. Irwing, T. Booth, and D.J. Hughes (Eds.), The Wiley Handbook of Psychometric Testing (pp. 211–252). Wiley.

See Also

rotations, Harman, Thurstone, CCAI

Examples

  data(GriffithMulaik, package = "GPArotation")

  ## Scree plot comparing raw versus reduced correlation matrix
  ## (Mulaik, 2018, recommends reduced matrix for ML factor analysis)
  fa.un <- factanal(factors = 6, covmat = GriffithMulaik,
                    n.obs = 523, rotation = "none")
  res.un <- quartimax(fa.un)

  ## Six factor oblimin rotation
  res.ob <- oblimin(fa.un, randomStarts = 100)
  summary(res.ob)

  ## Compare simple structure across rotation methods
  res.vm <- Varimax(fa.un)
  res.gm <- geominQ(fa.un, randomStarts = 100)

  cat(sprintf("%-12s  AUC = %.3f  Hyperplane = %.1f pct\n",
              "Varimax",  GPArotation:::calc_AUC(res.vm)$AUC_mean,
              GPArotation:::calc_hyperplane(res.vm)$HP_pct))
  cat(sprintf("%-12s  AUC = %.3f  Hyperplane = %.1f pct\n",
              "Oblimin",  GPArotation:::calc_AUC(res.ob)$AUC_mean,
              GPArotation:::calc_hyperplane(res.ob)$HP_pct))
  cat(sprintf("%-12s  AUC = %.3f  Hyperplane = %.1f pct\n",
              "GeominQ",  GPArotation:::calc_AUC(res.gm)$AUC_mean,
              GPArotation:::calc_hyperplane(res.gm)$HP_pct))

  ## --- Theory-guided rotation via partially specified target ---
  ## Mulaik (2018) Table 8.7 specifies a CFA hypothesis matrix:
  ## * = free parameter (NA), o = constrained to zero (0)
  GM_target <- matrix(c(
  ## F1  F2  F3  F4  F5  F6
    NA,  0,  0,  0, NA,  0,  ## ALLOW
    NA, NA,  0,  0, NA,  0,  ## OUTGO
     0, NA,  0,  0, NA,  0,  ## NONTALK
    NA,  0,  0,  0,  0,  0,  ## PERMISS
     0,  0, NA, NA,  0,  0,  ## FAIR
    NA, NA, NA, NA,  0,  0,  ## UNFRIEND
    NA, NA,  0, NA, NA,  0,  ## FORCELSS
    NA, NA,  0,  0, NA,  0,  ## RESTRAIN
    NA, NA,  0, NA, NA,  0,  ## NONDOMIN
     0,  0,  0, NA,  0,  0,  ## UNLOVING
     0,  0,  0,  0,  0, NA,  ## NONCNFRM
    NA,  0,  0,  0,  0, NA,  ## COMPLINT
     0, NA, NA,  0, NA,  0,  ## IMPARTL
     0,  0,  0,  0, NA,  0,  ## LEADER
     0, NA, NA, NA,  0,  0,  ## UNKIND
     0,  0,  0,  0,  0, NA,  ## UNCNVENT
     0,  0,  0, NA, NA,  0,  ## DIRECTNG
    NA,  0,  0,  0, NA, NA,  ## PROHIBTV
     0, NA,  0,  0,  0,  0,  ## UNGREGAR
     0,  0, NA, NA,  0,  0,  ## JUST
    NA, NA, NA, NA, NA,  0,  ## NEUTRAL
    NA, NA,  0, NA, NA, NA,  ## SOCIABLE
     0, NA,  0, NA,  0,  0,  ## UNAFFECT
     0, NA,  0, NA,  0, NA   ## REBELLOS
  ), nrow = 24, ncol = 6, byrow = TRUE)
  rownames(GM_target) <- rownames(GriffithMulaik)
  colnames(GM_target) <- paste0("F", 1:6)

  ## W: 1 = constrained to zero, 0 = free
  GM_target_pst <- ifelse(is.na(GM_target), 0, GM_target)
  GM_W          <- ifelse(is.na(GM_target), 0, 1)

  ## Oblique PST closely replicates Mulaik (2018) Table 8.8 CFA solution
  res.pstQ <- pstQ(fa.un, Target = GM_target_pst, W = GM_W)
  print(res.pstQ)

  ## Compare data-driven vs theory-guided rotation
  cat(sprintf("%-10s  AUC = %.3f  Hyperplane = %.1f pct\n",
              "PST-Q",   GPArotation:::calc_AUC(res.pstQ)$AUC_mean,
              GPArotation:::calc_hyperplane(res.pstQ)$HP_pct))
  cat(sprintf("%-10s  AUC = %.3f  Hyperplane = %.1f pct\n",
              "Oblimin", GPArotation:::calc_AUC(res.ob)$AUC_mean,
              GPArotation:::calc_hyperplane(res.ob)$HP_pct))

Harman's Eight Physical Variables

Description

Unrotated factor loading matrix extracted by the centroid method for Harman's 8 physical variables, a classic dataset in factor analysis.

Usage

  data(Harman)

Format

Harman8: a 8 \times 2 numeric matrix of unrotated centroid factor loadings. Row names are variable names; column names are CF1 and CF2.

Details

Harman8 is a 8 \times 2 matrix of unrotated factor loadings extracted by the centroid method (Table 14.3 in Harman, 1976) from the correlation matrix of 8 physical measurements taken on 305 girls (Table 2.3 in Harman, 1976). The centroid method was a standard extraction procedure before iterative methods became computationally feasible. The two columns represent unrotated centroid factors (CF1 and CF2). Row names are the 8 variable names.

The covariance matrix for these 8 physical variables is available in base R (see Harman23.cor).

Source

Harman, H.H. (1976). Modern Factor Analysis (3rd ed.). University of Chicago Press.

See Also

Harman23.cor, GPForth, rotations, Thurstone, WansbeekMeijer, CCAI, GriffithMulaik

Examples

  data(Harman, package = "GPArotation")

  # Centroid unrotated loadings from Harman (1976) Table 14.3
  print(Harman8)

  # Rotate the centroid solution
  quartimax(Harman8)
  oblimin(Harman8, randomStarts = 100)

Normalizing Weights for Factor Rotation

Description

Internal utility function that computes normalizing weights for factor loading matrices prior to rotation. Called by GPForth, GPFoblq, and the random-start wrappers.

Usage

    NormalizingWeight(A, normalize=FALSE)

Arguments

Details

NormalizingWeight is not exported from the NAMESPACE and is called internally by GPForth, GPFoblq, and the random-start wrapper functions. For a full description of the normalize argument and its options, see GPFRSorth.

The choice of normalization method can affect the rotation solution and its interpretation. For a detailed investigation of the effects of normalization on factor rotations, see Nguyen and Waller (2023).

Value

A numeric vector of normalizing weights. This function is not exported from the NAMESPACE and is only called internally by the gradient projection rotation functions. See GPFRSorth for details on the normalize argument.

References

Cureton, E.E. and Mulaik, S.A. (1975). The weighted varimax rotation and the promax rotation. Psychometrika, 40(2), 183–195. doi:10.1007/BF02291565

Nguyen, H.V. and Waller, N.G. (2023). Local minima and factor rotations in exploratory factor analysis. Psychological Methods, 28(5), 1122–1141. doi:10.1037/met0000467

See Also

GPFRSorth, GPForth

Examples

  data("CCAI", package = "GPArotation")
  # Kaiser normalization
  factanal(factors = 3, covmat = CCAI_R, n.obs = 461, rotation = "oblimin",
           control = list(rotate = list(normalize = TRUE)))

  data(Harman, package = "GPArotation")
  # two ways to do Kaiser normalization
  quartimin(Harman8, normalize = TRUE, randomStarts = 100)
  quartimin(Harman8, normalize = "Kaiser", randomStarts = 100)
  # Cureton-Mulaik normalization
  quartimin(Harman8, normalize = "CM", randomStarts = 100)

Random Starting Matrix for Factor Rotation

Description

Generates a random orthogonal matrix for use as an initial rotation matrix (Tmat) in gradient projection rotation functions.

Usage

    Random.Start(k = 2L)

Arguments

Details

The function generates a random orthogonal matrix using QR decomposition of a matrix of standard normal variates, with a sign correction applied to the diagonal of R to ensure uniform sampling from the Haar measure. This follows the approach of Stewart (1980) and Mezzadri (2007).

The naive approach of qr.Q(qr(matrix(rnorm(k*k), k))) does not guarantee uniform sampling from the Haar measure. The sign correction Q %*% diag(sign(diag(R))) is required to achieve this. This was updated in GPArotation 2024.2-1 following a suggestion by Yves Rosseel.

For oblique rotation a random starting transformation matrix can be generated by normalizing the columns of a random matrix: X %*% diag(1/sqrt(diag(crossprod(X)))) where X <- matrix(rnorm(k*k), k).

Value

A k \times k orthogonal matrix drawn uniformly from the Haar measure on the orthogonal group O(k). Columns have unit length and are mutually orthogonal.

Author(s)

Coen A. Bernaards and Robert I. Jennrich with some R modifications by Paul Gilbert. Updated following a suggestion by Yves Rosseel.

References

Stewart, G.W. (1980). The efficient generation of random orthogonal matrices with an application to condition estimators. SIAM Journal on Numerical Analysis, 17(3), 403–409. doi:10.1137/0717034

Mezzadri, F. (2007) How to generate random matrices from the classical compact groups. Notices of the American Mathematical Society, 54(5), 592–604. arXiv:math-ph/0609050

See Also

GPFRSorth, GPFRSoblq, GPForth, GPFoblq, rotations

Examples

  # Generate a 5 x 5 random orthogonal matrix
  Random.Start(5)

  # Generate a random orthogonal start matrix
  Q <- Random.Start(4)

  # Verify orthogonality: t(Q) %*% Q should equal the identity matrix
  stopifnot(isTRUE(all.equal(t(Q) %*% Q, diag(4), tolerance = 1e-10)))
  
  # Use as starting matrix for rotation
  data("Thurstone", package = "GPArotation")
  simplimax(box26, Tmat = Random.Start(3))

Thurstone Box Data

Description

Factor loading matrix derived from Thurstone's box data.

Usage

  data(Thurstone)

Format

A numeric matrix:

• box26: 26 x 3 matrix of unrotated factor loadings.

Details

Loading the data makes the following object available:

box26 is a 26 \times 3 unrotated factor loading matrix for 26 variables measured on a set of boxes. It appears frequently in the rotation literature as a dataset for comparing rotation criteria, particularly for assessing local minima.

The matrix is suitable as input to GPForth, GPFoblq, and all rotation wrapper functions in GPArotation. For a richer collection of factor analysis datasets, see the psych package.

Source

Thurstone, L.L. (1947). Multiple Factor Analysis. University of Chicago Press.

See Also

GPForth, rotations, Harman, CCAI, WansbeekMeijer, GriffithMulaik

Netherlands Television Viewership Data

Description

Correlation matrix for Netherlands television viewership, used as an example dataset for factor rotation. From Wansbeek and Meijer (2000), page 171.

Usage

  data(WansbeekMeijer)

Format

NetherlandsTV is a list with components:

• $cov: a 7 \times 7 numeric correlation matrix

• $n.obs: sample size (2154)

Details

NetherlandsTV is a list with components $cov (a 7 \times 7 correlation matrix for 7 television viewership variables measured in the Netherlands) and $n.obs (sample size, 2154). Use cov2cor(NetherlandsTV\$cov) to obtain the correlation matrix, or pass NetherlandsTV directly to factanal which handles the list structure automatically. It is used throughout the GPArotation documentation and vignettes as an example dataset for oblique rotation with 2 or 3 factors.

Source

Wansbeek, T. and Meijer, E. (2000). Measurement Error and Latent Variables in Econometrics. North-Holland.

See Also

GPForth, rotations, Thurstone, Harman, CCAI, GriffithMulaik

Examples

  data(WansbeekMeijer, package = "GPArotation")

  # Correlation matrix
  round(cov2cor(NetherlandsTV$cov), 2)

  # factanal picks up n.obs automatically from the list
  factanal(factors = 2, covmat = NetherlandsTV, rotation = "none")

  # Two-step oblique rotation
  fa.unrotated <- factanal(factors = 3, covmat = NetherlandsTV,
                           rotation = "none")
  oblimin(loadings(fa.unrotated), randomStarts = 100)

Echelon Rotation

Description

Rotates a factor loading matrix to an echelon parameterization.

Usage

    echelon(A, reference = seq(ncol(A)), ...)

Arguments

Details

The loading matrix is rotated so that the k rows indicated by reference form the Cholesky factorization given by t(chol(L[reference,] %*% t(L[reference,]))). This defines the rotation transformation, which is then applied to all rows to give the rotated loading matrix.

The optimization is not iterative and does not use the gradient projection algorithm. The function can be used directly or passed to factor analysis functions like factanal via the rotation argument.

This parameterization has several useful properties:

1. It can be useful for comparison with published results in this parameterization.

2. Standard errors are more straightforward to compute because the solution corresponds to an unconstrained optimization.

3. Models with k and k+1 factors are nested, making it straightforward to test the k-factor model versus the (k+1)-factor model. In particular, the Wald test and LM test can be used in addition to the LR test. The test of a k-factor model versus a (k+1)-factor model is a joint test of whether all free parameters (loadings) in the (k+1)st column are zero.

4. For some purposes, only the subspace spanned by the factors matters, not the specific parameterization within this subspace.

5. Back-predicted indicators (the explained portion of the indicators) do not depend on the rotation method. Combined with the greater ease of obtaining correct standard errors, this allows easier and more accurate prediction standard errors.

6. This parameterization and its standard errors can be used to detect identification problems (McDonald, 1999, pp. 181–182).

One use of echelon rotation is obtaining good starting values for subsequent rotation, though it may seem counterintuitive to rotate towards this solution afterwards.

Value

A GPArotation object which is a list with elements:

Author(s)

Erik Meijer and Paul Gilbert.

References

McDonald, R.P. (1999). Test Theory: A Unified Treatment. Erlbaum.

Wansbeek, T. and Meijer, E. (2000). Measurement Error and Latent Variables in Econometrics. North-Holland.

See Also

eiv, rotations, GPForth, GPFoblq

Examples

data("WansbeekMeijer", package = "GPArotation")
fa.un <- factanal(factors = 2, covmat = NetherlandsTV,
                  rotation = "none", n.obs = 2154)

# Echelon rotation — first k rows as reference (default)
fa.ech <- echelon(fa.un)
print(fa.ech)

# Equivalent via factanal rotation argument
fa.ech2 <- factanal(factors = 2, covmat = NetherlandsTV,
                    rotation = "echelon")

# Echelon rotation with a different reference set
fa.ech3 <- echelon(fa.un, reference = 6:7)
print(fa.ech3)

Errors-in-Variables Rotation

Description

Rotates a factor loading matrix to an errors-in-variables representation.

Usage

    eiv(A, identity = seq(ncol(A)), ...)

Arguments

Details

The loading matrix is rotated so that the k rows indicated by identity form an identity matrix, with the remaining M-k rows as free parameters. \Phi is also free.

The optimization is not iterative and does not use the gradient projection algorithm. The function can be used directly or passed to factor analysis functions like factanal via the rotation argument.

Viewed as a rotation method it is oblique, with an explicit solution. Given an initial loadings matrix L partitioned as L = (L_1^T, L_2^T)^T, the rotated loadings matrix is (I, (L_2 L_1^{-1})^T)^T and \Phi = L_1 L_1^T, where I is the k \times k identity matrix. It is assumed that \Phi = I for the initial loadings matrix.

Not all authors consider this representation to be a rotation in the strict sense.

This parameterization has several useful properties:

1. It can be useful for comparison with published results in this parameterization.

2. Standard errors are more straightforward to compute because the solution corresponds to an unconstrained optimization.

3. One may have prior knowledge about which reference variables load on only one factor without imposing restrictive constraints on other loadings — in this sense it has similarities to CFA.

4. For some purposes, only the subspace spanned by the factors matters, not the specific parameterization within this subspace.

5. Back-predicted indicators (the explained portion of the indicators) do not depend on the rotation method. Combined with the greater ease of obtaining correct standard errors, this allows easier and more accurate prediction standard errors.

One use of this parameterization is obtaining good starting values for subsequent rotation, though it may seem counterintuitive to rotate towards this solution afterwards.

Note that the rotated loadings for non-identity rows may be large in magnitude and should not be interpreted as standard factor loadings. \Phi is the factor covariance matrix L_1 L_1^T and is not a correlation matrix — diagonal values are not constrained to 1.

Value

A GPArotation object which is a list with elements:

Author(s)

Erik Meijer and Paul Gilbert.

References

Hägglund, G. (1982). Factor analysis by instrumental variables methods. Psychometrika, 47, 209–222. doi:10.1007/BF02296276

Lewin-Koh, S.C. and Amemiya, Y. (2003) Heteroscedastic factor analysis. Biometrika, 90(1), 85–97. doi:10.1093/biomet/90.1.85

Wansbeek, T. and Meijer, E. (2000). Measurement Error and Latent Variables in Econometrics. North-Holland.

See Also

echelon, rotations, GPForth, GPFoblq

Examples

data("WansbeekMeijer", package = "GPArotation")
fa.un <- factanal(factors = 2, covmat = NetherlandsTV,
                  rotation = "none")

# Eiv rotation: items 1 and 2 define the identity
fa.eiv <- eiv(fa.un)

# Equivalent via factanal rotation argument
fa.eiv2 <- factanal(factors = 2, covmat = NetherlandsTV,
                    rotation = "eiv")

# Rotated loadings
print(fa.eiv)

# Eiv rotation with a different identity set (items 6 and 7)
fa.eiv3 <- eiv(fa.un, identity = 6:7)
print(fa.eiv3)

L^p Rotation

Description

Performs L^p rotation to obtain sparse factor loadings.

Usage

GPForth.lp(A, Tmat = diag(ncol(A)), p = 1, normalize = FALSE, eps = 1e-05,
           maxit = 10000, gpaiter = 5)
GPFoblq.lp(A, Tmat = diag(ncol(A)), p = 1, normalize = FALSE, eps = 1e-05,
           maxit = 10000, gpaiter = 5)

Arguments

Details

These functions optimize an L^p rotation objective where 0 < p <= 1. A smaller value of p promotes greater sparsity in the loading matrix but increases computational difficulty. For guidance on choosing p, see the Concluding Remarks in Liu et al. (2023).

The user-facing wrapper functions lpT and lpQ provide random start functionality on top of GPForth.lp and GPFoblq.lp respectively, analogous to how GPFRSorth and GPFRSoblq wrap GPForth and GPFoblq for standard rotation criteria. For most users lpT and lpQ are the recommended entry points.

Since the L^p objective is nonsmooth, a different optimization method is required compared to smooth rotation criteria. GPForth.lp and GPFoblq.lp replace GPForth and GPFoblq for orthogonal and oblique L^p rotations, respectively.

The optimization uses an iterative reweighted least squares (IRLS) approach. The nonsmooth objective is approximated by a smooth weighted least squares function in the outer loop, which is then optimized using GPA in the inner loop (gpaiter controls the maximum inner iterations).

Normalization is not recommended for L^p rotation and may produce unexpected results.

Value

A GPArotation object, which is a list containing:

Author(s)

Xinyi Liu, with minor modifications for GPArotation by C. Bernaards.

References

Liu, X., Wallin, G., Chen, Y., and Moustaki, I. (2023). Rotation to sparse loadings using L^p losses and related inference problems. Psychometrika, 88(2), 527–553. doi:10.1007/s11336-023-09911-y

See Also

lpT, lpQ, vgQ.lp.wls

Examples

  data("WansbeekMeijer", package = "GPArotation")
  fa.unrotated <- factanal(factors = 2, covmat = NetherlandsTV, rotation = "none")
  options(warn = -1)

  # Orthogonal rotation — single start
  fa.lpT1 <- GPForth.lp(loadings(fa.unrotated), p = 1)

  # Orthogonal rotation — 10 random starts
  fa.lpT <- lpT(loadings(fa.unrotated), Tmat = Random.Start(2), p = 1,
                randomStarts = 10)
  print(fa.lpT, digits = 5, sortLoadings = FALSE, Table = TRUE, rotateMat = TRUE)

  # Oblique rotation — single start
  fa.lpQ1 <- GPFoblq.lp(loadings(fa.unrotated), p = 1)

  # Oblique rotation — 10 random starts
  fa.lpQ <- lpQ(loadings(fa.unrotated), p = 1, randomStarts = 10)
  summary(fa.lpQ, Structure = TRUE)

  # Compare Lp (p=1), Lp (p=0.5), and Geomin oblique rotations
  set.seed(1020)
  fa.lpQ1   <- lpQ(loadings(fa.unrotated), p = 1,   randomStarts = 10)
  fa.lpQ0.5 <- lpQ(loadings(fa.unrotated), p = 0.5, randomStarts = 10)
  fa.geo    <- geominQ(loadings(fa.unrotated),       randomStarts = 10)

  # With factor ordering using internal sortGPALoadings helper
  res <- round(cbind(GPArotation:::.sortGPALoadings(fa.lpQ1)$loadings,
                     GPArotation:::.sortGPALoadings(fa.lpQ0.5)$loadings,
                     GPArotation:::.sortGPALoadings(fa.geo)$loadings), 3)
  print(c("oblique --  Lp p=1        Lp p=0.5       Geomin"))
  print(res)

  # Without factor ordering
  res <- round(cbind(fa.lpQ1$loadings, fa.lpQ0.5$loadings, fa.geo$loadings), 3)
  print(c("oblique --  Lp p=1        Lp p=0.5       Geomin"))
  print(res)

  options(warn = 0)

Plot Method for GPArotation Objects

Description

Produces visual summaries of rotated factor matrices, including an APA-style salient table-style, a heatmap, a tripartite pairs plot of factor geometry, a macro-level profile bar chart, a target rotation schema, a rotation trajectory view, or a set of diagnostics.

Usage

  ## S3 method for class 'GPArotation'
plot(x, type = c("salient", "pairs", "profile", "vector", 
     "target", "heatmap", "trajectory", "diagnostics", "residuals"), 
     hide_hyperplane = FALSE, salient = 0.40, cutoff = .1, R = NULL, 
     factors = c(1,2), items = NULL, main = NULL, ...)
     
  

Arguments

Value

Invisibly returns the original x object.

Author(s)

Your Name

Examples

  ## Not run: 
    # Assuming 'res' is a fitted GPArotation object:
    plot(res, type = "salient")
    plot(res, type = "pairs")
    plot(res, type = "profile")
  
## End(Not run)

Algorithm Comparison for Two-Factor Orthogonal Rotation

Description

Produces a 2 \times 3 panel figure comparing the "legacy", "bb", and "cayley" optimization algorithms on a two-factor orthogonal rotation. The top row shows the criterion landscape and the algorithm trajectory for each algorithm; the bottom row shows the corresponding factor migration paths. An optional endgame inset zooms into the convergence region near the global minimum.

Usage

plot2fOrthComparison(A, method = "quartimax", items = NULL,
                     nang = 6000, maxit = 500, magnify = FALSE, 
                     inset_pos = "auto", randomStart = FALSE, ...)

Arguments

Value

Invisibly returns a list with components legacy, bb, and cayley, each a "GPArotation" object from the corresponding algorithm run.

Note

This function requires a two-factor loading matrix. The criterion landscape visualization is only defined for two-factor orthogonal rotations, where all rotation matrices are parameterized by a single angle \theta \in [0, 2\pi).

See Also

plot.GPArotation, GPForth, rotations

Examples

data("Harman", package = "GPArotation")
plot2fOrthComparison(Harman8, method = "quartimax", magnify = TRUE)

plot2fOrthComparison(Harman8, method = "tandemII", magnify = FALSE)

Print and Summary Methods for GPArotation Class Objects

Description

Print and summary methods for objects returned by GPFRSorth, GPFRSoblq, GPForth, or GPFoblq. Both print.GPArotation and summary.GPArotation apply consistent factor sorting and sign correction via the internal helper .sortGPALoadings when sortLoadings = TRUE. Factors are ordered by descending variance explained and signs are adjusted so that the sum of loadings per factor is positive. This convention matches that used by factanal. For oblique rotations, summary.GPArotation displays both the pattern matrix (regression coefficients of items on factors, controlling for factor intercorrelations) and the structure matrix (loadings %*% Phi, correlations between items and factors) when Structure = TRUE. The two matrices coincide for orthogonal rotations where \Phi = I. Output includes contributions of factors via SS loadings (sum of squared loadings); see Harman (1976), sections 2.4 and 12.4. For orthogonal rotations, Proportion Var and Cumulative Var are also shown. For oblique rotations these are suppressed as they are not meaningful when factors are correlated. Simple structure quality is reported via two criterion-free per-factor measures: AUC (Liu and Moustaki, 2024) and FSI (Lorenzo-Seva, 2003). Higher values indicate cleaner simple structure. These measures are shown for both orthogonal and oblique rotations and can be used to compare simple structure quality across rotation methods.

Usage

    ## S3 method for class 'GPArotation'
print(x, digits = 3, sortLoadings = TRUE, cutoff = 0.1,
                    rotateMat = FALSE, Table = FALSE, ...)
    ## S3 method for class 'GPArotation'
summary(object, digits = 3, Structure = TRUE, ...)
    ## S3 method for class 'summary.GPArotation'
print(x, cutoff = 0.1, ...)

Arguments

Details

Factor sorting and sign correction are applied consistently in both print and summary via the internal function .sortGPALoadings, adapted from factanal sorting conventions (R Core Team). This ensures that the pattern matrix shown by summary is consistent with the loadings shown by print.

The digits argument controls the number of decimal places shown in the loadings, structure matrix, and Phi.

Display options such as cutoff are controlled at print time: print(summary(x), cutoff = 0.4), consistent with standard R convention.

Accessing the rotated loading matrix

The rotated loading matrix can be accessed directly via x\$loadings for use in further computations. The formatted output with factor sorting, AUC, FSI, and SS loadings is produced by print(x) or simply by typing the object name at the console. The loadings function from the stats package can also be used but returns the unsorted matrix without simple structure measures.

Convergence and algorithm information

The print header reports convergence status, the algorithm used (algorithm), and the non-monotone line search window (fwindow). If convergence was not obtained, a diagnostic suggestion is provided based on the algorithm currently in use. See GPForth and GPFoblq for details on algorithm options.

For examples see GPFRSorth and the package vignettes: vignette("GPA1guide", package = "GPArotation").

Value

print.GPArotation returns the sorted GPArotation object invisibly when sortLoadings = TRUE, or the unsorted object when sortLoadings = FALSE. summary.GPArotation returns a summary.GPArotation object with sorted loadings and, for oblique rotations, the structure matrix. print.summary.GPArotation returns the object invisibly.

References

Harman, H.H. (1976) Modern Factor Analysis, 3rd ed. University of Chicago Press.

Liu, X., Wallin, G., Chen, Y., and Moustaki, I. (2023). Rotation to sparse loadings using L^p losses and related inference problems. Psychometrika, 88(2), 527–553. doi:10.1007/s11336-023-09911-y

Lorenzo-Seva, U. (2003) A factor simplicity index. Psychometrika, 68(1), 49–60. doi:10.1007/BF02296652

See Also

GPFRSorth, GPForth, factanal, summary

Examples

  data(Harman, package = "GPArotation")
  res <- oblimin(Harman8, normalize = TRUE, randomStarts = 100)

  # Print sorted loadings (default)
  print(res)

  # Print unsorted loadings
  print(res, sortLoadings = FALSE)

  # Summary with pattern and structure matrices
  summary(res, Structure = TRUE)

  # Summary without structure matrix
  summary(res, Structure = FALSE)

  # Print with iteration table
  print(res, Table = TRUE)

Rotations Functions Using Gradient Projection Algorithms

Description

Optimize factor loading rotation objective.

Usage

    oblimin(A, Tmat = NULL, gam=0, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    quartimin(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    targetT(A, Tmat = NULL, Target = NULL, normalize = FALSE, randomStarts = 0, 
        algorithm = "bb", fwindow = 10, ...) 
    targetQ(A, Tmat = NULL, Target = NULL, normalize = FALSE, randomStarts = 0, 
        algorithm = "bb", fwindow = 10, ...)
    pstT(A, Tmat = NULL, W = NULL, Target = NULL, normalize = FALSE,
        randomStarts = 0, algorithm = "bb", fwindow = 10, ...)
    pstQ(A, Tmat = NULL, W = NULL, Target = NULL, normalize = FALSE,
        randomStarts = 0, algorithm = "bb", fwindow = 10, ...)
    oblimax(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    entropy(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    quartimax(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    Varimax(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    simplimax(A, Tmat = NULL, k=NULL, normalize=FALSE, 
        randomStarts=0, algorithm = "bb", fwindow = 10,...)
    bentlerT(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    bentlerQ(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    tandemI(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    tandemII(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    geominT(A, Tmat = NULL, delta=0.01, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    geominQ(A, Tmat = NULL, delta=0.01, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    bigeominT(A, Tmat = NULL, delta=0.01, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    bigeominQ(A, Tmat = NULL, delta=0.01, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    cfT(A, Tmat = NULL, kappa=0, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    cfQ(A, Tmat = NULL, kappa=0, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    equamax(A, Tmat = NULL, kappa=NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    parsimax(A, Tmat = NULL, kappa=NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    infomaxT(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    infomaxQ(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    mccammon(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    varimin(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    bifactorT(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    bifactorQ(A, Tmat = NULL, normalize=FALSE, randomStarts=0, 
        algorithm = "bb", fwindow = 10, ...)
    lpT(A, Tmat=diag(ncol(A)), p=1, normalize=FALSE, eps=1e-05, maxit=2000, 
            randomStarts=0, gpaiter=5) 
    lpQ(A, Tmat=diag(ncol(A)), p=1, normalize=FALSE, eps=1e-05, maxit=2000, 
        randomStarts=0, gpaiter=5)
    

Arguments

Details

These functions optimize a rotation objective. They can be used directly or the function name can be passed to factor analysis functions like factanal. Several of the function names end in T or Q, which indicates if they are orthogonal or oblique rotations (using GPFRSorth or GPFRSoblq respectively). The gradient projection algorithms are described in Bernaards and Jennrich (2005).

The Varimax implementation in the list uses the gradient projection algorithm applied to vgQ.varimax. This implementation is different that the varimax rotation defined in the stats package. Additionally, varimax does Kaiser normalization by default whereas GPArotation::Varimax does not.

The argument kappa parameterizes the family for the Crawford-Ferguson method. If m is the number of factors and p is the number of indicators then kappa values having special names are 0=Quartimax, 1/p=Varimax, m/(2*p)=Equamax, (m-1)/(p+m-2)=Parsimax, 1=Factor parsimony.

Bifactor rotations, bifactorT and bifactorQ are called bifactor and biquartimin in Jennrich and Bentler (2011). For a comparison of exploratory bifactor analysis algorithms including those implemented here, see Garcia-Garzon, Abad and Garrido (2021).

The argument p is needed for L^p rotation. See Lp rotation for details on the rotation method.

update follows standard R convention — use it to modify arguments (e.g. randomStarts, normalize) not to change the rotation criterion. To re-rotate with a different criterion, use GPFRSorth or GPFRSoblq with attr(object, "A_unrotated").

Value

A GPArotation object which is a list with elements:

Author(s)

Coen A. Bernaards and Robert I. Jennrich with some R modifications by Paul Gilbert.

References

Bernaards, C.A. and Jennrich, R.I. (2005) Gradient Projection Algorithms and Software for Arbitrary Rotation Criteria in Factor Analysis. Educational and Psychological Measurement, 65, 676–696. doi:10.1177/0013164404272507

Bi, Y. and Barchard, K.A. (2024). Purchasing choices that reduce climate change: An exploratory factor analysis. Spectra Undergraduate Research Journal, 3(2), 8–14. doi:10.9741/2766-7227.1028

Fischer, R., & Fontaine, J. (2010). Methods for investigating structural equivalence. In D. Matsumoto & F. van de Vijver (Eds.), Cross-Cultural Research Methods in Psychology (pp. 179–215). Cambridge University Press. doi:10.1017/CBO9780511779381.010

Garcia-Garzon, E., Abad, F.J. and Garrido, L.E. (2021). On omega hierarchical estimation: A comparison of exploratory bi-factor analysis algorithms. Multivariate Behavioral Research, 56(1), 101–119. doi:10.1080/00273171.2020.1736977

Jennrich, R.I. and Bentler, P.M. (2011). Exploratory bi-factor analysis. Psychometrika, 76(4), 537–549. doi:10.1007/s11336-011-9218-4

For references to individual rotation criteria see vignette("GPA1guide", package = "GPArotation").

See Also

factanal, GPFRSorth, GPFRSoblq, vgQ, Harman8, NetherlandsTV, CCAI, box26

Examples

  # For extended examples see the vignettes:
  # vignette("GPA1guide",    package = "GPArotation")
  # vignette("GPA2local",    package = "GPArotation")
  # vignette("GPA3bifactor", package = "GPArotation")

  # --- Accessing rotated loadings ---
  data("Harman", package = "GPArotation") # 8 physical variables
  qHarman <- quartimax(Harman8)
  loadings(qHarman)                              # via extractor (recommended)
  qHarman$loadings                               # via direct list access
  all.equal(loadings(qHarman), qHarman$loadings) # identical

  # --- Rotating factanal loadings ---
  data("WansbeekMeijer", package = "GPArotation") # Netherlands TV viewership
  fa.unrotated <- factanal(factors = 2, covmat = NetherlandsTV,
                           normalize = TRUE, rotation = "none")
  quartimax(fa.unrotated, normalize = TRUE)
  geominQ(fa.unrotated, normalize = TRUE, randomStarts = 100)

  # --- Passing rotation to factanal ---
  # CCAI:Climate-Friendly Purchasing Choices domain of the Climate Change Action Inventory
  data("CCAI", package = "GPArotation") 
  factanal(factors = 3, covmat = CCAI_R, n.obs = 461, rotation = "infomaxT")
  factanal(factors = 3, covmat = CCAI_R, n.obs = 461, rotation = "infomaxT",
           control = list(rotate = list(normalize = TRUE, eps = 1e-6)))
           
  # --- Target rotation ---
  # Orthogonal target rotation of two varimax rotated matrices 
  # towards each other. Data from Fischer and Fontaine (2010).
  # See vignette("GPA1guide", package = "GPArotation") for further analyses.
  trBritain <- matrix(c(.783, -.163, .811, .202, .724, .209, .850, .064,
                       -.031, .592, -.028, .723, .388, .434, .141, .808,
                       .215, .709), byrow = TRUE, ncol = 2)

  trGermany <- matrix(c(.778, -.066, .875, .081, .751, .079, .739, .092,
                       .195, .574, -.030, .807, -.135, .717, .125, .738,
                       .060, .691), byrow = TRUE, ncol = 2)
  trx <- targetT(trGermany, Target = trBritain)
  print(trx$loadings - trBritain, cutoff = 0, digits = 3)  # difference from target

  # Plot discrepancy from target --- steelblue = below target, firebrick = above
  plot(trx, type = "target", Target = trBritain)

  # --- Partially specified target rotation ---
  # See vignette("GPA1guide", package = "GPArotation") for full context.
  # Unrotated loadings matrix A and partially specified target SPA.
  # NA entries in SPA are unspecified --- rotation is free there.
  # Numeric entries are the target values the rotation aims towards.
  A <- matrix(c(.664, .688, .492, .837, .705, .82, .661, .457, .765, .322,
                .248, .304, -0.291, -0.314, -0.377, .397, .294, .428,
                -0.075, .192, .224, .037, .155, -.104, .077, -.488, .009), ncol = 3)
  SPA <- matrix(c(rep(NA, 6), .7, .0, .7, rep(0, 3), rep(NA, 7),
                  0, 0, NA, 0, rep(NA, 4)), ncol = 3)
  comparison <- cbind(round(A, 3), rep(NA, nrow(A)), SPA)
  colnames(comparison) <- c("A.F1", "A.F2", "A.F3", "|", "T.F1", "T.F2", "T.F3")
  cat("Unrotated loadings (A) and partially specified target (SPA):\n")
  print(comparison, na.print = "NA")

  res.t <- targetT(A, Target = SPA)
  print(res.t)

  # Discrepancy from target --- specified elements only
  plot(res.t, type = "target", Target = SPA)
  
  # --- Random starts and Update calls ---
  # CCAI Climate-Friendly Purchasing Choices domain, 14 items, 3 oblique factors.
  # High factor intercorrelations make oblimin the natural choice.
  # Note: factanal uses MLE extraction; results differ somewhat from
  # PCA-based extraction used in Bi and Barchard (2024).
  data("CCAI", package = "GPArotation")
  fa.unrotated <- factanal(factors = 3, covmat = CCAI_R, n.obs = 461, rotation = "none")
  # Single random start via Tmat
  res.o <- oblimin(fa.unrotated, Tmat = Random.Start(3))
  res.o <- oblimin(fa.unrotated, randomStarts = 1)  # equivalent using randomStarts = 1
  res.o <- oblimin(fa.unrotated, randomStarts = 100) # multiple random starts 
  # Update arguments without re-specifying the full call
  res.o2 <- update(res.o, randomStarts = 250)  # randomStarts = 250
  res.o3 <- update(res.o, gam = -0.5) # gam = -0.5, randomStarts = 100
  res.o4 <- update(res.o, normalize = TRUE) # normalize = TRUE, randomStarts = 100
   
  # To change rotation criterion, use A_unrotated directly
  A <- attr(res.o, "A_unrotated")
  res.qt <- GPFRSoblq(A, method = "quartimin")

  # Directly via factanal call 
  factanal(factors = 3, covmat = CCAI_R, n.obs = 461, rotation = "oblimin",
      control = list(rotate = list(normalize = TRUE, gam = -0.1, randomStarts = 100)))

  # --- Assessing local minima ---
  # For detailed investigation of local minima across all random starts
  # see vignette("GPA2local", package = "GPArotation").
  data(Thurstone, package = "GPArotation")
  infomaxQ(box26, normalize = TRUE, randomStarts = 150)
  geominQ(box26,  normalize = TRUE, randomStarts = 150)

Simple Structure Quality Measures

Description

Functions for assessing the quality of simple structure in a rotated factor solution. Four criterion-free measures are available:

calc_AUC

Area Under the Curve measure of factor loading simplicity, computed per factor and overall (Liu et al., 2023).

calc_FSI

Factor Simplicity Index, computed per factor and overall (Lorenzo-Seva, 2003).

calc_simplicity

Three overall solution measures: Hoffman (1978), Gini (Kaiser, 1974), and Bentler (1977).

calc_hyperplane

Hyperplane count per factor and overall, measuring coordinate sparsity in the rotated solution.

calc_fitstats

Chi-square, SMRS, RMSEA with confidence interval, valid for MLE extracted loadings only.

All measures are criterion-free — they do not depend on how the rotation was obtained — making them useful for comparing simple structure quality across different rotation methods applied to the same unrotated solution. The table below summarizes the measures, their scope, range, and interpretation:

Note: No universally accepted cutoffs exist for these measures. The ranges above are approximate guidelines based on empirical experience. Always interpret in conjunction with substantive theory and model fit indices such as RMSEA and SRMR.

Note on Bentler: the Bentler index tends to be low when factor intercorrelations are high, even when simple structure is otherwise clean. A low Bentler value should therefore be interpreted in conjunction with AUC, FSI, and hyperplane counts rather than in isolation. High AUC and hyperplane counts alongside a low Bentler value typically indicate correlated factors rather than poor simple structure.

These functions are used internally by print.GPArotation and summary.GPArotation and are displayed automatically when printing or summarising a GPArotation object. They can also be called directly using the triple-colon operator, for example GPArotation:::calc_AUC(x).

Usage

## Not exported -- access via GPArotation:::calc_AUC(x) etc.
## calc_AUC(x, digits = 3L)
## calc_FSI(x, digits = 3L)
## calc_simplicity(x, digits = 3L)
## calc_hyperplane(x, cutoff = 0.1, digits = 3L)
## calc_fitstats(x)

Arguments

Details

All four measures are criterion-free — they do not depend on how the rotation was obtained — making them useful for comparing simple structure quality across different rotation methods applied to the same unrotated solution.

AUC (Liu et al., 2023)

For each factor j the squared loadings are sorted in descending order and their cumulative proportions of total factor variance are plotted against their rank. AUC is the area under this curve:

AUC_j = \frac{1}{p} \sum_{i=1}^{p} \frac{\sum_{r=1}^{i} l_{(r)j}^2}{\sum_{r=1}^{p} l_{(r)j}^2}

AUC = 1 indicates perfect simple structure (one non-zero loading); AUC = 0.5 indicates no simple structure (equal loadings). Adjusted AUC subtracts 0.5 so that 0 corresponds to no simple structure.

FSI (Lorenzo-Seva, 2003)

Based on the kurtosis of the squared loadings within each factor:

FSI_j = \frac{p \sum l_{ij}^4 - \left(\sum l_{ij}^2\right)^2} {(p-1)\left(\sum l_{ij}^2\right)^2}

FSI = 0 indicates no simple structure; FSI = 1 indicates perfect simple structure. More sensitive to dominant loadings than AUC.

Overall solution measures

Hoffman (1978):

H = 1 - \frac{\sum_j \sum_i l_{ij}^2 (1 - l_{ij}^2)}{p \cdot k \cdot 0.25}

Measures how far loadings are from 0.5. H = 1 when all loadings are 0 or 1; H = 0.25 when all loadings equal 0.5.

Gini (Kaiser, 1974): Based on the Gini coefficient of the absolute loadings within each factor, averaged across factors. Higher values indicate greater concentration of variance on fewer indicators per factor.

Bentler (1977): Ratio of between-factor to within-factor variance of squared loadings. Values near 1 indicate clean simple structure.

Hyperplane count (Mair, 2018)

For each factor j:

HP_j = \sum_{i=1}^{p} \mathbf{1}(|l_{ij}| < c)

where c is the hyperplane cutoff. The theoretical maximum is p(k-1), achieved when each indicator loads on exactly one factor. The percentage of the theoretical maximum provides a scale-free summary of hyperplane density.

AUC and FSI are per-factor measures shown in print.GPArotation. All four measures plus their summaries are shown in summary.GPArotation.

Author(s)

Coen A. Bernaards.

References

Bentler, P.M. (1977) Factor simplicity index and transformations. Psychometrika, 42(2), 277–295. doi:10.1007/BF02294054

Hoffman, P.J. (1978) The paramorphic representation of clinical judgment. Psychological Bulletin, 55(2), 116–131. doi:10.1037/h0047807

Kaiser, H.F. (1974) An index of factorial simplicity. Psychometrika, 39(1), 31–36. doi:10.1007/BF02291575

Liu, X., Wallin, G., Chen, Y., and Moustaki, I. (2023). Rotation to sparse loadings using L^p losses and related inference problems. Psychometrika, 88(2), 527–553. doi:10.1007/s11336-023-09911-y

Lorenzo-Seva, U. (2003) A factor simplicity index. Psychometrika, 68(1), 49–60. doi:10.1007/BF02296652

Mair, P. (2018). Modern Psychometrics with R. Springer. doi:10.1007/978-3-319-93177-7

See Also

print.GPArotation, summary.GPArotation

Examples

data("CCAI", package = "GPArotation")
fa.un <- factanal(factors = 3, covmat = CCAI_R, n.obs = 461,
                  rotation = "none")

res.tan <- tandemI(fa.un)
res.ob  <- oblimin(fa.un)

# Simple structure measures shown automatically in print
# Orthogonal: SS loadings, Proportion Var, Cumulative Var, AUC, FSI
print(res.tan)

# Oblique: SS loadings, AUC, FSI, Phi
print(res.ob)

# Summary adds overall measures: Hoffman, Gini, Bentler, hyperplane
summary(res.tan)
summary(res.ob)

# Direct access for custom analyses
GPArotation:::calc_AUC(res.ob)
GPArotation:::calc_FSI(res.ob)
GPArotation:::calc_simplicity(res.ob)
GPArotation:::calc_hyperplane(res.ob)

# Stricter hyperplane cutoff
GPArotation:::calc_hyperplane(res.ob, cutoff = 0.05)

# Compare simple structure across rotation methods
auc.tan <- GPArotation:::calc_AUC(res.tan)$AUC_mean
auc.ob  <- GPArotation:::calc_AUC(res.ob)$AUC_mean
fsi.tan <- GPArotation:::calc_FSI(res.tan)$FSI_mean
fsi.ob  <- GPArotation:::calc_FSI(res.ob)$FSI_mean

cat(paste0(formatC("Tandem I", width = 10, flag = "-"),
           "  AUC: ", round(auc.tan, 3),
           "  FSI: ", round(fsi.tan, 3), "\n"))
cat(paste0(formatC("Oblimin", width = 10, flag = "-"),
           "  AUC: ", round(auc.ob, 3),
           "  FSI: ", round(fsi.ob, 3), "\n"))           

Rotation Criteria: Value and Gradient Functions

Description

Rotation criterion functions that compute the value and gradient of the rotation objective Q. Not exported from the package NAMESPACE.

Usage

    vgQ.oblimin(L, gam=0)
    vgQ.quartimin(L)
    vgQ.target(L, Target=NULL)
    vgQ.pst(L, W=NULL, Target=NULL)
    vgQ.oblimax(L)
    vgQ.entropy(L)
    vgQ.quartimax(L)
    vgQ.varimax(L)
    vgQ.simplimax(L, k=nrow(L))
    vgQ.bentler(L)
    vgQ.tandemI(L)
    vgQ.tandemII(L)
    vgQ.geomin(L, delta=0.01)
    vgQ.bigeomin(L, delta=0.01)
    vgQ.cf(L, kappa=0)
    vgQ.infomax(L)
    vgQ.mccammon(L)
    vgQ.varimin(L)
    vgQ.bifactor(L)
    vgQ.lp.wls(L, W)

Arguments

Details

The vgQ.* functions are called internally by the gradient projection optimization routines and would typically not be used directly. They can be inspected using :::, for example: GPArotation:::vgQ.oblimin. The gradient projection algorithms are described in Bernaards and Jennrich (2005).

The T or Q suffix on rotation function names should be omitted for the corresponding vgQ.* function. For example, GPArotation:::vgQ.target is the criterion used by both targetT and targetQ.

See rotations for details on rotation arguments.

New rotation criteria can be added by defining a function named vgQ.newmethod. The function takes L as its first argument, plus any additional arguments, and must return a list with elements f, Gq, and Method.

Gradient projection without derivatives can be performed using the GPArotateDF package; type vignette("GPArotateDF", package = "GPArotation") at the command line.

Value

A list with:

Author(s)

Coen A. Bernaards and Robert I. Jennrich with some R modifications by Paul Gilbert.

References

Bernaards, C.A. and Jennrich, R.I. (2005) Gradient Projection Algorithms and Software for Arbitrary Rotation Criteria in Factor Analysis. Educational and Psychological Measurement, 65, 676–696. doi:10.1177/0013164404272507

See Also

rotations, GPFRSorth

Examples

  GPArotation:::vgQ.oblimin