数学建模常用算法汇总及python，MATLAB实现（七) —— sklearn和SPSS实现主成分分析

主成分分析算是比较重要的算法, 可以细看一下, 操作很简单, 理论部分主要看一下4, 5节, 不然不会用

主成分分析用python和SPSS实现都很方便

python相较于SPSS可以控制的更加细节

SPSS学起来则相对容易, 操作有图形界面

求求点赞关注呜呜

文章目录

数学建模常用算法汇总及python，MATLAB实现（七) —— sklearn和SPSS实现主成分分析
- 求求点赞关注呜呜
- 1. 何时需要主成分分析
- - 1. 降维
  - 2. 分析变量内部关系
- 2. 什么是主成分分析
- - (1) 基本思想
  - - 举例
    - 补充: 信息量和方差的关系
    - - 为什么方差越大携带信息越多
  - (2) 数学原理
  - - 数学原理一句话解释
    - 这样的超平面有以下两个性质*
    - 用二维来简单理解
    - 主成分的性质
- 3. 哪些变量适合做主成分分析
- - 1. Bartlett's球状检验
  - - python代码
  - 2. KMO检验
  - - python代码
- 4. 需要几个主成分
- - 1. 策略1:无脑2
  - 2. 策略2: 画碎石图
  - - 文字叙述: 怎么做?什么流程?
    - 举例
    - python实现代码
  - 策略3: 更进阶的用法
- 5. 如何评估分析结果
- 6. 实现
- - SPSS
  - pyhton
  - - 1. 导包
    - 2. 读数据
    - 3.基本信息
    - 4. 特征选择
    - 5. 开始用PCA降维

1. 何时需要主成分分析

主成分分析的作用:解决数据冗余问题, 给数据降维度, 对数据进行浓缩

1. 降维

我们在建模或者数据挖掘过程中可能会遇到:

数据的维度过大(人话就是:有很多属性, 比如一个样本有身高资产体重年龄一大堆不同属性的数据)
数据中含有噪声(错值)

等情况

这些情况下最适合用主成分分析先对数据降维再做后续处理

2. 分析变量内部关系

主成分也可以用作分析各个属性间的内部关系

因为主成分分析是将原来多个具有相关性的指标(属性)，重新组合成一组新的互相无关的 综合指标来代替原来的指标

tips

主成分分析主要的作用是降维, 分析变量内部关系一般是一个附属作用, 在分析之后可以得到一些系数, 从而顺便的观察到变量内部的相关关系

2. 什么是主成分分析

将相关性高的一组变量用一个变量代替, 新的变量尽量包含之前数据的信息, 这个新的变量就叫主成分, 主成分分析就是用已有的数据来得到一些主成分从而达到降维的目的

(1) 基本思想

将相关性高的一组变量用一个变量代替, 新的变量尽量包含之前数据的信息

通俗一点解释就是: 每个人都有身高和体重, 身高和体重这两个东西相关性就很高, 可以尝试用一个变量来代替身高和体重(当然会损失一些信息, 并且有偏差, 但是只要相关性足够高, 一般损失的信息不多, 偏差不大)

举例

如图所示, x1, x2 在散点图中的状态可以看出, 这俩玩意相关性很高1, 这俩玩意几乎就是一个函数关系(直线y1), 所以可以用一个变量代替这俩

先对x1, x2线性变换, 也就是把坐标轴改成y1, y2, 那么y1和y2都是x1, x2的线性组合
y1=k11x1+k12x2y2=k21x2+k22x2y_1 = k_{11}x_1 + k_{12}x_2\\ y_2 = k_{21}x_2 + k_{22}x_2 y1=k11x1+k12x2y2=k21x2+k22x2

这里就可以把y2丢弃了, 只保留y1作为x1, x2的替代变量

为什么可以这样做? 看下面

补充: 信息量和方差的关系

信息量和方差有关

看上图

在图中

看y1, y2坐标系, y1跨度比较大, y2变化非常小, 说明y1上方差比较大, 携带的信息多
这里就把y2甩掉

为什么方差越大携带信息越多

(例子是理解的关键)

举个简单例子理解

比如语文大家都是70分, 数学最低25, 最高90, 要看哪个学生优秀, 那就可以不看语文成绩了, 只看数学筛选学生优不优秀

(2) 数学原理

这部分要求不高, 看不看得懂都行, 有兴趣就简单了解即可

数学原理一句话解释

对于正交属性空间(所有属性也就是维度的空间)中的样本点, 用一个超平面(直线的高维推广, 比如在三维就是一个面)来表达所有的样本点(也就是把所有样本点映射到超平面上)

这样的超平面有以下两个性质*

样本点到这个超平面的距离都足够近
不同样本的投影在超平面上尽可能分开

用二维来简单理解

也就是样本点到直线的垂线长度加起来, 尽可能短

投影上去的方差尽可能大

因此这里红线的效果就比黑线好

我们的目标就是找到上面两个性质的最优值

主成分的性质

我们做完主成分分析之后会发现分析出来的主成分具有以下性质

分析出来的主成分有时不止一个, 主成分方差贡献率(就是上面图写的方差, 一个意思)都是逐级递减的, 第一个分析出来的主成分一定是方差贡献最多的
各个主成分之间是相互独立的

3. 哪些变量适合做主成分分析

我们遇到一堆数据, 不能因为我现在想要降维, 于是上来就无脑的主成分分析

在主成分分析之前有以下方法进行评估

下面的评估我只会用python评估, SPSS不知道, 可以自己搜一下怎么用, 其实大部分场景都不用评估

1. Bartlett’s球状检验

检验总体变量的相关矩阵是否是单位阵（相关系数矩阵对角线的所有元素均为1,所有非对角线上的元素均为零）；即检验各个变量是否各自独立。

python代码

(没有的库记得用pip安装在cmd(win)或者终端(linux/mac)输入pip install xx库)

from factor_analyzer.factor_analyzer import calculate_bartlett_sphericitychi_square_value, p_value = calculate_bartlett_sphericity(df)
print(chi_square_value, p_value)

2. KMO检验

python代码

检查变量间的相关性和偏相关性，取值在0-1之间；KOM统计量越接近1，变量间的相关性越强，偏相关性越弱，因子分析的效果越好。

from factor_analyzer.factor_analyzer import calculate_kmokmo_all, kmo_model = calculate_kmo(df)
print(kmo_all)

4. 需要几个主成分

所有主成分方差贡献率加起来应该大于等于80%

第五部分会讲怎么看方差贡献率

SPSS会自动帮你决定有几个主成分

如果使用python的sklearn包的话就需要自己判断一下了, 相对个性化一些

这里只讲python代码, 用SPSS的了解即可

1. 策略1:无脑2

先提供一种捷径, 大部分情况下2个左右的主成分是最好的

2. 策略2: 画碎石图

文字叙述: 怎么做?什么流程?

在真正进入主成分分析之前

我们先预先进行一次特殊的主成分分析

这次主成分分析我们设定主成分个数与我们输入数据的维度一样

然后再把每个主成分的方差贡献率输出或者可视化

看看需要几个主成分加起来才能达到80%(或者你想要的阈值, 一般阈值80%以上)的方差贡献率

举例

上图就是碎石图

比如输入数据有7维

那么分析之后的主成分也有7个

第一个主成分贡献率70%, 第二个20%, 后面的都差不多只有一点点, 这样判断下来, 就只需要两个主成分就可以获得原来数据90%的信息量了

那么后面的主成分分析, 我们只需要两个主成分来做主成分分析

python实现代码

用到了python的sklearn库

假设我们要分析的数据已经读出了, 是data

读数据这样读

from pandas import pd# 别的格式数据就改成read_xxx, 比如读excel就是read_excel
data = pd.read_csv('数据路径')

降维

from sklearn import decomposition
import seaborn as sns# 定义一个主成分分析(PCA)模型, 本来有个参数n_component
# 不填就是, 默认这次主成分分析我们设定主成分个数与我们输入数据的维度一样
pca_model = decomposition.PCA()# 对data进行降维, 确定参数(写在定义的模型内部了, 这个时候data还没有变)
pca_model.fit(data)# 改变data, 用上面一步解得的参数对data降维
data_test = pca_model.transform(data)# 得到贡献率
contribute_rate = pca_model.explained_variance_ratio_print(contribute_rate)# 画个图看看
sns.lineplot(x=np.arange(1, 8), y=contribute_rate)

策略3: 更进阶的用法

直接指定我需要多少的贡献率

比如我要97%的贡献率, 让sklearn自己来算如果要97%的贡献率需要几个主成分

那就在定义模型的时候直接指定 n_component=97%, 后面的操作一样

pca_model = decomposition.PCA(n_component=.97)

或者甚至可以用最大似然的方法决定有几个主成分, 篇幅有限这里不讲了, 可以自己搜一下

5. 如何评估分析结果

碎石图
把贡献率都输出出来看看

用SPSS分析之后就会有分析结果, 贡献率什么都有, 直接看就行, 或者再点一下碎石图画一下碎石图

6. 实现

SPSS

非常简单的操作

分析->降维->因子, 把想要降维的属性都点进去

然后降维就行了, 结果大概长这样

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-sIcXm7bp-1657561369934)(C:\Users\Lenovo\Desktop\SPSS.png)]

pyhton

1. 导包

In [4]:

import pandas as pd
import seaborn as sns
import numpy as np
import torch.nn as nn

2. 读数据

In [5]:

data = pd.read_excel('data/Corporation_evaluation.xlsx')

3.基本信息

In [6]:

data.head(3)

Out[6]:

	企业序号	净利润率/%	固定资产利润率/%	总产值利润率/%	销售收入利润率/%	产品成本利润率/%	物耗利润率/%	人均利润/(千元/人)	流动资金利润率/%
0	1	40.4	24.7	7.2	6.1	8.3	8.7	2.442	20.0
1	2	25.0	12.7	11.2	11.0	12.9	20.2	3.542	9.1
2	3	13.2	3.3	3.9	4.3	4.4	5.5	0.578	3.6

In [7]:

data.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 15 entries, 0 to 14
Data columns (total 9 columns):#   Column       Non-Null Count  Dtype
---  ------       --------------  -----  0   企业序号         15 non-null     int64  1   净利润率/%       15 non-null     float642   固定资产利润率/%    15 non-null     float643   总产值利润率/%     15 non-null     float644   销售收入利润率/%    15 non-null     float645   产品成本利润率/%    15 non-null     float646   物耗利润率/%      15 non-null     float647   人均利润/(千元/人)  15 non-null     float648   流动资金利润率/%    15 non-null     float64
dtypes: float64(8), int64(1)
memory usage: 1.2 KB

4. 特征选择

序号明显是无关数据, 去除掉, 复制到一个新的DataFrame里

In [31]:

data_without_serial = data.iloc[:, 1:8]

看一下协方差矩阵, 检查一下变量相关性

In [15]:

data_without_serial.corr()

Out[15]:

	净利润率/%	固定资产利润率/%	总产值利润率/%	销售收入利润率/%	产品成本利润率/%	物耗利润率/%	人均利润/(千元/人)
净利润率/%	1.000000	0.762986	0.701700	0.586838	0.595883	0.489560	0.597346
固定资产利润率/%	0.762986	1.000000	0.550389	0.466708	0.515790	0.419594	0.704630
总产值利润率/%	0.701700	0.550389	1.000000	0.840700	0.975962	0.816132	0.694095
销售收入利润率/%	0.586838	0.466708	0.840700	1.000000	0.866654	0.982280	0.492640
产品成本利润率/%	0.595883	0.515790	0.975962	0.866654	1.000000	0.866731	0.626012
物耗利润率/%	0.489560	0.419594	0.816132	0.982280	0.866731	1.000000	0.421641
人均利润/(千元/人)	0.597346	0.704630	0.694095	0.492640	0.626012	0.421641	1.000000

分析kmo系数, 检查特征是否适合降维(这套数据好像不太适合, 那就硬着头皮降)

In [17]:

import factor_analyzer.factor_analyzer as fa kmo_all, i = fa.calculate_kmo(data_without_serial)pd.DataFrame(zip(data_without_serial.columns, kmo_all), columns=['属性', 'kmo系数'])

Out[17]:

	属性	kmo系数
0	净利润率/%	0.417954
1	固定资产利润率/%	0.409777
2	总产值利润率/%	0.566945
3	销售收入利润率/%	0.594288
4	产品成本利润率/%	0.627004
5	物耗利润率/%	0.584813
6	人均利润/(千元/人)	0.455360

5. 开始用PCA降维

先用同等个数的特征降维, 观察贡献率来确定主成分个数

In [43]:

from sklearn import decompositionpca_model = decomposition.PCA()pca_model.fit(data_without_serial)data_test = pca_model.transform(data_without_serial)contribute_rate = pca_model.explained_variance_ratio_print(contribute_rate)sns.lineplot(x=np.arange(1, 8), y=contribute_rate)
[7.33866498e-01 2.02135729e-01 3.18938740e-02 2.93350166e-021.78296005e-03 7.18734974e-04 2.67187391e-04]

Out[43]:

<AxesSubplot:>

这里明显只需一两个主成分就可以占到80%到90%的信息量

这里设置主成分个数为2

In [53]:

pca_real = decomposition.PCA(n_components=2)pca_real.fit(data_without_serial)data_dir = pca_real.transform(data_without_serial)

评估降维结果

In [54]:

pd.DataFrame(zip([0, 1], pca_real.explained_variance_ratio_), columns=['特征序号', '信息量'])

Out[54]:

	特征序号	信息量
0	0	0.733866
1	1	0.202136

In [55]:

print('降维后总信息量', pca_real.explained_variance_ratio_.sum())
降维后总信息量 0.9360022269422483

In [56]:

data_dir

Out[56]:

array([[  5.5455228 ,  17.62345596],[  2.99790536,  -5.32484781],[-21.28863103,  -2.668911  ],[-12.90077775,   2.78150035],[ -1.25142989,   9.91612153],[ 20.83418515,  -7.94500664],[ -2.20824286,  -6.7415482 ],[ 13.38017708,  15.15368765],[ 34.97242118, -11.60230746],[ -1.76276938,  -3.48727368],[-19.49549195,  -2.14097427],[-34.52357753,  -6.42555221],[  1.89530172,   5.55388463],[  7.87905102,   5.58812265],[  5.92635609, -10.2803515 ]])