PCA（主成分分析）

引言

为什么会讲到PCA？

前一阵在做一些featureselection的工作，想到PCA可以用来降维（在某些方面来说Feature selection也要降维），因此学习了一下PCA才发现这是完全不同的过程（PCA中的降维（计算主成分的目的是将高维数据投影到较低维空间）与特征选择中的降维（选取特征）是两个概念）。PCA应该算是feature extraction的过程，feature extraction是通过变换把原始特征从最初的度量空间中变换到另一个维度减少了的度量空间（the original features in the measurement space areinitially transformed into a new dimension-reduced space via some specifiedtransformation）的过程，虽然通过featureextraction的方法是维度减少很多，但是转换后的feature几乎包含原始数据的全部信息。而feature selection则是选择一个subset，使得所选取的特征能够去掉原数据的冗余信息，是分类更好，有利于创建好的model等。扯远了，言归正传。

PCA（Principal Component Analysis）即主成分分析，是一种掌握事物主要矛盾的统计分析方法，是一种经典的用于数据降维的统计学方法。在本人所接触到的应用中，PCA被应用于面部识别，是面部识别领域经典的方法。下面便总结一下ＰＣＡ的算法与技术

PCA(统计学上)

奇异值分解

在讲解PCA之前，首先我们介绍一下什么是奇异值分解，对于一个方阵来说，可以对其进行特征值分解，而对于非方阵来说，我们没有办法进行特征值分解，因此引入了奇异值分解。

$$M = U \Sigma V^T$$

也就是说给定一个矩阵$M$， 他可以分解为由三个矩阵相乘。其中这三个矩阵满足下面的条件：

其中$U$的列是$MM^T$的单位正交特征向量 $U$的列是一组单位正交基

其中$\Sigma$是对角矩阵，对角线上的值为矩阵$M$的奇异值

其中$V$的列是$M^TM$的单位正交特征向量 $V$的列是一组单位正交基

PCA算法就用到了矩阵的这个性质

PCA算法的矩阵解释

给定我们一个机器学习问题，当我们提取完特征后，一般会将特征保存为一个矩阵进行存储，这个矩阵的行表示样本，列表示特征。假设我们的这个矩阵为$X$，根据上面介绍的奇异值分解，我们可以有：

$$X = W \Sigma V^T$$

其中$W$的列是$XX^T$的单位正交特征向量 $W$的列是一组正交基

由于$W$是一组单位正交基，那么我们可以将矩阵$X$（矩阵实际上是数据在一组默认的正交基空间（就是单位矩阵）中的映射）映射到$W$的单位正交基上，可以得到矩阵$X$在新的正交基上的表示$Y$，即:

$$Y^T=X^TW$$

将$X = W \Sigma V^T$带入有：

$$Y^T=X^TW = V \Sigma^T W^TW = V\Sigma^T$$

因为$W$是单位正交矩阵，所以$W^TW = WW^T = I$

因此有$Y=\Sigma V^T$， 我们的得到了将$X$变换为$Y$的方式。

其中$\Sigma$是对角矩阵，对角线上的值为矩阵$M$的奇异值，因此我们要将$X$映射到K为空间的话，那么只需要取前K个奇异值组成奇异值矩阵即可

$$Y_k =\Sigma_k V^T = I_{k*m}\Sigma V^T$$

那么还剩下一个问题，$V^T$是什么鬼。根据奇异值分解，我们知道：

矩阵$V$的列是$X^TX$的单位正交特征向量 $V$的列是一组单位正交基。

因此，给定了矩阵$X$，我们只要计算得到$X$的奇异值以及$XX^T$的特征向量，那么就可以计算得到矩阵的主成分矩阵$Y$

还有一个疑问点，$X^TX$是什么呢？ 这个就是$X$的协方差矩阵

$$Cov(x_i,x_j) = E((x_i - E(x_i))(x_j - E(x_j))) = E(XX^T) - \mu\mu^T$$

其中$\mu$是矩阵的列的均值向量。如果我们对矩阵首先按列归一化，那么就有$\mu = 0$

因此有$Cov(x_i,x_j) = E((x_i - E(x_i))(x_j - E(x_j))) = E(XX^T)$，即：

$$\Sigma = XX^T$$

以上就是PCA的数学表达。

没看懂？ 没看懂就对了。因为这里面省略了太多的细节，并且也隐藏了太多的数学上的意义。因此下面我们来更加详细的看下，为什么要这样做。

PCA算法的更加详细的数学解释

向量的内积

首先我们都学过向量的内积，那么向量的内积有什么意义呢？

我们知道向量的内积有如下的公式：

$$A \cdot B = |A||B|cos(\theta)$$

如果$|B| = 1$，即$B$为单位向量，那么实际上$A \cdot B$就是在$B$上的映射。

因此我们可以得到，向量的内积本质上就是一个向量在另一个向量上的映射投影。这个很关键，要理解。

对于我们的笛卡尔坐标系中的一个点$(x,y)$, 实际上可以表示为$x \cdot (1,0)^T + y \cdot (0,1)^T$ ，其中$(1,0)$和$(0,1)$就是一组正交基。

因此如果我们将$(x,y)$映射到一组新的基空间中，只需要用$(x,y)$与新的基向量做内积，就可以得到$(x,y)$在新的空间中的每个基向量上的坐标分量。

将其推广到多维向量中，并且利用矩阵来表示：

假设有M个N维向量，将其变换到R个N维向量表示的新空间中，那么可以将R个基组成矩阵$A_{N*R}$，向量按照行组成矩阵$B_{M*N}$,因此$BA$就可以得到向量在新空间中的向量表示。$BA$的第m行就是$B_{M*N}$中第m行变换的结果。

两个矩阵相乘，本质上就是左边矩阵的每一个行向量变换到右边矩阵的每个列向量为基所表示的空间中。这是也很关键，要理解。

那么现在有一个问题啊，怎么选择基呢。一个原则是“投影后的投影值尽可能的分散”

（因为如果投影后的投影值都挤在了一起，那么我们就没法分别不同的点了，也就没意义了，因此要尽可能的分散）

而方差是用来表示数据的离散程度的量，因此如果是一个二维的矩阵，如果要将其映射为一维，那么问题就可以转化为寻找一个一维基，使得所有的数据变换后，在这个基上的坐标，方差值最大。

推广到多维向量呢？那么问题就会变为：

将一组N维向量降到k维，目标就是选择k个单位正交基，使得变换后的不同特征的数据两两之间协方差为0，而特征内，方差尽可能的大。
（协方差为0保证了变换后的数据是线性不相关的，也就是特征之间没有冗余信息）

假设有2个样本的矩阵：

$$X = \begin{pmatrix} a_1 & b_1 \\ a_2 & b_2& \\ ... & ... \\ a_n & b_n \end{pmatrix}$$

那么$\frac{1}{m}XX^T = \begin{pmatrix} \frac{1}{m} \sum a_i^2 & \frac{1}{m} \sum a_i b_i \\ \frac{1}{m} \sum a_i b_i & \frac{1}{m} \sum b_i^2 \end{pmatrix}$

我们可以看到主对角线上是方差，其他的值为协方差。那么我们只要令主对角线上的值最大化，其他值为0，即可达到我们的目标。

假设原始数据为X， 其对应的协方差矩阵为C，P为基变换矩阵，那么有：

$$Y = PX$$

Y为变换后的矩阵，其协方差矩阵为D，那么我们的目的就是变换后的矩阵的协方差矩阵的主对角线的值最大化，其他值为0。因此有：

$$D = \frac{1}{m}YY^T = \frac{1}{m}(PX)(PX)^T = P(\frac{1}{m} XX^T )P^T = PCP^T$$

因此转换矩阵P就是能让原数据的协方差矩阵对角化的矩阵。

而P矩阵按特征值的大小排列的前k行主成的矩阵$P_k$乘以矩阵$X$就是X从n维降到k维的变换，$P_kX$就是前k个主成分组成的矩阵。

PCA的本质

PCA的本质就是将方差最大的方向作为主要特征，并且在各个正交方向上数据没有相关性。

主成分分析定理

关于主成分的求法有下面定理：

设$X$为$p$维随机变量，且$\Sigma = D(X)$存在，则$X$的第$i$个主成分$Y_i$与方差$D(Y_i)$分别为$Y_i=t_i^TX, D(Y_i)= \lambda_i$
其中$D(Y_2)$为的特征值，$t_i^T$为对应的单位特征向量。

看到这里如果你还没晕就证明你的统计学知识还是不错的，但当我们在machine learning或patternrecognition中用到PCA时，如果用上面的思想脉络总是感觉怪怪的，统计学家们一直在说做machine learning的人是在reinventingthe wheel，但我感觉做机器学习的人在一些方法的应用上更好理解，更能让人明白。

PCA（machine learning）

从上面我们可知，降维就是把高维的数据投影到低维的空间中，而线性方法的本质就在这里，所以上面我们用线性组合还表示数据。在machine leaning中，PCA的目的是找到在最小均方意义下最能代表原始数据的方法（在machine leaning中这个“最小均方”无处不在啊！）。

有$n$个$d$维的样本$X_1,X_2,...,X_n$，我们希望找到一个$X_0$来表示，即使平方误差准则函数$J_0 = \sum_{i=1}^n ||X_0-X_k||^2$达到最小，容易得到，这个值就是样本均值$X_0=m= \frac{1}{n} \sum_{i=1}^n X_k$。

通过把全部样本向通过样本均值的一条直线做投影，我们能够得到样本均值的直线上的一维向量

这条直线是$X = m+ae$， （e是通过样本均值的直线上的单位向量，a是一个实数）如果我们用$m+a_ke$来表示$X_k$，那么通过最小化均方误差准则函数，我们可以得到一组最优的$a_k$的集合，如下：

$$J_1(a_1,...,a_n,e) = \sum_{i=1}^n ||(m+a_ke)-X_k||^2 = \sum_{i=1}^n ||a_ke-(X_k-m)||^2$$ $$= \sum_{i=1}^n a_k^2||e||^2 -2\sum_{i=1}^n a_ke^T(x_k-m) + \sum_{i=1}^n ||m-X_k||^2$$

对$a_k$求导，可以得到$a_k =e^T(x_k-m)$

其实就是用一组数$a_k (k=1,2,...,n)$来表示原来的$n$个$d$维的样本$X_1,X_2,...,X_n$。

下面还有一个问题就是如何找到$e$的方向

引入样本离差矩阵$S = \sum_{i=1}^n (X_i-\overline{X})(X_i-\overline{X})^T$

将$a_k =e^T(x_k-m)$带入到$J_1(a_1,...,a_n,e) = \sum_{i=1}^n a_k^2||e||^2 -2\sum_{i=1}^n a_ke^T(x_k-m) + \sum_{i=1}^n ||m-X_k||^2$公式中，可以得到$J_1(e) = -e^TSe + \sum_{i=1}^n ||X_k - m||^2$, 显然是$J_1$最小的那个向量$e$，能使$-e^TSe$最大。在约束条件$||e||^2 = 1$下，利用拉格朗日乘子法，有：

$\mu = e^TSe - \lambda (e^Te-1)$,对$e$求偏导，有$\frac{\partial \mu}{\partial e} = 2Se-2\lambda e$,另其为0，有$Se =\lambda e$,因此，矩阵$S$最大的特征值对应的特征向量就是投影e的方向。

从几何上来讲，样本点$X_1,X_2,...,X_n$在$d$维空间上形成一个$d$维椭球形状的云团，样本离差矩阵$S$的特征向量就是这个云团的主轴，主成分分析就是提取云团散布最大的那些方向的方法，达到对特征空间进行降维的目的。

举个小例子直观说明一下PCA

对于一个训练集，100个对象模板，特征是10维，那么它可以建立一个10010的矩阵，作为样本。求这个样本的协方差矩阵，得到一个1010的协方差矩阵，然后求出这个协方差矩阵的特征值和特征向量，应该有10个特征值和特征向量，我们根据特征值的大小，取前四个特征值所对应的特征向量，构成一个104的矩阵，这个矩阵就是我们要求的特征矩阵，10010的样本矩阵乘以这个104的特征矩阵，就得到了一个1004的新的降维之后的样本矩阵，每个特征的维数下降了。

具体的应用（面部识别）

训练阶段

• 假设训练集有两百个样本，由灰度图组成每个样本大小M*N，则样本训练矩阵为$x =(x_1,x_2,...,x_{200})$,$x_i$ 为灰度图像沿纵轴的拉伸
• 计算平均脸$\Psi = \frac{1}{200} \sum_{i=1}^{200} x_i$
• 计算每一张脸与平均脸的差值$d_i = x_i -\Psi$
• 构造协方差矩阵$C = \frac{1}{200} \sum_{i=1}^{200}d_id_i^T$
• 求协方差矩阵的特征值与特征向量，构造特征脸空间
• 根据贡献率选取前p个最大特征值及其对应的特征向量,其中贡献率为$\varphi =\frac{\sum_{i=1}^{p} \lambda_i}{\sum_{i=1}^{p} \lambda_i} \geq a$，其中贡献率$a$ 取99%， 即训练样本的前p个特征向量集上的投影有99%的能量，特征脸空间为$w = (\mu_1,...,\mu_p)$
• 将每一幅人脸与平均脸的差值矢量投影到“特征脸”空间，即$\Omega =w^Td_i$

识别阶段

• 将待识别人脸与平均脸的差值投影到特征空间，得到其特征向量表示$\Omega^\Gamma =w^T( \Gamma -\Psi )$
• 定义阈值 $\theta =\frac{1}{2}max\left \{ ||\Omega_i - \Omega_j|| \right \},i,j=1,2,...,200$
• 采用欧式距离计算$\Omega^\Gamma$与每个人脸的距离$\varepsilon^2=||\Omega_i - \Omega^\Gamma ||^2 ,i=1,2,...,200$
• 为了区分人脸和非人脸，还需要计算原始图像$\Gamma$与由特征脸空间重建的图像$\Gamma_f$之间的距离$\varepsilon$，其中$\varepsilon^2 = ||\Gamma -\Gamma_f ||^2,\Gamma_f = w\Omega^\Gamma +\Psi$；
• 根据以下规则对人脸进行分类

  • 若$\varepsilon \geq \theta$则输入图像不是人脸图像
  • 若$\varepsilon \geq \theta$且$\forall i,\varepsilon_i \geq 0$则输入图像包含未知人脸
  • 若$\varepsilon < \theta$且$\forall i,\varepsilon_i < 0$则输入图像是库中第$k$个人的人脸