机器学习与高维信息检索 - Note 7 - 核主成分分析（Kernel Principal Component Analysis，K-PCA）

Note 7 - 核主成分分析（Kernel Principal Component Analysis）

标准PCA通过将观察到的数据投射到一个线性子空间来降低其维度。选择投影的方式是使以平方的标准欧氏准则衡量的误差最小，这也可以解释为减少白高斯噪声的一种方式。一个非常重要的应用是将PCA作为分类的预处理，因为分类器在减少噪声的特征空间中表现更好。

标准PCA的主要缺点是，它严重依赖数据的近似线性结构。在许多应用中，这是一个过于严格的假设。核PCA（K-PCA）是标准PCA的一个扩展，它没有这些缺点。K-PCA的关键思想是，它隐含地假设存在一个非线性映射

$$\varphi :{\mathbb{R}}^p\to \mathcal{H},\text{\hspace{1em}(7.1)}$$
其中$\mathcal{H}$是一个非常高维的向量空间（甚至可能是无限维的），其内积为${}^1⟨\cdot ,\cdot ⟩$。我们无妨把$\mathcal{H}$看作是一些${\mathbb{R}}^P$，有非常大的$P$。作为一个初步的步骤，我们重新表述众所周知的标准PCA，使其只涉及我们数据的内积。

${}^1$ 形式上，希尔伯特空间$\mathcal{H}$是一个实数或复数内积空间，就内积所引导的距离函数而言，它也是一个完整的公制空间。

7.1 用内积表示的线性PCA(Linear PCA expressed with inner products)

让$\mathbf{X}\in {\mathbb{R}}^{p\times n}$为居中的数据矩阵（这将是下面推导的一个关键假设），让$\mathbf{K}:={\mathbf{X}}^{\top }\mathbf{X}$为$(n\times n)$矩阵，由数据的所有内积组成。更确切地说，$(i,j)$项是$\mathbf{K}$第$i$和第$j$观测值的内积${\mathbf{x}}_i^{\top }{\mathbf{x}}_j$。

回顾一下，如果$\mathbf{X}=\mathbf{U}\mathbf{\Sigma }{\mathbf{V}}^{\top }$是数据矩阵的SVD。那么数据向量$\mathbf{y}\in {\mathbb{R}}^p$的前$k$主分量是由${\mathbf{U}}_k^{\top }$给出的。在目前的情况下，只给出了内积，因此不可能直接得到${\mathbf{U}}_k$。然而，$\mathbf{K}$的特征值分解允许对这个问题有一个优雅的解决方案。记住，我们有 $\mathbf{K}:={\mathbf{X}}^{\top }\mathbf{X}$ 并用它的SVD代替$\mathbf{X}$。我们得到

$$\mathbf{K}=\mathbf{V}{\mathbf{\Sigma }}^{\top }\mathbf{\Sigma }{\mathbf{V}}^{\top }\text{. }\text{\hspace{1em}(7.2)}$$

记住，$\mathbf{V}$是正交的，${\mathbf{\Sigma }}^{\top }\mathbf{\Sigma }$是对角的。因此，方程（7.2）是$\mathbf{K}$的特征值分解(EVD)，由于EVD的唯一性，通过计算$\mathbf{K}$的$k$最大特征值与它们各自的特征向量，我们得到${\sigma }_1^2,\cdots ,{\sigma }_k^2$和${\mathbf{V}}_k$。我们将假设${\sigma }_k>0$，否则我们可以减少目标子空间的维度而不损失任何数据信息。为了简单起见，我们定义对角矩阵${\Sigma }_k=\mathrm{diag}\left({\sigma }_1,\dots ,{\sigma }_k\right)$。

这使我们能够对于一个给定的观测值 $\mathbf{y}$只需使用内积就可以计算出 ${\mathbf{U}}_k^{\top }\mathbf{y}$的主成分。对于一个新的测量值$\mathbf{y}$，我们有

$$\begin{array}{rl}{\mathbf{V}}_k^{\top }{\mathbf{X}}^{\top }\mathbf{y}& ={\mathbf{V}}_k^{\top }\mathbf{V}\mathbf{\Sigma }{\mathbf{U}}^{\top }\mathbf{y}\\ & =\left[I_k\mid 0\right]\mathbf{\Sigma }{\mathbf{U}}^{\top }\mathbf{y}\\ & =\left[{\mathbf{\Sigma }}_k\mid 0\right]{\mathbf{U}}^{\top }\mathbf{y}\\ & ={\mathbf{\Sigma }}_k{\mathbf{U}}_k^{\top }\mathbf{y}.\end{array}\text{\hspace{1em}(7.3)}$$

将此方程与${\Sigma }_k^{-1}$相乘，我们可以得到

$${\mathbf{U}}_k^{\top }\mathbf{y}={\mathbf{\Sigma }}_k^{-1}{\mathbf{V}}_k^{\top }{\mathbf{X}}^{\top }\mathbf{y}={\mathbf{\Sigma }}_k^{-1}{\mathbf{V}}_k^{\top }\left[\begin{array}{c}{\mathbf{x}}_1^{\top }\mathbf{y}\\ ⋮\\ {\mathbf{x}}_n^{\top }\mathbf{y}\end{array}\right].\text{\hspace{1em}(7.4)}$$

注意，这个方程的右边可以通过只涉及数据的内积的数据来计算，即Gram-matrix $\mathbf{K}$和内积${\mathbf{x}}_n^{\top }\mathbf{y}$。

使数据居中
到目前为止，我们已经假定数据是居中的，也就是说，Gram-matrix $\mathbf{K}$来自居中化的数据。这个假设很关键，因为否则方程（7.4）的推导将不成立。现在让我们假设，我们没有可用的居中数据，而Gram-matrix $\mathbf{K}$是由非居中数据建立的。好消息是，可以通过以下公式直接从$\mathbf{K}$推导出对应于中心数据的格拉姆矩阵，例如$\tilde{\mathbf{K}}$，即

$$\tilde{\mathbf{K}}=\mathbf{H}\mathbf{K}\mathbf{H},\text{\hspace{1em}(7.5)}$$

其中$\mathbf{H}=\left({\mathbf{I}}_n-\frac{1}{n}{1}_n{1}_n^{\top }\right)$.

---

证明：
从方程（1.3）中可以看出，如果$\mathbf{X}$表示非居中的数据矩阵，那么居中的矩阵由以下公式给出

$$\overline{\mathbf{X}}=\mathbf{X}-\hat{\mu }{1}_n^{\top },\text{\hspace{1em}(7.6)}$$

有$\hat{\mu }=\frac{1}{n}\mathbf{X}{1}_n$。 由此，可以看出

$$\overline{\mathbf{X}}=\mathbf{X}-\frac{1}{n}\mathbf{X}{1}_n{1}_n^{\top }=\mathbf{X}\left({\mathbf{I}}_n-\frac{1}{n}{1}_n{1}_n^{\top }\right).\text{\hspace{1em}(7.7)}$$
用对称矩阵的速记符号来表示 $\mathbf{H}:={\mathbf{I}}_n-\frac{1}{n}{1}_n{1}_n^{\top }$。我们可以得到

$$\tilde{\mathbf{K}}={\overline{\mathbf{X}}}^{\top }\overline{\mathbf{X}}={\mathbf{H}\mathbf{X}}^{\top }\mathbf{X}\mathbf{H}=\mathbf{H}\mathbf{K}\mathbf{H}.\text{\hspace{1em}(7.8)}$$

---

因此，如果我们有来自非中心数据的Gram-matrix，我们可以很容易地计算出对应于中心数据的Gram-matrix（无需明确地将$\mathbf{X}$居中）。由此，我们可以–如上所述–计算${\mathbf{V}}_k$和${\mathbf{\Sigma }}_k$。为了计算新数据样本$\mathbf{y}$的主成分，我们首先要根据训练样本把$\mathbf{y}$居中处理。具体来说，这意味着我们必须从$\mathbf{y}$减去训练数据的经验平均值$\hat{\mu }=\frac{1}{n}{\sum }_i{\mathbf{x}}_i=\frac{1}{n}\mathbf{X}{1}_n$。然后我们必须用$\overline{\mathbf{X}}=\mathbf{X}\mathbf{H}$替换公式（7.3）和（7.4）中的$\mathbf{X}$。这就得到了

$${\mathbf{U}}_k^{\top }\left(\mathbf{y}-\frac{1}{n}\mathbf{X}{1}_n\right)={\mathbf{\Sigma }}_k^{-1}{\mathbf{V}}_k^{\top }(\mathbf{X}\mathbf{H}{)}^{\top }\left(\mathbf{y}-\frac{1}{n}\mathbf{X}{1}_n\right)={\mathbf{\Sigma }}_k^{-1}{\mathbf{V}}_k^{\top }{\mathbf{k}}_y\text{\hspace{1em}(7.9)}$$

其中

$${\mathbf{k}}_y=\mathbf{H}\left[\begin{array}{c}{\mathbf{x}}_1^{\top }\mathbf{y}\\ ⋮\\ {\mathbf{x}}_n^{\top }\mathbf{y}\end{array}\right]-\frac{1}{n}\mathbf{H}\mathbf{K}{1}_n$$

7.2 向核PCA过渡 (Transition to Kernel PCA)

现在，通过简单地将内积 ${\mathbf{x}}^{\top }\mathbf{y}$替换为$⟨\varphi (\mathbf{x}),\varphi (\mathbf{y})⟩$来扩展经典的PCA是很简单的。K-PCA的实际成功是由于计算$⟨\varphi (\mathbf{x}),\varphi (\mathbf{y})⟩$时，既不需要$\varphi$也不需要$⟨\cdot ,\cdot ⟩$ 的内积。相反，只需要有一个函数

$$\kappa :{\mathbb{R}}^p\times {\mathbb{R}}^p\to \mathbb{R},(\mathbf{x},\mathbf{y})\mapsto \kappa (\mathbf{x},\mathbf{y})\text{\hspace{1em}(7.10)}$$

这反映了$⟨\varphi (\mathbf{x}),\varphi (\mathbf{y})⟩$的属性，即对于所有$\mathbf{x}\in {\mathbb{R}}^p$来说是对称的，并且满足$\kappa (\mathbf{x},\mathbf{x})\ge 0$ 的正定属性。用$\kappa (\mathbf{x},\mathbf{y})$代替$⟨\varphi (\mathbf{x}),\varphi (\mathbf{y})⟩$，因此不需要知道特征映射$\varphi$，称为核技巧（Kernel trick）。

Definition 7.1 正定核 (Positive Definite Kernel)

1. 令 S : = { x 1 , … , x n } ⊂ R p S:=\left\{\mathbf{x}_{1}, \ldots, \mathbf{x}_{n}\right\} \subset \mathbb{R}^{p} S:={x1​,…,xn​}⊂Rp 与$\kappa (\cdot )$如上定义. 带有$(i,j)$ 项$k_{ij}=\kappa \left({\mathbf{x}}_i,{\mathbf{x}}_j\right)$的 $(n\times n)$-矩阵$\mathbf{K}$ 被称为$S$的Gram-或Kernel-矩阵，相对于$\kappa$。

2. 如果对于所有有限的、非空的集合 $S\subset {\mathbb{R}}^p$来说，相应的Gram-matrix是正半定的（正定的），那么这个函数$\kappa (\cdot )$被称为核（正定的核）。

在上一小节中，我们讨论了在潜在的希尔伯特空间中的数据居中问题。我们在这里也面临同样的问题。当一个新的数据点出现时，直接计算新数据点的非中心化内积向量与训练数据的关系是很简单的。我们将把它称为

$${\mathbf{k}}^{\text{new }}=\left[\begin{array}{c}⟨\varphi \left({\mathbf{x}}_1\right),\varphi (\mathbf{y})⟩\\ ⋮\\ ⟨\varphi \left({\mathbf{x}}_n\right),\varphi (\mathbf{y})⟩\end{array}\right].$$

那么，中心化数据的第$j$条是

$${\left({\mathbf{k}}_{cent}^{new}\right)}_j=⟨\varphi \left({\mathbf{x}}_j\right)-\frac{1}{n}\sum _i\varphi \left({\mathbf{x}}_i\right),\varphi (\mathbf{y})-\frac{1}{n}\sum _i\varphi \left({\mathbf{x}}_i\right)⟩$$

为了找到一个更简洁的表达方式，我们利用标量乘积的线性来得到

$$\begin{array}{c}⟨\varphi \left({\mathbf{x}}_j\right)-\frac{1}{n}\sum _i\varphi \left({\mathbf{x}}_i\right),\varphi (\mathbf{y})-\frac{1}{n}\sum _i\varphi \left({\mathbf{x}}_i\right)⟩\\ =⟨\varphi \left({\mathbf{x}}_j\right),\varphi (\mathbf{y})⟩-\frac{1}{n}⟨\sum _i\varphi \left({\mathbf{x}}_i\right),\varphi (\mathbf{y})⟩\\ -\frac{1}{n}⟨\varphi \left({\mathbf{x}}_j\right),\sum _i\varphi \left({\mathbf{x}}_i\right)⟩+\frac{1}{n^2}⟨\sum _i\varphi \left({\mathbf{x}}_i\right),\sum _i\varphi \left({\mathbf{x}}_i\right)⟩,\end{array}$$

而且很容易看出，

$$\begin{array}{c}{\mathbf{k}}_{\text{cent }}^{\text{new }}={\mathbf{k}}^{\text{new }}-\frac{1}{n}{1}_n{1}_n^{\top }{\mathbf{k}}^{\text{new }}-\frac{1}{n}\mathbf{K}{1}_n+\frac{1}{n^2}{1}_n{1}_n^{\top }\mathbf{K}{1}_n\\ ={\mathbf{H}\mathbf{k}}^{\text{new }}-\frac{1}{n}\mathbf{H}\mathbf{K}{1}_n.\end{array}$$

总之，对于一个训练集$\mathbf{X}=\left[{\mathbf{x}}_1,\dots {\mathbf{x}}_n\right],{\mathbf{x}}_i\in {\mathbb{R}}^p$ 的K-PCA包括以下步骤。

1. 找到一个合适的内核函数 $\kappa (\cdot )$ 并计算Gram-matrix

$$\mathbf{K}=\left[\begin{array}{cccc}\kappa \left({\mathbf{x}}_1,{\mathbf{x}}_1\right)& \kappa \left({\mathbf{x}}_1,{\mathbf{x}}_2\right)& \dots & \kappa \left({\mathbf{x}}_1,{\mathbf{x}}_n\right)\\ \kappa \left({\mathbf{x}}_2,{\mathbf{x}}_1\right)& \kappa \left({\mathbf{x}}_2,{\mathbf{x}}_2\right)& & \kappa \left({\mathbf{x}}_2,{\mathbf{x}}_n\right)\\ ⋮& ⋮& \ddots & ⋮\\ \kappa \left({\mathbf{x}}_n,{\mathbf{x}}_1\right)& \kappa \left({\mathbf{x}}_n,{\mathbf{x}}_2\right)& \dots & \kappa \left({\mathbf{x}}_n,{\mathbf{x}}_n\right)\end{array}\right]$$

2. 计算Gram-matrix $\tilde{\mathbf{K}}=\mathbf{H}\mathbf{K}\mathbf{H}$，其中$\mathbf{H}={\mathbf{I}}_n-\frac{1}{n}{1}_n{1}_n^{\top }$是居中化的数据。

3. 计算特征值分解$\tilde{\mathbf{K}}=\mathbf{V}\mathbf{\Lambda }{\mathbf{V}}^{\top }$。由于核函数的定义，矩阵$\mathbf{K}$是半正定的，因此$\mathbf{\Lambda }$的对角线项是非负的。因此，我们可以写成 $\mathbf{\Lambda }={\mathbf{\Sigma }}^2=\mathrm{diag}\left({\sigma }_1^2,\dots ,{\sigma }_n^2\right)$.

4. 定义缩减矩阵 ${\mathbf{\Sigma }}_k=\mathrm{diag}\left({\sigma }_1,\dots ,{\sigma }_k\right)\in {\mathbb{R}}^{k\times k}$ 与${\mathbf{V}}_k\in {\mathbb{R}}^{n\times k}$.

5. 缩减后的训练数据由以下公式得出

$$\mathbf{S}={\mathbf{\Sigma }}_k{\mathbf{V}}_k^{\top }\text{\hspace{1em}(7.11)}$$

6. 对于一个新的数据点$\mathbf{y}\in {\mathbb{R}}^p$ ，$k$核主成分的计算方法是

$${\mathbf{s}}_{\text{new }}:={\mathbf{\Sigma }}_k^{-1}{\mathbf{V}}_k^{\top }{\mathbf{k}}_{\text{cent }}^{\text{new }}\text{\hspace{1em}(7.12)}$$

其中

$$\begin{array}{rl}{\mathbf{k}}_{\text{cent }}^{\text{new }}& =\mathbf{H}\left({\mathbf{k}}^{\text{new }}-\frac{1}{n}\mathbf{K}{1}_n\right)\\ \text{ where }{\mathbf{k}}^{\text{new }}& ={\left[\kappa \left({\mathbf{x}}_1,\mathbf{y}\right),\dots ,\kappa \left({\mathbf{x}}_n,\mathbf{y}\right)\right]}^{\top }.\end{array}$$