I have to perform a factor analysis. My data is a questionnaire with same answers possible for each question, so no scale problem.
I tried to use factanal()
from base R and fa()
from the psych
package.
Even with having Maximum likelihood set for both methods, using the "promax" rotation gives different factor loadings.
The question has already been asked there : https://stackoverflow.com/questions/50573764/differences-between-fa-and-factanal-functions-in-r#:~:text=factanal%20performs%20a%20maximum%2Dlikelihood,least%20square%20regressions%20(OLR).
but the answer is not completely true, the article they link to say "same results" but it is same results only when we use "varimax" rotation.
So i wonder why there is a difference in factor loadings (psych
package offers more things, but why is it different when it is supposed to be the same?)
This is a reproducible example:
# Load the psych package for
# data analysis and visualization
library(psych)
library(GPArotation)
# Load the mtcars dataset
data(mtcars)
# Perform factor analysis on the mtcars dataset
factanal <- factanal(mtcars, factors=3, rotation="promax")
fa <- fa(mtcars, nfactors = 3, fm="ml", rotate="promax")
# Print the results
head(round(factanal$loadings, 2))
Factor1 Factor2 Factor3
mpg 0.54 0.16 -0.43
cyl -0.53 -0.58 0.08
disp -0.65 -0.32 0.19
hp -0.11 -0.48 0.47
drat 0.83 0.12 0.10
wt -0.71 0.18 0.54
head(round(fa$loadings, 2))
ML2 ML1 ML3
mpg 0.52 0.17 -0.43
cyl -0.51 -0.59 0.08
disp -0.63 -0.33 0.19
hp -0.09 -0.49 0.48
drat 0.83 0.13 0.10
wt -0.70 0.16 0.54
In R documentation, they recommend to use "oblique.Scores=TRUE", but that does not change the factor loadings:
fa2 <- fa(mtcars, nfactors = 3, fm="ml", rotate="promax", oblique.scores=T, scores="tenBerge")
head(round(fa2$loadings, 2))
ML2 ML1 ML3
mpg 0.52 0.17 -0.43
cyl -0.51 -0.59 0.08
disp -0.63 -0.33 0.19
hp -0.09 -0.49 0.48
drat 0.83 0.13 0.10
wt -0.70 0.16 0.54
You can get the same rotation as stats::promax
, used by factanal
, by using Promax
(note the upper case P). The help page of ?Promax
explains:
Promax is a very direct adaptation of the stats::promax function. The addition is that it will return the interfactor correlations as well as the loadings and rotation matrix.
I have not looked in the docs to see what promax
, from the psych
package, is doing (I'm sure it will be in there as they are very good) but a quick look at the source code for fa
reveals that it returns a rotation of the normalised loadings:
Promax = {pro <- Promax(loadings,...) #Promax without Kaiser normalization
loadings <- pro$loadings
Phi <- pro$Phi
rot.mat <- pro$rotmat},
promax = {#pro <- stats::promax(loadings,...) #from stats
pro <- kaiser(loadings,rotate="Promax",...) #calling promax will now do the Kaiser normalization before doing Promax rotation
Code showing the equivalence:
library(psych)
data(mtcars)
m1 <- factanal(mtcars, factors=3, rotation="promax")
m2 <- fa(mtcars, nfactors = 3, fm="ml", rotate="Promax")
loadings(m2)[] - loadings(m1)[]
# ML2 ML1 ML3
# mpg -1.705875e-06 1.473799e-06 -3.059439e-07
# cyl 1.594591e-07 -8.884557e-07 -1.818443e-07
# disp 8.168498e-07 4.712348e-07 2.310251e-06
# hp 6.522325e-07 -1.931424e-06 -1.722668e-06
# drat 8.599102e-07 1.509280e-06 8.609054e-07
# wt 2.121750e-06 8.515195e-07 3.925057e-06
# qsec 1.092731e-06 1.119313e-06 2.633606e-06
# vs 4.688156e-07 -8.692953e-07 -6.737374e-07
# am -1.884347e-06 1.230842e-06 -1.514240e-06
# gear -1.691097e-06 -2.111156e-06 -5.330317e-06
# carb 1.288899e-06 -2.468629e-06 -2.641025e-06