K-median 算法

1 Inroduction

1.1 背景和符号

k-中值问题几十年来一直受到重视，主要是计算机科学和运筹学研究。 这个问题通常形式化如下：

给定一个度量空间X，一组客户C∈X和一组设施F∈X，使得| X | = | F∪C| = n，在F中开放k个设施，以尽量减少每个客户到最近开放设施的距离总和。

我们将距离度量定义为对于i∈{1，...，n}，j∈{1，...，n}的，使得是度量空间X中点i和j之间的距离。Kariv和Hakim [1]证明在网络中发现这样的k个中值是一个NP难题，通过减少对它的支配集问题。一个简单的brute-force算法将检查F中的每个可能的size-k子集，为每个客户计算该集合中最接近的设施，并返回最佳设置。这个蛮力算法将在时间运行，其中。因此，对这个问题的学术研究主要集中在产生好的近似算法上。对于给定的算法，近似比率被定义为算法输出成本与最优成本之间可能的最差可能比率。然而，对于大多数问题实例，我们并不知道实际的最优成本，因此我们计算近似比率为算法返回的总成本除以1.2节讨论的宽松线性程序的最优值。 Jain等人[2]证明k-median问题是1+ 2/e≈“1.736” - 难 以在度量空间内逼近。我们注意到，在本文中，我们所有的距离度量都满足度量空间的属性。

1.2 与其他NP难问题的关系

k-median问题与设施位置问题（FLP）有许多相似之处。 在这个问题中，我们给出包含有客户C，设施F，使用设施i的客户j的成本dij以及与开放设施i相关联的成本fi的度量空间X. 有两个设施位置问题的变体：无容量和有容量。 在无容量设施位置问题（UFLP）中，任何数量的客户都可以使用任何给定的设施。 在容量化设施选址问题（CFLP）中，我们定义了一组变量V（i），使得使用设施i的客户数量必须小于或等于V（i）。 对于本文的其余部分，我们只关心自己的无容量设备位置问题。 在1966年，Balinski [3]为UFLP定义了以下整数程序。 这里，yi表示设施i是否打开，xij表示客户端j是否连接到设施i。

通过将整数约束xij∈{0,1}和yi∈{0,1}分别转换为xij≥0和yi≥0，可以放宽整数程序。UFLP的另一个变体是将开放设施的数量限制在k（UKFLP）。 我们可以通过添加yi值的总和必须小于或等于k的约束条件来构造整数（以及最终的线性松弛）程序。 通过设置fi = 0，我们可以很容易地将UKFLP连接到一个无容量限制的k-中值问题上。所得到的整数程序定义如下：

通过将整数约束xij∈{0,1}和yi∈{0,1}分别转换为xij≥0和yi≥0，可以放宽整数程序。 就本文而言，我们将无容量限制的k中介问题简称为k中值问题。 如果X = F = C，则k-中值算法简单地对度量空间中的点进行聚类

1.3 现有的理论研究

现有关于k-中值近似算法的大量研究。 前向贪婪算法通过选择下一个未打开的设施来迭代打开设施，以使上面定义的成本函数最小化。 Chrobak等人。 [4]证明了反向贪心算法（一种以所有设备开始打开和关闭最少增加成本函数的设施）的算法产生Ω（logn / loglogn）和O（logn）之间的近似比。 另一个经常使用的算法，围绕中心点划分进行分区（PAM），迭代计算交换  medoids（(网络)中心点） 的代价并选择交换，最大限度地降低总成本（如果有这样的交换仍然存在）。这是一种本地搜索算法。 在性能保证方面，2003年Arya等人 [5]证明了这种单次交换的局部搜索具有最坏情况下的近似比为5.在一般情况下，其中至多p设施可以一次交换，近似比为3 + 2 / p。

Charikar等人[6]为k-中值问题提供了第一个常数因子近似算法，实现了20/3的比率。这依赖于上一节中介绍的整数程序的LP松弛的舍入方案。此后的几篇论文使用了线性规划技术。我们在这里简单地考虑一些论文。 Jain和Vazirani [7]为UFLP获得了一个3-近似的主要算法。通过将UFL问题表达为k-中值问题的拉格朗日松弛和使用3-近似UFL算法作为子程序，这为k-中值问题提供了6-近似算法。在他现有的论文的基础上，Jain等人。 [8]将UFLP问题的这个比例提高到2，因此对于k-中值问题，这个比率提高到4。 2012年，Charikar和Li [9]开发了一种新的基于LP松弛的算法，该算法的逼近比率为3.25（1 +δ）。虽然比上述的3 + 2 / p局部搜索近似更差，但运行时间从减少到。 2013年，Li和Svensson [10]开发了一种近似算法，实现了的近似保证，最终改进了Arya等人的比率。

1.4现有实验研究

计算机科学和运筹学研究都对k-中值问题进行了各种实证研究。然而，就我们的知识（以及Charikar和Li的知识）而言，Charikar 2012 LP算法以前从未经验地实施或研究;这是我们论文的主要动机。下面，我们简要讨论最相关的先前的经验研究。 2011年，Nagarajan和Williamson [11]对各种k-中值近似算法进行了实证研究，包括Arya等人的局部搜索算法，Charikar的1999 LP舍入算法和Jain等人的4-近似算法。然而，这个实验研究在Charikar的2012 LP算法的开发之前，这是我们论文的主要动机。此外，由于多交换版本的运行时间长，Nagarajan和Williamson仅分析了单次交换的局部搜索（保证近似比率为3 + 2 / p = 5）。关于数据，他们在OR-Library的40个数据集和另外三个大小为100,150和316个节点的数据集上测试了他们的算法。在所有这些数据集中，F = C。Nagarajan和Williamson得出结论：局部搜索和Charikar的1999年LP舍入算法比Jain-Vazirani算法表现更好，并且Charikar算法通常表现比本地搜索更好。但是，对于大多数使用的数据集，LP解决方案通常是积分或非常接近积分。从所呈现的结果来看，这些算法如何在LP不返回整数解的数据集上执行并不清楚。

2013年，Bhattacharya等人。 [12]经验地分析了鲁棒k-中值问题，该问题试图将所有客户的最坏连接成本（而不是总成本）降到最小。虽然他们的论文专门针对强大的k-中值问题（与原来的k-中值问题相对），但其用于生成用于测试的设施和客户的随机度量空间的方法可以扩展到k-中值问题的任何变体。这些方法在4.2中进一步讨论。关于LP分析，2014年，Awasthi等人[13]对k均值和k中值问题进行了LP松弛完整性的理论和实证研究。他们发现，对于k-中值问题，LP解决方案比k-均值问题更频繁地积分。事实上，对于至少一些常数c和任意k的聚类分离，如果n足够大（尽管“足够大”没有根据问题的参数进行精确定义），则k中值LP解将是整数。这些结论激发了用基于LP松弛的近似算法研究大量k-中值问题实例。

我们的研究基于Nagarajan和Williamson的分析，通过对使用随机生成的数据集的各种k-median算法进行实证研究，如[12]中所述，它提供了一个更全面的算法性能图。 此外，我们将大部分分析重点放在Charikar的2012算法上，既未实施也未经过实证研究。 考虑到Charikar的2012算法具有至多3.25的近似比（理想情况下低于3并具有适当的参数调整），该算法与Arya等人的局部搜索相竞争，其提供了最知名的近似保证，直到Li-Svensson 因此确定Charikar算法在实践中的表现方式是非常重要的。

2 实现的算法

在本节中，我们将探索我们更详细地实现的算法。 一般来说，设C为城市集合，F为可能的设施位置集合，，。 我们让成本（S）为每个城市到S中最近设施的成本总和。

2.1 贪婪算法

2.1.1 向前贪婪

对于k中值问题，前向贪心算法在步骤t维持一个集合，并设置以使成本最小化，直到| S | = k。 这是一个在实践中快速执行的简单算法，但它至少是一个针对该问题的n-近似算法[4]。

2.1.2 反向贪婪

反向贪心算法首先在每个位置放置一个设施，反复从该组中点去除设施，使得成本增加最少，直到剩下k个设施。 [4]证明反向贪心算法在Ω（logn / loglogn）和O（logn）之间具有近似比率。

2.2 本地搜索

以最简单的形式，给定一组中位数S，局部搜索通过反复交换S和F-S中的一个中位数使得新的耗费更小从而得到一个k-median的解。在更一般的情况下，本地搜索具有可调整的参数p，允许一次更换多达p个中间值的集合。 初始化阶段通常通过将前向贪婪算法应用于问题来完成。 在这种情况下，该算法是完全确定性的。 如[5]中所示的局部搜索如下，其中B（S）是通过S进行p交换得到的可能集合的集合。

[5]证明了带p交换的局部搜索是一种3 + 2 / p近似算法，它在简单单交换情况下给出近似比率5。 算法中的每个交换都以时间运行。 他们还证明，只要能够提供至少一定数量的增益的交换被接受，步骤的数量也将是多项式的。 在实践中有很多实现技巧和数据结构来加速本地搜索，例如缓存哪些城市与每个设施相关联，但我们不会详细探讨这些技术。 值得注意的是，在统计和操作研究（OR）社区中，使用单一互换的局部搜索是一种长期存在的技术，它被称为分区中心（PAM）[14]。

2.3 Jain Vazirani

Jain和Vazirani [7]利用primal-dual 模式为无容量限制设施的位置（UFL）问题开发了一种3近似算法。 此外，他们表明UFL问题是k-中值问题的拉格朗日松弛。 因此，Jain和Vazirani能够将他们的primal - dual 方法用作k-中值问题的6-近似算法的子程序。 此外，他们的k-median算法具有非常快的运行时间：O（mlogm（L + log（n））），其中n和m是城市和设施的二分图中的顶点和边的数量。

该算法的直觉如下：每个城市都有一个贡献和一个连接状态。 随着时间t从0增加，每个城市都愿意考虑连接到其位置的半径t内的设施。 考虑距离贡献城市的距离r≤t的设施c。 如果f打开，则c连接到f并在未来所有时间切换到连接状态。 如果f未打开，则c有助于打开t-r。 当贡献总额大于设施开放成本时，设施开放。 当一个城市链接时，其现有的所有贡献仍然存在，但不再增加。 这是该算法的第一阶段。

第二阶段的重点是选择暂时开放城市的一个子集作为开放。在[8]中，Jain et al简化算法以消除清理阶段的需要。在第二阶段，如果dist（i，j）小于 i 链接到设施的时间，那么城市i和设施j之间的任何边i，j都被视为“特殊”边。该算法从这些边形成一个图T，找到T2，它有一个边i，j，如果在T中i和j之间存在长度为1或2的路径，并且找到H作为T2的子图，T2包括暂时开放设施（即仅保留对应于临时开放设施的顶点从阶段1开始）。这意味着对于任何一个城市c而言，所有临时开放的设施如果从c收到正面贡献，将在H中形成一个派系。最后，从H中选择一个最大独立子集并返回它。这意味着在最终的一套开放式设施中，每个城市最多只能贡献一个正数。虽然发现最大独立子集一般是一个NP难题，但在这种情况下，贪婪地添加不与当前独立集合冲突的设施在这种情况下产生最大独立集合[15]。

Jain和Vazirani在设施开放成本上使用二分查找来解决k-中值问题。 如果每个设施的开放成本为0，则上述两个阶段将开放nf设施。 如果设施开放成本任意高，则将打开一个设施。 该算法在这两个开放成本之间执行二分搜索以查找打开k个设施的成本。 如果不存在这样的成本，则算法采用两个大小为k1，k2的设施子集，使得k1 <k <k2。 该算法提供了一个四舍五入方案，以找到k个设施打开这两个子集。 舍入方案将成本最多增加2倍。出于分析的目的，我们没有实施舍入阶段;相反，我们仅比较了通过二分搜索成功找到 整数解 的实例的性能。

2.4 Charikar 2012 Rounding

Charikar和Li在[9]中开发了一种依赖LP-rounding方法，在运行时间中获得了保证的3.25（1+δ）近似，其中n = |F∪C|。 该算法依赖于杨氏近似的k-中值LP [16]，这产生了（1 +δ）k近似。 在我们的实现中，我们不使用Young的近似算法，所以运行时间稍慢，但保证界限不依赖于δ

2.4.1 预备步骤和符号

首先，该算法初始化1.3节所述的线性程序并生成分数解。 在下面的解释中，我们将yi表示为在LP解决方案中开放的设施i的部分组件，xij表示设施i和客户j之间开放的连接的“数量”。 对于C中的每个客户j，定义Fj为xij> 0的设施集合。这些是j在分数解中连接的设施。 因此分数解中j的连接代价为。 这是LP放松中客户j的成本。 定义为每个客户j的成本。 另外，定义表示具有严格小于r到j的距离的一组设施的函数B（j，r）。 所有yi = 0的设施都从F中移除，因为它们并不完全开放。 舍入计划分四步进行： 

2.4.2 过滤阶段

该算法通过选择子集C'⊆C来减少被考虑的客户端的数量,C'具有两个属性：（1）C'中的客户端彼此远离;（2）C \ C'中的客户端靠近C'中的某个客户端。 更具体地说，对于所有的客户端j, j'在C'中 并且并且对于C \ C'中的每个j，在C'中存在客户端j'，使得

2.4.3 绑定阶段

捆绑阶段尝试将可能的设施位置分成一组给C'中的客户端。 由于这些客户彼此之间距离很远，因此每个捆绑包在最终解决方案中都有很大的可能性产生开放式设施。 把对客户端j的一个捆绑组叫做Uj并将Rj定义为到下一个最接近的C'中的客户端的距离的一半。 对于每个客户j定义每个至少在中出现过一次的设施都被放入一个捆绑中。 出现在多个中的设施则分配给最近的一个客户端。

2.4.4 匹配阶段

匹配阶段将两个捆绑(包)Uj和Uj'分组，以便它们可以从联合分布中采样。 匹配的捆绑包应该靠近在一起，以便使用贪婪算法，其中客户j和j'按照距离的递增顺序匹配。 j中的每个客户端只匹配一次。 如果C'中的客户端数量为奇数，则与其他所有客户端最远的客户端不参与匹配。

2.4.5 取样阶段

定义函数采样阶段总结为如下算法4

采样阶段结束后，返回k个开放设施。 由于该算法仅返回期望中的k个设施，所以我们一直运行采样阶段，直到有k个设施开放。 或者，可以实现一个复杂的树结构，它依赖于对分布进行采样。

2.4.6实际改进和变体
多次迭代采样
该算法的唯一随机部分发生在采样阶段。考虑线性程序不返回整数解的情况。有多种可能的舍入方案可以将分数解转化为整数解。这些舍入计划可能会有不同的成本。在本文中，我们没有接受第一种方法，而是提出了一种算法来运行舍入方案一千次，并返回产生最低成本的一组设施。尽管多次采样不能提高保证的近似比率，但实际上它的结果至少与传统的Charikar算法一样好。我们在5.1节中分析这种算法的变体。在本文期间，我们将上面记录的传统Charikar和Li算法称为“Charikar 2012”，将多重迭代采样算法称为“Charikar多重采样”。
更改参数
在2.4.3节中，介绍了常数1.5。 1.5的常数提供3.25的近似比。在文章中简要地提到，如果使用1.0的常数，分析更简单。然而，保证的近似比降低到4倍。在与Charikar教授的讨论中，他建议∞的值会产生最佳的近似比率，可能会跳动（当时）十倍长的屏障。如果该常数等于∞，那么束将完全由xij> 0的设施组成。在5.2节中，我们探讨改变这个常数的结果;具体而言，我们研究了值1,1.5,10,100和∞。