TIT-2014《Randomized Dimensionality Reduction for $k$-means Clustering》

推荐深蓝学院的《深度神经网络加速：cuDNN 与 TensorRT》，课程面向就业，细致讲解CUDA运算的理论支撑与实践，学完可以系统化掌握CUDA基础编程知识以及TensorRT实战，并且能够利用GPU开发高性能、高并发的软件系统，感兴趣可以直接看看链接：深蓝学院《深度神经网络加速：cuDNN 与 TensorRT》

---

核心思想

论文的核心思想是研究 $k$-means 聚类的降维方法，旨在通过特征选择和特征提取技术降低高维数据的计算复杂性，同时保持聚类质量的理论保证。论文提出了一种新的特征选择方法和两种改进的特征提取方法，均基于随机化算法，提供了常数因子近似保证。核心观点是将 $k$-means 聚类问题从线性代数的视角重新表述，利用矩阵低秩近似技术（如奇异值分解 SVD 和随机投影）来实现降维，同时保持聚类目标函数的近似精度。

具体而言，论文通过以下方式实现降维：

1. 特征选择：从原始特征中选择一个小的子集，基于随机采样技术，结合近似 SVD 的思想，提出首个具有理论保证的 $k$-means 特征选择算法。
2. 特征提取：通过随机投影和快速近似 SVD 构建新的低维特征，优化时间复杂度和特征数量，同时保持聚类质量。

---

目标函数

$k$-means 聚类的目标是给定 $m$ 个欧几里得点 P = { p 1 , p 2 , … , p m } ⊆ R n \mathcal{P} = \{\mathbf{p}_1, \mathbf{p}_2, \ldots, \mathbf{p}_m\} \subseteq \mathbb{R}^n P={p1​,p2​,…,pm​}⊆Rn 和聚类数量 $k$，将这些点划分为 $k$ 个簇，使得每个点到其最近簇中心的平方欧几里得距离之和最小。目标函数定义为：

$$\mathcal{F}(\mathcal{P},\mathcal{S})=\sum _{i=1}^m\Vert {\mathbf{p}}_i-\mathbf{\mu }({\mathbf{p}}_i){\Vert }_2^2$$

其中：

• S = { S 1 , S 2 , … , S k } \mathcal{S} = \{\mathcal{S}_1, \mathcal{S}_2, \ldots, \mathcal{S}_k\} S={S1​,S2​,…,Sk​} 是 $\mathcal{P}$ 的一个 $k$ 分区，${\mathcal{S}}_j$ 表示第 $j$ 个簇，$s_j=|{\mathcal{S}}_j|$ 是簇的大小。
• ${\mathbf{\mu }}_j=\frac{{\sum }_{{\mathbf{p}}_i\in {\mathcal{S}}_j}{\mathbf{p}}_i}{s_j}$ 是第 $j$ 个簇的质心，$\mathbf{\mu }({\mathbf{p}}_i)$ 是点 ${\mathbf{p}}_i$ 所属簇的质心。

在线性代数视角下，数据点组成矩阵 $\mathbf{A}\in {\mathbb{R}}^{m\times n}$，每行表示一个数据点 ${\mathbf{p}}_i^{\top }$。聚类由簇指示矩阵 $\mathbf{X}\in {\mathbb{R}}^{m\times k}$ 表示，其中 ${\mathbf{X}}_{ij}=1/\sqrt{s_j}$ 如果点 ${\mathbf{p}}_i$ 属于簇 ${\mathcal{S}}_j$，否则为 0。目标函数可重写为：

$$\mathcal{F}(\mathbf{A},\mathbf{X})=\Vert \mathbf{A}-\mathbf{X}{\mathbf{X}}^{\top }\mathbf{A}{\Vert }_F^2$$

其中 $\Vert \cdot {\Vert }_F$ 表示 Frobenius 范数。目标是找到最优的簇指示矩阵 ${\mathbf{X}}_{\text{opt}}$，使得：

$${\mathbf{X}}_{\text{opt}}=\underset{\mathbf{X}\in \mathcal{X}}{\mathrm{argmin}}\Vert \mathbf{A}-\mathbf{X}{\mathbf{X}}^{\top }\mathbf{A}{\Vert }_F^2$$

其中 $\mathcal{X}$ 是所有 $m\times k$ 簇指示矩阵的集合，最优目标函数值为 $\mathcal{F}(\mathbf{A},{\mathbf{X}}_{\text{opt}})={\mathbf{F}}_{\text{opt}}$。

降维的目标是构建低维点集 P ~ = { p ~ 1 , p ~ 2 , … , p ~ m } ⊆ R r \tilde{\mathcal{P}} = \{\tilde{\mathbf{p}}_1, \tilde{\mathbf{p}}_2, \ldots, \tilde{\mathbf{p}}_m\} \subseteq \mathbb{R}^r P~={p~​1​,p~​2​,…,p~​m​}⊆Rr（$r\ll n$），使得在低维空间计算的最优 $k$-means 分区 ${\tilde{\mathcal{S}}}_{\text{opt}}$ 在原始空间的目标函数值满足：

$$\mathcal{F}(\mathcal{P},{\tilde{\mathcal{S}}}_{\text{opt}})\le \gamma \cdot \mathcal{F}(\mathcal{P},{\mathcal{S}}_{\text{opt}})$$

其中 $\gamma >0$ 是近似比率。

---

目标函数的优化过程

$k$-means 聚类的优化是一个 NP 难问题，传统方法如 Lloyd 算法通过迭代优化目标函数，但计算复杂性随维度 $n$ 增加而显著上升。论文通过降维降低计算复杂性，具体优化过程如下：

1. 特征选择（Theorem 11）：

   • 步骤：
     1. 给定数据矩阵 $\mathbf{A}\in {\mathbb{R}}^{m\times n}$ 和簇数 $k$，计算近似右奇异向量矩阵 $\mathbf{Z}\in {\mathbb{R}}^{n\times k}$，使用快速 Frobenius 范数 SVD 算法（Lemma 4），满足 ${\mathbf{Z}}^{\top }\mathbf{Z}={\mathbf{I}}_k$ 且 $\mathbb{E}\Vert \mathbf{A}-\mathbf{AZ}{\mathbf{Z}}^{\top }{\Vert }_F^2\le (1+\epsilon )\Vert \mathbf{A}-{\mathbf{A}}_k{\Vert }_F^2$。
     2. 使用随机采样方法（Definition 5），根据 $\mathbf{Z}$ 的行范数计算采样概率 $p_i=\frac{\Vert {\mathbf{Z}}_{(i)}{\Vert }_2^2}{\Vert \mathbf{Z}{\Vert }_F^2}$，从中选择 $r=O(k\log (k)/{\epsilon }^2)$ 个特征。
     3. 构建采样矩阵 $\mathbf{\Omega }\in {\mathbb{R}}^{n\times r}$ 和重缩放矩阵 $\mathbf{S}\in {\mathbb{R}}^{r\times r}$，得到低维矩阵 $\mathbf{C}=\mathbf{A}\mathbf{\Omega }\mathbf{S}$。
     4. 在 $\mathbf{C}$ 上运行 $k$-means 算法，得到低维最优分区 ${\tilde{\mathbf{X}}}_{\text{opt}}$，并将其应用于原始数据 $\mathbf{A}$。
   • 优化原理：通过近似 SVD 和随机采样，保留数据的主要结构，降低维度，同时保证目标函数的 $(3+\epsilon )$ 近似。

2. 特征提取 - 随机投影（Theorem 12）：

   • 步骤：
     1. 构造随机投影矩阵 $\mathbf{R}\in {\mathbb{R}}^{n\times r}$，其中 $r=O(k/{\epsilon }^2)$，元素为标准高斯分布或重缩放的随机符号。
     2. 计算低维矩阵 $\mathbf{C}=\mathbf{AR}\in {\mathbb{R}}^{m\times r}$，表示低维点集 $\tilde{\mathcal{P}}$。
     3. 在 $\mathbf{C}$ 上运行 $k$-means 算法，得到 ${\tilde{\mathbf{X}}}_{\text{opt}}$，并计算原始空间的目标函数值。
   • 优化原理：利用 Johnson-Lindenstrauss 引理，随机投影保留点之间的欧几里得距离，目标函数近似误差为 $(2+\epsilon )$。

3. 特征提取 - 近似 SVD（Theorem 13）：

   • 步骤：
     1. 使用快速近似 SVD 算法（Lemma 4）计算 $\mathbf{Z}\in {\mathbb{R}}^{n\times k}$，近似 $\mathbf{A}$ 的前 $k$ 个右奇异向量。
     2. 构建低维矩阵 $\mathbf{C}=\mathbf{AZ}\in {\mathbb{R}}^{m\times k}$。
     3. 在 $\mathbf{C}$ 上运行 $k$-means 算法，得到 ${\tilde{\mathbf{X}}}_{\text{opt}}$。
   • 优化原理：近似 SVD 提供接近最优的低秩近似，目标函数近似误差为 $(2+\epsilon )$，时间复杂度显著低于精确 SVD。

优化过程的关键是通过矩阵 $\mathbf{C}=\mathbf{AD}$（其中 $\mathbf{D}$ 为特征选择或提取矩阵）构造低维表示，确保 $\mathbf{C}\cdot \mathbf{H}$（适当的 $\mathbf{H}$）在 Frobenius 范数下接近 ${\mathbf{A}}_k$（$\mathbf{A}$ 的最佳秩 $k$ 近似）。利用矩阵 Pythagorean 定理和 SVD 的正交性，证明近似误差受控。

---

主要贡献点

1. 首个理论保证的特征选择算法：
   • 提出了一种随机化特征选择算法（Theorem 11），以 $O(mnk/\epsilon +k\log (k)/{\epsilon }^2\log (k\log (k)/\epsilon ))$ 时间复杂度选择 $r=O(k\log (k)/{\epsilon }^2)$ 个特征，达到 $(3+\epsilon )$ 近似误差。这是首个具有理论保证的 $k$-means 特征选择方法。
2. 改进的随机投影特征提取：
   • Theorem 12 提出了一种随机投影方法，所需维度从 $O(\log (m)/{\epsilon }^2)$ 减少到 $O(k/{\epsilon }^2)$，时间复杂度为 $O(mn[{\epsilon }^{-2}k/\log (n)])$，近似误差为 $(2+\epsilon )$，优于传统随机投影结果。
3. 快速近似 SVD 特征提取：
   • Theorem 13 利用快速近似 SVD，仅需 $r=k$ 个特征和 $O(mnk/\epsilon )$ 时间复杂度，达到 $(2+\epsilon )$ 近似误差，显著优于精确 SVD 的 $O(mn\min (m,n))$ 时间复杂度。
4. 线性代数视角：
   • 将 $k$-means 聚类问题重构为矩阵低秩近似问题，利用 SVD 和随机投影的理论工具，提供了统一的分析框架。
5. 实验验证：
   • 在合成数据集和真实数据集（如 USPS、COIL20、LIGHT、PIE、ORL）上验证了算法的有效性，表明小维度（如 $r=20$ 或 $30$）即可实现接近最优的聚类效果。

---

实验结果

实验在合成数据集和多个真实数据集（USPS、COIL20、LIGHT、PIE、ORL）上进行，评估了运行时间、目标函数值（归一化形式 $\mathcal{F}/\Vert \mathbf{A}{\Vert }_F^2$）和聚类准确率（基于标签的误分类率）。主要结果如下：

1. 合成数据集：
   • 降维方法显著优于朴素 $k$-means，运行时间大幅降低。
   • 当维度 $r\approx 20$ 时，聚类准确率接近最优，表明降维在分离良好的数据上非常有效。
2. 真实数据集：
   • 随着维度 $r$ 增加，归一化目标函数值逐渐接近朴素 $k$-means 的值。
   • 在 USPS、LIGHT 和 ORL 数据集上，提出的降维方法在准确率和目标函数值上优于 Laplacian Scores 方法。
   • 在 PIE 和 COIL20 数据集上，Laplacian Scores 在准确率上更优，但朴素 $k$-means 表现较差，表明这些数据可能分离性较差。
3. 运行时间：
   • 降维后 $k$-means 的每次迭代时间复杂度从 $O(kmn)$ 降至 $O(m{k}^2/{\epsilon }^2)$，显著提高效率。
   • 运行时间不随维度单调增加，可能是因为降维后 Lloyd 算法迭代次数变化。
4. 结论：
   • 实验表明 $r=20$ 或 $30$ 即可实现接近最优的聚类效果，验证了算法在理论和实践中的有效性。

---

算法实现过程

以下详细解释三种算法的实现过程，结合数学公式和伪代码。

1. 特征选择算法（Theorem 11）

算法描述：基于随机采样的特征选择，结合近似 SVD。
输入：数据矩阵 $\mathbf{A}\in {\mathbb{R}}^{m\times n}$，簇数 $k$，误差参数 $\epsilon$，失败概率 $\delta$。
输出：低维矩阵 $\mathbf{C}\in {\mathbb{R}}^{m\times r}$，簇指示矩阵 ${\tilde{\mathbf{X}}}_{\text{opt}}$。
步骤：

1. 计算近似右奇异向量：
   • 使用快速 Frobenius 范数 SVD 算法（Lemma 4）：
     $$\mathbf{Z}=\text{FastFrobeniusSVD}(\mathbf{A},k,\epsilon )$$
     其中 $\mathbf{Z}\in {\mathbb{R}}^{n\times k}$，满足 ${\mathbf{Z}}^{\top }\mathbf{Z}={\mathbf{I}}_k$，且 $\mathbb{E}\Vert \mathbf{A}-\mathbf{AZ}{\mathbf{Z}}^{\top }{\Vert }_F^2\le (1+\epsilon )\Vert \mathbf{A}-{\mathbf{A}}_k{\Vert }_F^2$。
   • 时间复杂度：$O(mnk/\epsilon )$。
2. 计算采样概率：
   • 对 $\mathbf{Z}$ 的每一行 ${\mathbf{Z}}_{(i)}$，计算概率：
     $$p_i=\frac{\Vert {\mathbf{Z}}_{(i)}{\Vert }_2^2}{\Vert \mathbf{Z}{\Vert }_F^2}$$
     其中 $\Vert \mathbf{Z}{\Vert }_F^2=k$（因为 ${\mathbf{Z}}^{\top }\mathbf{Z}={\mathbf{I}}_k$）。
   • 时间复杂度：$O(nk)$。
3. 随机采样特征：
   • 选择 $r=O(k\log (k)/{\epsilon }^2)$，调用随机采样算法（Definition 5）：
     $$[\mathbf{\Omega },\mathbf{S}]=\text{RandomizedSampling}(\mathbf{Z},r)$$
     其中 $\mathbf{\Omega }\in {\mathbb{R}}^{n\times r}$ 是采样矩阵，$\mathbf{S}\in {\mathbb{R}}^{r\times r}$ 是重缩放矩阵。
   • 时间复杂度：$O(n+r)$。
4. 构建低维矩阵：
   • 计算 $\mathbf{C}=\mathbf{A}\mathbf{\Omega }\mathbf{S}\in {\mathbb{R}}^{m\times r}$，表示选择的 $r$ 个重缩放特征。
   • 时间复杂度：$O(mr)$。
5. 运行 $k$-means：
   • 在 $\mathbf{C}$ 上运行 Lloyd 算法，得到 ${\tilde{\mathbf{X}}}_{\text{opt}}$。
   • 将 ${\tilde{\mathbf{X}}}_{\text{opt}}$ 应用于 $\mathbf{A}$，计算目标函数值 $\mathcal{F}(\mathbf{A},{\tilde{\mathbf{X}}}_{\text{opt}})$。
     总时间复杂度：$O(mnk/\epsilon +k\log (k)/{\epsilon }^2\log (k\log (k)/\epsilon ))$。
     理论保证：以至少 $1-3\delta$ 的概率，$\mathcal{F}(\mathbf{A},{\tilde{\mathbf{X}}}_{\text{opt}})\le (3+\epsilon )\mathcal{F}(\mathbf{A},{\mathbf{X}}_{\text{opt}})$。

2. 随机投影特征提取（Theorem 12）

算法描述：基于随机投影的特征提取，利用 Johnson-Lindenstrauss 变换。
输入：$\mathbf{A}\in {\mathbb{R}}^{m\times n}$，$k$，$\epsilon$，$\delta$。
输出：$\mathbf{C}\in {\mathbb{R}}^{m\times r}$，${\tilde{\mathbf{X}}}_{\text{opt}}$。
步骤：

1. 生成随机投影矩阵：
   • 构造 $\mathbf{R}\in {\mathbb{R}}^{n\times r}$，其中 $r=O(k/{\epsilon }^2)$，元素为 i.i.d. $\mathcal{N}(0,1)$ 或重缩放随机符号 $\pm 1/\sqrt{r}$。
2. 投影数据：
   • 计算 $\mathbf{C}=\mathbf{AR}\in {\mathbb{R}}^{m\times r}$。
   • 时间复杂度：$O(mnr)=O(mn[{\epsilon }^{-2}k/\log (n)])$。
3. 运行 $k$-means：
   • 在 $\mathbf{C}$ 上运行 Lloyd 算法，得到 ${\tilde{\mathbf{X}}}_{\text{opt}}$。
   • 计算 $\mathcal{F}(\mathbf{A},{\tilde{\mathbf{X}}}_{\text{opt}})$。
     总时间复杂度：$O(mn[{\epsilon }^{-2}k/\log (n)])$。
     理论保证：以至少 0.97 的概率，$\mathcal{F}(\mathbf{A},{\tilde{\mathbf{X}}}_{\text{opt}})\le (2+\epsilon )\mathcal{F}(\mathbf{A},{\mathbf{X}}_{\text{opt}})$。

3. 近似 SVD 特征提取（Theorem 13）

算法描述：基于快速近似 SVD 的特征提取。
输入：$\mathbf{A}\in {\mathbb{R}}^{m\times n}$，$k$，$\epsilon$，$\delta$。
输出：$\mathbf{C}\in {\mathbb{R}}^{m\times k}$，${\tilde{\mathbf{X}}}_{\text{opt}}$。
步骤：

1. 计算近似 SVD：
   • 使用 Lemma 4 计算 $\mathbf{Z}=\text{FastFrobeniusSVD}(\mathbf{A},k,\epsilon )$。
   • 时间复杂度：$O(mnk/\epsilon )$。
2. 构建低维矩阵：
   • 计算 $\mathbf{C}=\mathbf{AZ}\in {\mathbb{R}}^{m\times k}$。
   • 时间复杂度：$O(mnk)$。
3. 运行 $k$-means：
   • 在 $\mathbf{C}$ 上运行 Lloyd 算法，得到 ${\tilde{\mathbf{X}}}_{\text{opt}}$。
   • 计算 $\mathcal{F}(\mathbf{A},{\tilde{\mathbf{X}}}_{\text{opt}})$。
     总时间复杂度：$O(mnk/\epsilon )$。
     理论保证：以至少 0.99 的概率，$\mathcal{F}(\mathbf{A},{\tilde{\mathbf{X}}}_{\text{opt}})\le (2+\epsilon )\mathcal{F}(\mathbf{A},{\mathbf{X}}_{\text{opt}})$。

---

总结

该论文通过线性代数视角和随机化技术，为 $k$-means 聚类提供了高效的降维方法，显著降低了计算复杂性，同时保持理论上的近似保证。特征选择和特征提取算法的实现过程清晰，结合了现代矩阵分解和随机投影技术。实验结果进一步验证了算法在实际数据集上的有效性，为未来研究 $(1+\epsilon )$ 近似误差的降维方法提供了方向。