PCA数学推导及原理（转）

原文：https://zhuanlan.zhihu.com/p/26951643

在多元统计分析中，主成分分析（Principal components analysis，PCA）是一种分析、简化数据集的技术。主成分分析经常用于减少数据集的维数，同时保持数据集中的对方差贡献最大的特征。这是通过保留低阶主成分，忽略高阶主成分做到的。这样低阶成分往往能够保留住数据的最重要方面。

PCA在机器学习中经常被用到，是数据预处理的重要步骤。它主要基于以下考虑：

高维特征中很多特征之间存在相关性，含有冗余信息
相比于低维数据，高维数据计算更复杂

PCA的数学原理

如下图，平面上有很多二维空间的特征点，如果想对这些特征点做特征降维（变为一维），应该怎么做呢？大家应该都知道需要进行投影，但还要考虑在哪个方向上进行投影，例如图中需要投影到长箭头方向即可，但考虑为什么不在短箭头上投影？

PCA本质上是一个有损的特征压缩过程，但是我们期望损失的精度尽可能地少，也就是希望压缩的过程中保留最多的原始信息。要达到这种目的，我们希望降维（投影）后的数据点尽可能地分散。如图，相比于长箭头，如果在短箭头上进行投影，那么重叠的点会更多，也就意味着信息丢失的更多，因而选择长箭头方向。

基于这种思想，我们希望投影后的数据点尽可能地分散。而这种分散程度在数学上可以利用方差来表示。设降维后的特征为 $A$ ，也就是希望 $var(A)=\frac{1}{m}\sum_i^m(a_i-\mu_a)^2$ 尽可能地大（ $a_i$ 为特征 $A$ 中的值， $\mu_a$ 为均值），而由于在PCA降维前，一般已经做了特征零均值化处理，为了方便，记 $var(A)=\frac{1}{m}\sum_i^ma_i^2$ 。

同样，为了减少特征的冗余信息，我们希望降维后的各特征之间互不相关。而不相关性可以用协方差来衡量。设降维后的两个特征为 $A$ 、 $B$ ，则希望 $Cov(A,B)=\frac{1}{m}\sum_i^ma_ib_i$ 为0。

现假设我们的数据为

$\begin{align} X = \left[ \begin{matrix} a_1&b_1\\ a_2&b_2\\ \vdots&\vdots\\ a_m&b_m \end{matrix} \right] \end{align}$

构造出协方差矩阵，并乘以系数 $\frac{1}{m}$ ，则

$\begin{align} \frac{1}{m}X^TX = \left[ \begin{matrix} \frac{1}{m}\sum_i^ma_i^2&\frac{1}{m}\sum_i^ma_ib_i\\ \frac{1}{m}\sum_i^ma_ib_i&\frac{1}{m}\sum_i^mb_i^2\\ \end{matrix} \right] \end{align}$

可以看出 $\frac{1}{m}X^TX$ 的对角线元素就是各特征的方差，其他各位置的元素就是各特征之间的协方差。因而只需要降维后的数据协方差矩阵满足对角矩阵的条件即可。

设 $Y$ 为原始数据 $X$ 做完PCA降维后的数据，满足 $Y=XP$ （矩阵乘法相当于映射，若 $P$ 为的列向量为基向量，那么就相当于映射到新的坐标系）， $Y_c$ ， $X_c$ 分别为对应的协方差矩阵，那么

$\begin{align} &Y_c=\frac{1}{m}Y^TY\\ &=\frac{1}{m}(XP)^TXP\\ &=\frac{1}{m}P^TX^TXP\\ &=P^T(\frac{1}{m}X^TX)P\\ &=P^TX_cP \end{align}$

因而，我们只需要计算出 $P$ ，使 $Y_c=P^TX_cP$ 满足对角矩阵的条件即可。而 $X_c$ 为实对称矩阵，我们只需要对它做矩阵对角化即可。

PCA的原理基本就是这样，还是挺简单的。

PCA的推导证明

PCA的构建：PCA需要构建一个编码器 $f$ ，由输入 $x\in R^n$ 得到一个最优编码 $c\in R^l$ （若 $l<n$ ，则做了降维编码）；同时有一个解码器 $g$ ，解码后的输出 $g(c)$ 尽可能地与 $x$ 相近。

PCA由我们所选择的解码器决定，在数学上，我们使用矩阵将 $c$ 映射回 $R^n$ ，即 $g(c)=Dc$ ，其中 $D\in R^{n\times l}$ 定义解码的矩阵。

为了限制PCA的唯一性，我们限制 $D$ 中所有列向量彼此正交且均有单位范数（否则 $D$ 、 $c$ 同比例增加、减少会产生无数个解）。

在数学上，为了满足PCA构建中的条件，我们利用 $L_2$ 范数来衡量 $g(c)$ 与 $x$ 的相近程度。即 $c^*=argmin_c||x-g(c)||_2$ ，也就是 $c^*=argmin_c||x-g(c)||_2^2$

该最小化函数可以简化为

$\begin{align} &(x-g(c))^T(x-g(c))\\ &=x^Tx-x^Tg(c)-g(c)^Tx+g(c)^Tg(c)\\ &=x^Tx-2x^Tg(c)+g(c)^Tg(c) \end{align}$

因而，优化目标变为 $c^*=argmin_c-2x^Tg(c)+g(c)^Tg(c)$ ，再带入 $g(c)=Dc$ ，

$\begin{align} &c^*=argmin_c-2x^TDc+c^TD^TDc\\ &=argmin_c-2x^TDc+c^Tc(D^TD=I_l) \end{align}$

再求偏导

$\begin{align} &\nabla_c(-2x^TDc+c^Tc)=0\\ &-2D^Tx+2c=0\\ &c=D^Tx \end{align}$

于是我们可以得到编码函数 $f(x)=D^Tx$ ，PCA的重构操作也就可以定义为 $r(x)=g(c)=g(f(x))=DD^Tx$ 。问题接着就转化成如何求编码矩阵 $D$ 。由于PCA算法是在整个数据矩阵上进行编码，因而也要用 $D$ 对所有数据进行解码，所以需要最小化所有维上的误差矩阵的Frobenius范数：

$D^*=argmin_D\sqrt{\sum_{i,j}(x^{(i)}-r(x^{(i)}))_j^2}~~subject~to~D^TD=I_l$

我们考虑 $l=1$ 的情况，则 $D$ 是一个单一向量 $d$ ，则上式可以转化为

$d^*=argmin_d\sum_{i}||(x^{(i)}-dd^Tx^{(i)}||_2^2~~subject~to~||d||_2=1$

而 $d^Tx^{(i)}$ 为标量，转置与自身相等，上式通常写作

$d^*=argmin_d\sum_{i}||(x^{(i)}-x^{(i)T}dd||_2^2~~subject~to~||d||_2=1$

再将每一个输入点叠加起来，我们得到

$d^*=argmin_d\sum_{i}||X-X^{T}dd||_F^2~~subject~to~d^Td=1$

Frobenius范数简化成（考虑约束条件 $d^Td=1$ ）

$\begin{align} &argmin_d\sum_{i}||X-X^{T}dd||_F^2\\ &=argmin_dTr((X-X^{T}dd)^T(X-X^{T}dd))\\ &=argmin_dTr(X^TX-X^TXd^T-dd^TX^TX+dd^TX^TXdd^T)\\ &=argmin_d-Tr(X^TXd^T)+Tr(dd^TX^TX)+Tr(dd^TX^TXdd^T)\\ &=argmin_d-2Tr(X^TXdd^T)+Tr(dd^TX^TXdd^T)\\ &=argmin_d-2Tr(X^TXdd^T)+Tr(X^TXdd^Tdd^T)\\ &=argmin_d-Tr(X^TXdd^T)\\ &=argmax_dTr(X^TXdd^T)\\ &=argmax_dTr(d^TX^TXd)~~subject~to~d^Td=1 \end{align}$

最后的优化目标可以利用 $\frac{\partial Tr(ABA^TC)}{\partial A}=CAB+C^TAB^T$ 以及拉格朗日乘数法来求解，可得最优的 $d$ 是 $X^TX$ 的最大特征值对应的特征向量。

上面的推导特定于 $l=1$ 的情况，仅有一个主成分。一般来说，矩阵 $D$ 由 $X^TX$ 的前 $l$ 个最大的特征值对应的特征向量组成（利用归纳法，将 $D_{l+1}$ 表示为 $D_l$ 的函数即可，需要两个辅助矩阵：单位对角矩阵 $R^{(l+1)\times l}$ 以及 $(0,0\cdots,0,1)^T\in R^{l+1}$ ，省去证明过程）。

参考

主成分分析
CodingLabs - PCA的数学原理
《Deep Learning》 Ian Goodfellow et al.