Should I reorder principal components after rotation?

Question

I recently noticed that psych::principal reorders principal components on (automatic) rotation, according to their Eigenvalues (from highest to lowest). (Recall that rotation matrix-multiplies the loadings, and so the order of their squared column sums (fka. the Eigenvalues) can change, too).

Here's an example:

library(qmethod)
library(psych)
data("lipset")  # this dataset is used because it causes a re-ordering of components
Lipset <- cor(x = lipset[[1]], method = "pearson")  # must calculate cor matrix first

  # calculate unrotated loadings:
principal.unrotated <- principal(r = Lipset, nfactors = 4, rotate = "none")$loa  
  # calculate varimax rotated loadings:
principal.varimax   <- principal(r = Lipset, nfactors = 4, rotate = "varimax")$loa  
  # manually calculate varimax rotmat on unrotated loadings:
rot.mat.varimax     <- varimax(x = principal.unrotated)$rotmat  

  # should manually reproduce the varimax rotation:
repr.varimax <- unclass(principal.unrotated) %*% rot.mat.varimax  
repr.varimax
#>             [,1]        [,2]        [,3]        [,4]
#> US1 -0.229427530  0.15096123  0.81283285  0.06465534
#> US2  0.002334831 -0.11383236  0.89114173 -0.06389700
#> US3 -0.009194167  0.79325633  0.03603989  0.18951141
#> US4  0.255174168  0.76681887  0.26319367 -0.08631099
#> JP5  0.003227787 -0.87396713  0.21605277  0.05615580
#> CA6  0.922371369  0.08409883 -0.01191349 -0.08215986
#> UK7  0.823285358  0.07913797 -0.17255592  0.03003107
#> US8 -0.447664930  0.02677531  0.37686206 -0.60826777
#> FR9 -0.158333934  0.06509079  0.11083125  0.87837316
  # notice how eigenvalue order is out of whack:
apply(X = repr.varimax, MARGIN = 2, FUN = function(x) sum(x^2))  
#> [1] 1.871892 2.035122 1.756305 1.203962

unclass(principal.varimax)  # notice how cols 1 and 2 have changed
#>             PC2          PC1         PC3         PC4
#> US1  0.15096123 -0.229427530  0.81283285  0.06465534
#> US2 -0.11383236  0.002334831  0.89114173 -0.06389700
#> US3  0.79325633 -0.009194167  0.03603989  0.18951141
#> US4  0.76681887  0.255174168  0.26319367 -0.08631099
#> JP5 -0.87396713  0.003227787  0.21605277  0.05615580
#> CA6  0.08409883  0.922371369 -0.01191349 -0.08215986
#> UK7  0.07913797  0.823285358 -0.17255592  0.03003107
#> US8  0.02677531 -0.447664930  0.37686206 -0.60826777
#> FR9  0.06509079 -0.158333934  0.11083125  0.87837316
  # eigenvalue order is fine again:
apply(X = principal.varimax, MARGIN = 2, FUN = function(x) sum(x^2))  
#>      PC2      PC1      PC3      PC4 
#> 2.035122 1.871892 1.756305 1.203962

I get how and why this works.

My question is simply: is there a reason – other than convention or convenience – why this reordering would be necessary?

(This causes a bunch of problems for my use case in another function, and I'd like to avoid if that is statistically sound).

The way I see it, rotated principal components are no longer principal components anyway, so you might as well leave them in any order they come in.

Answer 1

As the author of the psych package I will try to explain what is going on and offer a solution.

First, as @ttnhpns correctly points out, rotated principal components are no longer principal components, they are merely components. @ttnhpns is also correct that we should not call the sums of squares of these rotated components "eigen values" but rather "after rotation Sums of squares". (Indeed, that is how they are labeled in the output.) However, common usage, unfortunately, tends to call them eigen values.

Then, when rotating a solution, the order of the sums of squares of the rotated components will not necessarily be monotonically decreasing. A conventional solution to this problem is to order by Sums of Squares and then just name the components RC1 ... RCn and forget about it. I prefer to keep the names of the unrotated components to highlight what rotation does. (Which is to change the amount of importance of each component.)

This has been asked a lot and in fact has led to the final comment on the fa help page (and probably should be added the principal help page as well):

A frequently asked question is why are the factor names of the rotated solution not in ascending order? That is, for example, if factoring the 25 items of the bfi, the factor names are MR2 MR3 MR5 MR1 MR4, rather than the seemingly more logical "MR1" "MR2" "MR3" "MR4" "MR5". This is for pedagogical reasons, in that factors as extracted are orthogonal and are in order of amount of variance accounted for. But when rotated (orthogonally) or transformed (obliquely) the simple structure solution does not preserve that order. The factor names are, of course, arbitrary, and are kept with the original names to show the effect of rotation/transformation.

To make the order of the rotated / transformed components match that of the unrotated, you can just reorder the components using the fa.organize function.

library(qmethod)
library(psych)
data("lipset")
lip   <- lipset$ldata
p4   <- principal(lip,4)               # the default is to rotate using varimax
p4.n <- principal(lip,4,rotate="none") # extract the unrotated components
v4   <- varimax(p4.n$loadings)         # r otate them using varimax
fa.congruence(p4,v4)
      PC1  PC2   PC3   PC4
PC2  0.14 1.00  0.03  0.06
PC1  1.00 0.14 -0.25  0.03
PC3 -0.25 0.03  1.00 -0.10
PC4  0.03 0.06 -0.10  1.00
p4.o <- fa.organize(p4,paste0("PC",1:4))
fa.congruence(v4,p4.o)  # just show the rotated and organized solutions

      PC1  PC2   PC3   PC4
PC1  1.00 0.14 -0.25  0.03
PC2  0.14 1.00  0.03  0.06
PC3 -0.25 0.03  1.00 -0.10
PC4  0.03 0.06 -0.10  1.00

If enough people request it, the next psych release (which won't be for several months) can include a do not reorder option, but I prefer the fa.reorganize option.