主成分分析(PCA)-Python代码

主成分分析是利用降维的思想，在损失很少信息的前提下把多个指标转化为几个综合指标的多元统计方法。
算法的具体步骤如下：
1）对向量X进行去中心化。
2）计算向量X的协方差矩阵，自由度可以选择0或者1。
3）计算协方差矩阵的特征值和特征向量。
4）选取最大的k个特征值及其特征向量。
5）用X与特征向量相乘。

代码如下：

一、导入库

# 数据处理
import pandas as pd
import numpy as np
 
# 绘图
import seaborn as sns
import matplotlib.pyplot as plt

二、读取数据集

df = pd.read_csv(r"C:\Users\1\Desktop\mydata.csv", encoding='gbk', index_col=0).reset_index(drop=True)
## 数据来源：国家统计局  2020年31个城市空气质量情况  http://www.stats.gov.cn/tjsj/ndsj/2021/indexch.htm
df

三、球状检验

# Bartlett's球状检验
from factor_analyzer.factor_analyzer import calculate_bartlett_sphericity
 
chi_square_value, p_value = calculate_bartlett_sphericity(df)
print(chi_square_value, p_value)

四、KMO检验

# 检查变量间的相关性和偏相关性，取值在0-1之间；KOM统计量越接近1，变量间的相关性越强，偏相关性越弱，因子分析的效果越好。
# 通常取值从0.6开始进行因子分析
from factor_analyzer.factor_analyzer import calculate_kmo
 
kmo_all, kmo_model = calculate_kmo(df)
print(kmo_all)

五、相关矩阵

### step-1  标准化
from sklearn import preprocessing
df = preprocessing.scale(df)
df

### step-2 相关系数矩阵
covX = np.around(np.corrcoef(df.T),decimals=3)
covX

### step-3 特征值和特征向量
featValue, featVec=  np.linalg.eig(covX.T)  #求解系数相关矩阵的特征值和特征向量
featValue, featVec

### step-4 对特征值进行排序并输出  降序
featValue = sorted(featValue)[::-1]
featValue

六、绘制散点图和折线图

# 同样的数据绘制散点图和折线图
plt.scatter(range(1, df.shape[1] + 1), featValue)
plt.plot(range(1, df.shape[1] + 1), featValue)
 
# 显示图的标题和xy轴的名字
# 最好使用英文，中文可能乱码
plt.title("Scree Plot")  
plt.xlabel("Factors")
plt.ylabel("Eigenvalue")
 
plt.grid()  # 显示网格
plt.show()  # 显示图形

七、特征值贡献度

gx = featValue/np.sum(featValue)
gx

八、特征值的累计贡献率

lg = np.cumsum(gx)
lg

九、排序选出主成分

#选出主成分
k=[i for i in range(len(lg)) if lg[i]<0.85]
k = list(k)
print(k)

十、选出主成分对应的特征向量矩阵

selectVec = np.matrix(featVec.T[k]).T
selectVe=selectVec*(-1)
selectVec

十一、主成分得分

finalData = np.dot(df,selectVec)
finalData

十二、热力图

# 绘图
 
plt.figure(figsize = (14,14))
ax = sns.heatmap(selectVec, annot=True, cmap="BuPu")
 
# 设置y轴字体大小
ax.yaxis.set_tick_params(labelsize=15)
plt.title("Factor Analysis", fontsize="xx-large")
 
# 设置y轴标签
plt.ylabel("Sepal Width", fontsize="xx-large")
# 显示图片
plt.show()
 
# 保存图片
# plt.savefig("factorAnalysis", dpi=500)