ML Note 3.4 - 数据降维算法 PCA / t-SNE

维数灾难的提出背景

1957年，数学家 Richard Bellman 在研究动态规划问题时首次提出了维数灾难（curse of dimensionality）的概念，用于描述高维空间中的一系列数学现象。例如，高维数据会带来指数级增长的计算复杂度，在 GPU 广泛应用之前严重限制了算法的可应用性。

维数灾难的内涵

维数灾难描述了与高维数据相关的一系列数学现象。特别地，在机器学习领域，维数灾难通常指数据集维数和数据集大小之间的指数级关系。机器学习算法将一组特征作为输入，这些特征可能具有许多不同类别的变量。以临床医学机器学习为例，特征可能包括患者人口统计数据（例如年龄、性别、体重）、临床特征（例如 CRP、心率、体温）或来自医学成像的信息（例如每个单独像素的灰度值）。一个数据集的维数，指的是数据集中属性或特征的数量。当数据集的维数达到一定规模，就被称为高维数据。一般来说，随着特征数量的增长，机器模型训练算法所需的样本数量呈指数级增长。由于高维数据而导致训练机器学习模型的困难被称为“维度灾难”。 下面介绍维数灾难的两个重要方面：数据稀疏性（data sparsity）和距离聚集性（distance concentration）。

数据稀疏性

经过训练的机器学习模型可以预测给定输入数据样本的标签。在训练模型时，使用部分数据用于训练模型，其余数据用于评估模型的泛化性能。泛化性能是指模型对于训练集以外的数据样本的预测准确率。泛化性能强的模型对测试数据的预测准确度应该非常接近其在训练数据上的准确度。

建立泛化模型的一种有效方法是捕获预测变量值和相应目标的不同可能组合。例如，一个数据集包含两个特征：性别和年龄段。理想情况下，我们应该捕获两个属性值的所有可能组合的目标。对于两种性别和四个年龄段，我们需要 8 个训练样本来覆盖所有可能组合。如果目标取决于第三个属性，比如体型，覆盖所有组合所需的训练样本数量会显着增加。这个例子表明，随着特征数量的增加，泛化模型所需的训练样本数量也显着增加。

实际应用中，可用的训练样本可能没有观察到所有属性组合的目标。 这是因为某些组合比其他组合出现的概率更高。训练样本没有捕获所有组合的性质在高维数据中被称为“数据稀疏性”或简称为“稀疏性”。数据稀疏性是维度灾难的一个方面。用稀疏数据训练模型可能会导致高方差或过拟合情况。这是因为在训练模型时，模型已经从频繁出现的属性组合中学习，并且可以准确地预测结果。当将不太频繁出现的组合实时输入模型时，它可能无法准确预测结果。

距离聚集性

维数灾难的另一个方面是距离聚集性。随着数据维数的增加，任意两个样本之间的距离会趋于相同。因此，高维空间中的距离度量可能无法反应样本之间的相似性。对于一些基于距离度量的机器学习模型，例如聚类算法，在低维空间中验证有效的距离度量可能不是用于高维空间。

数据降维

数据降维是解决维数灾难的有效方法之一，下面介绍数据降维中常用的 PCA 和 t-SNE 方法。数据降维的目的是将高维空间的一系列数据点 ${\mathbf{x}}_1,{\mathbf{x}}_2,\cdots ,{\mathbf{x}}_m$ 映射到低维空间 ${\mathbf{y}}_1,{\mathbf{y}}_2,\cdots ,{\mathbf{y}}_m$ ，同时最大程度保留数据集信息。一般来说，低维空间的维数远小于原高维空间的维数 $n$ 。

Principal Components Analysis

PCA 尝试寻找 $n$ 维原空间的一个 $M$ 维子空间，使样本集在其上的投影保留了尽量多的信息。在 PCA 之前，我们首先需要对数据集归一
$$x_i\leftarrow x_i-\bar{x}$$

使用设计矩阵的定义 $X=[x_1,x_2,\dots ,x_m{]}^T\in {\mathbb{R}}^{m\times n}$ 我们可以完成从原空间到子空间的映射
$$Y=XP$$

其中 $P=[u_1,u_2,\dots ,u_M]\in {\mathbb{R}}^{n\times M}$ 是由一组单位正交向量构成的变换矩阵。对于保留尽量多的原始信息这一目标，存在着下面两种理解

1. $Y$ 的样本方差尽可能大
2. 从 $Y$ 恢复到原空间后偏差尽可能小

下面将会证明这两种理解是等价的。

最大方差思想

首先考虑 $M=1$ 的情况，也就是样本集被映射为一维空间的一组数。因为变换矩阵 $P=[u_1]\in {\mathbb{R}}^{n\times 1}$ ，所以算法只需找到一个向量 $u=u_1$ 满足
$$\begin{array}{rl}\underset{u}{\max }& ||Xu|{|}^2/m\\ s.t.& u^{T}u=1\end{array}$$

样本集的协方差矩阵表示为
$$\Sigma =\frac{1}{m}{X}^{T}X$$

则优化目标可以化为 $u^T\Sigma u$ 。应用拉格朗日乘数法可以得到
$$\frac{\partial }{\partial u}(u^T\Sigma u+\lambda (u^{T}u-1))=\Sigma u+\lambda u=0$$

因此 $u$ 和 $\lambda$ 分别是 $\Sigma$ 的特征向量和对应特征值。重新考察目标函数
$$u^T\Sigma u=\lambda u^{T}u=\lambda$$

也就是说，样本在特征向量 $u$ 方向上的投影，所构成的样本集的方差即为其特征值。因此 $u$ 应该是 $\Sigma$ 的最大特征值对应的特征向量，称为 first principal components 。一般地对于 $M\ne 1$ 的情况，可以循环地执行以下步骤

1. 选取原空间中的第一主成分 $u$
2. 取原空间对于向量 $u$ 的补空间作为下一次循环的原空间

最小均方误差思想

假设存在原空间的某个坐标系 $[u_1,u_2,\dots ,u_n]$ ，取其中前 $M$ 个坐标轴 $[u_1,u_2,\dots ,u_M]$ 构成子空间坐标变换矩阵 $P$ 。不妨设某个样本 $x$ 在这样的坐标系下具有坐标 $(x_1,x_2,\dots ,x_n)$ ，则容易证明其对应的 $y=(x_1,x_2,\dots ,x_M)$ 。

现在我们希望找到一种方法使 $y$ 恢复到原空间，并与 $x$ 尽可能一致。但是由于缺少 $x_{M+1},x_{M+2},\dots ,x_n$ 的信息，我们只能用其均值 $0$ 来进行估计。也就是说重构得到的向量
$$\tilde{x}=(y_1,y_2,\dots ,y_M,0,\dots ,0)$$

这时就会引入重构误差
$$J=||x-\tilde{x}|{|}^2=\sum _{i=M+1}^n{x}_i^2=\sum _{i=M+1}^n{u}_i^{T}x{x}^T{u}_i$$

对于所有样本，均方误差的定义如下
$$J=\sum _{i=M+1}^n{u}_i^T\Sigma u_i$$

接下来的讨论与最大方差思想类似。不同点在于这里需要选出 $n-M$ 个特征向量使特征值最小。这些特征向量将作为补空间的一组基，而子空间的基则是剩余的 $M$ 个主成分。因此在结论上，最小均方误差思想和最大方差思想一致。

使用 SVD 处理高维数据

模型求解可以通过 SVD 完成，但是对于高维数据 $n\gg m$ 计算量会大一些。对于这种情况我们可以计算矩阵 $X{X}^T$ 的特征值 $\lambda$ 和特征向量 $u$ 。这个矩阵的维数是 $m$ 因此比较容易计算。算出来之后可以证明
$$X^{T}X(X^{T}u)=X^T\cdot \lambda u$$

因此 $X^{T}u$ 即为 $\Sigma$ 的特征向量，同时特征值 $\lambda$ 并没有改变。因为
$$\begin{array}{rcl}||A{v}_i|{|}^2& =& (A{v}_i{)}^{T}A{v}_i\\ & =& v_i^T{A}^{T}A{v}_i\\ & =& v_i^T({\lambda }_i{v}_i)\\ & =& {\lambda }_i\end{array}$$

所以有 $||X^{T}u|{|}^2=\lambda$ 。归一化可得特征向量
$$v=\frac{X^{T}u}{\sqrt{\lambda }}$$

相关系数矩阵

如果原空间中数据的各个特征量纲不同，可以使用数据的相关系数矩阵代替原数据。但是计算相关系数时，需要保证特征的标准差非零。当且仅当某个特征的取值为常数，其标准差才会为零。因此可以将标准差为零的特征提前删除，而不会造成误差。

结果分析

设 $F_j=X{u}_j$ 表示第 $j$ 个主成分的取值。不难证明
$$F_i^T{F}_j=u_i^T{X}^{T}X{u}_j=\{\begin{array}{ll}m{\lambda }_j& i=j\\ 0& i\ne j\end{array}$$

因为 $X$ 经过了归一化，因此
$$\begin{array}{rcl}Var(F_j)& =& {\lambda }_j\\ Cov(F_i,F_j)& =& 0\end{array}$$

贡献率

从最大方差的角度思考，我们的目标是最大化 ${\sum }_{j=1}^{M}Var(F_j)$ 。因此某个主成分对整体的贡献可以被 $\lambda$ 所度量。定义第 $j$ 主成分贡献率
$$c_j=\frac{{\lambda }_j}{{\sum }_{j=1}^n{\lambda }_j}$$

不难看出随着 $j$ 的增大，贡献率在逐渐减小，但累积贡献率 ${\sum }_{i=1}^j{c}_i$ 却在增大。实际应用中，一般对累积贡献率设定某个阈值 $\xi$ 来确定主成分的个数。

被提取率

对于原数据中的每个特征，其信息可能分散到多个主成分中
$$x=[u_1,u_2,\dots ,u_n]y$$

我们同样使用方差作为信息的度量，设
$$Var(x_i)={\sigma }_i^2=Var\left(\sum _{j=1}^n(u_j{)}_i{F}_j\right)=\sum _{j=1}^n(u_j{)}_i^2{\lambda }_j$$

其中 $(u_j{)}_i$ 表示主成分 $u_j$ 的第 $i$ 个元素。定义
$${\Omega }_i=\sum _{j=1}^M\frac{(u_j{)}_i^2{\lambda }_j}{{\sigma }_i^2}$$

为特征 $x_i$ 的被提取率。事实上相关系数
$$\rho (x_i,F_j)=\frac{(u_j{)}_i\sqrt{{\lambda }_j}}{{\sigma }_i}$$

因此被提取率还可以被定义为
$${\Omega }_i=\sum _{j=1}^M{\rho }^2(x_i,F_j)$$

Stochastic Neighbor Embedding

SNE 方法由 Hinton 于 2002 年提出，其基本思想是在低维空间中保持样本在原高维空间中的距离分布。具体来说，对于样本对 $(i,j)$ ，使用 $d_{ij}$ 表示其距离。 SNE 使用类似 Softmax 的方式激活距离度量

$$p_{ij}=\frac{\exp \left(-d_{ij}^2\right)}{{\sum }_{k\ne i}\exp \left(-d_{ik}^2\right)},$$

类似地，低维空间中的距离度量表示为

$$q_{ij}=\frac{\exp \left(-{\tilde{d}}_{ij}^2\right)}{{\sum }_{k\ne i}\exp \left(-{\tilde{d}}_{ik}^2\right)},$$

其中 ${\tilde{d}}_{ij}$ 表示样本对 $(i,j)$ 在低维空间中的距离。 为了保持降维前后，样本之间的距离关系，SNE 最小化 $p_{ij}$ 和 $q_{ij}$ 之间的 KL 散度

$$C=\sum _{i}KL(P_i\Vert Q_i)=\sum _{ij}{p}_{ij}\log \frac{p_{ij}}{q_{ij}}.$$

距离度量 $d_{ij}$ 可以作为输入，也可以通过高维空间的样本计算出来

$$d_{ij}^2=\frac{{\Vert {\mathbf{x}}_i-{\mathbf{x}}_j\Vert }^2}{2{\sigma }_i^2},$$

其中 $x_i$ 和 $x_j$ 分别表示样本 $i$ 和 $j$ 在高维空间中的特征表示。${\sigma }_i$ 可以人工设置，也可以通过求解 $\log k=H(P_i)=-{\sum }_j{p}_{ij}\log p_{ij}$ 得到，其中 $k$ 是一个手工设置的超参数，称为困惑度（perplexity）。类似地，低维空间中的距离定义为

$${\tilde{d}}_{ij}^2={\Vert {\mathbf{y}}_i-{\mathbf{y}}_j\Vert }^2,$$

其中 ${\mathbf{y}}_i$ 和 ${\mathbf{y}}_j$ 分别表示样本 $i$ 和 $j$ 在低维空间中的特征表示。

Visualizing Data using t-SNE

针对 SNE 方法的两个问题

1. 损失函数难以优化
2. 拥挤问题（crowding problem）

t-SNE 做出了两点改进

1. 将高维空间的距离度量改为对称形式，以简化 SNE 损失函数的梯度
2. 将低维空间的距离度量改为学生 t 分布

对称 SNE

SNE 方法中，样本对 $(i,j)$ 的距离度量 $p_{ij}$ 和样本对 $(j,i)$ 的距离度量 $p_{ji}$ 不相等，这是因为 ${\sigma }_i$ 和 ${\sigma }_j$ 不相等导致的。一个直观的想法是对所有样本对使用相同的 $\sigma$ 。但是当数据集中存在离群点时，这样的修改会带来问题。设 ${\mathbf{x}}_i$ 是一个离群点，则对于任意一个样本 ${\mathbf{x}}_j$ ， $p_{ij}$ 的值很小。因此 ${\mathbf{y}}_i$ 的位置基本不影响 SNE 的损失函数，导致 ${\mathbf{y}}_i$ 无法反应 ${\mathbf{x}}_i$ 所包含的信息。

为了解决这个问题，对称 SNE 重新定义邻居概率为

$$p_{ij}\leftarrow \frac{p_{ij}+p_{ji}}{2}.$$

需要注意，这一定义只对高维空间中的距离度量生效，因为低维空间的距离度量已经是对称形式了。对称 SNE 的损失函数具有简单的梯度函数

$$\frac{\delta C}{\delta {\mathbf{y}}_i}=4\sum _j\left({\mathbf{y}}_i-{\mathbf{y}}_j\right)\left(p_{ij}-q_{ij}\right).$$

作为对比， SNE 的损失函数如下

$$\frac{\partial C}{\partial {\mathbf{y}}_i}=2\sum _j\left({\mathbf{y}}_i-{\mathbf{y}}_j\right)\left(p_{ij}-q_{ij}+p_{ji}-q_{ji}\right).$$

t-SNE

t-SNE 方法是 SNE 方法的改进，同样由 Hinton 提出。想象二维曲线流形上的一系列点，将其映射到高维空间中。下图中的“瑞士卷”数据集就是这样的一个例子。显然，我们可以在二维空间中重建这些点之间的距离关系。

下面想象一个十维空间中的流形（例如 MNIST 手写数字数据集），将其表示为高维空间中的一系列样本点。有许多可能的原因使得这些样本被映射到二维空间后不能保持距离关系。首先，十维空间中可能有至多 11 个相互等距的样本点。想要在二维空间中找到 11 个相互等距的样本点却是不可能的。另一个相关的问题是拥挤问题。如果要在二维空间中准确建模相近样本的距离，那么距离稍远的样本就需要被放在二维空间中非常远的位置。在 SNE 中使用样本间的弹簧连接来避免样本在低维空间中过于分散。

t-SNE 使用学生 t 分布建模低维空间的样本距离。相比 SNE 中使用的正态分布， t 分布的尾部概率更大，从而在样本间提供适当的引力。具体来说， t-SNE 使用下式计算 $q_{ij}$

$$q_{ij}=\frac{{\left(1+{\Vert {\mathbf{y}}_i-{\mathbf{y}}_j\Vert }^2\right)}^{-1}}{{\sum }_{k\ne l}{\left(1+{\Vert {\mathbf{y}}_k-{\mathbf{y}}_l\Vert }^2\right)}^{-1}}.$$

t-SNE 的梯度函数如下

$$\frac{\delta C}{\delta {\mathbf{y}}_i}=4\sum _j\left({\mathbf{y}}_i-{\mathbf{y}}_j\right)\left(p_{ij}-q_{ij}\right){\left(1+{\Vert y_i-y_j\Vert }^2\right)}^{-1}.$$

参考资料

1. Bellman, R. (1957). Dynamic Programming.
2. Chmiel, E., & Moore, C. Curse of dimensionality. Reference article, Radiopaedia.org. (accessed on 27 Mar 2022) https://doi.org/10.53347/rID-69156
3. Understanding Curse of Dimensionality
4. Shlens, J. (2003). A Tutorial on Principal Component Analysis. http://arxiv.org/abs/1404.1100
5. Hinton, G., & Roweis, S. (2002). Stochastic Neighbor Embedding. Annual Conference on Neural Information Processing Systems.
6. van der Maaten, L., & Hinton, G. (2008). Visualizing Data using t-SNE. Journal of Machine Learning Research, 9, 2579–2605.