徐航帆，劉 叢，唐堅剛，彭敦陸

(上海理工大學 光電信息與計算機工程學院,上海 200082)

現(xiàn)今世界數(shù)據(jù)量爆炸增長，對大數(shù)據(jù)進行高效地分析將在各行各業(yè)起到關(guān)鍵性作用。聚類分析[1]是一種無監(jiān)督學習方法，也是數(shù)據(jù)挖掘和機器學習中一個重要的研究方向。聚類是指將數(shù)據(jù)在不需要先驗知識的情況下根據(jù)某種相似度劃分為多個簇的過程。傳統(tǒng)的聚類算法，例如K-means算法[2]具有簡單高效的特點，在球形結(jié)構(gòu)的數(shù)據(jù)上有著較高的聚類精確度。但是，其缺乏處理復雜結(jié)構(gòu)數(shù)據(jù)的能力[3]。當樣本空間為非凸時，K-means算法往往達不到理想的效果[4]。

為了能在任意形狀的樣本空間上聚類，且收斂于全局最優(yōu)解，研究人員開始研究新型的聚類算法，也被稱為譜聚類[5]算法(Spectral Clustering Algorithm，SC)。譜圖理論[6]是該算法的理論基礎(chǔ)，其本質(zhì)是將聚類問題轉(zhuǎn)化為圖的最優(yōu)劃分問題。近幾十年來，研究人員提出了多種譜聚類算法，例如Ratio Cut算法[7]和Normalized Cut算法[8]。這兩種算法將數(shù)據(jù)轉(zhuǎn)換為基于相似性的加權(quán)無向圖， 然后通過圖優(yōu)化算法進行求解。文獻[9]提出了NJW算法，該算法是一種應(yīng)用比較廣泛的譜聚類算法。然而譜聚類算法在面對大規(guī)模數(shù)據(jù)集時，其構(gòu)建相似度矩陣所需的空間復雜度(O(n2))和拉普拉斯矩陣特征分解所需的時間復雜度(O(n3))十分龐大，產(chǎn)生的計算成本將難以承受。

為了克服譜聚類的在大數(shù)據(jù)集上的可擴展性問題，一些加速譜聚類的算法被提出來，其中最自然的想法就是降低拉普拉斯矩陣的特征分解時間。2004年，文獻[10]采用經(jīng)典的Nystr?m[11]方法有效地計算了特征分解的近似解，從原始數(shù)據(jù)集中隨機選取p個樣本，然后計算相似度矩陣W∈p×p，根據(jù)這個矩陣計算特征分解來近似表示原始相似度矩陣的特征分解。文獻[12]提出了可擴展的Nystr?m方法，通過隨機低秩矩陣近似算法大大減少了算法的運行時間。文獻[13]提出了一種自適應(yīng)采樣的方法，改進了Nystr?m譜聚類算法的聚類效果。

后來，文獻[14]提出了一種通過地標點表示的加速譜聚類算法(Landmark-based Spectral Clustering，LSC)，其效果優(yōu)于前述算法聚類。該算法中，需首先對原始數(shù)據(jù)進行采樣，選取p個地標點，通過p個點與原始n個數(shù)據(jù)點成對相似度來構(gòu)建相似度矩陣Z∈n×p。然后，利用稀疏編碼技術(shù)[15]調(diào)整矩陣Z，使其成為稀疏相似度矩陣[16]，從而將計算特征分解的時間復雜度降低為O(p3+p2n)。文獻[14]提出了兩種采樣算法：(1)隨機采樣選取地標點；(2)K-means算法采樣選取聚類中心作為地標點。經(jīng)實驗證明，K-means采樣在大部分數(shù)據(jù)集上的最終聚類精確度比隨機采樣要高，但K-means采樣選取的地標點需重復讀取數(shù)據(jù)，時間復雜度大。隨機采樣選取地標點隨機性較大，地標點有時往往選取的不夠均勻，效果較差，文獻[17]提出了一種基于Pagerank算法選取地標點的方法，通過構(gòu)建數(shù)據(jù)點之間成對相似度矩陣W∈p×p選取地標點，并將該方法通過并行計算框架實現(xiàn)，有效減少了運行時間，獲得了更好的聚類精度。但是該算法的空間復雜度(O(n2))較大，面對大型數(shù)據(jù)集時需要更多的計算節(jié)點，成本較大。文獻[18]提出了一種快速地標點選取算法，有效保留了數(shù)據(jù)的原始信息，在聚類精度和時間復雜度之間取得較好的平衡，但該方法仍會受到離群地標點的影響。

綜上所述，基于地標點的加速譜聚類算法存在以下兩個問題：(1)地標點的選取難以在時間空間成本和最終聚類精度之間保持有效平衡；(2)聚類結(jié)果易受分布不均勻地標點和離群地標點的影響。

針對上述問題，本文提出了改進地標點采樣(Improved Landmark Selection，ILS)的加速譜聚類算法(LSC-ILS)。該方法通過隨機多組待選地標點集，根據(jù)地標點之間相似度標準差信息，選擇分布最均勻的地標點集，結(jié)合地標點周圍原始數(shù)據(jù)點局部密度分布信息，去除周圍原始點局部密度分布較小的離群地標點。實驗證明此算法能用較少的時間成本獲得較高的聚類精確度,在時間空間成本和聚類精度上取得較好的平衡。

1 相關(guān)工作

1.1 譜聚類

給定數(shù)據(jù)集X={x1,x2,…,xn}∈m×n，其中n表示樣本個數(shù)，m表示數(shù)據(jù)維度。譜聚類首先根據(jù)數(shù)據(jù)點之間成對相似度矩陣W∈n×n構(gòu)造無向G=(V,E)，其中W中第i行第j個元素wij≥ 0表示點xi和點xj之間的相似度。X中的樣本點對應(yīng)V的頂點，任意兩點之間權(quán)重由W給出。度矩陣D是由相似度矩陣W的行和組成的對角線矩陣,如式(1)所示。

(1)

L=D-W被稱為拉普拉斯矩陣。然后譜聚類對拉普拉斯矩陣進行特征分解，并對前k個最小特征值對應(yīng)的特征向量組成的點進行聚類(一般利用K-means聚類算法)，作為最終聚類結(jié)果?；A(chǔ)譜聚類算法流程如下文所示。

輸入：n個數(shù)據(jù)點x1,x2,x3,…,xn∈m；聚類數(shù)目k。

輸出：聚類結(jié)果標簽。

步驟1根據(jù)數(shù)據(jù)點之間相似度構(gòu)造相似度矩陣W∈n×n，根據(jù)式(2)計算度矩陣D∈n×n；

步驟2根據(jù)拉普拉斯矩陣L=D-1/2(D-W)D-1/2計算前k個最小特征值對應(yīng)的特征向量，組成點集Q={q1,q2,…,qk}；

步驟3Q的每一行視作一個點，通過K-means算法得到最終聚類結(jié)果。

由于傳統(tǒng)譜聚類算法空間復雜度(O(n2))和時間復雜度(O(n3))，很難將其擴展到大規(guī)模數(shù)據(jù)集的應(yīng)用。因此本文提出了擴展到大規(guī)模數(shù)據(jù)集上的加速譜聚類的算法。

1.2 基于地標點的加速譜聚類算法

基于地標點的加速譜聚類算法是一種在大規(guī)模數(shù)據(jù)集上加速譜聚類的方法。該算法通過矩陣分解技術(shù)找到一組點，通過這組點與原始點之間的關(guān)系來近似表示原始數(shù)據(jù)。X={x1,x2,…,xn}∈m×n是原始數(shù)據(jù)矩陣，矩陣分解技術(shù)就是找到p個m維的數(shù)據(jù)點集U∈m×p和與原始點的關(guān)系矩陣Z∈p×n，通過式(2)近似表示原始數(shù)據(jù)矩陣X。

X≈UZ

(2)

式中，矩陣U為地標點集，每一列為每一個地標點，這些地標點為原始數(shù)據(jù)點的代表點。對于任意一個數(shù)據(jù)點xi∈X，它的近似點φi可以通過式(3)表示。

(3)

式中，uj是U中第j列向量，表示第j個地標點；zji是矩陣Z中第j行第i列的元素。如果uj不在點xi最近的r(≤p)個鄰域中，則zji置為0，因此Z就變成了稀疏表示矩陣。通過U(i)∈m×r表示U的子矩陣，由xi的r個最近的地標點組成。zji可由式(4)計算。

(4)

K(·)是和函數(shù)，通常使用高斯核函數(shù)

(5)

根據(jù)矩陣Z計算相似度矩陣W，如式(6)

W=(D-1/2Z)T(D-1/2Z)

(6)

其中，D是一個p×p對角矩陣，由Z的行和構(gòu)成。D-1/2Z的奇異值分解(Singular Value Decomposition，SVD)如下

D-1/2Z=AΣBT

(7)

其中，Σ=diag(?1,?2,…,?p)，?1≥?2≥…≥?p≥0是D-1/2Z的奇異值；A={a1,a2,…,ap}∈p×p是左奇異向量矩陣；B={b1,b2,…,bp}∈n×p是右奇異向量矩陣；?i2為特征值。很明顯，B的列向量為W的特征向量，A的列向量為(D-1/2Z)(D-1/2Z)T的特征向量。

由于矩陣(D-1/2Z)(D-1/2Z)T的維度是p×p，可以先計算A，只需要O(p3)，B因此可以由下式(8)計算。

BT=∑-1AT(D-1/2Z)

(8)

由于p?n，所以拉普拉斯矩陣的分解時間就從O(n3)減少到O(p3+p2n)，存儲空間由O(n2)減少到O(np)，大幅減少了譜聚類所需的時間和空間。

2 改進地標點采樣的加速譜聚類算法

地標點的選取作為基于地標點的加速譜聚類算法的最關(guān)鍵部分，在很大程度上決定了最后聚類結(jié)果。隨機采樣的方法在大規(guī)模數(shù)據(jù)上選取的地標點易集中于某一局部區(qū)域，分布不均勻。K-means采樣算法存在著反復讀取數(shù)據(jù)，時間消耗較大的問題；Pagerank采樣算法需先構(gòu)造數(shù)據(jù)點之間成對相似度矩陣，空間復雜度較大。此外，聚類結(jié)果也易受離群地標點的影響。針對以上問題，本文提出了改進地標點采樣的加速譜聚類算法。該算法首先多次隨機采樣選取多組地標點集，計算各組地標點之間相似度的標準差(Standard Deviation，SD)，通過標準差來衡量地標點的均勻程度；然后計算地標點周圍原始數(shù)據(jù)點密度分布來衡量地標點的代表性，去除代表性較差的離群地標點，進一步在不影響聚類精度的情況下減少時間復雜度。

2.1 改進的地標點采樣算法

給定數(shù)據(jù)集X={x1,x2,…,xn}∈m×n，其中n表示樣本個數(shù)，m表示數(shù)據(jù)維度。首先通過均勻隨機采樣，從原始數(shù)據(jù)集中進行c次地標點采樣，每次采樣p個地標點，形成c個地標點集，其中第i個地標點集表示為Ui={u1,u2,…,up}∈m×p。為了選擇最佳的地標點集，通過地標點集Ui中點之間成對相似度矩陣Si∈p×p的行標準差之和來衡量地標點的分布均勻程度。其中，點uj與點ul之間相似度計算如式(9)所示。

(9)

其中，K(·)是核函數(shù)，通常由式(5)計算。

第j行的行標準差計算如下

(10)

(11)

各個地標點集的行標準差之和即為SD={SD1,SD2,…,SDc},其中max(SD)對應(yīng)的地標點集Ui={u1,u2,…,up}為最佳地標點集。

接下來，計算各個地標點的局部密度權(quán)重，通過局部密度權(quán)重反映地標點周圍原始數(shù)據(jù)點的分布情況。局部密度權(quán)重越大表示地標點周圍數(shù)據(jù)點分布越密集，地標點就越重要，局部密度權(quán)重越小表示地標點周圍數(shù)據(jù)點分布越稀疏，地標點就越不重要；然后去除密度權(quán)重小于閾值γ的地標點(Landmarks Reduction，LR)。具體步驟為：首先構(gòu)造地標點與原始數(shù)據(jù)點之間的相似度矩陣Z′∈p×n,z′ji表示第j個地標點與原始數(shù)據(jù)點xi之間的相似度，計算式為

(12)

其中，K(·)是核函數(shù)，通常由式(5)計算；然后計算每個地標點的密度權(quán)重。地標點uj權(quán)重dj計算如下

(13)

其中，count(z′ji>s,z′ji∈z′j)表示原始樣本中與地標點uj相似度大于s(0～1)的數(shù)量；最后根據(jù)閾值γ(0～1)，選出權(quán)重dj>γ的p′個地標點。

2.2 改進地標點采樣的加速譜聚類算法流程

圖1 LSC-ILS流程圖

圖1為改進地標點采樣的加速譜聚類算法流程，其具體流程如下文所示。

輸入：原始數(shù)據(jù)集X={x1,x2,…,xn}∈m×n，地標點數(shù)量p；系數(shù)相似度矩陣近鄰個數(shù)r；地標點集個數(shù)c；閾值s和γ；聚類數(shù)目k。

輸出：聚類結(jié)果標簽。

步驟1對原始數(shù)據(jù)集進行c次隨機采樣選取c個地標點集U={U1,U2,U3,…,Ui,…,Uc}；

步驟2根據(jù)式(9)計算每個地標點集的相似度矩陣Si；

步驟3根據(jù)式(10)和式(11)計算出每個地標點集相似度矩陣的標準差SDi；

步驟4計算max(SD)，輸出其對應(yīng)的地標點集；

步驟5根據(jù)式(12)計算地標點密度權(quán)重dj；

步驟6選出閾值dj>γ的p′個地標點；

步驟7根據(jù)式(5)計算地標點與原始數(shù)據(jù)點之間的稀疏相似度矩陣Z；

步驟8計算(D-1/2Z)(D-1/2Z)T的前k個最小特征值對應(yīng)的特征向量A={a1,a2,…,ap′}；

步驟9根據(jù)式(7)計算B={b1,b2,…,bp′}；

步驟10B的每一行視作一個點，通過K-means算法得到最終聚類結(jié)果標簽ri∈R。

2.3 時間復雜度和空間復雜度分析

本節(jié)是對本文算法的時間復雜度和空間復雜度分析。數(shù)據(jù)樣本數(shù)量為n，采樣次數(shù)為c次，每次采樣地標點數(shù)量為p，最后輸出地標點數(shù)量為d，其余算法采樣地標點數(shù)量為p。

采用隨機采樣選取c組待選地標點的時間復雜度可以忽略不記，構(gòu)建地標點相似度矩陣需要O(cp2)的時間復雜度，需要O(cp2)空間復雜度。根據(jù)地標點與原始點之間的相似度矩陣計算地標點密度權(quán)重以及構(gòu)建稀疏相似度矩陣需要O(pn)時間復雜度。計算稀疏相似度矩陣需要O(pn2)，根據(jù)稀疏相似度矩陣計算左奇異向量需要O(d3)的時間復雜度，計算右奇異向量需要O(d2n)的時間復雜度。

表1 時間復雜度分析 (p?n,d

表1是4種算法的時間復雜度分析，由表1可以看出本文算法的時間復雜度總和最小。

3 實驗驗證

為了驗證本文算法的有效性，進行實驗分析，從算法的聚類精確度和運行時間兩指標進行評估。本文選取原始譜聚類[5](記為SC)、Nystr?m 近似譜聚類[11](記為Nystr?m)、基于隨機采樣的地標點譜聚類[13](記為LSC-R)和基于K均值中心的地標點譜聚類[13](記為 LSC-K)與本文算法(記為LSC-ILS)進行對比。為了證明本文提出的去除地標點算法(LR)有效性，本文將LSC-R和LSC-K分別結(jié)合LR算法(記為LSC-R-LR和LSC-K-LR)與原始LSC-R和LSC-K算法進行實驗比較。

3.1 數(shù)據(jù)集

為了驗證本文算法面對大規(guī)模數(shù)據(jù)集的性能，選擇兩個大規(guī)模UCI數(shù)據(jù)集進行實驗。第一個為MINIST，該數(shù)據(jù)集是一個手寫數(shù)字的數(shù)據(jù)集，共有70 000個樣本，共10類，每個樣本被視為784維數(shù)據(jù)。第二個為Pendigits，該數(shù)據(jù)集也是一個手寫數(shù)字的數(shù)據(jù)集，包含了44個作者寫的250個手寫樣例，經(jīng)過處理得到10 992個樣本，每個對象有16個屬性，共被分成10類。

3.2 評價指標

為了評價聚類有效性(聚類算法對數(shù)據(jù)進行劃分的正確程度)，本節(jié)采取了兩個指標，即聚類精確性(Accuracy，ACC)[19]和標準化互信息(Normalized Mutual Information，NMI)[20]進行數(shù)值分析。實驗在英特爾Core i3-8100@3.6 GHz，8 GB內(nèi)存的計算機上運行，代碼在MATLAB環(huán)境下編寫。

給定數(shù)據(jù)點xi，令ri∈R和fi∈F分別為聚類結(jié)果標簽和預(yù)定義標簽，ACC定義如下

(14)

式中，n是樣本數(shù)量；δ(x,y)表示若x=y，為1，否則為0。map(ri)是聚類結(jié)果標簽ri映射到數(shù)據(jù)語料庫中的真正標簽，采用Kuhn-Munkres算法[21]可以找到最佳映射。

MI(R，F(xiàn))互信息表示為

(15)

式中，p(ri)和p(fj)是從數(shù)據(jù)集中任意選擇的樣本屬于標簽ri和標簽fj的概率；p(ri,fj)是從數(shù)據(jù)集中任意選擇的樣本屬于簇ri和fj的聯(lián)合概率。NMI標準化互信息計算如下

(16)

H(R)=-∑ri∈Rp(ri)logp(ri)

(17)

式中，H(R)是聚類結(jié)果R的熵，表示R中變量的混亂程度，R中變量越混亂，H(R)的值越接近1；反之H(R)的值越接近0。很明顯ACC指標和NMI指標取值范圍都是[0,1]，越接近1表示聚類結(jié)果越精確。

3.3 實驗結(jié)果

7種聚類算法在MINIST和Pendigits上的聚類結(jié)果如表2和表3所示，聚類指標分別為聚類精確度(ACC)和標準化互信息(NMI)。為了比較的公平性，6種加速譜聚類算法的抽樣數(shù)均為1 000，基于地標點的加速譜聚類算法稀疏相似度矩陣的最近鄰個數(shù)r為5。其中，本文算法中s參數(shù)，γ參數(shù)和c參數(shù)均為人工設(shè)定，選取較優(yōu)的值。

表2和表3分別表示MINIST數(shù)據(jù)集和Pendigits數(shù)據(jù)集熵算法的各項指標對比。

表2 MINIST數(shù)據(jù)集上算法性能對比

表3 Pendigits數(shù)據(jù)集上算法性能對比

由表2和表3中可以看出，本文算法LSC-ILS聚類精度比Nystr?m和LSC-R算法高，比LSC-K算法低，運行時間短于LSC-K算法和Nystr?m算法，與LSC-R算法運行時間幾乎持平。LSC-R-LR算法和LSC-K-LR算法的聚類精度分別與LSC-R算法和LSC-K算法幾乎持平，但運行時間有所減少。因此，綜合聚類精確度與運行時間本文算法做到了更好的平衡，并且本文中LR算法可擴展到其它算法中，在保持聚類精確度的情況下減少運行時間。

為了驗證地標點個數(shù)對加速聚類算法的影響，將本文算法與其余3種加速譜算法在數(shù)據(jù)集MINIST上針對不同地標點數(shù)量在3個指標上進行對比。地標點數(shù)量從300開始，增量為300，到2 100結(jié)束，結(jié)果如圖2～圖4所示。

圖2 MINIST數(shù)據(jù)集上ACC指標隨地標點增長變化

由圖2～圖4可以看出，除了Nystr?m算法，其余算法的聚類精確性均隨著地標點個數(shù)增加而增加。本文算法聚類精確度始終保持在LSC-K與LSC-R之間，但運行時間幾乎與LSC-R算法持平，遠小于LSC-K算法，運行時間增長慢于LSC-K算法，并且隨著地標點個數(shù)的增加，本文算法精確度增長加快。

圖3 MINIST數(shù)據(jù)集上NMI指標隨地標點數(shù)量增長變化

圖4 MINIST數(shù)據(jù)集上時間指標隨地標點數(shù)量增長變化

為了進一步驗證本文中LR算法對LSC-R和LSC-K算法隨地標點數(shù)量增長的變化的影響，將LSC-R和LSC-K分別與LSC-R-LR和LSC-K-LR在數(shù)據(jù)集MINIST上針對不同地標點數(shù)量在3個指標上進行對比。地標點數(shù)量同樣從300開始，增量為300，到2 100結(jié)束。結(jié)果如圖5～圖10所示。

由圖5～圖10可以看出，LSC-R-LR算法和LSC-K-LR算法在不同的地標點數(shù)量情況下的聚類精度幾乎與LSC-R算法和LSC-K算法持平，運行時間有所減少，進一步證明了本文LR算法的穩(wěn)定性和有效性。

圖5 ACC指標隨地標點數(shù)量變化

圖6 NMI指標隨地標點數(shù)量變化

圖7 Time指標隨地標點數(shù)量變化

圖8 ACC指標隨地標點數(shù)量變化

圖9 NMI指標隨地標點數(shù)量變化

圖10 Time指標隨地標點數(shù)量變化

4 結(jié)束語

本文針對傳統(tǒng)加速譜聚類算法存在的一些問題提出了改進地標采樣的加速譜聚類算法。通過多次隨機采樣選取多組地標點集，利用各個地標點集內(nèi)地標點之間相似度矩陣的行標準差和衡量地標點的分布情況，選取最佳地標點集，并根據(jù)地標點的密度權(quán)重去除一些離群地標點。本文從理論和實驗分析證明了該算法的有效性。下一步計劃將本文算法擴展到分布式計算框架上，進一步提高算法效率。