主成分分析

主成分能用来做什么

聚类、回归、寻找变量之间的共线性关系！！！

主成分分析的一般目的由两点组成：（参考王学民老师）

（1）将多个有相关关系的变量压缩成少数几个不相关的主成分（综合变量），并保留绝大部分信息；

（2）给出各主成分的具有实际背景和意义的解释

雷点在于！！！

多数论文以及博客使用主成分进行综合评价，王学民老师认为这是不对的，个人参考了一些这方面的研究后觉得王老师说的对。。。

因为个人理解能力有限，一知半解，只好放老师的文章来观摩学习

下面是复现一篇论文时写的python代码，除了最后评价分析这里不科学，其他的地方还行吧，权当存档留念，以后再也不傻乎乎的用主成分进行综合评价了。。。

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import seaborn as sns
from sklearn.preprocessing import scale
from sklearn.decomposition import PCA
from factor_analyzer.factor_analyzer import calculate_bartlett_sphericity
from factor_analyzer.factor_analyzer import calculate_kmo

# 设置中文字符及负号
plt.rcParams['font.sans-serif']=['Kaiti']
plt.rcParams['axes.unicode_minus']=False

# 读取文件
MVdata=pd.read_excel(r"C:\Users\Jr.hu\Desktop\指标正向化后.xlsx",index_col=0)

# 展示数据
print("原始数据展示:\n", MVdata)

# 首先对数据进行检验分析
# 1.Bartlett's球状检验
chi_square_value, p_value = calculate_bartlett_sphericity(MVdata)
print("bartlett球状检验参数：\n卡方值为：{}，p值为：{}".format(chi_square_value, p_value))

# 2.KMO检验
kmo_all, kmo_model = calculate_kmo(MVdata)
print("KMO检验参数：\n", kmo_model)
# 参数检验通过


# 这里需要先进行一些分析，为主成分分析做铺垫
# 对原始数据进行标准差标准化  这里采用课件中的方法
S=(MVdata-MVdata.mean())/MVdata.std()
print("标准化处理后的数据:\n", round(S,4))

# 求相关系数矩阵
CovX=np.around(np.corrcoef(MVdata.T),decimals=3)
print("相关系数矩阵为：\n",CovX)

# 求解特征值和特征向量
featValue,featVec=np.linalg.eig(CovX.T)
print("特征值为：\n",featValue)
print("特征向量为：\n",featVec)

# 对特征值进行排序输出
featValue=sorted(featValue,reverse=True)
print("特征值从高到低依次为：\n",featValue)

# 绘制特征值碎石图
plt.scatter(range(1,MVdata.shape[1]+1),featValue)
plt.plot(range(1,MVdata.shape[1]+1),featValue)
plt.hlines(y=1,xmin=0,xmax=14,colors='red',linestyles='--')
plt.title("Scree Plot")
plt.xlabel("Factors")
plt.ylabel("Eigenvalue")
plt.grid()  # 显示网格
plt.show()  # 显示图形

# 绘制累计贡献度折线图
gx=featValue/np.sum(featValue)   # 贡献度
lg=np.cumsum(gx)    # 累计贡献度
plt.scatter(range(1,MVdata.shape[1]+1),lg)
plt.plot(range(1,MVdata.shape[1]+1),lg)
plt.hlines(y=0.85,xmin=0,xmax=14,colors='red',linestyles='--')
plt.text(1,0.86,"85%水平线")
plt.title("累计贡献度")
plt.grid()  # 显示网格
plt.show()  # 显示图形
# 由碎石图及折线图综合分析可知，主成分分析应分为四类，k=4


# 第一步，计算主成分对象
pcal =PCA(n_components=4).fit(S)

# 第二步，确定主成分
V=pcal.explained_variance_   # 各个主成分所解释的数据方差
W=pcal.explained_variance_ratio_  # 各个主成分的方差解释比例
W.sum()    # 所有主成分解释的总方差
print("各个主成分所解释的数据方差为:\n", V)
print("各个主成分的方差贡献率为:\n", W)
print("累计方差贡献率为:\n", W.sum())

# 第三步，主成分负荷
Comp_P=pd.DataFrame(pcal.components_.T)   # 主成分载荷矩阵
print("主成分载荷矩阵为:\n", Comp_P)

# 这里可以绘制一张热力图表示原先的13个因子与降维后的4个主成分之间的关系
ax = sns.heatmap(Comp_P, annot=True, cmap="coolwarm")
ax.set_xticks([0.5, 1.5, 2.5, 3.5])
ax.set_xticklabels(['Comp1', 'Comp2', 'Comp3', 'Comp4'])  # 设置x轴的标签
ax.set_yticks([0.5,1.5,2.5,3.5,4.5,5.5,6.5,7.5,8.5,9.5,10.5,11.5,12.5])
ax.set_yticklabels(['x1','x2','x3','x4','x5','x6','x7','x8','x9','x10','x11','x12','x13'])  # 设置y轴的标签
plt.title("Factor Analysis")
plt.show()

# 第四步，先算各个主成分得分
Si=pcal.fit_transform(S)
SR=pd.DataFrame(Si,columns=['Comp1','Comp2','Comp3','Comp4'],index=MVdata.index)
print("主成分得分:\n", SR)

# 第五步，再算综合评价得分
SR['Comp']=Si.dot(W)
SR['Rank']=SR.Comp.rank(ascending=False)
print("主成分综合评价得分:\n", SR)


# 数据可视化一：箱线图
# 箱线图1.0(最终得分的箱线图)
sns.boxplot(y=SR['Comp'])   # y轴为各地最终的评价得分
plt.ylabel('综合评价最终得分')
plt.title('最终得分箱线图')
plt.tight_layout()
plt.show()

# 箱线图2.0(四个主成分得分的箱线图)
fig, axes = plt.subplots(2, 2, figsize=(10, 8))   # 创建包含四个子图的画布
sns.boxplot(y=SR['Comp1'], ax=axes[0, 0])
sns.boxplot(y=SR['Comp2'], ax=axes[0, 1])
sns.boxplot(y=SR['Comp3'], ax=axes[1, 0])
sns.boxplot(y=SR['Comp4'], ax=axes[1, 1])
plt.tight_layout()
plt.show()

# 数据可视化二：小提琴图
# 小提琴图1.0(最终得分小提琴图)
sns.violinplot(x='Comp',data=SR)
plt.title('最终得分小提琴图')
plt.xlabel('综合评价最终得分')
plt.tight_layout()
plt.show()

# 小提琴图2.0(四个主成分得分的小提琴图)
fig, axes = plt.subplots(2, 2, figsize=(10, 8))   # 创建包含四个子图的画布
sns.violinplot(y=SR['Comp1'], ax=axes[0, 0])
sns.violinplot(y=SR['Comp2'], ax=axes[0, 1])
sns.violinplot(y=SR['Comp3'], ax=axes[1, 0])
sns.violinplot(y=SR['Comp4'], ax=axes[1, 1])
plt.tight_layout()
plt.show()

# 可视化三：条图
# 这个没想到合适的数据，先放一放

# 可视化四：概率分布图
# 概率分布图1.0---数据标准化处理之后的概率分布图
sns.histplot(S,kde=True,color='blue',bins=30)   # 这里可以单独看S["x1"]...
plt.title('数据标准化处理后的概率分布图')  # 标准化处理后满足正态分布
plt.ylabel('概率密度')
plt.xlabel('数值')
plt.show()

# 概率分布图2.0---这里我们选取了最终评价得分Comp为y轴,主成分一得分Comp1为x轴,可以呈现出较好的效果
sns.jointplot(x=SR['Comp1'],y=SR['Comp'],data=SR)
plt.show()

# 概率分布图3.0---多个变量,四个主成分的得分
sns.pairplot(SR, vars=['Comp1','Comp2','Comp3','Comp4'])
plt.show()

# 可视化五：柱状图
# 柱状图1.0---各个主成分的方差贡献率
x = ['Comp1','Comp2','Comp3','Comp4']
y = ['Comp1','Comp1、2','Comp1、2、3','Comp1、2、3、4']
plt.grid(ls="--", alpha=0.5)
plt.title("各个主成分的方差贡献率")
bars = plt.bar(x, W)
for i, bar in enumerate(bars):   # 添加百分比标签
    plt.text(bar.get_x() + bar.get_width() / 2, bar.get_height(), f"{W[i]*100:.2f}%",
             ha='center', va='bottom')
plt.show()

# 柱状图2.0---累计方差贡献率柱状图
cumulative_variance = [sum(W[:i+1]) for i in range(len(W))]    # 计算累计方差贡献率
plt.figure()
plt.grid(ls="--",alpha=0.5)
plt.title("累计方差贡献率")
bars = plt.bar(y, cumulative_variance)
for i, bar in enumerate(bars):      # 添加百分比标签
    plt.text(bar.get_x() + bar.get_width() / 2, bar.get_height(), f"{cumulative_variance[i]*100:.2f}%",
             ha='center', va='bottom')
plt.show()

# 柱状图3.0---综合评级得分柱状图
plt.figure()
plt.grid(ls="--",alpha=0.5)
plt.title("各省份综合评价得分条形图")
bars = plt.bar(SR['Rank'], SR['Comp'])
text= [
    "北京市", "天津市", "河北省", "山西省", "内蒙古自治区",
    "辽宁省", "吉林省", "黑龙江省", "上海市", "江苏省",
    "浙江省", "安徽省", "福建省", "江西省", "山东省",
    "河南省", "湖北省", "湖南省", "广东省", "广西壮族自治区",
    "海南省", "重庆市", "四川省", "贵州省", "云南省",
    "西藏自治区", "陕西省", "甘肃省", "青海省", "宁夏回族自治区",
    "新疆维吾尔自治区"]
for bar, province in zip(bars, text):
    yval = bar.get_height()
    plt.text(bar.get_x() + bar.get_width()/2, yval, province, ha='center', va='bottom')
plt.show()

# 主成分图（取前两个主成分）
plt.plot(SR.Comp1,SR.Comp2,'.')
A1=SR.Comp1
A2=SR.Comp2
AI=SR.index
plt.hlines(0,A1.min(),A1.max(),linestyles="dotted")
plt.vlines(0,A2.min(),A2.max(),linestyles="dotted")
for i in range(SR.index.shape[0]):
    plt.text(A1[i],A2[i],AI[i])
plt.show()

关于相关系数矩阵与协方差矩阵

主成分分析（PCA）既可以使用协方差矩阵，也可以使用相关系数矩阵，具体选择哪种取决于数据的特性和分析的目的。以下是对两种矩阵在PCA中使用的考量因素：

协方差矩阵： - 适用情况：当原始数据的各个变量具有不同的量纲（如长度、重量、时间等）或者变量的尺度差异较大时，使用协方差矩阵更为合适。协方差矩阵直接反映了变量间因数值大小变化引起的线性关系，保留了原始变量的量纲信息。 - 优点：能够捕捉变量间原始的、未经标准化的波动关系，对于那些量纲差异显著且这种差异对分析有重要意义的情况，使用协方差矩阵能够保持这些差异对主成分的影响。

相关系数矩阵： - 适用情况：当原始数据的各个变量在同一量纲下或者经过标准化处理（即所有变量均转化为无量纲的、均值为0、标准差为1的状态）时，使用相关系数矩阵更为合适。相关系数矩阵衡量的是变量间关系的方向和强度，不受变量尺度的影响。 - 优点： - 消除量纲影响：由于相关系数矩阵中的元素是标准化数据的协方差，它已经消除了变量各自尺度的影响，只反映变量间相对变化的一致性。这使得不同量纲或不同尺度的变量在PCA中得以公平比较，避免了量纲较大或波动范围较大的变量对分析结果的主导。 - 易于解释：相关系数的取值范围在-1到1之间，可以直接解读为变量间的线性相关程度，更便于理解和解释变量间的关系。

选择原则： - 如果数据集中的变量具有相同的量纲或已经过标准化处理，通常建议使用相关系数矩阵，因为它能够更纯粹地反映变量间的相对关系，不受变量自身尺度的干扰。 - 如果变量间存在显著的量纲差异，并且这些差异对分析结果有重要含义

sklearn库默认使用协方差矩阵

sklearn 库在执行主成分分析（PCA）时，默认使用的是协方差矩阵。这是因为 sklearn.decomposition.PCA 类在计算主成分时，默认计算的是数据的协方差矩阵。然而，用户可以根据实际需求选择是否对输入数据进行预处理（如标准化），从而转换为使用相关系数矩阵进行PCA。

以下是使用 sklearn 进行PCA时与矩阵选择相关的几点说明：

1. 默认行为：如果不进行额外的预处理，直接将原始数据传递给 sklearn.decomposition.PCA 类的实例，并进行拟合（fit），那么PCA将基于数据的协方差矩阵来确定主成分。

2. 标准化选项：如果希望使用相关系数矩阵进行PCA，通常需要先对数据进行标准化，即将数据转换为无量纲的、均值为0、标准差为1的形式。在 sklearn 中，可以使用 preprocessing 模块中的 StandardScaler 类来完成这一操作。对数据标准化后再进行PCA，此时虽然PCA内部依然计算的是标准化后数据的协方差矩阵，但由于标准化数据的均值和方差已调整为一致，此时的协方差矩阵实际上等价于相关系数矩阵，因为标准化消除了量纲和尺度的影响。

例如，一个完整的流程可能如下：

from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA

# 假设 data 是未标准化的原始数据
scaler = StandardScaler()
scaled_data = scaler.fit_transform(data)

# 创建并拟合PCA模型，此时使用的是等价于相关系数矩阵的标准化数据的协方差矩阵
pca = PCA()
pca.fit(scaled_data)

3. 参数设置：sklearn 的 PCA 类提供了 covariance 参数（从 scikit-learn 1.1 及更高版本开始），允许用户直接指定是否基于数据的协方差矩阵还是基于相关系数矩阵来进行PCA计算。不过，需要注意的是，设置 covariance='corr' 时，sklearn 会自动对数据进行标准化，然后计算相关系数矩阵。因此，对于用户来说，直接对数据进行标准化仍然是推荐的做法，而不需要显式设置 covariance 参数。

总结起来，在使用 sklearn 库进行主成分分析时，默认是基于协方差矩阵。若要基于相关系数矩阵进行分析，通常做法是先对数据进行标准化处理，确保变量间的关系不受量纲和尺度差异的影响。尽管 sklearn 提供了参数来直接指定基于相关系数矩阵计算，但实践中更常见的是通过预处理步骤标准化数据，间接实现基于相关系数矩阵的PCA。

错上加错（载荷矩阵）

关于载荷矩阵的一点思考与说明

在上面的处理过程中，使用到了载荷矩阵，但是按照同样的思路使用spss来分析，得到的主成分系数矩阵不一样，结果相差很大

# 第三步，主成分负荷
Comp_P=pd.DataFrame(pcal.components_.T)   # 主成分载荷矩阵
print("主成分载荷矩阵为:\n", Comp_P)

上面的代码中关于载荷矩阵的计算是没有问题的，载荷矩阵描述了主变量与原始变量之间的相关系数，各个主成分的得分系数矩阵需要使用载荷矩阵除以对应特征值的平方根，这个系数才能表示各个变量与主变量之间的线性关系

在上面的代码中，使用sklrean库进行分析直接用载荷矩阵描述各个变量与主变量之间的线性关系，缺少了一步转换

# 主成分得分系数矩阵  
for i in range(4): 
    Comp_P[i]=Comp_P[i]/math.sqrt(V[i])    # 计算公式为载荷矩阵除以对应主成分的特征值，这个操作是参考这篇文章：https://zhuanlan.zhihu.com/p/49481213
print("主成分得分矩阵:\n", Comp_P)

至于使用spss分析时，做主成分分析不要选择旋转,否则就变成了因子分析,结果当然不一样

关于因子分析与主成分分析的区别和联系，下面放了一篇很好的文章

但是换句话说，反正都错了，错上加错也无所谓了。。。。 ## 参考借鉴

做综合评价：

http://t.csdnimg.cn/VvKCZ

http://t.csdnimg.cn/ak1oY

http://t.csdnimg.cn/V8Cak

做降维、聚类：

https://www.cnblogs.com/RHadoop-Hive/p/9441505.html

做回归：

主成分分析 - Sophie的文章 - 知乎 https://zhuanlan.zhihu.com/p/653310400

关于载荷矩阵的说明：

主成分分析PCA之Python与SPSS展示的区别 - canhui87的文章 - 知乎

https://zhuanlan.zhihu.com/p/76622543

如何做主成分分析和因子分析？它们的区别与联系在哪里？ - 胡保强的文章 - 知乎

https://zhuanlan.zhihu.com/p/49481213

一起进步！