因子分析factor_analyzer工具包介绍
#此代码是我基于因子分析的专用工具包(factor_analyzer)搭建的
#此代码能解决需要输出多个(n_factors=1,2,3,...m)潜变量的因子分析问题
#如需将此代码放在自己个人社交平台，请给个引用"B站 耿大哥讲算法"
import numpy as np
from sklearn import preprocessing
from factor_analyzer import calculate_kmo
from factor_analyzer import FactorAnalyzer
from factor_analyzer import calculate_bartlett_sphericity
#读取数据&校准化
# 指标1 指标2 指标3 指标4指标5 指标6 指标7 指标8
data=np.array([[40.4,24.7,7.20,6.10,8.30,8.70,2.442,20.0], #样本1
[25.0,12.7,11.2,11.0,12.9,20.2,3.542,9.10], #样本2
[13.2,3.30,3.90,4.30,4.40,5.50,0.578,3.60], #样本3
[22.3,6.70,5.60,3.70,6.00,7.40,0.176,7.30], #样本4
[34.3,11.8,7.10,7.10,8.00,8.90,1.726,27.5], #样本5
[35.6,12.5,16.4,16.7,22.8,29.3,3.017,26.6], #样本6
[22.0,7.80,9.90,10.2,12.6,17.6,0.847,10.6], #样本7
[48.4,13.4,10.9,9.90,10.9,13.9,1.772,17.8], #样本8
[40.6,19.1,19.8,19.0,29.7,39.6,2.449,35.8], #样本9
[24.8,8.00,9.80,8.90,11.9,16.2,0.789,13.7], #样本10
[12.5,9.70,4.20,4.20,4.60,6.50,0.874,3.90], #样本11
[1.80,0.60,0.70,0.70,0.80,1.10,0.056,1.00], #样本13
[32.3,13.9,9.40,8.30,9.80,13.3,2.126,17.1], #样本13
[38.5,9.10,11.3,9.50,12.2,16.4,1.327,11.6]])#样本14
# 对样本Z-score标准化
X=preprocessing.scale(data,axis=0)
#KMO和barlett检验
bartlett,KMO=calculate_bartlett_sphericity(data),calculate_kmo(data)[1]
print('-'*54);print('KMO和bartlett检验')
print('KMO='+str('%.3f'%calculate_kmo(data)[1]),end='\t')
if KMO<=1 and KMO>=0.9:print('非常好')
elif KMO<0.9 and KMO>=0.8:print('好')
elif KMO<0.8 and KMO>=0.7:print('一般')
elif KMO<0.7 and KMO>=0.6:print('差')
elif KMO<0.6 and KMO>=0.5:print('很差')
else:print('不能接受')
print('bartlett 统计量='+str('%.3f'%bartlett[0]))
print('bartlett 自由度='+str(int(data.shape[1]*(data.shape[1]-1)/2)))
print('bartlett P值='+str('%.3f'%bartlett[1]))
#(常用)method='minres'迭代主因子,'principal'主成分
#(不常用)method='ml'极大似然,'mle'极大似然估计,'uls'无条件最小二乘,
#(常用)rotation=varimax标准化最大方差法
#(不常用)rotation=promax,oblimin,oblimax,quartimin,quartimax,equamax
FA=FactorAnalyzer(n_factors=2,method='minres',rotation='varimax')
FA.fit(X)
#计算载荷矩阵
A=FA.loadings_
#计算样本得分矩阵
F=FA.transform(X)
#计算得分系数矩阵
S=np.mat(X).I*F
#计算共同度
h=FA.get_communalities()
#计算方差贡献
g=FA.get_factor_variance()
#计算旋转矩阵
Q=FA.rotation_matrix_
# 输出使载荷间的总方差函数达到最大的旋转矩阵Q
print('-'*50)
print('旋转矩阵Q')
for i in range(Q.shape[0]):
for j in range(Q.shape[1]):print('%.3f'%Q[i,j],end='\t')
print()
# 输出载荷矩阵的旋转解
print('-'*50);print('载荷矩阵的旋转解')
print('指标|组件',end='\t')
for i in range(A.shape[1]):print('潜变量f'+str(i+1),end='\t')
print('共同度')
for i in range(A.shape[0]):
print('指标'+str(i+1),end='\t')
for j in range(A.shape[1]):print('%.3f'%A[i,j],end='\t')
print('%.3f'%h[i],end='\t')
print()
print('总贡献',end='\t')
for i in g[0]:print('%.3f'%i,end='\t')
print('%.3f'%sum(g[0]))
print('贡献率',end='\t')
for i in g[1]:print('%.1f'%(i*100)+'%',end='\t')
print('-----')
print('贡献累积',end='\t')
for i in g[2]:print('%.1f'%(i*100)+'%',end='\t')
print('-----')
print('-'*50)
#输出旋转后的成分得分系数矩阵
print('成分得分系数矩阵(旋转后)')
print('指标|组件',end='\t')
for i in range(A.shape[1]):print('潜变量f'+str(i+1),end='\t')
print()
for i in range(S.shape[0]):
print('指标'+str(i+1),end='\t')
for j in range(S.shape[1]):print('%.3f'%S[i,j],end='\t')
print()
#计算各样本的综合得分并排序
print('-'*50)
print('各样本的得分及排序');print('样本|fi',end='\t')
for i in range(A.shape[1]):print('潜变量f'+str(i+1),end='\t')
print('综合变量f',end='\t')
for i in range(A.shape[1]):print('按f'+str(i+1)+'排序',end='\t')
print('按f排序')
def xvhao(x):
list2=[i for i in x]
list2.sort(reverse=True)
return [list2.index(i)+1 for i in x]
list3=[i/sum(g[1]) for i in g[1]]
list4=[sum(list3[j]*F[i,j] for j in range(F.shape[1])) for i in range(F.shape[0])]
for i in range(F.shape[0]):
print('样本'+str(i+1),end='\t')
for j in range(F.shape[1]):print('%.3f'%F[i,j],end='\t')
print('%.3f'%list4[i],' ',end='\t')
for j in range(F.shape[1]):print(xvhao(F[:,j])[i],' ',end='\t')
print(xvhao(list4)[i],' ',end='\t')
print()
print('-'*50)
