Optimize factor loading rotation objective.

oblimin(A, Tmat=diag(ncol(A)), gam=0, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0)
    quartimin(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0)
    targetT(A, Tmat=diag(ncol(A)), Target=NULL, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0, L=NULL)
    targetQ(A, Tmat=diag(ncol(A)), Target=NULL, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0, L=NULL)
    pstT(A, Tmat=diag(ncol(A)), W=NULL, Target=NULL, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0, L=NULL)
    pstQ(A, Tmat=diag(ncol(A)), W=NULL, Target=NULL, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0, L=NULL)
    oblimax(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    entropy(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    quartimax(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5,maxit=1000,randomStarts=0)
    Varimax(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    simplimax(A, Tmat=diag(ncol(A)), k=nrow(A), normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0)
    bentlerT(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    bentlerQ(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    tandemI(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    tandemII(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    geominT(A, Tmat=diag(ncol(A)), delta=.01, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0)
    geominQ(A, Tmat=diag(ncol(A)), delta=.01, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0)
    bigeominT(A, Tmat=diag(ncol(A)), delta=.01, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0)
    bigeominQ(A, Tmat=diag(ncol(A)), delta=.01, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0)
    cfT(A, Tmat=diag(ncol(A)), kappa=0, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0)
    cfQ(A, Tmat=diag(ncol(A)), kappa=0, normalize=FALSE, eps=1e-5, 
        maxit=1000, randomStarts=0)
    equamax(A, Tmat=diag(ncol(A)), kappa=ncol(A)/(2*nrow(A)), normalize=FALSE,
        eps=1e-5, maxit=1000, randomStarts = 0)
    parsimax(A, Tmat=diag(ncol(A)), kappa=(ncol(A)-1)/(ncol(A)+nrow(A)-2), 
        normalize=FALSE, eps=1e-5, maxit=1000, randomStarts = 0)
    infomaxT(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    infomaxQ(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    mccammon(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    varimin(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000, randomStarts=0)
    bifactorT(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000,randomStarts=0)
    bifactorQ(A, Tmat=diag(ncol(A)), normalize=FALSE, eps=1e-5, maxit=1000,randomStarts=0)

Arguments

A

an initial loadings matrix to be rotated.

Tmat

initial rotation matrix.

gam

0=Quartimin, .5=Biquartimin, 1=Covarimin.

Target

rotation target for objective calculation.

W

weighting of each element in target.

k

number of close to zero loadings.

delta

constant added to Lambda^2 in objective calculation.

kappa

see details.

normalize

parameter passed to optimization routine (GPForth or GPFoblq).

eps

parameter passed to optimization routine (GPForth or GPFoblq).

maxit

parameter passed to optimization routine (GPForth or GPFoblq).

randomStarts

parameter passed to optimization routine (GPFRSorth or GPFRSoblq).

L

provided for backward compatibility in target rotations only. Use A going forward.

Value

A list (which includes elements used by factanal) with:

loadings

Lh from GPFRSorth or GPFRSoblq.

Th

Th from GPFRSorth or GPFRSoblq.

Table

Table from GPForth or GPFoblq.

method

A string indicating the rotation objective function.

orthogonal

A logical indicating if the rotation is orthogonal.

convergence

Convergence indicator from GPFRSorth or GPFRSoblq.

Phi

t(Th) %*% Th. The covariance matrix of the rotated factors. This will be the identity matrix for orthogonal rotations so is omitted (NULL) for the result from GPFRSorth and GPForth.

randStartChar

Vector indicating results from random starts from GPFRSorth or GPFRSoblq

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).

Rotations which are available are

obliminobliqueoblimin family
quartiminoblique
targetTorthogonaltarget rotation
targetQobliquetarget rotation
pstTorthogonalpartially specified target rotation
pstQobliquepartially specified target rotation
oblimaxoblique
entropyorthogonalminimum entropy
quartimaxorthogonal
varimaxorthogonal
simplimaxoblique
bentlerTorthogonalBentler's invariant pattern simplicity criterion
bentlerQobliqueBentler's invariant pattern simplicity criterion
tandemIorthogonalTandem principle I criterion
tandemIIorthogonalTandem principle II criterion
geominTorthogonal
geominQoblique
bigeominTorthogonal
bigeominQoblique
cfTorthogonalCrawford-Ferguson family
cfQobliqueCrawford-Ferguson family
equamaxorthogonalCrawford-Ferguson family
parsimaxorthogonalCrawford-Ferguson family
infomaxTorthogonal
infomaxQoblique
mccammonorthogonalMcCammon minimum entropy ratio
variminorthogonal
bifactorTorthogonalJennrich and Bentler bifactor rotation
bifactorQobliqueJennrich and Bentler biquartimin rotation

Note that Varimax defined here uses vgQ.varimax and is not varimax defined in the stats package. stats:::varimax does Kaiser normalization by default whereas Varimax defined here 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.

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.

Bifactor rotation, bifactorT and bifactorQ are called bifactor and biquartimin in Jennrich, R.I. and Bentler, P.M. (2011) Exploratory bi-factor analysis. Psychometrika, 76.

Author

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

Examples

  # see GPFRSorth and GPFRSoblq for more examples
  
  # getting loadings matrices
  data("Harman", package="GPArotation")
  qHarman  <- GPFRSorth(Harman8, Tmat=diag(2), method="quartimax")
  qHarman <- quartimax(Harman8) 
  loadings(qHarman) - qHarman$loadings   #2 ways to get the loadings
#>      [,1] [,2]
#> [1,]    0    0
#> [2,]    0    0
#> [3,]    0    0
#> [4,]    0    0
#> [5,]    0    0
#> [6,]    0    0
#> [7,]    0    0
#> [8,]    0    0

  # factanal loadings used in GPArotation
  data("WansbeekMeijer", package="GPArotation")
  fa.unrotated  <- factanal(factors = 2, covmat=NetherlandsTV, normalize=TRUE, rotation="none")
  quartimax(loadings(fa.unrotated), normalize=TRUE)
#> Orthogonal rotation method Quartimax converged.
#> Loadings:
#>          Factor1 Factor2
#> NL1        0.327   0.720
#> TV2        0.424   0.725
#> NL3        0.323   0.703
#> RTL4       0.669   0.214
#> RTL5       0.738   0.183
#> Veronica   0.811   0.109
#> SBS6       0.719   0.106
#> 
#>                Factor1 Factor2
#> SS loadings      2.556   1.641
#> Proportion Var   0.365   0.234
#> Cumulative Var   0.365   0.600
  geominQ(loadings(fa.unrotated), normalize=TRUE, randomStarts=100)
#> Oblique rotation method Geomin converged at lowest minimum.
#> Of 100 random starts 100% converged, 100% at the same lowest minimum.
#> Loadings at lowest minimum:
#>          Factor1 Factor2
#> NL1      -0.0159  0.8002
#> TV2       0.0897  0.7848
#> NL3      -0.0111  0.7796
#> RTL4      0.6306  0.1126
#> RTL5      0.7235  0.0601
#> Veronica  0.8437 -0.0459
#> SBS6      0.7435 -0.0295
#> 
#>                Factor1 Factor2
#> SS loadings      2.254   1.943
#> Proportion Var   0.322   0.278
#> Cumulative Var   0.322   0.600
#> 
#> Phi:
#>         Factor1 Factor2
#> Factor1   1.000   0.582
#> Factor2   0.582   1.000

  # passing arguments to factanal (See vignette for a caution)
  # vignette("GPAguide", package = "GPArotation")
  data(ability.cov)
  factanal(factors = 2, covmat = ability.cov, rotation="infomaxT")
#> 
#> Call:
#> factanal(factors = 2, covmat = ability.cov, rotation = "infomaxT")
#> 
#> Uniquenesses:
#> general picture  blocks    maze reading   vocab 
#>   0.455   0.589   0.218   0.769   0.052   0.334 
#> 
#> Loadings:
#>         Factor1 Factor2
#> general 0.495   0.547  
#> picture 0.151   0.623  
#> blocks  0.199   0.861  
#> maze    0.105   0.469  
#> reading 0.955   0.189  
#> vocab   0.783   0.231  
#> 
#>                Factor1 Factor2
#> SS loadings      1.844   1.738
#> Proportion Var   0.307   0.290
#> Cumulative Var   0.307   0.597
#> 
#> Test of the hypothesis that 2 factors are sufficient.
#> The chi square statistic is 6.11 on 4 degrees of freedom.
#> The p-value is 0.191 
  factanal(factors = 2, covmat = ability.cov, rotation="infomaxT", 
    control=list(rotate=list(normalize = TRUE, eps = 1e-6)))
#> 
#> Call:
#> factanal(factors = 2, covmat = ability.cov, rotation = "infomaxT",     control = list(rotate = list(normalize = TRUE, eps = 1e-06)))
#> 
#> Uniquenesses:
#> general picture  blocks    maze reading   vocab 
#>   0.455   0.589   0.218   0.769   0.052   0.334 
#> 
#> Loadings:
#>         Factor1 Factor2
#> general 0.568   0.471  
#> picture 0.629   0.124  
#> blocks  0.869   0.161  
#> maze    0.473          
#> reading 0.231   0.946  
#> vocab   0.265   0.772  
#> 
#>                Factor1 Factor2
#> SS loadings      1.821   1.761
#> Proportion Var   0.304   0.293
#> Cumulative Var   0.304   0.597
#> 
#> Test of the hypothesis that 2 factors are sufficient.
#> The chi square statistic is 6.11 on 4 degrees of freedom.
#> The p-value is 0.191 
  # when using factanal for oblique rotation it is best to use the rotation command directly
  # instead of including it in the factanal command (see Vignette).  
  fa.unrotated  <- factanal(factors = 3, covmat=NetherlandsTV, normalize=TRUE, rotation="none")
  quartimin(loadings(fa.unrotated), normalize=TRUE)
#> Oblique rotation method Quartimin converged.
#> Loadings:
#>          Factor1  Factor2 Factor3
#> NL1       0.8085 -0.03829  0.0166
#> TV2       0.7866  0.07540  0.0116
#> NL3       0.7801 -0.00504 -0.0113
#> RTL4      0.0533  0.74602 -0.0768
#> RTL5      0.0483  0.65405  0.0903
#> Veronica -0.0721  0.81431  0.0553
#> SBS6      0.0205  0.02007  0.9755
#> 
#>                Factor1 Factor2 Factor3
#> SS loadings      1.926   1.724   1.026
#> Proportion Var   0.275   0.246   0.147
#> Cumulative Var   0.275   0.521   0.668
#> 
#> Phi:
#>         Factor1 Factor2 Factor3
#> Factor1   1.000   0.615   0.381
#> Factor2   0.615   1.000   0.687
#> Factor3   0.381   0.687   1.000

  # oblique target rotation of 2 varimax rotated matrices towards each other
  # See vignette for additional context and computation,
  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 <- targetQ(trGermany, Target = trBritain)
  # Difference between rotated loadings matrix and target matrix 
  y <- trx$loadings - trBritain
  
  # partially specified target; See vignette for additional method
  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)
  targetT(A, Target=SPA)
#> Orthogonal rotation method Target rotation converged.
#> Loadings:
#>        [,1]    [,2]   [,3]
#>  [1,] 0.611  0.0231 0.4203
#>  [2,] 0.734  0.0554 0.1731
#>  [3,] 0.607 -0.0861 0.0926
#>  [4,] 0.608  0.6218 0.1756
#>  [5,] 0.549  0.5642 0.0172
#>  [6,] 0.499  0.7130 0.2614
#>  [7,] 0.705 -0.0691 0.3139
#>  [8,] 0.235  0.0234 0.6910
#>  [9,] 0.768 -0.0392 0.4216
#> 
#>                 [,1]  [,2]  [,3]
#> SS loadings    3.341 1.231 1.068
#> Proportion Var 0.371 0.137 0.119
#> Cumulative Var 0.371 0.508 0.627

  # using random starts
  data("WansbeekMeijer", package="GPArotation")
  fa.unrotated  <- factanal(factors = 3, covmat=NetherlandsTV, normalize=TRUE, rotation="none")
  # single rotation with a random start
  oblimin(loadings(fa.unrotated), Tmat=Random.Start(3))
#> Oblique rotation method Oblimin Quartimin converged.
#> Loadings:
#>           Factor1  Factor2  Factor3
#> NL1       0.81156 -0.04402  0.01607
#> TV2       0.79138  0.06961  0.00987
#> NL3       0.78396 -0.01027 -0.01260
#> RTL4      0.06646  0.74487 -0.08584
#> RTL5      0.05754  0.65088  0.08480
#> Veronica -0.06034  0.81206  0.04782
#> SBS6      0.00682  0.00651  0.99012
#> 
#>                Factor1 Factor2 Factor3
#> SS loadings      1.947   1.701   1.029
#> Proportion Var   0.278   0.243   0.147
#> Cumulative Var   0.278   0.521   0.668
#> 
#> Phi:
#>         Factor1 Factor2 Factor3
#> Factor1   1.000   0.610   0.402
#> Factor2   0.610   1.000   0.704
#> Factor3   0.402   0.704   1.000
  oblimin(loadings(fa.unrotated), randomStarts=1)
#> Oblique rotation method Oblimin Quartimin converged.
#> Loadings:
#>           Factor1  Factor2  Factor3
#> NL1       0.81155 -0.04402  0.01607
#> TV2       0.79137  0.06961  0.00987
#> NL3       0.78396 -0.01027 -0.01260
#> RTL4      0.06646  0.74487 -0.08584
#> RTL5      0.05754  0.65088  0.08480
#> Veronica -0.06034  0.81206  0.04782
#> SBS6      0.00683  0.00652  0.99011
#> 
#>                Factor1 Factor2 Factor3
#> SS loadings      1.947   1.701   1.029
#> Proportion Var   0.278   0.243   0.147
#> Cumulative Var   0.278   0.521   0.668
#> 
#> Phi:
#>         Factor1 Factor2 Factor3
#> Factor1   1.000   0.610   0.402
#> Factor2   0.610   1.000   0.704
#> Factor3   0.402   0.704   1.000
  # multiple random starts
  oblimin(loadings(fa.unrotated), 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
#> NL1       0.81155 -0.04402  0.01607
#> TV2       0.79137  0.06961  0.00987
#> NL3       0.78395 -0.01027 -0.01260
#> RTL4      0.06646  0.74487 -0.08584
#> RTL5      0.05754  0.65088  0.08480
#> Veronica -0.06034  0.81206  0.04782
#> SBS6      0.00682  0.00652  0.99012
#> 
#>                Factor1 Factor2 Factor3
#> SS loadings      1.947   1.701   1.029
#> Proportion Var   0.278   0.243   0.147
#> Cumulative Var   0.278   0.521   0.668
#> 
#> Phi:
#>         Factor1 Factor2 Factor3
#> Factor1   1.000   0.610   0.402
#> Factor2   0.610   1.000   0.704
#> Factor3   0.402   0.704   1.000

  # assessing local minima for box26 data
  data(Thurstone, package = "GPArotation")
  infomaxQ(box26, normalize = TRUE, randomStarts = 150)
#> Oblique rotation method Infomax converged at lowest minimum.
#> Of 150 random starts 100% converged, 33% at the same lowest minimum.
#> Random starts converged to 9 different local minima
#> Loadings at lowest minimum:
#>          [,1]      [,2]    [,3]
#>  [1,] -0.6979  0.733141  0.7412
#>  [2,]  0.6320 -0.445887  0.6788
#>  [3,]  0.7337  0.664606 -0.4890
#>  [4,] -0.0120  0.102265  0.9371
#>  [5,]  0.0554  0.906382  0.0942
#>  [6,]  0.8940  0.107167  0.0572
#>  [7,] -0.3082  0.407151  0.8923
#>  [8,]  0.3103 -0.144537  0.8407
#>  [9,] -0.2544  0.895283  0.3483
#> [10,]  0.3288  0.933975 -0.1188
#> [11,]  0.8560 -0.112853  0.2744
#> [12,]  0.8807  0.301270 -0.1431
#> [13,] -1.0675  1.018748 -0.0251
#> [14,]  1.0675 -1.018748  0.0251
#> [15,] -1.1473  0.000244  1.0721
#> [16,]  1.1473 -0.000244 -1.0721
#> [17,]  0.0080 -1.057386  1.0279
#> [18,] -0.0080  1.057386 -1.0279
#> [19,]  0.0904 -0.005924  0.9215
#> [20,]  0.0962  0.922324  0.0229
#> [21,]  0.8755  0.108684  0.0816
#> [22,]  0.1000  0.000594  0.8995
#> [23,]  0.1219  0.886070  0.0131
#> [24,]  0.8371  0.101842  0.1155
#> [25,]  0.3478  0.405151  0.4190
#> [26,]  0.4730  0.309916  0.3582
#> 
#>                 [,1]  [,2]  [,3]
#> SS loadings    8.668 8.559 8.183
#> Proportion Var 0.333 0.329 0.315
#> Cumulative Var 0.333 0.663 0.977
#> 
#> Phi:
#>       [,1]  [,2]  [,3]
#> [1,] 1.000 0.548 0.607
#> [2,] 0.548 1.000 0.543
#> [3,] 0.607 0.543 1.000
  geominQ(box26, normalize = TRUE, randomStarts = 150)
#> Oblique rotation method Geomin converged at lowest minimum.
#> Of 150 random starts 100% converged, 30% at the same lowest minimum.
#> Random starts converged to 4 different local minima
#> Loadings at lowest minimum:
#>           [,1]     [,2]      [,3]
#>  [1,] -0.01480 -0.00966  0.992934
#>  [2,]  0.94194  0.04837  0.059132
#>  [3,]  0.06094  0.96469 -0.000941
#>  [4,]  0.64216 -0.01368  0.637319
#>  [5,]  0.00200  0.64560  0.596311
#>  [6,]  0.61311  0.64401 -0.024703
#>  [7,]  0.38294  0.00698  0.834795
#>  [8,]  0.81242  0.03353  0.384702
#>  [9,] -0.02006  0.41809  0.788315
#> [10,]  0.02744  0.85762  0.444477
#> [11,]  0.76608  0.45213 -0.018941
#> [12,]  0.44155  0.78388 -0.028514
#> [13,] -0.82910  0.00876  0.747107
#> [14,]  0.82910 -0.00876 -0.747107
#> [15,]  0.00630 -0.82506  0.816393
#> [16,] -0.00630  0.82506 -0.816393
#> [17,]  0.84823 -0.79818 -0.008815
#> [18,] -0.84823  0.79818  0.008815
#> [19,]  0.71018 -0.02017  0.548317
#> [20,] -0.02338  0.68851  0.556578
#> [21,]  0.61810  0.63120 -0.006276
#> [22,]  0.70019 -0.00779  0.537441
#> [23,] -0.00946  0.68113  0.525157
#> [24,]  0.61769  0.59906  0.015742
#> [25,]  0.47824  0.46620  0.452515
#> [26,]  0.52777  0.48684  0.340444
#> 
#>                 [,1]  [,2]  [,3]
#> SS loadings    8.824 8.795 7.790
#> Proportion Var 0.339 0.338 0.300
#> Cumulative Var 0.339 0.678 0.977
#> 
#> Phi:
#>       [,1]  [,2]  [,3]
#> [1,] 1.000 0.268 0.206
#> [2,] 0.268 1.000 0.278
#> [3,] 0.206 0.278 1.000
  # for detailed investigation of local minima, consult package 'fungible' 
  # library(fungible)
  # faMain(urLoadings=box26, rotate="geominQ", rotateControl=list(numberStarts=150))
  # library(psych) # package 'psych' with random starts:
  # faRotations(box26, rotate = "geominQ", hyper = 0.15, n.rotations = 150)