LNCS-2009《Adaptive Sampling for $k$-Means Clustering》

---

核心思想

论文的核心思想是通过自适应采样（adaptive sampling）改进$k$-means聚类的初始化过程，提出一种高效的算法，生成$O(k)$个中心，以常数概率获得$k$-means问题的常数因子双标准（bi-criteria）近似解，并在这些中心中通过线性规划（LP）技术提取$k$个中心，获得常数因子近似解。传统$k$-means++算法通过$D^2$采样获得期望上$O(\log k)$近似，但其时间复杂度为$\Theta (nkd)$，对大规模数据集计算成本高昂。论文通过采样$O(k)$个中心（而非$k$个），结合LP技术，绕过了$k$-means++的归纳分析，实现了时间复杂度为$\mathcal{O}(nkd)$的常数因子近似算法。该方法不仅提升了效率，还通过理论分析和简化证明，揭示了自适应采样的潜力，适用于$k$-means及相关聚类问题（如$k$-median）。

---

目标函数

$k$-means聚类的目标是最小化量化误差（平方误差和，Sum of Squared Error, SSE），将数据集$X\subseteq {\mathbb{R}}^d$（包含$n$个点）划分为$k$个簇。给定中心集合$C\subseteq {\mathbb{R}}^d$（$|C|=k$），目标函数定义为：

$$\varphi (C)=\sum _{x\in X}\underset{c\in C}{\min }\Vert x-c{\Vert }^2$$

其中，$\Vert x-c{\Vert }^2$表示点$x$到中心$c$的平方欧几里得距离，${\min }_{c\in C}\Vert x-c{\Vert }^2$表示点$x$到最近中心的平方距离。对于子集$A\subseteq X$，定义其贡献为：

$${\varphi }_A(C)=\sum _{x\in A}\underset{c\in C}{\min }\Vert x-c{\Vert }^2$$

令$C_{OPT}$为最优$k$中心集，满足$\varphi (C_{OPT})={\varphi }_{OPT}$为最小量化误差。一个解$C$为$\alpha$近似解，若：

$$\varphi (C)\le \alpha \cdot \varphi (C_{OPT})$$

论文的目标是通过自适应采样生成中心集$S$（$|S|=O(k)$），使$\varphi (S)\le 20\varphi (C_{OPT})$（双标准近似，中心数超过$k$），并从中提取$k$个中心$C$，使$\varphi (C)\le \mathcal{O}(1)\cdot \varphi (C_{OPT})$（标准常数近似）。

此外，论文引入加权$k$-means问题，用于从$S$中选择$k$个中心。给定点集 X ′ = { μ i } X' = \{\mu_i\} X′={μi​}（每个${\mu }_i$为簇$X_i$的均值）和权重$w_i$，加权目标函数为：

$${\varphi }^{\prime }(C)=\sum _{x_i\in X^{\prime }}\underset{c\in C}{\min }{w}_i\Vert x_i-c{\Vert }^2$$

其中，${\varphi }_A^{\prime }(C)={\sum }_{x_i\in A}{\min }_{c\in C}{w}_i\Vert x_i-c{\Vert }^2$为子集$A\subseteq X^{\prime }$的贡献。

---

目标函数的优化过程

论文提出的算法通过两阶段优化目标函数$\varphi (C)$：首先通过自适应采样生成$O(k)$个中心（双标准近似），然后通过LP技术从这些中心中选择$k$个中心（标准近似）。优化过程如下：

1. 双标准近似（自适应采样）

• 输入：数据集$X\subseteq {\mathbb{R}}^d$（$n$个点），簇数$k$。
• 输出：中心集$S\subseteq X$，大小$t=\lceil 16(k+\sqrt{k})\rceil =O(k)$。
• 算法步骤：
  1. 初始化$S_0=\varnothing$。
  2. 对于$i=1$到$t$：
     • 按概率分布采样点$x$，概率与当前成本成正比：
       Pr ⁡ ( x ) ∝ ϕ { x } ( S i − 1 ) = min ⁡ c ∈ S i − 1 ∥ x − c ∥ 2 \operatorname{Pr}(x) \propto \phi_{\{x\}}(S_{i-1}) = \min_{c \in S_{i-1}} \|x - c\|^2 Pr(x)∝ϕ{x}​(Si−1​)=c∈Si−1​min​∥x−c∥2
       （当$i=1$时，均匀采样）。
     • 更新 S i ← S i − 1 ∪ { x } S_i \leftarrow S_{i-1} \cup \{x\} Si​←Si−1​∪{x}。
  3. 返回$S=S_t$。
• 时间复杂度：每次采样需计算所有点到$S_{i-1}$的距离，复杂度$\mathcal{O}(nd)$，共$t=O(k)$步，总复杂度$\mathcal{O}(nkd)$。
• 理论分析：
  • 定义“好”簇和“坏”簇：
    Good ⁡ i = { A j : ϕ A j ( S i − 1 ) ≤ 10 ϕ A j ( C O P T ) } , Bad ⁡ i = { A 1 , … , A k } ∖ Good ⁡ i \operatorname{Good}_i = \{A_j : \phi_{A_j}(S_{i-1}) \leq 10 \phi_{A_j}(C_{OPT})\}, \quad \operatorname{Bad}_i = \{A_1, \ldots, A_k\} \setminus \operatorname{Good}_i Goodi​={Aj​:ϕAj​​(Si−1​)≤10ϕAj​​(COPT​)},Badi​={A1​,…,Ak​}∖Goodi​
  • 引理1：若$\varphi (S_{i-1})>20\varphi (C_{OPT})$，则采样到${\mathrm{Bad}}_i$中某簇的概率$\ge 1/2$。
  • 引理2：对于坏簇$A\in {\mathrm{Bad}}_i$，若采样点 b ∈ B ( α ) = { x ∈ A : ∥ x − μ ∥ ≤ α r } b \in B(\alpha) = \{x \in A : \|x - \mu\| \leq \alpha r\} b∈B(α)={x∈A:∥x−μ∥≤αr}（$\mu$为$A$的均值，$r=\sqrt{{\varphi }_A(C_{OPT})/|A|}$，$\alpha \le 3$），则：
    ϕ A ( S i − 1 ∪ { b } ) ≤ 10 ϕ A ( C O P T ) \phi_A(S_{i-1} \cup \{b\}) \leq 10 \phi_A(C_{OPT}) ϕA​(Si−1​∪{b})≤10ϕA​(COPT​)
  • 引理3：$|B(\alpha )|\ge |A|(1-1/{\alpha }^2)$，表明$B(\alpha )$包含大部分点。
  • 引理4：条件概率$\mathrm{Pr}(x\in B(\alpha )\mid x\in A,A\in {\mathrm{Bad}}_i)\ge \frac{(3-\alpha {)}^2}{10}(1-1/{\alpha }^2)$。
  • 引理5：若$x\in A\in {\mathrm{Bad}}_i$，则$\mathrm{Pr}({\varphi }_A(S_i)\le 10{\varphi }_A(C_{OPT}))\ge 0.126$（取$\alpha =1.44225$）。
  • 定理1：采样$t=\lceil 16(k+\sqrt{k})\rceil$个点后，$\varphi (S)\le 20\varphi (C_{OPT})$的概率$\ge 0.03$。通过重复$\mathcal{O}(\log (1/\delta ))$次，可将成功概率提升至$1-\delta$。
  • 定理3：推广到$(4+ϵ)$近似，需采样$t=\mathcal{O}(k/ϵ\cdot \log (1/ϵ))$个点。
• 超鞅分析：
  • 定义指示变量$X_i$：若${\mathrm{Bad}}_{i+1}={\mathrm{Bad}}_i$，则$X_i=1$，否则$X_i=0$。
  • 构造超鞅$J_i={\sum }_{j=1}^i(X_j-(1-p))$，其中$p=0.063$。
  • 应用Azuma-Hoeffding不等式，证明在$t=(k+\sqrt{k})/p$步后，坏簇数降为0的概率$\ge 0.03$。

2. 标准近似（加权$k$-means）

• 输入：双标准解$S$，诱导数据集$X$的分区$X_1\cup \cdots \cup X_t$，每个簇$X_i$的均值${\mu }_i$和大小$n_i$。
• 加权$k$-means问题：
  • 点集 X ′ = { μ 1 , … , μ t } X' = \{\mu_1, \ldots, \mu_t\} X′={μ1​,…,μt​}，权重$w_i=n_i$。
  • 目标：最小化${\varphi }^{\prime }(C)={\sum }_{i=1}^t{n}_i{\min }_{c\in C}\Vert {\mu }_i-c{\Vert }^2$，$|C|=k$。
• 引理6：
  • 加权问题的最优解$C_{OPT}^{\prime }$满足：
    $${\varphi }^{\prime }(C_{OPT}^{\prime })\le 2\varphi (S)+2\varphi (C_{OPT})$$
  • 证明基于三角不等式，分析${\mu }_i$到$C_{OPT}$的距离。
• 定理4：
  • 若$C$为加权$k$-means问题的$\beta$近似解（即${\varphi }^{\prime }(C)\le \beta {\varphi }^{\prime }(C_{OPT}^{\prime })$），则：
    $$\varphi (C)\le (2\beta +1)\varphi (S)+2\beta \varphi (C_{OPT})$$
  • 结合$\varphi (S)\le 20\varphi (C_{OPT})$，若$\beta$为常数，则$\varphi (C)=\mathcal{O}(1)\cdot \varphi (C_{OPT})$。
• LP技术：
  • 使用Jain和Vazirani [JV01]或Charikar等人 [CGTS02]的LP方法解决加权$k$-means问题。
  • 由于点数为$O(k)$，LP求解时间为$\text{poly}(k,\log n)$，总时间仍为$\mathcal{O}(nkd)$。
• 理论保证：
  • 双标准解$S$以常数概率满足$\varphi (S)\le 20\varphi (C_{OPT})$。
  • 从$S$提取的$k$个中心$C$以常数概率满足$\varphi (C)=\mathcal{O}(1)\cdot \varphi (C_{OPT})$。

3. 后续$k$-means迭代（隐含）：

• 论文未明确讨论Lloyd算法，但生成的$k$个中心可作为Lloyd算法的初始中心，进一步优化$\varphi (C)$。
• 由于初始中心已接近常数近似，Lloyd迭代通常快速收敛至局部最优。

---

主要的贡献点

1. 常数因子双标准近似：
   • 提出自适应采样算法，生成$O(k)$个中心，以常数概率（$\ge 0.03$）达到$20$倍近似，可通过重复提升概率至$1-\delta$。
2. 常数因子标准近似：
   • 从$O(k)$个中心中通过LP技术提取$k$个中心，获得常数因子近似，时间复杂度为$\mathcal{O}(nkd)$，优于此前超线性时间的常数近似算法。
3. 绕过复杂归纳：
   • 相较于$k$-means++的归纳分析，论文使用超鞅和简单概率分析，证明更简洁，适用于其他问题（如$k$-median）。
4. 加权$k$-means框架：
   • 提出将双标准解转化为加权$k$-means问题，点数降至$O(k)$，显著降低LP求解复杂度。
5. 简化下界证明：
   • 提供Arthur和Vassilvitskii [AV07]中$\Omega (\log k)$期望误差下界的简化证明，揭示其误导性：尽管期望误差为$\Omega (\log k)$，自适应采样以高概率获得常数近似。
6. 广泛适用性：
   • 证明技术可推广至$k$-median和${\ell }_p$聚类问题（最小化距离的$p$次方和），因$p$次方欧几里得距离满足弱三角不等式。

---

实验结果

论文为理论分析工作，未提供实验结果或数据集的实证评估。作者专注于算法的理论保证，包括：

• 双标准近似的概率保证：采样$t=\lceil 16(k+\sqrt{k})\rceil$个中心后，$\varphi (S)\le 20\varphi (C_{OPT})$的概率$\ge 0.03$，可通过重复提升。
• 标准近似的常数因子：从$S$提取$k$个中心，获得$\mathcal{O}(1)$近似，时间复杂度$\mathcal{O}(nkd)$。
• 下界分析：通过简化证明，展示自适应采样在特定构造数据集（$k$个正则单纯形）上期望误差为$\Omega (\log k)\varphi (C_{OPT})$，但以高概率（$1-\Theta (\frac{{\delta }^2}{{\Delta }^2}k\log k)$）获得常数近似。

由于缺乏实验，论文未比较算法在真实数据集上的性能（如量化误差、运行时间）或与其他方法（如$k$-means++、$k$-means||）的对比。这可能是因为工作重点在于理论突破，而非实证验证。

---

算法的实现过程

论文提出了两个主要算法：双标准近似算法（自适应采样生成$O(k)$个中心）和标准近似算法（通过加权$k$-means选择$k$个中心）。以下详细解释实现过程，结合伪代码和说明。

1. 双标准近似算法（自适应采样）

Algorithm Bi-criteria Approximation by Adaptive Sampling
输入: 数据集 X ⊆ ℝ^d (n 个点)，簇数 k
输出: 中心集 S ⊆ X，大小 t = ⌈16(k + √k)⌉1. 初始化:S_0 ← ∅2. 循环采样 t 个中心:for i = 1 to t do// 按当前成本采样计算每个点 x 的当前成本 φ_{{x}}(S_{i-1}) = min_{c ∈ S_{i-1}} ||x - c||^2if i = 1 then均匀随机选择 x ∈ Xelse按概率 Pr(x) ∝ φ_{{x}}(S_{i-1}) 采样 xS_i ← S_{i-1} ∪ {x}end3. 返回 S ← S_t

实现细节：

• 初始化：
  • 设置空中心集$S_0$，准备存储$t=\lceil 16(k+\sqrt{k})\rceil$个中心。
• 概率采样：
  • 对于$i=1$，由于$S_0=\varnothing$，无法计算距离，采用均匀采样（$\mathrm{Pr}(x)=1/n$）。
  • 对于$i\ge 2$，计算每个点$x$到当前中心集$S_{i-1}$的最近距离平方 ϕ { x } ( S i − 1 ) = min ⁡ c ∈ S i − 1 ∥ x − c ∥ 2 \phi_{\{x\}}(S_{i-1}) = \min_{c \in S_{i-1}} \|x - c\|^2 ϕ{x}​(Si−1​)=minc∈Si−1​​∥x−c∥2。
  • 归一化概率：
    Pr ⁡ ( x ) = ϕ { x } ( S i − 1 ) ∑ x ′ ∈ X ϕ { x ′ } ( S i − 1 ) \operatorname{Pr}(x) = \frac{\phi_{\{x\}}(S_{i-1})}{\sum_{x' \in X} \phi_{\{x'\}}(S_{i-1})} Pr(x)=∑x′∈X​ϕ{x′}​(Si−1​)ϕ{x}​(Si−1​)​
  • 使用累积分布函数（CDF）或轮盘赌选择法实现概率采样：
    • 计算所有点的成本总和 ϕ ( S i − 1 ) = ∑ x ∈ X ϕ { x } ( S i − 1 ) \phi(S_{i-1}) = \sum_{x \in X} \phi_{\{x\}}(S_{i-1}) ϕ(Si−1​)=∑x∈X​ϕ{x}​(Si−1​)。
    • 生成均匀随机数$u\in [0,1]$，选择满足${\sum }_{j=1}^{m-1}\mathrm{Pr}(x_j)<u\le {\sum }_{j=1}^m\mathrm{Pr}(x_j)$的点$x_m$。
• 距离计算：
  • 每次迭代需计算$n$个点到$S_{i-1}$（最多$i-1$个中心）的距离，单点距离计算为$\mathcal{O}(d)$，总复杂度$\mathcal{O}(nd)$。
  • 共$t=O(k)$步，总复杂度$\mathcal{O}(nkd)$。
• 存储与更新：
  • 存储中心集$S_i$（最多$O(k)$个点，每个点$d$维），空间复杂度$\mathcal{O}(kd)$。
  • 更新$S_i$仅需添加新点$x$，复杂度$\mathcal{O}(1)$。
• 随机性：
  • 使用伪随机数生成器确保采样随机性，实验中可设置不同种子以重复运行。
• 数值稳定性：
  • 距离平方可能导致大数值，需使用双精度浮点数避免溢出。
  • 归一化概率时，成本总和$\varphi (S_{i-1})$可能很大，需小心处理除法。

优化技巧：

• 增量距离计算：在$S_{i-1}$到$S_i$时，仅更新每个点到新中心$x$的距离，保留到$S_{i-1}$中最近中心的距离，取最小值，降低单步复杂度。
• 并行化：计算$n$个点的距离可并行分配到多核CPU，加速采样。
• 数据结构：使用k-d树或ball树加速最近邻查询，降低距离计算复杂度（对高维数据效果有限）。

理论保证：

• 采样$t=O(k)$个点后，$S$以概率$\ge 0.03$满足$\varphi (S)\le 20\varphi (C_{OPT})$。
• 重复$\mathcal{O}(\log (1/\delta ))$次，选择$\varphi (S)$最小的解，提升成功概率至$1-\delta$。

2. 标准近似算法（加权$k$-means）

Algorithm Standard Approximation by Weighted k-Means
输入: 双标准中心集 S (|S|=t)，数据集 X，簇数 k
输出: 中心集 C，|C|=k1. 分区数据集:对每个点 x ∈ X，分配到最近的中心 s_i ∈ S，得到分区 X_1 ∪ ... ∪ X_t计算每个簇 X_i 的均值 μ_i 和大小 n_i2. 构造加权 k-means 问题:点集 X' ← {μ_1, ..., μ_t}权重 w_i ← n_i for i=1 to t3. 解决加权 k-means:使用 LP 方法（[JV01], [CGTS02]）求解：最小化 φ'(C) = ∑_{i=1}^t n_i min_{c ∈ C} ||μ_i - c||^2，|C|=k得到中心集 C4. 返回 C

实现细节：

• 分区与均值计算：
  • 对每个点$x\in X$，计算到$S$中最近中心的距离，分配到对应簇$X_i$，复杂度$\mathcal{O}(ntd)$。
  • 计算每个簇$X_i$的均值：
    $${\mu }_i=\frac{1}{n_i}\sum _{x\in X_i}x$$
    复杂度$\mathcal{O}(nd)$，共$t$个簇，总复杂度$\mathcal{O}(ntd)$。
  • 记录簇大小$n_i=|X_i|$。
• 加权$k$-means问题：
  • 构造点集 X ′ = { μ 1 , … , μ t } X' = \{\mu_1, \ldots, \mu_t\} X′={μ1​,…,μt​}（$t=O(k)$），权重$w_i=n_i$。
  • 目标函数：
    $${\varphi }^{\prime }(C)=\sum _{i=1}^t{n}_i\underset{c\in C}{\min }\Vert {\mu }_i-c{\Vert }^2$$
• LP求解：
  • 应用Jain和Vazirani [JV01]或Charikar等人 [CGTS02]的LP算法，针对平方欧几里得距离（满足弱三角不等式）。
  • LP问题规模为$O(k)$个点，变量数和约束数为$\text{poly}(k)$，求解时间为$\text{poly}(k,\log n)$。
  • LP输出$k$个中心$C$，为加权问题的$\beta$近似解（$\beta$为常数）。
• 复杂度分析：
  • 分区和均值计算：$\mathcal{O}(ntd)=\mathcal{O}(nkd)$（$t=O(k)$）。
  • LP求解：$\text{poly}(k,\log n)$，通常远小于$\mathcal{O}(nkd)$。
  • 总时间：$\mathcal{O}(nkd+\text{poly}(k,\log n))\approx \mathcal{O}(nkd)$。
• 数值稳定性：
  • 均值计算需累加高维向量，建议使用Kahan求和算法减少误差。
  • LP求解需高质量线性规划库（如Gurobi、CPLEX），确保数值稳定性。

优化技巧：

• 快速分区：使用k-d树或近似最近邻算法加速点到$S$的分配。
• LP并行化：LP求解可利用并行线性代数库（如Eigen、BLAS）。
• 预处理：若$S$质量较高，可尝试启发式方法（如贪心选择$k$个中心）代替LP，降低计算成本。

理论保证：

• 若$\varphi (S)\le 20\varphi (C_{OPT})$，则${\varphi }^{\prime }(C_{OPT}^{\prime })\le 2\varphi (S)+2\varphi (C_{OPT})$。
• 若$C$为$\beta$近似解，则：
  $$\varphi (C)\le (2\beta +1)\varphi (S)+2\beta \varphi (C_{OPT})=\mathcal{O}(1)\cdot \varphi (C_{OPT})$$

3. 与$k$-means++的对比

• $k$-means++：
  • 采样$k$个中心，每次按$D^2$采样，概率$\mathrm{Pr}(x)\propto {\min }_{c\in S_{i-1}}\Vert x-c{\Vert }^2$。
  • 期望上获得$O(\log k)$近似，时间复杂度$\Theta (nkd)$。
• 本文算法：
  • 采样$O(k)$个中心，生成双标准解，概率上获得$20$倍近似。
  • 通过加权$k$-means提取$k$个中心，获得常数近似，时间复杂度仍为$\mathcal{O}(nkd)$。
  • 避免$k$-means++的复杂归纳分析，使用超鞅和LP技术，理论更简洁。

4. 潜在扩展

• Lloyd迭代：将$C$作为Lloyd算法的初始中心，进一步优化$\varphi (C)$，通常只需少量迭代。
• 并行实现：采样和LP求解可分布式实现，适用于大规模数据集。
• 流式数据：自适应采样可改编为在线算法，处理数据流聚类。

---

总结

论文通过自适应采样和加权$k$-means框架，提出了一种高效的$k$-means聚类算法，以常数概率获得常数因子近似，时间复杂度为$\mathcal{O}(nkd)$。其核心创新在于采样$O(k)$个中心生成双标准解，并通过LP技术提取$k$个中心，绕过了$k$-means++的复杂分析。理论证明简洁，适用于多种聚类问题，但缺乏实验验证。算法实现清晰，结合概率采样和LP求解，适合理论研究和潜在的大规模应用。