详解主成分分析PCA

2019-04-08 本文已影响433人霞客环肥

主成分分析（ Principal components analysis），简称PCA，是最主要的数据降维方法之一。本文从PCA的思想开始，一步一步推导PCA。

1.0 PCA的最大可分性的思想

对于 $X = \begin{bmatrix} x_1 \\ x_2 \\ ... \\ x_n \end{bmatrix}$ , $X \in R^n$ 。我们希望 $X$ 从 $n$ 维降到 $n^{'}$ 维，同时希望信息损失最少。比如，从 $n = 2$ 维降到 $n^{'} = 1$ ：

image.png

我们既可以降维到第一主成分轴，也可以降维到第二主成分轴。那么如何找到这这些主成分轴并且选择最优成分轴呢？

直观上，第一主成分轴优于第二主成分轴，即具有最大可分性。
下面解决一些基本概念。

2.0 基变换

欲获得原始数据新的表示空间，最简单的方法是对原始数据进行线性变换（基变换）：

$Y = PX$

其中 $X$ 是原始样本， $P$ 是基向量， $Y$ 是新表达。

数学表达：
$\begin{bmatrix} p_1 \\ p_2 \\ \vdots \\ p_R \end{bmatrix}_{R \times N} \begin{bmatrix} x_1 & x_2 & \cdots & x_M \end{bmatrix}_{N \times M} = \begin{bmatrix} p_1 x_1 & p_1 x_2 & \cdots & p_1 x_M \\ p_2 x_1 & p_2 x_2 & \cdots & p_2 x_M \\ \vdots & \vdots & \ddots & \vdots \\ p_R x_1 & p_R x_2 & \cdots & p_R x_M\end{bmatrix}_{R\times M}$

其中 $p_i$ 是行向量，表示第 $i$ 个基， $x_j$ 是一个列向量，表示第 $j$ 个原始数据记录.
当 $R < N$ 时即基的维度 < 数据维度时，可达到降维的目的。即：
$X \in R^{N \times M} \rightarrow Y \in R^{R \times M}$

以直角坐标系下的点(3,2)为例，欲将点(3,2)变换为新基上的坐标，就是用(3,2)与第一个基做内积运算，作为第一个新的坐标分量，然后用(3,2)与第二个基做内积运算，作为第二个新坐标的分量。

image.png
实际上，我们可以用矩阵相乘的形式简洁的表示这个变换：

可以稍微推广一下，如果我们有m个二维向量，只要将二维向量按列排成一个两行m列矩阵，然后用“基矩阵”乘以这个矩阵，就得到了所有这些向量在新基下的值。例如(1,1)，(2,2)，(3,3)，想变换到刚才那组基上，则可以这样表示：
$\begin{bmatrix}\frac{1}{\sqrt 2} & \frac{1}{\sqrt 2} \\ -\frac{1}{\sqrt 2} & \frac{1}{\sqrt 2} \end{bmatrix} \begin{bmatrix} 1 & 2 & 3 \\ 1 & 2 & 3\end{bmatrix} = \begin{bmatrix} 2\sqrt 2 & 4\sqrt2 & 6\sqrt2 \\ 0 & 0 & 0 \end{bmatrix}$

3.0 方差

回顾一下，我们的目的是希望在降维过程中损失最少，换言之，我们希望投影后的数据尽可能分散开。这种分散程度可以用方差来表达，方差越大，数据越分散。

定义方差 $Var$ ：对于单一随机变量 $a$ ，
$Var(a) = \frac{1}{m} \sum_{i = 1}^m (a_i - \mu)^2$
对数据做去中心化（方便后面操作）：
$Var(a) = \frac{1}{m} \sum_{i = 1}^m a_i ^2$

随机变量 $a$ 表达了 $a$ 的取值与其数学期望之间的偏离程度。若 $Var(a)$ 较小，意味着 $a$ 的取值主要集中在期望 $\mu$ 也就是 $E(a)$ 的附近，反之，若 $Var(a)$ 较大，意味着 $a$ 的取值比较分散。

为了避免过于抽象，我们以一个具体的例子展开。假设我们5个样本数据，分别是 $x_1 = \begin{bmatrix} 1 \\ 1\end{bmatrix}, x_2 = \begin{bmatrix} 1 \\ 3\end{bmatrix}, x_3 = \begin{bmatrix} 2 \\ 3\end{bmatrix}, x_4 = \begin{bmatrix} 4 \\ 4\end{bmatrix} ,x_5 = \begin{bmatrix} 2 \\ 4 \end{bmatrix}$ ，将它们表示成矩阵形式：
$X = \begin{bmatrix} 1 & 1 & 2 & 4 & 2 \\ 1 & 3 & 3 & 4 & 4 \end{bmatrix}$
为了后续处理方便，我们首先将每个字段内所有值都减去字段均值，其结果是将每个字段都变为均值为0.

我们看上面的数据，设第一个特征为 $a$ ，第二个特征为 $b$ , 此时某一个样本可以写作： $x_i = \begin{bmatrix} a \\ b \end{bmatrix}$
且特征 $a$ 的均值为2, 特征 $b$ 的均值为3，所以变换后：
$X = \begin{bmatrix} -1 & -1 & 0 & 2 & 0 \\ -2 & 0 & 0 & 1 & 1 \end{bmatrix}$

$Var(a ) = \frac{\sqrt 6} {5}$ $Var(b ) = \frac{\sqrt 6} {5}$

4.0 协方差

协方差（Covariance）在概率论和统计学中用于衡量两个变量的总体误差。

比如对于二维随机变量 $x_i = \begin{bmatrix} a \\ b \end{bmatrix}$ ，特征 $a,b$ 除了自身的数学期望和方差，还需要讨论 $a,b$ 之间互相关系的数学特征。

定义协方差 $Cov$ ：
$Cov(a, b) = \frac{1}{m}\sum_{i = 1}^m a_i b_i$

当 $Cov(a, b) = 0$ 时，变量 $a,b$ 完全独立，这也是我们希望达到的优化目标。

方差是协方差的一种特殊情况，即当两个变量是相同的情况:
$Cov(a, a) = Var(a)$

5.0 协方差矩阵

对于二维随机变量 $x_i = \begin{bmatrix} a \\ b \end{bmatrix}$ ,

定义协方差矩阵 $C$ :
$C = \begin{bmatrix} Var(a) & Cov(a, b) \\ Cov(b, a) &Var(b) \end{bmatrix}$

对于n维随机变量 $x_i = \begin{bmatrix} x_1 \\ x_2 \\ \vdots \ x_n \end{bmatrix}$ ,

$C = \begin{bmatrix} Var(x_1) & Cov(x_1, x_2) &\cdots & Cov(x_1, x_n)\\ Cov(x_2, x_1)& Var(x_2) & \cdots & Cov(x_1, x_n)\\ \vdots & \vdots & \ddots & \vdots \\ Cov(x_n, x_1) & Cov(x_n, x_2) & \cdots & Var(x_n)\\ \end{bmatrix}$

可见，协方差矩阵是 $n$ 行 $n$ 列的对称矩阵，主对角线上是方差，而协对角线上是协方差。

依然我们以一个具体的例子展开，还是这5个样本数据， $x_1 = \begin{bmatrix} 1 \\ 1\end{bmatrix}$ , $x_2 = \begin{bmatrix} 1 \\ 3\end{bmatrix}, x_3 = \begin{bmatrix} 2 \\ 3\end{bmatrix}, x_4 = \begin{bmatrix} 4 \\ 4\end{bmatrix} ,x_5 = \begin{bmatrix} 2 \\ 4 \end{bmatrix}$ ，将它们去中心化后表示成矩阵形式：
$X = \begin{bmatrix} -1 & -1 & 0 & 2 & 0 \\ -2 & 0 & 0 & 1 & 1 \end{bmatrix}$
那如果有 $m$ 个样本的话，
$X =\begin{bmatrix} a_1 & a_2 & \cdots &a_m \\ b_1 & b_2 & \cdots & b_m \end{bmatrix}$
对 $X$ 做一些变换，用 $X$ 乘以 $X$ 的转置，并乘上系数1/m：
$\frac{1}{m}XX^T = \frac{1}{m}\begin{bmatrix} a_1 & a_2 & \cdots &a_m \\ b_1 & b_2 & \cdots & b_m \end{bmatrix} \begin{bmatrix} a_1 & b_1 \\ a_2 & b_2 \\ \vdots & \vdots \\ a_m &b_m \end{bmatrix}$ $= \begin{bmatrix} \frac{1}{m} \sum_{i = 1}^m a_i ^2 & \frac{1}{m}\sum_{i = 1}^m a_i b_i \\ \frac{1}{m}\sum_{i = 1}^m a_i b_i& \frac{1}{m} \sum_{i = 1}^m b_i^2 \end{bmatrix}$

这不正是协方差矩阵嘛！

现在我们可以说：

设我们有m个n维数据记录，将其按列排成n乘m的矩阵X，设 $C = \frac{1}{m}XX^T$ ，则 $C$ 是一个对称矩阵，其对角线分别个各个特征的方差，而第i行j列和j行i列元素相同，表示i和j两个特征之间的协方差。

6.0 协方差矩阵对角化

回顾一下：

现在我们有 $m$ 个样本数据，每个样本有 $n$ 个特征，那么设这些原始数据为 $X$ ， $X$ 为 $n$ 行 $m$ 列的矩阵。
想要找到一个基 $P$ ，使 $Y_{r \times m} = P_{r \times n}X_{n \times m}$ ，其中 $r<n$ ，达到降维的目的。

设 $X$ 的协方差矩阵为 $C$ ， $Y$ 的协方差矩阵为 $D$ ，且 $Y = PX$ 。

我们的目的变为：对原始数据 $X$ 做PCA后，得到的 $Y$ 的协方差矩阵 $D$ 的各个方向方差最大，协方差为0。
那么 $C$ 与 $D$ 是什么关系呢？

$D = \frac{1}{m}YY^T$
$= \frac{1}{m} (PX)(PX)^T$
$= \frac{1}{m}PXX^TP^T$
$= \frac{1}{m}P(XX^T)P^T$
$= PCP^T$
$= P \begin{bmatrix} \frac{1}{m} \sum_{i = 1}^m a_i ^2 & \frac{1}{m}\sum_{i = 1}^m a_i b_i \\ \frac{1}{m}\sum_{i = 1}^m a_i b_i& \frac{1}{m} \sum_{i = 1}^m b_i^2 \end{bmatrix} P^T$

我们要找的 $P$ 不是别的，而是能让原始协方差矩阵对角化的 $P$ 。

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

现在所有焦点都聚焦在了协方差矩阵对角化问题上。

由上文知道，协方差矩阵 $C$ 是一个是对称矩阵，在线性代数上，实对称矩阵有一系列非常好的性质：

1）实对称矩阵不同特征值对应的特征向量必然正交。

2）设特征向量 $\lambda$ 重数为 $r$ ，则必然存在 $r$ 个线性无关的特征向量对应于 $\lambda$ ，因此可以将这 $r$ 个特征向量单位正交化。

由上面两条可知，一个 $n$ 行 $n$ 列的实对称矩阵一定可以找到 $n$ 个单位正交特征向量，设这 $n$ 个特征向量为 $e_1,e_2,⋯,e_n$ ，我们将其按列组成矩阵：
$E = \begin{bmatrix} e_1 & e_2 & \cdots \ e_n \end{bmatrix}$

则对协方差矩阵 $C$ 有如下结论：
$E^T C E = \Lambda = \begin{bmatrix} \lambda_1 \\ & \lambda_2 \\ &&\ddots \\ &&&\lambda_n\end{bmatrix}$

其中 $\Lambda$ 为对角矩阵，其对角元素为各特征向量对应的特征值（可能有重复）。

结合上面的公式：
$D = PCP^T$
其中， $D$ 为对角矩阵，我们可以得到：
$P = E^T$
$P$ 是协方差矩阵 $C$ 的特征向量单位化后按行排列出的矩阵，其中每一行都是 $C$ 的一个特征向量。如果设 $P$ 按照 $\Lambda$ 中特征值的从大到小，将特征向量从上到下排列，则用 $P$ 的前 $K$ 行组成的矩阵乘以原始数据矩阵 $X$ ，就得到了我们需要的降维后的数据矩阵 $Y$ 。

7.0 PCA算法

总结一下PCA的算法步骤：

设有 $m$ 条 $n$ 维数据。

1）将原始数据按列组成 $n$ 行 $m$ 列矩阵X

2）将 $X$ 的每一行（代表一个特征）进行零均值化，即减去这一行的均值

3）求出协方差矩阵 $C=\frac{1}{m}XX^𝖳$

4）求出协方差矩阵 $C$ 的特征值及对应的特征向量

5）将特征向量按对应特征值大小从上到下按行排列成矩阵，取前 $k$ 行组成矩阵 $P$

6） $Y=PX$ 即为降维到 $k$ 维后的数据

8.0 实例

这里以上文提到的：
$x_1 = \begin{bmatrix} 1 \\ 1\end{bmatrix}, x_2 = \begin{bmatrix} 1 \\ 3\end{bmatrix}, x_3 = \begin{bmatrix} 2 \\ 3\end{bmatrix}, x_4 = \begin{bmatrix} 4 \\ 4\end{bmatrix} ,x_5 = \begin{bmatrix} 2 \\ 4 \end{bmatrix}$ ，将它们表示成矩阵形式：
$X = \begin{bmatrix} 1 & 1 & 2 & 4 & 2 \\ 1 & 3 & 3 & 4 & 4 \end{bmatrix}$

我们用PCA方法将这组二维数据其降到一维。

为了后续处理方便，我们首先将每个特征内所有值都减去字段均值，其结果是将每个字段都变为均值为0.
$X = \begin{bmatrix} -1 & -1 & 0 & 2 & 0 \\ -2 & 0 & 0 & 1 & 1 \end{bmatrix}$
因为这个矩阵的每行已经是零均值，这里我们直接求协方差矩阵：
$C = \frac{1}{5} \begin{bmatrix} -1 & -1 & 0 & 2 & 0 \\ -2 & 0 & 0 & 1 & 1 \end{bmatrix}\begin{bmatrix} -1 & -2 \\ -1 & 0 \\ 0 & 0 \\ 2 & 1\\ 0 & 1 \end{bmatrix} = \begin{bmatrix} \frac{6}{5} & \frac{4}{5} \\ \frac{4}{5} & \frac{6}{5} \end{bmatrix}$
对于矩阵 $C$ :
$C = \begin{bmatrix} \frac{6}{5} & \frac{4}{5} \\ \frac{4}{5} & \frac{6}{5} \end{bmatrix}$
$\lambda$ 和 $v$ 分别是特征值和特征向量，
$\because Cv = \lambda v$ ，则：
$(C - \lambda I)v = 0$
为了使这个方程式有非零解，矩阵 $(C - \lambda I)$ 的行列式必须是0：
$det(C - \lambda I) = 0$
即：
$det(\begin{bmatrix} \frac{6}{5}-\lambda & \frac{4}{5} \\ \frac{4}{5} & \frac{6}{5}-\lambda \end{bmatrix}) = 0$
则：
$(\frac{6}{5}-\lambda) ^2 -\frac{16}{25} = 0$
分解得：
$(\lambda -2)(5\lambda -2) = 0$
找到2个特征值， $\lambda = 2$ , $\lambda = \frac{2}{5}$ ,

when $\lambda = 2$ :
$(C - \lambda I)v = 0$
即：
$\begin{bmatrix} -\frac{4}{5} & \frac{4}{5} \\ \frac{4}{5} & - \frac{4} {5} \end{bmatrix} \begin{bmatrix} v_1 \\ v_2 \\ \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\ \end{bmatrix}$
则：
$v_1 - v_2 = 0$
$v_1$ 和 $v_2$ 可以取任意值，我们取归一化的 $v_1$ 和 $v_2$ ，即： $v_1^2 + v_2^2 = 1$ ,
此时 $v_1 = \frac{\sqrt{2} } {2}$ 和 $v_2 = \frac{\sqrt{2} } {2}$
$v = \begin{bmatrix}\frac{\sqrt{2} } {2} \\ \sqrt{2} \over 2 \end{bmatrix}$

when $\lambda = \frac{2}{5}$ :
$(C - \lambda I)v = 0$
即：
$\begin{bmatrix} \frac{4}{5} & \frac{4}{5} \\ \frac{4}{5} & \frac{4} {5} \end{bmatrix} \begin{bmatrix} v_1 \\ v_2 \\ \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\ \end{bmatrix}$
则：
$v_1 + v_2 = 0$
$v_1$ 和 $v_2$ 可以取任意值，我们取归一化的 $v_1$ 和 $v_2$ ，即： $v_1^2 + v_2^2 = 1$
此时 $v_1 = \frac{\sqrt{2} } {2}$ 和 $v_2 = -\frac{\sqrt{2} } {2}$
$v = \begin{bmatrix} -\frac{\sqrt{2} } {2} \\ \sqrt{2} \over 2 \end{bmatrix}$

所以：
$P = \begin{bmatrix} \sqrt{2} \over 2 & \sqrt{2} \over 2 \\ -\sqrt{2} \over 2 & \sqrt{2} \over 2 \\ \end{bmatrix}$

可以验证协方差矩阵C的对角化：
$PCP^T = \begin{bmatrix} \sqrt{2} \over 2 & \sqrt{2} \over 2 \\ -\sqrt{2} \over 2 & \sqrt{2} \over 2 \\ \end{bmatrix} \begin{bmatrix} \frac{6}{5} & \frac{4}{5} \\ \frac{4}{5} & \frac{6}{5} \end{bmatrix} \begin{bmatrix} \sqrt{2} \over 2 & -\sqrt{2} \over 2 \\ \sqrt{2} \over 2 & \sqrt{2} \over 2 \\ \end{bmatrix} = \begin{bmatrix} 2 & 0 \\ 0 & \frac{2}{5}\end{bmatrix}$
最后我们用 $P$ 的第一行乘以数据矩阵，就得到了降维后的表示：
$Y = PX = \begin{bmatrix} \sqrt{2} \over 2 & \sqrt{2} \over 2 \end{bmatrix} \begin{bmatrix} -1 & -1 & 0 & 2 & 0 \\ -2 & 0 & 0 & 1 & 1 \end{bmatrix} = \begin{bmatrix} -\frac{3}{2} \sqrt 2 & -\frac{\sqrt 2} {2} & 0 & \frac{3}{2} \sqrt 2 & \frac{\sqrt 2} {2} \end{bmatrix}$

降维投影结果如下图：