Core Algorithms and Random-Start Wrappers
GPA.RdGradient 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=1000,
method="varimax", methodArgs=NULL)
GPFoblq(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000,
method="quartimin", methodArgs=NULL)
GPFRSorth(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000,
method="varimax", methodArgs=NULL, randomStarts=0, ...)
GPFRSoblq(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000,
method="quartimin", methodArgs=NULL, randomStarts=0, ...)Arguments
- A
initial factor loadings matrix for which the rotation criterion is to be optimized.
- Tmat
initial rotation matrix.
- normalize
see details.
- eps
convergence is assumed when the norm of the gradient is smaller than
eps.- maxit
maximum number of iterations allowed in the main loop.
- method
rotation objective criterion.
- methodArgs
a list of additional arguments passed to the rotation objective.
- randomStarts
number of random starts (
GPFRSorthandGPFRSoblqonly).- ...
additional arguments passed to
GPForthorGPFoblq, such asepsandmaxitwhen calling viaGPFRSorthorGPFRSoblq.
Value
A GPArotation object which is a list with elements:
- loadings
The rotated loadings matrix, one column per factor. If random starts were requested, this is the solution with the lowest criterion value.
- Th
The rotation matrix, satisfying
loadings %*% t(Th) = Afor orthogonal rotation andloadings = A %*% solve(t(Th))for oblique rotation.- Table
A matrix recording the iteration history: iteration number, criterion value, log10 of the gradient norm, and step size (alpha).
- method
A string indicating the rotation criterion.
- orthogonal
A logical indicating if the rotation is orthogonal.
- convergence
A logical indicating if convergence was obtained.
- Phi
t(Th) %*% Th, the covariance matrix of the rotated factors. Omitted (NULL) for orthogonal rotations.- Gq
The gradient of the criterion at the rotated loadings.
- randStartChar
A named vector summarising random start results:
randomStarts,Converged,atMinimum,localMins. Only present whenrandomStarts > 1.
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 matrixAsum 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 ofAare 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
Aas an argument and return a vector, which is then applied in the same manner as the vector input described above. SeeNormalizingWeightfor 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 matrixTmat. The initial rotation matrixTmatcan also be set by the user.Setting
randomStarts = 1initializesTmatwith a single random matrix.Setting
randomStarts > 1attempts multiple random starts and returns the rotated loadings matrix that achieved the lowest criterion valuefacross 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)).
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. You can see a direct comparison of this in the examples section.
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
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., & 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
Examples
# --- Basic rotation calls ---
data(Harman, package = "GPArotation") # 8 physical variables
quartimax(Harman8) # direct rotation call
#> Orthogonal rotation method Quartimax converged.
#> Loadings:
#> CF1 CF2
#> height 0.899 0.195
#> arm.span 0.934 0.130
#> forearm 0.902 0.104
#> lower.leg 0.877 0.171
#> weight 0.316 0.876
#> bitro.diameter 0.251 0.773
#> chest.girth 0.198 0.715
#> chest.width 0.308 0.659
#>
#> CF1 CF2
#> SS loadings 3.559 2.407
#> Proportion Var 0.445 0.301
#> Cumulative Var 0.445 0.746
GPFRSorth(Harman8, method = "quartimax") # equivalent via wrapper
#> Orthogonal rotation method Quartimax converged.
#> Loadings:
#> CF1 CF2
#> height 0.899 0.195
#> arm.span 0.934 0.130
#> forearm 0.902 0.104
#> lower.leg 0.877 0.171
#> weight 0.316 0.876
#> bitro.diameter 0.251 0.773
#> chest.girth 0.198 0.715
#> chest.width 0.308 0.659
#>
#> CF1 CF2
#> SS loadings 3.559 2.407
#> Proportion Var 0.445 0.301
#> Cumulative Var 0.445 0.746
GPFRSoblq(Harman8, method = "quartimin", normalize = TRUE)
#> Oblique rotation method Quartimin converged.
#> Loadings:
#> CF1 CF2
#> height 0.892 0.056
#> arm.span 0.954 -0.023
#> forearm 0.929 -0.046
#> lower.leg 0.877 0.034
#> weight 0.014 0.925
#> bitro.diameter -0.017 0.821
#> chest.girth -0.052 0.765
#> chest.width 0.086 0.683
#>
#> CF1 CF2
#> SS loadings 3.362 2.604
#> Proportion Var 0.420 0.325
#> Cumulative Var 0.420 0.746
#>
#> Phi:
#> CF1 CF2
#> CF1 1.000 0.473
#> CF2 0.473 1.000
loadings(quartimin(Harman8, normalize = TRUE)) # extract loadings directly
#> CF1 CF2
#> height 0.89176370 0.05613350
#> arm.span 0.95362150 -0.02312061
#> forearm 0.92909390 -0.04638106
#> lower.leg 0.87662670 0.03377464
#> weight 0.01365796 0.92502211
#> bitro.diameter -0.01730126 0.82126806
#> chest.girth -0.05247672 0.76496236
#> chest.width 0.08585315 0.68314140
# --- 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
#> Oblique rotation method Equamax converged.
#> Loadings:
#> [,1] [,2] [,3]
#> [1,] 0.627 0.681 -0.380
#> [2,] 0.711 -0.203 0.568
#> [3,] -0.145 0.630 0.692
#> [4,] 0.876 0.241 0.137
#> [5,] 0.257 0.844 0.222
#> [6,] 0.302 0.236 0.810
#> [7,] 0.810 0.467 -0.075
#> [8,] 0.828 0.049 0.362
#> [9,] 0.408 0.829 -0.017
#> [10,] 0.137 0.875 0.436
#> [11,] 0.454 0.066 0.767
#> [12,] 0.148 0.378 0.805
#> [13,] -0.124 0.763 -0.758
#> [14,] 0.124 -0.763 0.758
#> [15,] 0.689 0.005 -0.857
#> [16,] -0.689 -0.005 0.857
#> [17,] 0.765 -0.792 -0.038
#> [18,] -0.765 0.792 0.038
#> [19,] 0.869 0.155 0.205
#> [20,] 0.205 0.853 0.250
#> [21,] 0.320 0.239 0.797
#> [22,] 0.853 0.159 0.211
#> [23,] 0.196 0.822 0.265
#> [24,] 0.340 0.233 0.767
#> [25,] 0.545 0.485 0.434
#> [26,] 0.505 0.407 0.518
#>
#> [,1] [,2] [,3]
#> SS loadings 8.685 8.417 8.307
#> Proportion Var 0.334 0.324 0.320
#> Cumulative Var 0.334 0.658 0.977
#>
#> Phi:
#> [,1] [,2] [,3]
#> [1,] 1.000 0.179 0.200
#> [2,] 0.179 1.000 0.161
#> [3,] 0.200 0.161 1.000
GPFRSoblq(box26, method = "cf", methodArgs = list(kappa = 2/27)) # Parsimax
#> Oblique rotation method Parsimax converged.
#> Loadings:
#> [,1] [,2] [,3]
#> [1,] 0.629 0.680 -0.376
#> [2,] 0.705 -0.204 0.574
#> [3,] -0.150 0.630 0.692
#> [4,] 0.874 0.240 0.142
#> [5,] 0.255 0.843 0.223
#> [6,] 0.296 0.236 0.813
#> [7,] 0.810 0.466 -0.070
#> [8,] 0.825 0.048 0.368
#> [9,] 0.408 0.829 -0.015
#> [10,] 0.133 0.874 0.437
#> [11,] 0.448 0.065 0.772
#> [12,] 0.142 0.378 0.807
#> [13,] -0.118 0.764 -0.761
#> [14,] 0.118 -0.764 0.761
#> [15,] 0.695 0.004 -0.854
#> [16,] -0.695 -0.004 0.854
#> [17,] 0.764 -0.793 -0.033
#> [18,] -0.764 0.793 0.033
#> [19,] 0.867 0.154 0.211
#> [20,] 0.203 0.852 0.251
#> [21,] 0.314 0.238 0.800
#> [22,] 0.850 0.158 0.217
#> [23,] 0.194 0.822 0.266
#> [24,] 0.334 0.232 0.770
#> [25,] 0.541 0.484 0.438
#> [26,] 0.501 0.406 0.522
#>
#> [,1] [,2] [,3]
#> SS loadings 8.641 8.411 8.357
#> Proportion Var 0.332 0.324 0.321
#> Cumulative Var 0.332 0.656 0.977
#>
#> Phi:
#> [,1] [,2] [,3]
#> [1,] 1.000 0.179 0.202
#> [2,] 0.179 1.000 0.163
#> [3,] 0.202 0.163 1.000
# --- 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 reported
# by Bernaards and others and 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(loadings(fa.unrotated), kappa = 0.3, normalize = TRUE,
randomStarts = 100)
fa.cf
#> Oblique rotation method Crawford-Ferguson: kappa=0.3 converged at lowest minimum.
#> Of 100 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> Factor1 Factor2 Factor3
#> NL1 0.713 0.111 0.073
#> TV2 0.700 0.131 0.162
#> NL3 0.692 0.085 0.099
#> RTL4 0.089 0.107 0.618
#> RTL5 0.066 0.265 0.525
#> Veronica -0.030 0.251 0.642
#> SBS6 -0.066 1.066 -0.085
#>
#> Factor1 Factor2 Factor3
#> SS loadings 1.672 1.569 1.435
#> Proportion Var 0.239 0.224 0.205
#> Cumulative Var 0.239 0.463 0.668
#>
#> Phi:
#> Factor1 Factor2 Factor3
#> Factor1 1.000 0.313 0.392
#> Factor2 0.313 1.000 0.623
#> Factor3 0.392 0.623 1.000
# 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")
}
#> Oblique rotation method Crawford-Ferguson: kappa=0.3 converged at lowest minimum.
#> Of 100 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> Factor1 Factor2 Factor3
#> NL1 0.713 0.111 0.073
#> TV2 0.700 0.131 0.162
#> NL3 0.692 0.085 0.099
#> RTL4 0.089 0.107 0.618
#> RTL5 0.066 0.265 0.525
#> Veronica -0.030 0.251 0.642
#> SBS6 -0.066 1.066 -0.085
#>
#> Factor1 Factor2 Factor3
#> SS loadings 1.672 1.569 1.435
#> Proportion Var 0.239 0.224 0.205
#> Cumulative Var 0.239 0.463 0.668
#>
#> Phi:
#> Factor1 Factor2 Factor3
#> Factor1 1.000 0.313 0.392
#> Factor2 0.313 1.000 0.623
#> Factor3 0.392 0.623 1.000
#> Maximum difference in loadings between two-step and single-step:
#> [1] 0
# --- 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(loadings(fa.unrotated), gam = -0.5, randomStarts = 20)
# gam = -0.5: more orthogonal than quartimin
res # default print
#> Oblique rotation method Oblimin g=-0.5 converged at lowest minimum.
#> Of 20 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> Factor1 Factor2 Factor3
#> CCAI8 0.912 0.039 0.001
#> CCAI6 0.801 0.045 0.046
#> CCAI7 0.759 -0.011 0.133
#> CCAI11 0.533 0.272 0.218
#> CCAI12 0.487 0.360 0.166
#> CCAI10 0.383 0.295 0.288
#> CCAI14 0.046 0.946 0.002
#> CCAI13 0.062 0.899 0.049
#> CCAI5 0.063 0.613 0.291
#> CCAI2 0.088 -0.042 0.591
#> CCAI4 0.042 0.083 0.751
#> CCAI1 0.086 0.041 0.661
#> CCAI3 0.155 0.143 0.606
#> CCAI9 0.229 0.203 0.498
#>
#> Factor1 Factor2 Factor3
#> SS loadings 3.662 3.176 3.061
#> Proportion Var 0.262 0.227 0.219
#> Cumulative Var 0.262 0.488 0.707
#>
#> Phi:
#> Factor1 Factor2 Factor3
#> Factor1 1.000 0.523 0.578
#> Factor2 0.523 1.000 0.528
#> Factor3 0.578 0.528 1.000
print(res) # equivalent to above
#> Oblique rotation method Oblimin g=-0.5 converged at lowest minimum.
#> Of 20 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> Factor1 Factor2 Factor3
#> CCAI8 0.912 0.039 0.001
#> CCAI6 0.801 0.045 0.046
#> CCAI7 0.759 -0.011 0.133
#> CCAI11 0.533 0.272 0.218
#> CCAI12 0.487 0.360 0.166
#> CCAI10 0.383 0.295 0.288
#> CCAI14 0.046 0.946 0.002
#> CCAI13 0.062 0.899 0.049
#> CCAI5 0.063 0.613 0.291
#> CCAI2 0.088 -0.042 0.591
#> CCAI4 0.042 0.083 0.751
#> CCAI1 0.086 0.041 0.661
#> CCAI3 0.155 0.143 0.606
#> CCAI9 0.229 0.203 0.498
#>
#> Factor1 Factor2 Factor3
#> SS loadings 3.662 3.176 3.061
#> Proportion Var 0.262 0.227 0.219
#> Cumulative Var 0.262 0.488 0.707
#>
#> Phi:
#> Factor1 Factor2 Factor3
#> Factor1 1.000 0.523 0.578
#> Factor2 0.523 1.000 0.528
#> Factor3 0.578 0.528 1.000
print(res, Table = TRUE) # include iteration table
#> Oblique rotation method Oblimin g=-0.5 converged at lowest minimum.
#> Of 20 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> Factor1 Factor2 Factor3
#> CCAI8 0.912 0.039 0.001
#> CCAI6 0.801 0.045 0.046
#> CCAI7 0.759 -0.011 0.133
#> CCAI11 0.533 0.272 0.218
#> CCAI12 0.487 0.360 0.166
#> CCAI10 0.383 0.295 0.288
#> CCAI14 0.046 0.946 0.002
#> CCAI13 0.062 0.899 0.049
#> CCAI5 0.063 0.613 0.291
#> CCAI2 0.088 -0.042 0.591
#> CCAI4 0.042 0.083 0.751
#> CCAI1 0.086 0.041 0.661
#> CCAI3 0.155 0.143 0.606
#> CCAI9 0.229 0.203 0.498
#>
#> Factor1 Factor2 Factor3
#> SS loadings 3.662 3.176 3.061
#> Proportion Var 0.262 0.227 0.219
#> Cumulative Var 0.262 0.488 0.707
#>
#> Phi:
#> Factor1 Factor2 Factor3
#> Factor1 1.000 0.523 0.578
#> Factor2 0.523 1.000 0.528
#> Factor3 0.578 0.528 1.000
#>
#> Iteration table:
#> iter f log10(s) alpha
#> [1,] 0 1.401 0.601 1.0000
#> [2,] 1 0.803 0.130 0.0625
#> [3,] 2 0.644 -0.202 0.1250
#> [4,] 3 0.578 -0.392 0.2500
#> [5,] 4 0.511 -0.445 0.5000
#> [6,] 5 0.478 -0.366 0.5000
#> [7,] 6 0.464 -0.726 0.1250
#> [8,] 7 0.456 -0.819 0.2500
#> [9,] 8 0.447 -1.017 0.5000
#> [10,] 9 0.444 -1.118 0.5000
#> [11,] 10 0.444 -1.384 0.1250
#> [12,] 11 0.444 -1.557 0.2500
#> [13,] 12 0.443 -1.725 0.5000
#> [14,] 13 0.443 -1.989 0.1250
#> [15,] 14 0.443 -2.152 0.2500
#> [16,] 15 0.443 -2.314 0.5000
#> [17,] 16 0.443 -2.559 0.1250
#> [18,] 17 0.443 -2.703 0.2500
#> [19,] 18 0.443 -2.872 0.5000
#> [20,] 19 0.443 -3.055 0.2500
#> [21,] 20 0.443 -3.219 0.2500
#> [22,] 21 0.443 -3.369 0.2500
#> [23,] 22 0.443 -3.381 0.5000
#> [24,] 23 0.443 -3.720 0.1250
#> [25,] 24 0.443 -3.859 0.2500
#> [26,] 25 0.443 -3.933 0.5000
#> [27,] 26 0.443 -4.207 0.1250
#> [28,] 27 0.443 -4.340 0.2500
#> [29,] 28 0.443 -4.472 0.5000
#> [30,] 29 0.443 -4.654 0.2500
#> [31,] 30 0.443 -4.819 0.2500
#> [32,] 31 0.443 -4.969 0.2500
#> [33,] 32 0.443 -4.965 0.5000
#> [34,] 33 0.443 -5.315 0.1250
print(res, rotateMat = TRUE) # include rotating matrix
#> Oblique rotation method Oblimin g=-0.5 converged at lowest minimum.
#> Of 20 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> Factor1 Factor2 Factor3
#> CCAI8 0.912 0.039 0.001
#> CCAI6 0.801 0.045 0.046
#> CCAI7 0.759 -0.011 0.133
#> CCAI11 0.533 0.272 0.218
#> CCAI12 0.487 0.360 0.166
#> CCAI10 0.383 0.295 0.288
#> CCAI14 0.046 0.946 0.002
#> CCAI13 0.062 0.899 0.049
#> CCAI5 0.063 0.613 0.291
#> CCAI2 0.088 -0.042 0.591
#> CCAI4 0.042 0.083 0.751
#> CCAI1 0.086 0.041 0.661
#> CCAI3 0.155 0.143 0.606
#> CCAI9 0.229 0.203 0.498
#>
#> Factor1 Factor2 Factor3
#> SS loadings 3.662 3.176 3.061
#> Proportion Var 0.262 0.227 0.219
#> Cumulative Var 0.262 0.488 0.707
#>
#> Phi:
#> Factor1 Factor2 Factor3
#> Factor1 1.000 0.523 0.578
#> Factor2 0.523 1.000 0.528
#> Factor3 0.578 0.528 1.000
#>
#> Rotating matrix:
#> [,1] [,2] [,3]
#> [1,] 0.329 0.662 0.179
#> [2,] 0.936 -1.013 0.294
#> [3,] -0.826 -0.273 1.248
print(res, digits = 2) # rounded to 2 decimal places
#> Oblique rotation method Oblimin g=-0.5 converged at lowest minimum.
#> Of 20 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> Factor1 Factor2 Factor3
#> CCAI8 0.91 0.04 0.00
#> CCAI6 0.80 0.05 0.05
#> CCAI7 0.76 -0.01 0.13
#> CCAI11 0.53 0.27 0.22
#> CCAI12 0.49 0.36 0.17
#> CCAI10 0.38 0.30 0.29
#> CCAI14 0.05 0.95 0.00
#> CCAI13 0.06 0.90 0.05
#> CCAI5 0.06 0.61 0.29
#> CCAI2 0.09 -0.04 0.59
#> CCAI4 0.04 0.08 0.75
#> CCAI1 0.09 0.04 0.66
#> CCAI3 0.16 0.14 0.61
#> CCAI9 0.23 0.20 0.50
#>
#> Factor1 Factor2 Factor3
#> SS loadings 3.66 3.18 3.06
#> Proportion Var 0.26 0.23 0.22
#> Cumulative Var 0.26 0.49 0.71
#>
#> Phi:
#> Factor1 Factor2 Factor3
#> Factor1 1.00 0.52 0.58
#> Factor2 0.52 1.00 0.53
#> Factor3 0.58 0.53 1.00
summary(res) # pattern and structure matrices for oblique rotation
#> Oblique rotation method Oblimin g=-0.5 converged in 33 iterations.
#> Pattern (loadings):
#> Factor1 Factor2 Factor3
#> CCAI8 0.912 0.039 0.001
#> CCAI6 0.801 0.045 0.046
#> CCAI7 0.759 -0.011 0.133
#> CCAI11 0.533 0.272 0.218
#> CCAI12 0.487 0.360 0.166
#> CCAI10 0.383 0.295 0.288
#> CCAI14 0.046 0.946 0.002
#> CCAI13 0.062 0.899 0.049
#> CCAI5 0.063 0.613 0.291
#> CCAI2 0.088 -0.042 0.591
#> CCAI4 0.042 0.083 0.751
#> CCAI1 0.086 0.041 0.661
#> CCAI3 0.155 0.143 0.606
#> CCAI9 0.229 0.203 0.498
#>
#> Structure:
#> Factor1 Factor2 Factor3
#> CCAI8 0.933 0.517 0.549
#> CCAI6 0.851 0.489 0.533
#> CCAI7 0.830 0.457 0.566
#> CCAI11 0.801 0.666 0.669
#> CCAI12 0.771 0.702 0.637
#> CCAI10 0.704 0.648 0.665
#> CCAI14 0.542 0.971 0.527
#> CCAI13 0.561 0.957 0.559
#> CCAI5 0.552 0.799 0.650
#> CCAI2 0.407 0.315 0.619
#> CCAI4 0.520 0.501 0.819
#> CCAI1 0.490 0.435 0.733
#> CCAI3 0.581 0.545 0.772
#> CCAI9 0.623 0.585 0.737
summary(res, Structure = FALSE) # pattern matrix only
#> Oblique rotation method Oblimin g=-0.5 converged in 33 iterations.
#> Loadings:
#> Factor1 Factor2 Factor3
#> CCAI8 0.912 0.039 0.001
#> CCAI6 0.801 0.045 0.046
#> CCAI7 0.759 -0.011 0.133
#> CCAI11 0.533 0.272 0.218
#> CCAI12 0.487 0.360 0.166
#> CCAI10 0.383 0.295 0.288
#> CCAI14 0.046 0.946 0.002
#> CCAI13 0.062 0.899 0.049
#> CCAI5 0.063 0.613 0.291
#> CCAI2 0.088 -0.042 0.591
#> CCAI4 0.042 0.083 0.751
#> CCAI1 0.086 0.041 0.661
#> CCAI3 0.155 0.143 0.606
#> CCAI9 0.229 0.203 0.498
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
#> randomStarts Converged atMinimum localMins
#> 50 50 16 4
# --- 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
#> CF1 CF2
#> height -0.05613414 -0.89176336
#> arm.span 0.02311995 -0.95362121
#> forearm 0.04638043 -0.92909364
#> lower.leg -0.03377527 -0.87662638
#> weight -0.92502250 -0.01365713
#> bitro.diameter -0.82126839 0.01730199
#> chest.girth -0.76496264 0.05247738
#> chest.width -0.68314175 -0.08585252
max(abs(print(xusl)$loadings - xusl$loadings)) == 0 # FALSE: print() reorders
#> Oblique rotation method Quartimin converged at lowest minimum.
#> Of 100 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> CF1 CF2
#> height 0.892 0.056
#> arm.span 0.954 -0.023
#> forearm 0.929 -0.046
#> lower.leg 0.877 0.034
#> weight 0.014 0.925
#> bitro.diameter -0.017 0.821
#> chest.girth -0.052 0.765
#> chest.width 0.086 0.683
#>
#> CF1 CF2
#> SS loadings 3.362 2.604
#> Proportion Var 0.420 0.325
#> Cumulative Var 0.420 0.746
#>
#> Phi:
#> CF1 CF2
#> CF1 1.000 0.473
#> CF2 0.473 1.000
#> [1] FALSE
xsl <- print(xusl) # capture sorted result
#> Oblique rotation method Quartimin converged at lowest minimum.
#> Of 100 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> CF1 CF2
#> height 0.892 0.056
#> arm.span 0.954 -0.023
#> forearm 0.929 -0.046
#> lower.leg 0.877 0.034
#> weight 0.014 0.925
#> bitro.diameter -0.017 0.821
#> chest.girth -0.052 0.765
#> chest.width 0.086 0.683
#>
#> CF1 CF2
#> SS loadings 3.362 2.604
#> Proportion Var 0.420 0.325
#> Cumulative Var 0.420 0.746
#>
#> Phi:
#> CF1 CF2
#> CF1 1.000 0.473
#> CF2 0.473 1.000
max(abs(print(xsl)$loadings - xsl$loadings)) == 0 # TRUE: already sorted
#> Oblique rotation method Quartimin converged at lowest minimum.
#> Of 100 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> CF1 CF2
#> height 0.892 0.056
#> arm.span 0.954 -0.023
#> forearm 0.929 -0.046
#> lower.leg 0.877 0.034
#> weight 0.014 0.925
#> bitro.diameter -0.017 0.821
#> chest.girth -0.052 0.765
#> chest.width 0.086 0.683
#>
#> CF1 CF2
#> SS loadings 3.362 2.604
#> Proportion Var 0.420 0.325
#> Cumulative Var 0.420 0.746
#>
#> Phi:
#> CF1 CF2
#> CF1 1.000 0.473
#> CF2 0.473 1.000
#> [1] TRUE
# --- Normalization ---
# Kaiser normalization
data("CCAI", package = "GPArotation")
fa.unrotated <- factanal(factors = 3, covmat = CCAI_R, n.obs = 461, rotation = "none")
oblimin(loadings(fa.unrotated), normalize = TRUE, randomStarts = 100)
#> Oblique rotation method Oblimin Quartimin converged at lowest minimum.
#> Of 100 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#> Factor1 Factor2 Factor3
#> CCAI8 0.994 -0.074 -0.015
#> CCAI6 0.866 -0.015 -0.004
#> CCAI7 0.811 0.091 -0.066
#> CCAI11 0.541 0.186 0.244
#> CCAI12 0.494 0.126 0.343
#> CCAI10 0.368 0.275 0.274
#> CCAI14 0.011 -0.061 1.003
#> CCAI13 0.025 -0.007 0.949
#> CCAI5 0.005 0.283 0.632
#> CCAI2 0.020 0.659 -0.083
#> CCAI4 -0.057 0.835 0.044
#> CCAI1 0.005 0.733 0.002
#> CCAI3 0.084 0.660 0.109
#> CCAI9 0.176 0.528 0.174
#>
#> Factor1 Factor2 Factor3
#> SS loadings 3.628 3.136 3.136
#> Proportion Var 0.259 0.224 0.224
#> Cumulative Var 0.259 0.483 0.707
#>
#> Phi:
#> Factor1 Factor2 Factor3
#> Factor1 1.000 0.713 0.625
#> Factor2 0.713 1.000 0.638
#> Factor3 0.625 0.638 1.000
# --- Legacy engine ---
# Demonstrates numerical equivalence of current and legacy implementations.
# The maximum absolute difference should be zero or at machine epsilon.
data("Harman", package = "GPArotation")
# Rotate using current engine
res_new <- oblimin(Harman8)
# Rotate using legacy engine
res_legacy <- GPArotation:::GPFoblq.legacy(Harman8, method = "oblimin")
# Numerically identical
max(abs(loadings(res_new) - loadings(res_legacy)))
#> [1] 0