主成分分析(PCA)
好多小伙伴在学习数学建模的时候会遇到这么一个算法—主成分分析,顾名思义,是用来分析主要成分的。简单点说,就是对评价指标集进行分析,得到新的若干个成分,来代替原始的评价指标(听不懂没关系!继续往下看,保证可以听得懂!)。
它是一种普遍的降维算法,在机器学习、生物分析等各个领域都有着广泛的应用。接下来,跟着我的思路,带你掌握这个算法,并学会如何去应用!
1 例题引入
大家设想一下这个例子,为了衡量一个人的体育综合能力,我们设有如下评价指标集{BMI(身高体重),肺活量,400米跑成绩,800米跑成绩,跳高,跳远}。接着我们分析问题的时候,觉得这些指标太多了(似乎并不多,但是假设你觉得多),对我们分析问题比较麻烦,这个时候我需要对指标进行处理,这个时候我就进行主成分分析。很多同学对这有个误解,认为主成分分析出来的是从原来指标集中直接挑出来的,绝对不是这样!大家请看下图:
![](https://i-blog.csdnimg.cn/blog_migrate/6e26c010331c5359b43313124e367aac.png)
在这个图中,左边是我们刚刚设定的变量,这里我用1-6来表示,右边是我们求出来的主成分(PC),为什么会有6个主成分呢?实际上我们用matlab软件求出来的就是多少个指标对应着多少个主成分。有人会问,不是说好的降维吗,为啥还是六个主成分呢?别急,我们先保留这个问题。
大家还需要注意的一点就是,每一个主成分和原来的指标都是存在一个线性关系的,我这里图中只画出来PC1,其他的PC也是这样的,也就说,我们可以把PC1用下面公式表示:
P
C
1
=
a
X
1
+
b
X
2
+
.
.
.
+
f
X
6
P{C_1} = a{X_1} + b{X_2} + ... + f{X_6}
PC1=aX1+bX2+...+fX6
什么意思呢?也就是说,每一个主成分都和原来的变量有关系,并不是说主成分就是从原来变量中直接挑出来的,但是他是有一定意义的。比如说,假设这个主成分PC1在X2,X3,X4(对应着肺活量,400米跑,800米跑)前面的系数比较大且为正值,那我们可以分析得到这个主成分主要侧重于人的耐力方面(肺活量,400米,800米这三个指标本身就体现了一个人的耐力);同理如果某个主成分PC2,在X5,X6,对着的前面系数比较大,说明我们分析得到的这个主成分是衡量人的爆发力方面。但是注意,并不是每一个主成分都可以分析出它们的意义。
我们可以看到前面两个主成分是彼此独立的,衡量人体育成绩不同方面的,这也是符合我们寻找主成分的定义的。接下来我们看一下,主成分分析的具体原理又是怎么实现的呢?
2 基本思想
在展开PCA数学原理之前,我们先整体把握一下其思想:
PCA算法是一种降维算法,降维就不可避免的涉及到数据的损失,就好比在三维世界中你看到的一个球,在二维世界中你只能看到一个面。怎么让我们的数据经过降维之后,尽可能的保留原始信息,这是我们需要解决的。
我们期望找到主成分彼此之间尽可能的独立,就比如前面例题中找到的主成分,分别来侧重说明某一方面,相互独立。如果不是相互独立,说明我们找的主成分不具有代表性,是不符合要求的。
总结一下:
- 降维后数据尽可能损失小
- 得到主成分之间尽可能独立
以上就是我们PCA算法的基本思想,想深入继续理解的同学可以继续往下看,不想深入了解的,直接想调用代码的可以直接略过下小节。
3 数学原理
在这里,我们不多赘述复杂的过程与推导,我们只关心理解算法,了解算法,知道算法是怎么来的,怎么解决的,有大体的了解就ok了。
我们首先有数据集如下图所示:
![](https://i-blog.csdnimg.cn/blog_migrate/33093f6137b99c13577c9db0419dd806.png)
这里每一行代表一个特征(N),每一列代表一个样本(M),同时我们选择N个R维的正交基,如下图所示:
![](https://i-blog.csdnimg.cn/blog_migrate/4d8ba2cb5a15d96426efc107e4bc9e75.png)
接下来我们会有:
![](https://i-blog.csdnimg.cn/blog_migrate/bc8b91dfea80a5dfc4ce389df239bfd7.png)
![](https://i-blog.csdnimg.cn/blog_migrate/cee95db298c359598e2daf7ef8cf6b80.png)
大家对这一步不必过于纠结,这一步就是进行一个线性的变化,进行降维的过程,可以看到经过一个矩阵P进行处理后,我们得到的矩阵X’的维度发生了变化,每一行仍然代表是特征,但是此时的特征数变成R(R<N),样本的数目没有变,达到了我们降维的目的。同时对于变化后的每一个变量,有:
p
i
x
j
=
p
i
1
x
1
j
+
p
i
2
x
2
j
+
.
.
.
+
p
i
N
x
N
N
{p_i}{x_j} = {p_{i1}}{x_{1j}} + {p_{i2}}{x_{2j}} + ... + {p_{iN}}{x_{NN}}
pixj=pi1x1j+pi2x2j+...+piNxNN
是线性的关系,也就是为什么我们说PCA是一种线性降维。当然我们也有非线性降维,后续我们还会继续讨论这个问题。
我们回头再看我们的问题,我们可以找任何一个P矩阵对我们原始矩阵进行降维,但是怎么才能找到我们最好的,最优的P矩阵呢?此时问题转化为求P矩阵的最优问题,优化条件则是我们上文讲过的两个问题:一是降维后数据尽可能的损失小;二是降维得到的主成分尽可能的独立。
对于问题一,我们可以自然而然想到,损失尽可能的小,即让数据在降维后的点尽可能的分散,这比较好理解,比如原来三维空间中的点彼此分离,如果将它们投影到二维空间那么必定会有一些点相互重叠,造成数据的损失,而我们要做的即让经过P矩阵处理后,点尽可能的分散。大家注意,这里的点指的是样本(也就是每列,在这里强调是要和后面区分),我们采用的方法是用方差(刻画样本分散程度的统计量)来进行刻画:
var
(
x
j
)
=
1
M
∑
k
=
1
M
(
x
j
k
−
μ
)
2
\operatorname{var} ({x_j}) = \frac{1}{M}\sum\limits_{k = 1}^M {{{({x_{jk}} - \mu )}^2}}
var(xj)=M1k=1∑M(xjk−μ)2
对于问题二,主成分之间尽可能独立,我们采用的是协方差来刻画:
C
o
v
(
x
i
,
x
j
)
=
E
[
x
i
x
j
]
=
1
M
∑
k
=
1
M
x
i
k
x
j
k
Cov({x_i},{x_j}) = E[{x_i}{x_j}] = \frac{1}{M}\sum\limits_{k = 1}^M {{x_{ik}}{x_{jk}}}
Cov(xi,xj)=E[xixj]=M1k=1∑Mxikxjk
在这个公式中x和前面公式表述并不是一个意思,这里x代表我们的主成分,也就是降维后特征,我们还没有进行选择,但是要求特征两两之间尽可能的独立,用来刻画事务的某一方面。
那么如何去构造一个矩阵,能够包括以上的信息呢?首先我们对原始数据矩阵进行0均值处理得到矩阵A(方便计算,对我们来说无所谓,因为我们交给matlab计算,在这里处理方便观察解释公式):
![](https://i-blog.csdnimg.cn/blog_migrate/cf7d68acf3e0cf00f3bbe5aec7e052d5.png)
对矩阵A进行如下构造:
![](https://i-blog.csdnimg.cn/blog_migrate/a102b5a260903e3d4b7bc8b14f2fd9d1.png)
整理后可得结果:
![](https://i-blog.csdnimg.cn/blog_migrate/992cdb73b96e7c6fede89ac5bf730763.png)
到此为止,聪明的同学可以看到,对角线上的都是方差,除了对角线上的都是协方差,我们的优化目标是期望方差尽可能的大,协方差为零。到这里我们的主成分思想全部讲完了,它的精髓可能还需要你多看几遍才能明白,我们要求的是矩阵P,但是优化目标在后面的构造的矩阵中。
4 MATLAB 求解步骤
上面是我总结的主成分分析的步骤,实际上如果你仅仅是为了降维的话,第4步就已经结束了,接下来的5,6是进一步的扩展,是利用主成分分析进行一个综合评价的过程,它的思想有点类似于熵权法求综合评价值,后面我也会讨论这个专题。
最重要的环节!上代码!!
x=[1.0000 40.4000 24.7000 7.2000 6.1000 8.3000 8.7000 2.4420 20.0000;
2.0000 25.0000 12.7000 11.2000 11.0000 12.9000 20.2000 3.5420 9.1000;
3.0000 13.2000 3.3000 3.9000 4.3000 4.4000 5.5000 0.5780 3.6000;
4.0000 22.3000 6.7000 5.6000 3.7000 6.0000 7.4000 0.1760 7.3000;
5.0000 34.3000 11.8000 7.1000 7.1000 8.0000 8.9000 1.7260 27.5000;
6.0000 35.6000 12.5000 16.4000 16.7000 22.8000 29.3000 3.0170 26.6000;
7.0000 22.0000 7.8000 9.9000 10.2000 12.6000 17.6000 0.8470 10.6000;
8.0000 48.4000 13.4000 10.9000 9.9000 10.9000 13.9000 1.7720 17.8000;
9.0000 40.6000 19.1000 19.8000 19.0000 29.7000 39.6000 2.4490 35.8000;
10.0000 24.8000 8.0000 9.8000 8.9000 11.9000 16.2000 0.7890 13.7000;
11.0000 12.5000 9.7000 4.2000 4.2000 4.6000 6.5000 0.8740 3.9000;
12.0000 1.8000 0.6000 0.7000 0.7000 0.8000 1.1000 0.0560 1.0000;
13.0000 32.3000 13.9000 9.4000 8.3000 9.8000 13.3000 2.1260 17.1000;
14.0000 38.5000 9.1000 11.3000 9.5000 12.2000 16.4000 1.3270 11.6000;
15.0000 26.2000 10.1000 5.6000 15.6000 7.7000 30.1000 0.1260 25.9000];
g=zscore(x);
r=corrcoef(g);
[x,y,z]=pcacov(r);
num=3;%这个可以进行修改,根据计算得到的Z值进行判断
df=g*x(:,1:num);%行代表公司,所以全部选择
tf=df*z(1:num)/100;
[stf,ind]=sort(tf,'descend'); %stf是综合得分,ind是公司标签
stf=stf', ind=ind'
还剩一个小尾巴没给大家说,大家尝试运行代码后发现得到的主成分个数和你的变量数是一样的,在开头的时候给大家埋下了伏笔,怎么去选择主成分呢?
在pcacov这个函数返回三个值,x每一列代表每一个主成分PC与原始变量之间的系数,y返回的是特征值,z返回的是贡献率(其实是根据某一个特征值在所有特征值的占比得到的,CCA算法也有涉及),主成分都是按照贡献率大小排序的,我们选择累计贡献率大于85%的主成分即可,毕竟我们降维不可能避免数据损失,按照上面的例题,我们三个主成分贡献率达89%,也就说能够解释原始数据的89%,剩下的数据在降维过程丢失了(其实是人为舍去了)。
奈何本人实力有限,必有疏漏,欢迎批评指正!后续推出各种算法,欢迎大家订阅!