While performing an exploratory factor analysis, I want to see the explained variance by the factors. It is straightforward and I printed the output of the fa() method. One can get the same summary table as the output of the $loadings, as shown below. The problem is that the explained variance deviates in the 2 tables. What is the reason for this, is it a bug in the psych package? Which output is the correct one?
library(psych)
fa_results <- fa(df_sq1sq9, nfactors = 5, rotate = "oblimin")
print(fa_results)
fa_results$loadings
Please consider that I am not talking about the factor loadings but the Cumulative var row in the below tables. The discrepancies appear in the case of orthogonal and also in oblique rotations, and the deviations of the two outputs are greater in the case of oblique rotations. Therefore, the example is shown with oblimin rotation.
Oblimin Rotation:
>print(fa_results)
[...]
MR1 MR3 MR5 MR2 MR4
SS loadings 2.45 1.50 1.07 1.05 0.60
Proportion Var 0.27 0.17 0.12 0.12 0.07
Cumulative Var 0.27 0.44 0.56 0.68 0.74
Proportion Explained 0.37 0.22 0.16 0.16 0.09
Cumulative Proportion 0.37 0.59 0.75 0.91 1.00
With factor correlations of
MR1 MR3 MR5 MR2 MR4
MR1 1.00 0.59 0.64 -0.01 0.30
MR3 0.59 1.00 0.33 0.12 0.21
MR5 0.64 0.33 1.00 0.28 0.11
MR2 -0.01 0.12 0.28 1.00 -0.17
MR4 0.30 0.21 0.11 -0.17 1.00
Oblimin Rotation:
>fa_results$loadings
[...]
MR1 MR3 MR5 MR2 MR4
SS loadings 2.310 1.350 0.975 1.038 0.563
Proportion Var 0.257 0.150 0.108 0.115 0.063
Cumulative Var 0.257 0.407 0.515 0.630 0.693
The question I also raised as an update in https://stats.stackexchange.com/questions/657692/exploratory-factor-analysis-oblique-rotation-variance-explained
$endgroup$
The problem is in the print function in core R for an object of class loading. It automatically finds the column sums of the squared loadings. This ignores the fact that the factors are correlated.
Thus, this value will vary by rotation (which it should not, because the total amount of variance accounted for remains constant across rotations.
library(psych)
f5 <- fa(bfi[1:25],5) # 5 factors default to oblimin
f5l <- f5$loadings
tr(t(f5l) %*% f5l %*% f5$Phi) #this adjusts for the correlations
tr(t(f5l)%*% f5l) #this is the same sum of squares and is fine for orthogonal factors but give an incorrect answer for oblique factors.
Compare this to the unrotated solution
f5.none <- fa(bfi[1:25],5, rotate="none")
sum(f5.none$Vaccounted[1,]) #total variance accounted for matches the oblique above.
Then, although you did not ask, you might look at the communalities which in the orthogal solution are just the row sums of squares. In the oblique case, I find them by the diagonal of the reproduced correlation matrix:
h2 <- diag(f5l %*% f5$Phi %*% t(f5l))
sum(h2)
#is the same as the total variance accounted for from above.
$endgroup$
1