图文并茂的PCA教程

声明：

1. 参考：PCA数学原理、*

PCA——主成分分析

简介

PCA全称Principal Component Analysis，即主成分分析，是一种常用的数据降维方法。它可以通过线性变换将原始数据变换为一组各维度线性无关的表示，以此来提取数据的主要线性分量。

z=wTx $$z=w^{T}x$$

  其中，z为低维矩阵，x为高维矩阵，w为两者之间的映射关系。

划重点：

1. 线性无关是因为构建的特征轴是正交的
2. 主要线性分量，因为只选取了方差够大的特征（或者说是主成分）

相应的，PCA解释方差并对离群点很敏感：少量原远离中心的点对方差有很大的影响，从而也对特征向量有很大的影响。

线性变换

一个矩阵与一个列向量A相乘，等到一个新的列向量B，则称该矩阵为列向量A到列向量B的线性变换。

我们希望投影后投影值尽可能分散，而这种分散程度，可以用数学上的方差来表述。

Var(a)=1m∑i=1m(ai−μ)2 $$Var(a)=\frac{1}{m}\underset{i=1}{\overset{m}{\sum }}(a_i-\mu {)}^2$$

即寻找一个一维基，使得所有数据变换为这个基上的坐标表示后，方差值最大。

解释：方差越大，说明数据越分散。通常认为，数据的某个特征维度上数据越分散，该特征越重要。

对于更高维度，还有一个问题需要解决，考虑三维降到二维问题。与之前相同，首先我们希望找到一个方向使得投影后方差最大，这样就完成了第一个方向的选择，继而我们选择第二个投影方向。如果我们还是单纯只选择方差最大的方向，很明显，这个方向与第一个方向应该是“几乎重合在一起”，显然这样的维度是没有用的，因此，应该有其他约束条件——就是正交

解释：从直观上说，让两个字段尽可能表示更多的原始信息，我们是不希望它们之间存在（线性）相关性的，因为相关性意味着两个字段不是完全独立，必然存在重复表示的信息。

数学上可以用两个字段的协方差表示其相关性：

Cov(a,b)=1m∑i=1m(ai−μa)(bi−μb) $$Cov(a,b)=\frac{1}{m}\underset{i=1}{\overset{m}{\sum }}(a_i-{\mu }_a)(b_i-{\mu }_b)$$

当协方差为0时，表示两个字段线性不相关。

总结一下，PCA的优化目标是：
将一组N维向量降为K维（K大于0，小于N），其目标是选择K个单位正交基，使得原始数据变换到这组基上后，各字段两两间协方差为0，而字段的方差则尽可能大。

---

所以现在的重点是方差和协方差

协方差

在统计学上，协方差用来刻画两个随机变量之间的相关性，反映的是变量之间的二阶统计特性。考虑两个随机变量 Xi $X_i$和 Xj $X_j$，它们的协方差定义为

cov(Xi,Xj)=E[(Xi−E(Xi))(Xj−E(Xj))] $$cov(X_i,X_j)=E[(X_i-E(X_i))(X_j-E(X_j))]$$

tips：独立，不相关与协方差为零三者的关系
只讨论离散型随机变量的情形。
独立：随机变量 ξ,η $\xi ,\eta$独立是指对于任意的常数a,b，都有

P(ξ=a,η=b)=P(ξ=a)⋅P(η=b) $$P(\xi =a,\eta =b)=P(\xi =a)\cdot P(\eta =b)$$

.
相关性，相关系数

ρξη=cov(ξ,η)var(ξ)−−−−−√var(η)−−−−−√ $${\rho }_{\xi \eta }=\frac{cov(\xi ,\eta )}{\sqrt{var(\xi )}\sqrt{var(\eta )}}$$

相关系数其实是“线性相关系数”
相关系数和协方差在描述相关性方面是等价的，但独立与相关性的关系是：

独立=>不相关

协方差矩阵：
假设有两个变量a和b，特征维度为m，那么构成的数据集矩阵为：

X=(a1b1a2b2......ambm) $$X=\left(\begin{array}{cccc}{a}_1& a_2& ...& a_m\\ b_1& b_2& ...& b_m\end{array}\right)$$

再假设它们的均值都是0，对于有两个均值为0的m维向量组成的向量组，

1mXXT=(1m∑mi=1a2i1m∑mi=1aibi1m∑mi=1aibi1m∑mi=1b2i) $$\frac{1}{m}X{X}^T=\left(\begin{array}{cc}\frac{1}{m}\underset{i=1}{\overset{m}{\sum }}{a}_i^2& \frac{1}{m}\underset{i=1}{\overset{m}{\sum }}{a}_i{b}_i\\ \frac{1}{m}\underset{i=1}{\overset{m}{\sum }}{a}_i{b}_i& \frac{1}{m}\underset{i=1}{\overset{m}{\sum }}{b}_i^2\end{array}\right)$$

可以发现对角线上的元素是两个字段的方差，其他元素是两个字段的协方差，两者都被统一到了一个矩阵——协方差矩阵中。

回顾一下前面所说的PCA算法的目标：方差max，协方差min！！

要达到PCA降维目的，等价于将协方差矩阵对角化：即除对角线外的其他元素化为0，并且在对角线上将元素按大小从上到下排列，这样我们就达到了优化目的。

设原始数据矩阵X对应的协方差矩阵为C，而P是一组基按行组成的矩阵，设Y=PX，则Y为X对P做基变换后的数据。设Y的协方差矩阵为D，我们推导一下D与C的关系：

D=1mYYT=1m(PX)(PX)T=1mPXXTPT=P(1mXXT)PT=PCPT $$D=\frac{1}{m}Y{Y}^T=\frac{1}{m}(PX)(PX{)}^T=\frac{1}{m}PX{X}^T{P}^T=P(\frac{1}{m}X{X}^T)P^T=PC{P}^T$$

解释：想让原始数据集X =>pca成数据集Y，使得Y的协方差矩阵是个对角矩阵。
有上述推导可得，若有矩阵P能使X的协方差矩阵对角化，则P就是我们要找的PCA变换。

优化目标变成了寻找一个矩阵 P $P$，满足 PCPT $PC{P}^T$是一个对角矩阵，并且对角元素按从大到小依次排列，那么 P $P$的前 K $K$行就是要寻找的基，用 P $P$的前 K $K$行组成的矩阵乘以 X $X$就使得 X $X$从 N $N$维降到了 K $K$维并满足上述优化条件。

矩阵对角化

首先，原始数据矩阵X的协方差矩阵C是一个实对称矩阵，它有特殊的数学性质：

1. 实对称矩阵不同特征值对应的特征向量必然正交。
2. 设特征值 λ $\lambda$重数为r，则必然存在r个线性无关的特征向量对应于 λ $\lambda$，因此可以将这r个特征向量单位正交化。

一个n行n列的实对称矩阵一定可以找到n个单位正交特征向量，设这n个特征向量为 e1,e2,...,en $e_1,e_2,...,e_n$，我们将其按列组成矩阵：

E=(e1 e2 ... en) $$E=(e_1{e}_2...e_n)$$

则对协方差矩阵C有如下结论：

ETCE=Λ=⎛⎝⎜⎜⎜λ1λ2...λn⎞⎠⎟⎟⎟ $$E^{T}CE=\mathrm{\Lambda }=\left(\begin{array}{cccc}{\lambda }_1& & & \\ & {\lambda }_2& & \\ & & ...& \\ & & & {\lambda }_n\end{array}\right)$$

这里不懂的朋友可以查阅线性代数相关书籍。

P=ET $$P=E^T$$

P是协方差矩阵的特征向量单位化后按行排列出的矩阵，其中每一行都是C的一个特征向量。如果设P按照中特征值的从大到小，将特征向量从上到下排列，则用P的前K行组成的矩阵乘以原始数据矩阵X，就得到了我们需要的降维后的数据矩阵Y。

在解释一下，特征值 λ $\lambda$为什么要从大到小排列，为什么要选较大的 λ $\lambda$？？？
因为我们协方差矩阵的对角线元素是方差，我们想要找方差交大的特征维度，所以要选择较大的对角线元素。
而对角矩阵 Λ $\mathrm{\Lambda }$虽然是C经过线性变化后的矩阵，但它在对角线上元素的大小关系没变，特征维度 i $i$对应的特征值 λi ${\lambda }_i$越大，该维度上数据的方差越大。

另一种解释思路

该思路基于拉格朗日问题的求解方法。
回到一开始， z=wTx $z=w^{T}x$。其中，最主要的成分是这样的 w1 $w_1$，样本投影到 w1 $w_1$上之后最分散，使得样本点之间的差别变得最明显。为了得到唯一解且是该方向成为最重要的因素，我们要求 ||w1||=1 $||w_1||=1$. 如果 z1=w1Tx $z_1={w_1}^{T}x$且 Cov(x)=∑ $Cov(x)=\sum$，则

Var(z1)=E[(wTx−wTμ)2]=w1t∑w1 $$Var(z_1)=E[(w^{T}x-w^T\mu {)}^2]={w_1}^t\sum w_1$$

寻找 w1 $w_1$ ，使得 w1 $w_1$ 在约束下最大化。将这写成拉格朗日问题，则有：

maxw1w1T∑w1−α(wT1w1−1) $$\underset{w_1}{max}{w_1}^T\sum w_1-\alpha (w_1^T{w}_1-1)$$

关于 w1 $w_1$ 求导并让它等于0，有：

2∑w1−2αw1=0 $$2\sum w_1-2\alpha w_1=0$$

因此， ∑w1=αw1 $\sum w_1=\alpha w_1$
如果 w1 $w_1$ 是协方差矩阵 ∑ $\sum$ 的特征向量， a是对应的特征值，则上式成立。因为我们想最大化

Var(z1)=wT1∑w1=αwT1w1=α $$Var(z_1)=w_1^T\sum w_1=\alpha w_1^T{w}_1=\alpha$$

所以为了方差最大，我们选择具有最大特征值的特征向量。因此，主成分是输入样本的协方差矩阵的具有最大特征值 λ1=α ${\lambda }_1=\alpha$ 的特征向量。
第二个主成分 w2 $w_2$ 也应该是最大化方差，具有单位长度，并且与 w1 $w_1$ 正交。后一个要求是使得投影后 z2=wT2x $z_2=w_2^{T}x$ 与 z1 $z_1$ 不相关。对于第二个主成分，有

maxw2wT2∑w2−α(wT2w2−1)−β(wT2w1−0) $$\underset{w_2}{max}{w}_2^T\sum w_2-\alpha (w_2^T{w}_2-1)-\beta (w_2^T{w}_1-0)$$

最后，该式简化为 ∑w2=αw2 $\sum w_2=\alpha w_2$ ，这表明 w2 $w_2$ 应该是 ∑ $\sum$ 的具有第二大特征值 λ2=α ${\lambda }_2=\alpha$ 的特征向量。类次的，我们可以证明其他维被具有递减特征值的特征向量给出。

算法及实例

PCA算法

总结一下PCA的算法步骤：
设有n条m维数据。

1. 将原始数据按列组成m行n列矩阵X
2. 将X的每一行(代表一个属性字段）进行零均值化
3. 求出协方差矩阵 C=1mXXT $C=\frac{1}{m}X{X}^T$
4. 求出协方差矩阵的特征值及对应的特征向量
5. 将特征相浪按对应特征值大小从上到下按行排列成矩阵，取前k行组成矩阵P
6. Y=PX $Y=PX$即为降维到k维后的数据

关于PCA的python实现代码可以参考这里，不过ipynb文件可能在github上刷不出来，建议下载下来用jupyter notebook打开。

实例

原始数据集矩阵X：

(1113234424) $$\left(\begin{array}{ccccc}1& 1& 2& 4& 2\\ 1& 3& 3& 4& 4\end{array}\right)$$

求均值后：

(−1−1−10002101) $$\left(\begin{array}{ccccc}-1& -1& 0& 2& 0\\ -1& 0& 0& 1& 1\end{array}\right)$$

再求协方差矩阵

C=15(−1−1−10002101)⋅⎛⎝⎜⎜⎜⎜⎜⎜−1−1020−10011⎞⎠⎟⎟⎟⎟⎟⎟=(65454565) $$C=\frac{1}{5}\left(\begin{array}{ccccc}-1& -1& 0& 2& 0\\ -1& 0& 0& 1& 1\end{array}\right)\cdot \left(\begin{array}{cc}-1& -1\\ -1& 0\\ 0& 0\\ 2& 1\\ 0& 1\end{array}\right)=\left(\begin{array}{cc}\frac{6}{5}& \frac{4}{5}\\ \frac{4}{5}& \frac{6}{5}\end{array}\right)$$

特征值：

λ1=2,λ2=25 $${\lambda }_1=2,{\lambda }_2=\frac{2}{5}$$

对应的特征向量：

c1⎛⎝12√12√⎞⎠,c1⎛⎝−12√12√⎞⎠ $$c1\left(\begin{array}{c}\frac{1}{\sqrt{2}}\\ \frac{1}{\sqrt{2}}\end{array}\right),c1\left(\begin{array}{c}-\frac{1}{\sqrt{2}}\\ \frac{1}{\sqrt{2}}\end{array}\right)$$

标准化（其实不标准化也一样，只是稍显不专业）

P=⎛⎝12√−12√12√12√⎞⎠ $$P=\left(\begin{array}{cc}\frac{1}{\sqrt{2}}& \frac{1}{\sqrt{2}}\\ -\frac{1}{\sqrt{2}}& \frac{1}{\sqrt{2}}\end{array}\right)$$

选择较大特征值对应的特征向量：

(12√12√) $$\left(\begin{array}{cc}\frac{1}{\sqrt{2}}& \frac{1}{\sqrt{2}}\end{array}\right)$$

执行PCA变换：Y=PX，得到的Y就是PCA降维后的值数据集矩阵：

Y=(12√12√)⋅(−1−1−10002101)=(−32√−12√032√−12√) $$Y=\left(\begin{array}{cc}\frac{1}{\sqrt{2}}& \frac{1}{\sqrt{2}}\end{array}\right)\cdot \left(\begin{array}{ccccc}-1& -1& 0& 2& 0\\ -1& 0& 0& 1& 1\end{array}\right)=\left(\begin{array}{ccccc}-\frac{3}{\sqrt{2}}& -\frac{1}{\sqrt{2}}& 0& \frac{3}{\sqrt{2}}& -\frac{1}{\sqrt{2}}\end{array}\right)$$

降维过程的示意图

进一步讨论

根据上面对PCA的数学原理的解释，我们可以了解到一些PCA的能力和限制。PCA本质上是将方差最大的方向作为主要特征，并且在各个正交方向上将数据“离相关”，也就是让它们在不同正交方向上没有相关性。

因此，PCA也存在一些限制，例如它可以很好的解除线性相关，但是对于高阶相关性就没有办法了，对于存在高阶相关性的数据，可以考虑Kernel PCA，通过Kernel函数将非线性相关转为线性相关，关于这点就不展开讨论了。另外，PCA假设数据各主特征是分布在正交方向上，如果在非正交方向上存在几个方差较大的方向，PCA的效果就大打折扣了。

最后需要说明的是，PCA是一种无参数技术，也就是说面对同样的数据，如果不考虑清洗，谁来做结果都一样，没有主观参数的介入，所以PCA便于通用实现，但是本身无法个性化的优化。