基于最小聚類求解k-means問題算法

王守強(qiáng)，朱大銘

（1. 山東交通學(xué)院 信息工程系，山東 濟(jì)南 250023；2. 山東大學(xué) 計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，山東 濟(jì)南 250100）

1 引言

給定 d維空間中的點(diǎn)集P，k-means聚類問題要求在空間中選取k個中心點(diǎn)，使P中點(diǎn)與其距離最近的中心點(diǎn)的距離平方和最小。形式化描述為

實(shí)例：點(diǎn)集 P ∈Rd，正整數(shù)k∈Z+。

k-means問題相當(dāng)于在d維空間中計(jì)算k個中心點(diǎn)，以中心點(diǎn)為核心將給定點(diǎn)集P劃分為k個子集，優(yōu)化目標(biāo)為給定點(diǎn)到其所屬子集中心點(diǎn)的距離平方和最小。該問題是NP-Hard問題[1]。其教科書算法為Lloyd給出的啟發(fā)式算法[2,3]，Lloyd算法簡單而容易實(shí)現(xiàn)，但運(yùn)行結(jié)果依賴于初始值，算法無法保證一個確切的求解近似度。Kanungo[4]等采用局部搜索技術(shù)給出k-means問題（9+ε）近似度算法。在執(zhí)行算法前需要對點(diǎn)集空間結(jié)構(gòu)進(jìn)行劃分求得一個候選中心點(diǎn)集[5]，候選中心點(diǎn)集的求解十分復(fù)雜，算法顯得不夠?qū)嵱?。Song等[6]進(jìn)一步證明，如果將給定點(diǎn)集P作為中心點(diǎn)的候選點(diǎn)集，對候選點(diǎn)集執(zhí)行局部搜索，可使算法的近似度達(dá)到O(1)。2004年Kumar[7]等人給出求解k-means問題的(1+ε)隨機(jī)近似算法，算法的時間復(fù)雜度為

2) 對于m最小聚類問題，改進(jìn)了Kumar提出的 k-means問題的(1+ε)隨機(jī)近似算法，將 Kumar給出算法的時間復(fù)雜度由 O ( 2(k/ε)O(1)dn)改進(jìn)為證明改進(jìn)算法求到(1+ε)近似解的成功概率至少為如果該算法運(yùn)行次，則以接近于 1的概率求到該問題的(1+ε)近似解。

3) 設(shè)計(jì)出求解k-means問題的局部搜索隨機(jī)算法。從k個初始中心點(diǎn)的選取以及生成候選中心點(diǎn)集2個方面分別對Song的局部搜索算法提出新的隨機(jī)策略。Song局部搜索算法的時間復(fù)雜度為本文算法的時間復(fù)雜度為O(nk3dln(k)2)/α)。實(shí)驗(yàn)結(jié)果表明：新算法所求解的精度和算法的運(yùn)算時間均優(yōu)于 Song給出的局部搜索算法。

2 符號標(biāo)記與基本結(jié)論

定義 2對于點(diǎn)集 P的 m最小聚類問題，稱α=km/|P|為該問題的最小聚類參數(shù)。

定義3設(shè)是d維空間中的k個點(diǎn)，若存在 l個實(shí)數(shù)滿足：

d對于 d維空間中某點(diǎn) x ∈ R ，如果成立，則稱 x為的凸組合點(diǎn)。

引理1[8]點(diǎn)集P的質(zhì)心點(diǎn)為P的1-means問題最優(yōu)解的中心點(diǎn)。

Inaba等[7]指出在歐氏空間中，只需從P中隨機(jī)選取部分點(diǎn)，計(jì)算取樣點(diǎn)的質(zhì)心點(diǎn)，則該質(zhì)心點(diǎn)以較大的概率成為P的 ε 近似質(zhì)心點(diǎn)。該定理描述如下。

引理2[7]給定點(diǎn)集P，從P中均勻地隨機(jī)選取m個點(diǎn)，設(shè)S為取樣點(diǎn)的集合，m=|S|， c(S)為S的質(zhì)心點(diǎn)。存在δ (0＜δ＜1)，使下述不等式成立的概率至少為1-δ。

引理3[8]給定點(diǎn)集P，對于任一個點(diǎn)x ∈ Rd，

3 基于取樣技術(shù)k-means近似算法

3.1 近似度期望值為2的算法

記 ?k2(P )為給定實(shí)例點(diǎn)集 P的最優(yōu)解的解值。在給出該算法之前，首先探討當(dāng)中心點(diǎn)集包含每個最優(yōu)子集 Pi(1≤i≤k)中的一個點(diǎn)時，算法近似度的期望值。

定理 1設(shè) P所對應(yīng)的最優(yōu)聚類劃分子集為如果從每個 Pi中均勻地隨機(jī)選取一個點(diǎn)ci，以這k個點(diǎn)作為中心點(diǎn)求到的解值記為?，則

證明設(shè) P的最優(yōu)解所對應(yīng)的中心點(diǎn)為從每個 Pi中均勻地隨機(jī)選取一個點(diǎn) ci，由所有 ci構(gòu)成中心點(diǎn)集合記為不失一般性，假設(shè) ci∈Pi，根據(jù)已知條件，Pi中任一點(diǎn)被選作中心點(diǎn) ci的概率為針對點(diǎn) x ∈Pi,定義{distance(x,ci)}，D(x)實(shí)際表示為點(diǎn)x與集合 C - { ci}中最近點(diǎn)的距離，則由此可得：

(根據(jù)引理3)

給定正整數(shù)m，如果P為m最小聚類問題。下面證明只需從P中隨機(jī)選取一個樣本點(diǎn)集S，則S以較高的概率滿足它包含每個最優(yōu)子集 Pi中至少一個點(diǎn)，即

定理2如果從P中均勻地隨機(jī)選取k ln(2k)個點(diǎn)，記取樣點(diǎn)集為S，那么對于所有P，αi均成立的概率至少為1，即?iPr[|S∩P|2i

證明不失一般性，設(shè)|記定義則 nis的期望值根據(jù)Chernoff Bounds不等式得：

其中，0＜λ＜1。定義事件Ai為滿足足 nis小于λ|S |nni的事件，則：

對每個最優(yōu)子集Pi，如果均成立，枚舉S中所有k個點(diǎn)的子集，則至少存在一個子集滿足：k個點(diǎn)分別來自于最優(yōu)子集 P1,P2,… ,Pk中的一個點(diǎn)。以這k個點(diǎn)作為中心點(diǎn)，根據(jù)定理1知，所求解的近似性能比的期望值至多為 2。據(jù)此，給出下述算法。

算法1k-means問題近似算法

輸入：點(diǎn)集P，參數(shù)α。

2) 從S中任取k個點(diǎn)，以這k個點(diǎn)作為中心點(diǎn)，對P進(jìn)行一次劃分，計(jì)算劃分后的值?。

3) 重復(fù)2)的操作，枚舉所有可能k個點(diǎn)，選擇? 最小的值所對應(yīng)的k個中心點(diǎn)作為算法最終解。

4) 返回所求的中心點(diǎn)。

定理3算法1至少以1/2的概率求到近似度期望值為2的解。

證明根據(jù)定理1及定理2，該定理結(jié)論顯然成立。

算法在求得k個中心點(diǎn)時，通過枚舉取樣點(diǎn)集中所有k個點(diǎn)的子集找出最小解。枚舉子集的個數(shù)為 C|kS|，由于將S值代入該式可得枚舉子集個數(shù)為由于每次將實(shí)例點(diǎn)集P中的點(diǎn)分配給k個中心點(diǎn)時，時間復(fù)雜度為 O(nkd)。因此，算法 1的時間復(fù)雜度為

3.2 k-means聚類(1+ε)近似算法

定理4對滿足m最小聚類問題實(shí)例點(diǎn)集P，如果從P中均勻地隨機(jī)選取8 k ln(4k)個點(diǎn)，記取樣εα點(diǎn)集為S，則S包含每個最優(yōu)子集Pi中至少2/ε個點(diǎn)的概率大于等于3/4。

證明證明過程類似于定理 2。記 ni=|Pi|，定義則 nis的期望值根據(jù)Chernoff Bounds不等式得：

定義事件Ai為 nis小于λ|S |nni的事件，則：

由于|Pk|≥m，根據(jù)最小聚類參數(shù)定義α=km/|P|可得成立。

引理4給定點(diǎn)集P，設(shè)均屬于P中同一個最優(yōu)子集Pi中的點(diǎn)，如果x∈P并且 x是由所組成凸組合的點(diǎn)，則x∈Pi

證明由于 x是屬于中凸組合內(nèi)的點(diǎn)，根據(jù)凸組合定義，存在l個大于等于0的數(shù)值滿 足 ：其 中設(shè) ci為 Pi的質(zhì)心點(diǎn)，則：

算法2k-means問題(1+ε)近似算法

輸入：點(diǎn)集P，參數(shù)ε，每個子集最少聚類個數(shù)m。

1) 計(jì)算參數(shù)α＝(mk)/|P|，置C=φ。

3) 調(diào)用函數(shù) Cost=Online-k-means(S, C, ε)。

4) 返回代價(jià)Cost及集合C。

函數(shù) Online-k-means(·)的實(shí)現(xiàn)采用遞歸方法，從i＝1開始，首先從S中枚舉所有2/ε個點(diǎn)的子集，則至少存在一個子集 Si′滿足 Si′?Si，由 Si′求出 Si′′并進(jìn)一步計(jì)算 Si′∪ Si′的質(zhì)心點(diǎn) ci′，然后從S中刪除屬于集合 Si′∪ Si′中的所有點(diǎn)。刪除 Si′∪ Si′中的所有點(diǎn)的目的在于當(dāng)求解下一個子集時，可以減少枚舉子集的個數(shù)，從而提高算法的運(yùn)行效率。函數(shù)的實(shí)現(xiàn)過程描述如下。

函數(shù) Online-k-means(S, C, ε)

輸入：樣本點(diǎn)集S，已求中心點(diǎn)集C，參數(shù)ε 值。

1) If |C|=k then

2) 以C作為中心點(diǎn)，計(jì)算解值Sum。

3) if Sum＜MinCost then MinCost=Sum，保存C和MinCost。

4) If |C|＜k and |S|＜2/ε 返回。

5) Repeat。

6) 從S中取2/ε個點(diǎn),設(shè)點(diǎn)的集合為S′。

7) 從S-S′中找出滿足S′凸組合中的點(diǎn)，設(shè)點(diǎn)集為 S ''。

9) Online-k-means(S,C, ε)。

10) Until 窮舉完畢S中所有2/ε個點(diǎn)的子集。

11) 返回MinCost的解值。

引理 5如果對每個最優(yōu)子集 Pi，均成立，則函數(shù) Online-k-means(·)能夠從 S中求出一個子集 Si′，滿足成立。

證明記由上述第 6)步知，算法是從 S中枚舉所有 2/ε個點(diǎn)，則枚舉子集中必存在一個子集S′屬于Si。上述第7）步求出S中滿足S′的凸組合的點(diǎn)集 S′，根據(jù)引理 4，，因此并且成立。

引理 6對每個最優(yōu)子集 Pi，如果成立，則函數(shù) Online-k-means(·)求出 k-means問題(1+ε)近似解的成功概率不小于(1/2)k。

證明根據(jù)引理5，對每個最優(yōu)子集Pi，函數(shù)Online-k-means(·)從 S中求出一個子集 Si′，Si′滿足成立。以 Si′的質(zhì)心點(diǎn)作為 Pi的中心點(diǎn)，根據(jù)引理2，該質(zhì)心點(diǎn)至少以1/2的概率滿足Pi的 1-means問題的(1+ε)近似解。因此，算法求解 k-means問題(1+ε)近似解的成功概率至少為(1/2)k。

定理5對于滿足m最小聚類的點(diǎn)集P，算法2至少以2的概率求出該問題的(1+ε)近似算法，算

2k+2法的時間復(fù)雜度為

證明根據(jù)定理4，對任意最優(yōu)子集成立的概率不小于 3/4。當(dāng)條件成立時，由引理 6知，函數(shù) Online-k-means(·)求解k-means問題(1+ε)近似解的成功概率至少為(1/2)k。因此，算法2的成功概率為

記 r=2/ε，T(|S|)代表求函數(shù) Online-k-means(·)中以 S作為實(shí)例點(diǎn)集的枚舉個數(shù)，則：由此可得T(|S|)至多為

將|S|及r值代入上式可得：枚舉所有可能k個中心點(diǎn)的子集個數(shù)至多為因此整個算法的時間復(fù)雜度為

3.3 局部搜索隨機(jī)算法

基于文獻(xiàn)[6]的局部搜索算法，本文提出k-means問題的局部搜索的隨機(jī)算法，其隨機(jī)策略主要體現(xiàn)如下。

1) 文獻(xiàn)[6]候選中心點(diǎn)集選自給定實(shí)例點(diǎn)集P，新的隨機(jī)算法則以 P的一個取樣子集作為候選中心點(diǎn)集，S滿足以較高的概率包含每個最優(yōu)子集至少一個點(diǎn)。

2) 在 k個初始中心點(diǎn)的選取方面。文獻(xiàn)[6]是從候選中心點(diǎn)集中任取k個點(diǎn)，新的隨機(jī)算法則采用非均勻地隨機(jī)選取策略從P中選取k個點(diǎn)，使得這k個點(diǎn)盡可能地分屬于k個不同最優(yōu)子集中的點(diǎn)。

初始中心點(diǎn)選取算法實(shí)現(xiàn)描述如下。

算法3k個初始中心點(diǎn)選取算法

輸入：點(diǎn)集P

1) 從P中隨機(jī)選取2個點(diǎn) x1、 x2，選擇概率

2) While 選取點(diǎn)數(shù)≤k。

3) 從P中隨機(jī)選取一點(diǎn)xi(i≥3)，選取概率為

4) End While。

5) 返回選取點(diǎn){x1,…, xk}。

算法3首先從P中隨機(jī)選取2個點(diǎn) S = { x1,x2}，遵循兩點(diǎn)距離越大，被選取概率越大的原則。再依次隨機(jī)選取點(diǎn) x3,… ,xk，選取第i(i＞2)個點(diǎn)時，遵循一個點(diǎn)與已選擇的點(diǎn)距離平方和越大，則該點(diǎn)被選取概率越大的原則。因此算法規(guī)定第一次選擇2個點(diǎn){x1,x2}的概率表達(dá)式為設(shè)規(guī)定選擇第i個點(diǎn)xi的概率表達(dá)式為算法4k-means局部搜索隨機(jī)算法

輸入：點(diǎn)集P，最小聚類參數(shù)α。

3) Repeat。

4) 以C為中心點(diǎn)，對S作一次劃分，劃分子集為{S1,S2,…,Sk}。

5) For i=1 to k。

6) 對于Si中每一個點(diǎn)x，以C′為中心點(diǎn),計(jì)算給定點(diǎn)集P的k-means解值?′,如果?′＜?,置

7) Next i。

9) 返回C以及解值?。

與文獻(xiàn)[6]局部搜索算法相比，選取k個初始中心點(diǎn)時，算法4中的第2)步不是從給定實(shí)例點(diǎn)集P中任取k個點(diǎn)，而是調(diào)用算法3從P中非均勻地隨機(jī)選取k個初始中心點(diǎn)。算法4中的第3)～8)步局部搜索時以取樣子集S作為候選中心點(diǎn)集。根據(jù)定理2，該取樣子集能夠較好地代表P。在執(zhí)行中心點(diǎn)交換時，取樣子集能夠減少交換次數(shù)，從而達(dá)到降底算法時間復(fù)雜度的目的。

根據(jù)算法4第6)步，算法每次迭代時，所求P的值下降(1-ε/k)倍。記?k2(P)為給定實(shí)例點(diǎn)集 P的最優(yōu)解值，?2(P,C)為算法所求最終解值，?2(P,C0)為算法初始中心點(diǎn)所對應(yīng)的值，算法的迭代次數(shù)記為t。則：

由此可得：算法迭代次數(shù)

4 實(shí)驗(yàn)結(jié)果

本文選用 Iris、RuspIni、Spath Postal Zone Data、Cloud 和SPAM 數(shù)據(jù)集測試本文相關(guān)算法。選用數(shù)據(jù)集 Iris、RuspIni、Spath Postal Zone Data來驗(yàn)證算法1的運(yùn)算結(jié)果；選用UCI中2個高維數(shù)據(jù)集Cloud以及SPAM測試局部搜索隨機(jī)算法的運(yùn)算性能。表1列出5個數(shù)據(jù)集的基本屬性。

表1 選用數(shù)據(jù)集說明

4.1 近似度算法實(shí)驗(yàn)

根據(jù)算法 1，在執(zhí)行算法前，需要給出算法的最小聚類參數(shù)α。受篇幅所限、僅給出α=0.5作為參數(shù)值，對不同的k值進(jìn)行實(shí)驗(yàn)結(jié)果。表2中最優(yōu)值一列來自于文獻(xiàn)[10～12]。由于算法的隨機(jī)性，表2中實(shí)驗(yàn)結(jié)果取自算法運(yùn)算 20次后所求解值的最小值，實(shí)驗(yàn)結(jié)果參見表2。

表2 算法1實(shí)驗(yàn)結(jié)果(α＝0.5)

在表2中，近似度一列值等于實(shí)驗(yàn)結(jié)果與最優(yōu)值的比值。由表2近似度一列可知：本文實(shí)驗(yàn)所得到的算法近似度均小于2。通過表2實(shí)驗(yàn)結(jié)果可以看出：對隨機(jī)取樣點(diǎn)集S，枚舉S中所有可能k個點(diǎn)子集作為中心點(diǎn)，則以較高概率獲得近似度期望值為2的解。1的常數(shù)，因此，迭代此數(shù)可簡化為O(kln(k))。由算法4的第7)步知，每次迭代時，需要涉及|S|個中心點(diǎn)交換；每次交換后，算法重新計(jì)算k-means解值的時間復(fù)雜度為O(nkd)。所以，整個算法4的時

4.2 改進(jìn)局部搜索算法實(shí)驗(yàn)

為檢驗(yàn)改進(jìn)后局部搜索算法性能，本文將數(shù)據(jù)集Cloud以及SPAM分別應(yīng)用在文獻(xiàn)[6]局部搜索算法以及本文算法4上進(jìn)行測試。表3給出2種算法的實(shí)驗(yàn)結(jié)果，由于初始中心點(diǎn)選取隨機(jī)性，將程序執(zhí)行20次。表2中運(yùn)算值所對應(yīng)的兩列是指兩算法運(yùn)算20次后所求的最小值，而表2中運(yùn)算時間所對應(yīng)的2列值則指程序20次運(yùn)算后的平均時間。

為便于比較和描述，表中 2運(yùn)算值=算法實(shí)際運(yùn)行結(jié)果/點(diǎn)集個數(shù)。

表3 局部搜索算法實(shí)驗(yàn)結(jié)果

相對于文獻(xiàn)[6]，定義算法 4的改進(jìn)值=[1-(算法4)/(文獻(xiàn)[6]算法)]×100%。對于Cloud數(shù)據(jù)集，當(dāng)k=10、25、50時，算法2第2）步的運(yùn)算時間分別改進(jìn)了96.75%、93.32%以及89.21%，因此算法時間得到顯著提高。除運(yùn)算時間外，局部搜索隨機(jī)算法的運(yùn)算值也均小于文獻(xiàn)[6]的算法。對于SPAM數(shù)據(jù)，當(dāng) k=10、25、50時，算法所求的解值分別改進(jìn)了51.91%、87.38%和91.63%，而算法的運(yùn)算時間則分別改進(jìn)了98.98%、97.98%以及96.39%。由表3可以看出，針對給定的實(shí)例點(diǎn)集P，采用局部搜索隨機(jī)算法，無論是算法所求值的精度還是算法的運(yùn)算時間均優(yōu)于文[6]所給出的算法。

5 結(jié)束語

本文分析了基于最小聚類k-means問題的隨機(jī)近似算法，利用取樣技術(shù)，給出了該問題近似度期望值為2的隨機(jī)算法。同時探討了該子問題的(1+ε)近似算法的求解方案，將 Kumar所給出的(1+ε)近似方案的時間復(fù)雜度進(jìn)行了改進(jìn)，并分析了算法的成功概率。利用取樣技術(shù)，本文設(shè)計(jì)局部搜索隨機(jī)算法。最后，選取了部分實(shí)驗(yàn)數(shù)據(jù)，對算法近似度以及局部搜索隨機(jī)算法進(jìn)行了驗(yàn)證。但本文還有幾個問題需要進(jìn)行探討。1) 本文近似算法是通過枚舉樣本點(diǎn)集中部分點(diǎn)的求到的，顯然算法的時間復(fù)雜度高，能否不需要進(jìn)行組合，而是通過某種策略直解在取樣子集中找出滿足給定條件的某些點(diǎn)？2) 該算法能否進(jìn)一步提高成功概率？3) 如何找出k個初始中心點(diǎn)，使得這k個點(diǎn)以較高的概率分別來自于k個不同最優(yōu)子集中的一個點(diǎn)。

[1] DRINEAS P, FRIEZE A, KANNAN R, et al. Clustering large graphs via the singular value decomposition[J]. Machine Learning, 2004,56(1-3)∶9-33.

[2] MACQUEEN J B. Some methods for classification and analysis of multivariate observations[A]. Proceedings of the 5th Berkeley Symposium on Mathematical Statistics and Probability[C]. California, USA,1967. 281-297.

[3] LLOYD S P. Least squares quantization in PCM [J]. IEEE Transactions on Information Theory, 1982, 28(2)∶129-137.

[4] KANUNGO T, MOUNT D M, NETANYAHU N, et al. A local search approximation algorithm for k-means clustering[J]. Computational Geometry, 2004, 28∶ 89-112.

[5] MATOUSEK J. On approximate geometric k-clustering[J]. Discrete and Computational Geometry, 2000, 24∶ 61-84.

[6] SONG M J, RAJASEKARAN S. Fast k-means algorithms with constant approximation[A]. Proceedings of the 16th Annual International Symposium on Algorithms and Computation[C]. Sanya, Hainan,China, 2005,1029-1038.

[7] INABA M, KAOTH N, IMAI H. Application of weighted voronoi diagrams and randomization to variance-based k-clustering(extended abstract)[A]. Proceedings of the tenth annual symposium on Computational Geometry[C]. Stony Brook, New York, USA, 1994. 332-339.

[8] KUMAR A, SABHARWAL Y, SEN S. A sample linear time (1+ε)algorithm for k-means clustering in any dimensions[A]. Proceedings of the 45th IEEE Symposium on the Foundations of Computer science[C].Washington, DC, USA, 2004. 454-462.

[9] MOTAWANI R, RAGHAVAN P. Randomized Algorithms[M]. Cambridge University Press, Cambridge, UK, 1995.

[10] RUSPINI E H. Numerical methods for fuzzy clustering[J]. Inform Science,1970, 2(3)∶319-350.

[11] SPAETH H. Cluster Analysis Algorithms for Data Reduction and Classification of Objects[M]. John Wiley & Sones, 1980.

[12] PENG J M, XIA Y. A New Theoretical Framework for k-Means-Type Clustering[R]. McMaster University, Advanced Optimization Laboratory, Tech Rep∶ ADVOL2004/06, 2004.