CH20-Dimensionality Reduction

无监督学习的一个常见问题是降维，需要学习一个从高维空间 $x\in \mathbb{R}^D$ 到低维空间 $z\in\mathbb{R}^L$ 的映射。这个映射既可以是一个可以用于任意输入的参数模型 $z = f(x;\theta)$ ，又可以是一个非参数模型，计算输入数据集中数据点的embedding，但不能用于其他数据。后一种方法最常用来数据可视化，而前一种方法还可以用来作为其他学习算法的预处理步骤。例如，我们可以先降维，把 $x$ 映射到 $z$ ，再用 $z$ 训练一个线性分类器，把 $z$ 映射到 $y$ .

20.1 Principal components analysis#

最常用的降维方法就是主成分分析（PCA），其基本想法是为高维数据 $x\in \mathbb{R}^D$ 寻找一个低维空间 $z\in\mathbb{R}^L$ 的线性且正交的投影，我们认为这样的低维表示是原始数据的一个较好的近似，原因如下：

如果我们先投影 or 编码 $x$ 得到了 $z = \mathrm{W}^{\mathsf{T}}x$ ，再反投影 or 解码 $z$ 得到 $\hat{x} = \mathrm{W}z$ ，我们希望 $\hat{x}$ 和 $x$ 的 $\ell_2$ 距离尽可能近，特殊情况下，我们可以定义如下重建误差 or distortion：

\mathcal{L}(\mathrm{W})\triangleq \sum_{n=1}^N|| x_n - \mathrm{decode}(\mathrm{encode}(x_n;\mathrm{W});\mathrm{W}) ||^2_2

其中 encode 和 decode 函数均为线性映射，原因我们在下面会解释。

在 Section 20.1.2 中，我们会显示上面这个目标函数的最小值在 $\hat{\mathrm{W}}=\mathrm{U}_L$ 处取到， $\mathrm{U}_L$ 包含了经验协方差矩阵的特征值最大的 $L$ 个特征向量，经验协方差为：

\hat{\Sigma}=\frac{1}{N} \sum_{n=1}^N (x_n-\bar{x})(x_n-\bar{x})^{\mathsf{T}}=\frac{1}{N} \mathrm{X}_c^{\mathsf{T}} \mathrm{X}_c

其中 $\mathrm{X}_c$ 是 $N\times D$ 的设计矩阵的中心化版本。

在 Section 20.2.2 中，我们会显示这等价于最大化似然一个潜在线性高斯模型，又称为 probabilistic PCA.

20.1.1 Examples#

如果我们把 2d 的数据投影到 1d，这个 1d 的方向就是原数据方差最大的方向。

20.1.2 Derivation of the algorithm#

假设我们有一个数据集 $\mathcal{D} = \{x_n:n=1:N\}$ ，其中 $x_n\in \mathbb{R}^D$ . 我们可以用一个 $N\times D$ 的数据矩阵 $X$ 来代表这个数据集；我们假设 $\bar{x}=\frac{1}{N}\sum_{n=1}^Nx_n=0$ ，可以用中心化数据确保这个条件。

现在我们要用一个低维表示 $z_n\in \mathbb{R}^L$ 来近似 $x_n$ . 我们认为每个 $x_n$ 都可以用一系列基函数 $w_1,\dots,w_L$ 的加权组合来“解释”， $w_k\in\mathbb{R}^D$ , 而这些基函数的权重就是 $z_n\in \mathbb{R}^L$ ，也就是说，我们假设 $x_n\approx\sum_{k=1}^L z_{nk}w_k$ .

$z_n$ 向量是 $x_n$ 的低维表示，也称作潜在向量，因为它是由数据中未观察到的隐变量组成的。这些隐变量的集合被称作潜在因子。

损失函数：

\begin{aligned} \mathcal{L}(\mathrm{W},\mathrm{Z}) &=\frac{1}{N} ||\mathrm{X}-\mathrm{ZW}^{\mathsf{T}}||_F^2\\ &=\frac{1}{N} ||\mathrm{X}^{\mathsf{T}}-\mathrm{WZ}^{\mathsf{T}}||_F^2\\ &=\frac{1}{N} ||\mathrm{X}^{\mathsf{T}}-\mathrm{WZ}^{\mathsf{T}}||_F^2\\ &=\frac{1}{N} \sum_{n=1}^N||x_n-\mathrm{W}z_n||^2\\ \end{aligned}

这也被称作（平均）重建误差。

我们想在“ $W$ 是一个正交矩阵 ”这个限制条件下优化这个损失函数，下面我们证明最优解在 $\hat{\mathrm{W}}=\mathrm{U}_L$ 取到。

20.1.2.1 Base case#

从一维解的情况开始讨论， $w_1\in\mathbb{R}^D$ ，接下来我们再讨论 $w_2,w_3,\dots$ 的解。

记数据点关于第一个基向量的系数为 $\tilde{\mathrm{z}}_1 = [z_{11},\dots,z_{N1}]\in\mathbb{R}^N$ . 重建误差为：

\begin{aligned} \mathcal{L}(w_1,\tilde{\mathrm{z}}_1)&= \frac{1}{N} \sum_{n=1}^N || x_n-z_{n1}w_1||^2\\ &=\frac{1}{N} \sum_{n=1}^N (x_n-z_{n1}w_1)^{\mathsf{T}}(x_n-z_{n1}w_1)\\ &=\frac{1}{N} \sum_{n=1}^N [x_n^{\mathsf{T}}x_n-2z_{n1}w_1^{\mathsf{T}}x_n+z_{n1}^2w_1^{\mathsf{T}}w_1]\\ &=\frac{1}{N} \sum_{n=1}^N [x_n^{\mathsf{T}}x_n-2z_{n1}w_1^{\mathsf{T}}x_n+z_{n1}^2] \end{aligned}

其中由正交假设可知 $w_1^{\mathsf{T}}w_1=1$ . 对 $z_{n1}$ 求导可得：

\frac{\partial}{\partial{z_{n1}}}\mathcal{L}(w_1,\tilde{\mathrm{z}_1})= \frac{1}{N}[-2w_1^{\mathsf{T}}x_n+2z_{n1}]=0\\ \Longrightarrow z_{n1}=w_1^{\mathsf{T}}x_n

所以最佳的 embedding 就是把数据正交投影到 $w_1$ . 把这个解再代入损失函数就得到了该权重的损失：

\begin{aligned} \mathcal{L}(w_1)=&\mathcal{L}(w_1,\tilde{\mathrm{z}}_1^*(w_1))\\ &=\frac{1}{N}\sum_{n=1}^N[x_n^\mathsf{T}x_n-z_{n1}^2]\\ &=\mathrm{const}-\frac{1}{N}\sum_{n=1}^{N}z_{n1}^2 \end{aligned}

接下来对 $w_1$ 进行优化, 注意到

\mathcal{L}(w_1)=-\frac{1}{N}\sum_{n=1}^N z_{n1}^2 = -\frac{1}{N}\sum_{n=1}^N w_1^{\mathsf{T}} x_n x_n^{\mathsf{T}} w_1 =-w_1^{\mathsf{T}} \hat{\Sigma} w_1

其中 $\Sigma$ 是经验协方差矩阵（因为数据中心化了）。

对于这个问题，如果不限制 $||w_1||$ ，我们很难优化，所以我们加入约束条件： $||w_1||=1$ 并运用拉格朗日乘子法可得：

\hat{\mathcal{L}}(w_1)= w_1^\mathsf{T}\hat{\Sigma}w_1-\lambda_1(w_1^\mathsf{T}w_1-1)

对 $w_1$ 求导得：

\frac{\partial}{\partial{w_{1}}}\hat{\mathcal{L}}(w_1)=2\hat{\Sigma}w_1-2\lambda_1 w_1=0\\ \Longrightarrow \hat\Sigma w_1 = \lambda_1 w_1

因此最佳投影方向是协方差矩阵的特征值。两边乘上 $w_1^\mathsf{T}$ 并使用条件 $w_1^\mathsf{T}w_1 = 1$ 我们发现：

w_1^\mathsf{T}\hat\Sigma w_1 = \lambda_1

由于我们要最小化损失，也就是最大化这个量，我们需要选取最大的特征值对应的特征向量。

20.1.2.2 最优化的权重向量最大化了投影数据的方差#

由于数据被中心化了，我们有：

\mathbb{E}[z_{n1}]=\mathbb{E}[x_n^\mathsf{T}w_1]=\mathbb{E}[x_n]^\mathsf{T}w_1=0

因此投影数据的方差为：

\mathbb{V}[\tilde{\mathrm{z}}_1]=\mathbb{E}[\tilde{\mathrm{z}}_1^2]-(\mathbb{E}[\tilde{\mathrm{z}}_1])^2=\frac{1}{N}\sum_{n=1}^Nz_{n1}^2-0=-\mathcal{L}(w_1)+\mathrm{const}

因此最小化重建误差就等价于最大化投影数据的方差，这也是为什么经常说 PCA 是找到方差最大的方向。

20.1.2.3 Induction step#

现在我们找另一个方向 $w_2$ 来进一步最小化重建误差，约束条件： $w_1^\mathsf{T}w_2=0$ 和 $w_2^\mathsf{T}w_2 = 1$ . 误差为：

\mathcal{L}(w_1,\tilde{\mathrm{z}}_1,w_2,\tilde{\mathrm{z}}_2)= \frac{1}{N}\sum_{n=1}^N ||x_n-z_{n1}w_1-z_{n2}w_2 ||^2

对 $w_1$ 和 $z_1$ 优化得到和上一个部分一样的结果，对 $w_2$ 和 $z_2$ 优化得到：

\hat\Sigma w_2=\lambda_2w_2

也就是第二大的特征值对应的特征向量。

由此继续下去，可知 $\hat{\mathrm{W}}=\mathrm{U}_L$ .

20.1.3 Computational issues#

这一部分我们讨论一些使用 PCA 相关的实际问题。

20.1.3.1 协方差矩阵 vs 相关性矩阵#

我们一直在对协方差矩阵做特征值分解，但是更好的选择是相关性矩阵。

原因：PCA 容易被未标准化的数据误导，有时候方差大的方向仅仅是由于测量尺度不同导致的，下图展示了一个这样的例子。

Note：标准化后的样本数据的协方差矩阵就是原始样本数据的相关矩阵

pca_standardization

20.1.3.2 处理高维数据#

现在我们已经把 PCA 转换成了一个找 $D\times D$ 的协方差矩阵 $\mathrm{X}^\mathsf{T}\mathrm{X}$ 的特征向量的问题。如果 $D>N$ ，此时计算 $N\times N$ 的 Gram matrix $\mathrm{X}\mathrm{X}^\mathsf{T}$ 速度更快。

令 $\mathrm{U}$ 为一个包含 $\mathrm{X}\mathrm{X}^\mathsf{T}$ 特征向量的正交矩阵，对应的特征值在 $\Lambda$ 中。由定义有 $(\mathrm{X}\mathrm{X}^\mathsf{T})\mathrm{U}=\mathrm{U}\Lambda$ . 两边同乘 $\mathrm{X}^\mathsf{T}$ 得到： $(\mathrm{X}^\mathsf{T}\mathrm{X})(\mathrm{X}^\mathsf{T}\mathrm{U})=(\mathrm{X}^\mathsf{T}\mathrm{U})\Lambda$ .

从上式我们可以看到协方差矩阵 $\mathrm{X}^\mathsf{T}\mathrm{X}$ 的特征向量是 $\mathrm{V}=\mathrm{X}^\mathsf{T}\mathrm{U}$ ，特征值由 $\Lambda$ 给出。但是，这些特征向量是未标准化的，因为 $||v_j||^2=u_j^\mathsf{T}\mathrm{XX}^\mathsf{T}u_j=\lambda_ju_j^\mathsf{T}u_j=\lambda_j$ .

标准化后的特征向量由 $\mathrm{V}=\mathrm{X}^\mathsf{T}\mathrm{U}\Lambda^{-\frac{1}{2}}$ 给出。

上述方法提供了计算 PCA 基函数的另一种方法，同时让我们可以使用核技巧，见 Section 20.4.6.

20.1.3.3 用 SVD 计算 PCA#

这一部分我们显示用特征向量法计算 PCA 等价于截断 SVD.

假设 $\mathrm{X}$ 以被中心化，令协方差矩阵 $\Sigma \propto \mathrm{X}^{\mathsf{T}}\mathrm{X}$ 的 top L 特征分解矩阵为 $\mathrm{U}_{\Sigma}\Lambda_\Sigma\mathrm{U}^{\mathsf{T}}_\Sigma$ . 由 Section 20.1.2 的推导可知投影权重的最佳估计是 top L 特征值，所以 $\mathrm{W}=\mathrm{U}_\Sigma$ .

现在令 $\mathrm{U}_{\Sigma}\Lambda_\Sigma\mathrm{V}^{\mathsf{T}}_\Sigma \approx \mathrm{X}$ 为数据矩阵 $\mathrm{X}$ 的 L-truncated SVD 近似. 由第 7 章可知 $\mathrm{X}$ 的右奇异向量是 $\mathrm{X}^{\mathsf{T}}\mathrm{X}$ 的特征向量，所以 $\mathrm{V}_X=\mathrm{U}_\Sigma=\mathrm{W}$ . （另外，协方差矩阵的特征值和数据矩阵的奇异值之间差一个系数 $\lambda_k=s_k^2/N$ .）

现在假设我们对投影后的结果而不是投影矩阵感兴趣，我们有：

\mathrm{Z}=\mathrm{XW}=\mathrm{U}_X \mathrm{S}_X \mathrm{V}^{\mathsf{T}}_X \mathrm{V}_X = \mathrm{U}_X \mathrm{S}_X

最后，如果我们想近似重建数据，我们有：

\hat{X}=\mathrm{ZW}^{\mathsf{T}}=\mathrm{U}_X \mathrm{S}_X \mathrm{V}^{\mathsf{T}}_X

这和截断 SVD 近似是一模一样的。

因此我们发现对协方差矩阵 $\Sigma$ 做特征值分解和对数据矩阵 $\mathrm{X}$ 做奇异值分解都是做 PCA. 通常倾向于使用后一种方法，主要是计算上方便。对于非常高维的问题，我们可以用随机化的 SVD 算法来减少时间复杂度。

20.1.4 选择潜在维度的数目#

20.1.4.1 重建误差#

定义使用 L 个潜在维度的模型在数据集 $\mathcal{D}$ 的重建误差：

\mathcal{L}_L=\frac{1}{|\mathcal{D}|}\sum_{n\in\mathcal{D}}||x_n-\hat{x}_n||^2\\ where\quad \hat{x}_n=\mathrm{W}z_n+\mu ,\quad z_n= \mathrm{W}^{\mathsf{T}}(x_n-\mu)

$\mu$ 是经验误差， $\mathrm{W}$ 按如上方法估计。在 MNIST 数据集上画出 $\mathcal{L}_L$ 随 $L$ 变化图像，可以发现它下降很快，意味着我们可以用很小数目的因子来捕捉到 pixel 之间绝大部分的经验相关性。

20.1.4.2 Scree plots#

我们可以证明：

\mathcal{L}_L=\sum_{j=L+1}^{D}\lambda_j

当维数增加，特征值变小，重建误差也随之变小，可以画出对应的 scree plot 来选择合适的维数。

另一个常用量是方差解释比例，定义如下：

F_L = \frac{\sum_{j=1}^L \lambda_j}{\sum_{j^\prime=1}^{L^{max}}\lambda_{j^\prime}}

20.1.4.3 Profile likelihood#

尽管在 PCA 的损失函数图中无 U 形，但有时也会有一个 knee 或 elbow，即损失突然变小的点。这个想法可以理解为： $L^*$ 是“true”潜在维度，当 $L<L^*$ ，损失函数的下降速率会快，当 $L>L^*$ ，速率就会变慢。一种自动检测这种变化点的方法是计算 profile likelihood.