使用Promax旋转时R中fa()与factanal()的因子载荷差异原因

我尝试使用base R中的factanal()函数和psych包中的fa()函数，对响应量表一致的问卷数据进行因子分析。两者均设置为*极大似然（Maximum likelihood）*方法并使用Promax旋转，但得到的因子载荷结果却不相同。以下是可复现示例：

# 加载所需分析包
library(psych) 
library(GPArotation)

# 加载mtcars数据集
data(mtcars) 

# 分别用两个函数执行因子分析
factanal_result <- factanal(mtcars, factors=3, rotation="promax") 
fa_result <- fa(mtcars, nfactors = 3, fm="ml", rotate="promax")

# 输出factanal的载荷结果（取前6行，保留2位小数）
head(round(factanal_result$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

# 输出fa的载荷结果（取前6行，保留2位小数）
head(round(fa_result$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

按照R文档建议设置oblique.Scores=TRUE后，因子载荷结果并未改变：

fa2_result <- fa(mtcars, nfactors = 3, fm="ml", rotate="promax", oblique.scores=T, scores="tenBerge")
head(round(fa2_result$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

此前相关讨论中提到的“结果一致”结论并不准确，该结论仅在使用varimax旋转时成立。

内容的提问来源于stack exchange，提问作者BPeif