龔運鴻，周新志，雷印杰

(四川大學 電子信息學院，成都 610065 )

0 引言

隨著程序需要處理的數(shù)據(jù)量越來越龐大，現(xiàn)如今以GB或TB為單位的數(shù)據(jù)集已經(jīng)十分普遍了，數(shù)據(jù)挖掘中必須重視的一個問題就是如何高效得處理如此龐大的數(shù)據(jù)。即使算法的復雜程度是線性增長的，時間和空間的消耗也不容忽視。K均值聚類算法由Stuart Lloyd等人在1957年第一次提出[1]，K均值聚類算法源于信號處理的一種矢量量化方法，由于其概念簡單、收斂速度快、易于實現(xiàn)等特點，現(xiàn)如今在數(shù)據(jù)挖掘領域聚類分析中十分流行。然而K均值聚類算法的復雜度比較高，如何高效進行算法計算是一個急切的研究方向。目前，陶冶[2]等人證明并實現(xiàn)了并行K均值聚類算法。喻金平[3]等人提出一種改進的人工蜂群算法，解決了K均值算法的搜索能力差的問題?；衄撉颷4]等人提出分塊、并行的K均值聚算法，采用“合并訪問”、“多級規(guī)約求和”和“負載均衡”等優(yōu)化策略優(yōu)化并行算法，提高了算法的運行速度。對于K均值聚類并行計算的研究還存在許多缺陷，比如沒有針對CUDA并行計算平臺進行優(yōu)化，也沒有針對中心點更新效率問題提出解釋等。根據(jù)以上研究的缺陷，本文利用NVIDIA的CUDA并行平臺，在傳統(tǒng)K均值算法基礎上，采用了一種滑動門并行計算中心點算法優(yōu)化K均值算法更新中心點的耗時問題。通過與規(guī)約法計算中心點算法相比，獲得了很好的加速比。

1 CUDA并行計算平臺

隨著社會的發(fā)展，CPU逐漸達到了速度極限并且購置成本也在急速上升。通過對顯示圖像的優(yōu)化，GPU技術卻逐漸完善了起來，在通用計算領域，GPU的能力逐漸超過了傳統(tǒng)CPU。通過近幾年的發(fā)展，計算技術正在從單一CPU串行計算方式向GPGPU并行協(xié)同計算發(fā)展。為了讓開發(fā)者無需學習復雜的著色語言和圖像處理原語，方便更多的開發(fā)者能夠進行GPU編程，顯卡廠家NVIDIA在2007年提供了一個方便開發(fā)者使用的接口——Compute Unified Device Architecture,即CUDA。CUDA可以讓開發(fā)者直接訪問GPU的虛擬指令設定和并行計算元素。與在GPU上使用圖形API進行計算的傳統(tǒng)方式相比，CUDA具有十分巨大的優(yōu)勢：1)程序可以在內(nèi)存的任意位置進行讀??；2)在CUDA4.0以上的版本擁有統(tǒng)一虛擬內(nèi)存；3)在CUDA6.0以上的版本擁有統(tǒng)一內(nèi)存；4)CUDA為開發(fā)者開辟一個快速共享內(nèi)存區(qū)域，其數(shù)據(jù)可以在線程中共享。這可以幫助開發(fā)者管理緩存、建立更高的帶寬；5)可以更快地從GPU上下載和回寫數(shù)據(jù)。在科技研究方面，CUDA技術可謂炙手可熱，越來越多的研發(fā)人員投入到了CUDA并行技術研究當中，所以在科研領域，CUDA也擁有了許多研究成果：在醫(yī)療領域，可以將CUDA加速技術運用在CT或者MRI掃描圖像的VR技術上；在物理建模上，特別是流體動力學領域有很好的成果；在神經(jīng)網(wǎng)絡訓練領域中，CUDA技術可以很好地解決機器學習遇到的瓶頸。

CUDA將C語言進行擴展，這樣可以使開發(fā)者使用C語言等高級語言在GPU設備上進行并行程序編程，方便開發(fā)者使用。使用CUDA進行編寫的代碼，既可以在CPU上使用，也可以在GPU上使用。使用CUDA并行計算平臺時，主機端為CPU，設備端為GPU。在GPU運行時，基于NVIDIA自身底層硬件特點，采用了單程序多數(shù)據(jù)(SPMD)的并行計算方式來處理數(shù)據(jù)，屬于單指令多數(shù)據(jù)(SIMD)的一種變體。并行編程的核心是線程的概念，運行的程序稱為內(nèi)核函數(shù)(kernel)，并行計算時，每一個線程都會同時處理一個kernel，CUDA中數(shù)據(jù)通過線程-塊-網(wǎng)格計算方式進行分配。每個線程塊都有自己唯一的識別id，一個或者多個線程塊會由流多處理器SM來并行處理，每一個SM由很多個32位寄存器組成。每個SM中有8個流處理器(SP)構成。每一個流處理器都有一個共享存儲器，所有SP都可以訪問共享存儲器中的內(nèi)容。如果線程塊是一維的，那么blockId.x就可以確定線程塊的id了。如果線程塊是二維的，那么需要blockId.x和blockId.y才能確定線程塊的id。每一個線程也有自己唯一識別id，通過這個id查找需要處理數(shù)據(jù)的位置。如果線程是一維的，那么threadIdx.x就可以確定線程的id了。如果線程是二維的，那么需要threadIdx.x和threadIdx.y才能確定線程的id。如果線程是三維的，那么需要threadIdx.x，threadIdx.y，threadIdx.z這3個參數(shù)才能確定線程的id。同一線程塊中的線程如果需要進行數(shù)據(jù)傳輸，一般是在共享存儲器中進行的。

2 K均值聚類算法

2.1 K均值聚類算法介紹

聚類是一種無監(jiān)督的學習，它將相似的對象歸到同一簇中，簇內(nèi)的對象越相似，聚類的效果越好。K均值聚類算法可以發(fā)現(xiàn)k個不同的簇，且每個簇的中心采用簇中所含值的均值計算而成。

K均值聚類算法是將含有n個數(shù)據(jù)點劃分為k個簇Cj(j=1,2,...,k;k≤n)。首先在n個數(shù)據(jù)點中隨機選取k個數(shù)據(jù)點代表k個簇的質心，再根據(jù)剩余數(shù)據(jù)點與k個質心的距離，將剩余數(shù)據(jù)點劃分到與其距離最近的簇Cj中。然后重新計算每個簇中所有值的均值，并將這個均值作為新的質心。該過程不斷重復，直到誤差平方和函數(shù)SSE收斂，其定義如公式(1)所示：

(1)

(2)

式中,k代表選擇的簇的數(shù)量；Cj代表第j(j=1,2,...,k;k≤n)個簇；x代表簇Cj中的任意一個數(shù)據(jù)點；cj是簇Cj的均值；mj代表簇Cj中數(shù)據(jù)點的個數(shù)。

2.2 K均值聚類流程

1)初始化聚類中心。在n個數(shù)據(jù)點鐘隨機選取k個數(shù)據(jù)點作為初始聚類中心C1,C2,...,Ck。

2)計算剩余數(shù)據(jù)點到每個初始聚類中心的距離，將剩余數(shù)據(jù)點幾個劃分到與其最近的簇中心代表的簇中，根據(jù)公式(1)計算SSE的值。

3)分別算出各個簇中所有數(shù)據(jù)點的均值，用這些均值替換初始聚類中心，用新的聚類中心重復步驟2)。

4)將上一次SSE的值和本次SSE的值進行比較，差值絕對值大于閥值，代表SSE收斂，則進行步驟5)，否則進行步驟5)。

5)算法結束。

排除其他因素只考慮算法核心部分，把每一次的比較、乘法、加法操作都當做一次浮點運算，用Tflop來代表一次浮點運算花費的時間，則算法核心部分每次迭代所用的時間見表1。

表1 K均值聚類算法核心部分所用運算時間

3 數(shù)據(jù)準備階段并行設計

3.1 初始化數(shù)據(jù)的并行設計

選取初始化聚類中心的過程，不是完全隨機產(chǎn)生的，需要先確定數(shù)據(jù)點在所有維度的最大最小值，如圖1所示。

圖1 二維數(shù)據(jù)點分布圖

在選取初值時，求最大最小值是典型的規(guī)約運算，而規(guī)約運算是可以并行化的。對于n個輸入數(shù)據(jù)和操作⊕，規(guī)約可表示為：

⊕a1⊕a2⊕...⊕an

(3)

圖2展示了處理N個元素規(guī)約操作的實現(xiàn)。

圖2 3種不同規(guī)約操作的實現(xiàn)

由圖2可以發(fā)現(xiàn)，不同的實現(xiàn)方式，其時間復雜度也是不同的。其中，串行實現(xiàn)完成規(guī)約操作需要n-1步，兩種對數(shù)步長的規(guī)約操作只需要lgn步，大大減少了時間復雜度。但是，由于成對方式不能合并內(nèi)存事務，成對方式在CUDA實現(xiàn)中性能比較差。在CUDA中，交替方式在全局內(nèi)存和共享內(nèi)存都有很好的效果。在全局內(nèi)存中，將blockDim.x*gridDim.x的倍數(shù)作為作為交替因子，這樣可以獲得比較好的性能。在共享內(nèi)存中，需要保持線程塊相鄰的線程保持活躍狀態(tài)，同時，為了避免內(nèi)存的沖突，需按確定的交替因子來累計結果，這樣可以獲得良好的效果。

3.2 分配數(shù)據(jù)的并行設計

基于GPU的K均值聚類算法在進行數(shù)據(jù)對象分配這一步驟的時候有兩種策略。第一種策略是面向每個簇的中心，通過計算出該簇的中心與每個數(shù)據(jù)對象的距離，然后將每個數(shù)據(jù)對象歸并到距離簇中心最近的那個簇中。這種策略適應于GPU的處理核心數(shù)量比較少的情況，此時，GPU中的每個處理核心可以對應一個簇的中心，并且能夠連續(xù)的處理所有的數(shù)據(jù)對象。另一種策略是面向每個數(shù)據(jù)對象，每個數(shù)據(jù)對象計算與所有簇中心的距離，然后數(shù)據(jù)對象將會被分配到距離簇中心最近的那個簇當中。該策略適應于GPU的處理核心比較大的情況，目前主流的GPU一般都有超過100個處理核心，因此第二種處理策略比較合適。

圖3 數(shù)據(jù)分配過程不同策略對比圖

4 滑動門中心點并行算法設計

在目前的CUDA并行計算算法優(yōu)化中,一般使用帶狀劃分的并行劃分方法[5-8]，如圖4所示，K均值聚類算法的中心點存儲在共享內(nèi)存中，而樣本保存在全局存儲器中，每一行數(shù)據(jù)由一個線程進行處理。帶狀劃分方法能最大化地利用GPU的計算內(nèi)核。

采用帶狀劃方法的方法，線程1～m將在同一時間去讀取共享儲存器索引位置0的數(shù)據(jù)。在CUDA中，寄存器的訪問速度最快。所以為了減少重復訪問的時間，可以把共享存儲器中的數(shù)據(jù)轉化到寄存器當中進行存儲。

圖4 帶狀劃分的并行聚類算法

采用滑動門并行算法時，也是采用帶劃分的劃分方式，在每個塊上開啟m個線程，相同簇內(nèi)的樣本都會在這些線程中進行計算。每個塊中的樣本點值最后合并至臨時中心點存儲區(qū)s_cust中。m的最佳值由式(4)計算。

(4)

式中，L代表樣本向量的維度；wsize代表wrap的大??；Ci代表簇i中的樣本總數(shù)；tNum代表每個block中可以開啟的最大線程數(shù)量；SMsize代表共享內(nèi)存的大小。

在長度為L的s_cust上分配大小為m的滑動門，block內(nèi)的線程按照線程號分配其計算的數(shù)據(jù)。在同一時刻，block內(nèi)所有線程的存儲數(shù)據(jù)均或落在滑動門范圍內(nèi)，且每個線程會計算獨立維度的數(shù)據(jù)，這樣就避免了存儲數(shù)據(jù)時的線程同步。并行計算完m個數(shù)據(jù)之后，滑動門會向前移動，同時線程也會向前啟動，計算下一批m個數(shù)據(jù)，直到滑動門回到初始位置為止?；瑒娱T并行計算中心點的算法過程偽代碼如下：

圖4 滑動門并行計算中心點

輸入：所有樣本的聚類結果；總數(shù)為n的樣本集D，其中樣本為L維向量；簇數(shù)目k。

輸出：所有簇的中心點集合U

for每個簇do

1)根據(jù)(4)選取合適的m值。

2)for j=ni,do:

(1) j = j/m

(2) for i=0 to L

for 每個線程 paraplle do：

s=blockIdx.x × block_Size + (theadIdx.x +i)%L;

s_cust[i] += data[s];

end loop;

(3) s=blockIdx.x × block_Size + i;

data[blockIdx.x × L + i]=s_cust[i];

end loop;

3)or 每個線程 paralle do:

U[threadIdx.x] = data[0][threadIdx.x] / ;

end loop

根據(jù)實驗環(huán)境的GPU配置，由(1)可以計算出最合適的m值。通過輸入的總數(shù)為n的樣本集，可以計算出程序所需的迭代次數(shù)為L×logmn。其中，步驟(2)的迭代次數(shù)是L，m個維的值將會通過滑動門的方式被并行計算。同一個塊中所有

的線程都會對存儲器s_cust中的值進行訪問，為了降低訪問數(shù)據(jù)時的延遲，可以采用共享存儲器來存儲s_cust中的所有的值。通過步驟3)，共享存儲器中的所有的數(shù)據(jù)都會按照塊的編號轉移到全局存儲器當中，這樣就不用重新分配其他的存儲空間，可以降低運算的速度。當步驟2)運行結束之后，該簇內(nèi)的樣本的和將會存儲在全局存儲器的的第一行中，中心點存儲區(qū)將會保存根據(jù)步驟3)得到的均值。通過分析，滑動門中心點并行算法的時間復雜度為logmn，相比傳統(tǒng)規(guī)約法的計算中心點的時間復雜度相比，得到了明顯的提高。

4 實驗結果

為了驗證滑動門并行算法的有效性，本實驗對一個含有25789個數(shù)據(jù)點的樣本集進行K均值聚類算法并行計算，采用單個簇一次迭代中分別使用滑動門法與規(guī)約法在不同m值下并行計算中心點的運行時間比。

表2 實驗環(huán)境

表3 滑動門法與規(guī)約法計算中心點時間加速比

通過表(3)分析，可以看出當m值分配較小值時，傳統(tǒng)規(guī)約并行算法可以得到比較好的加速比，但是與滑動門法相比，優(yōu)勢并不十分明顯，只有0.05的優(yōu)勢。但是，當分配的m達到相當規(guī)模之后，滑動門中心點并行算法的優(yōu)勢逐漸顯示了出來，并且優(yōu)勢逐漸擴大，從m=64的0.09到m=512時的0.47。

5 結語

本文基于CUDA平臺，對K均值聚類算法的并行加速計算研究，討論了CUDA在實現(xiàn)并行計算的關鍵技術。本文介紹了CUDA平臺的基本內(nèi)容和K均值聚類算法的基本原理。在并行步驟中，基于以往研究的不足之處，提出了對初始值規(guī)約計算的并行方法進行的討論和選擇、數(shù)據(jù)分配時的并行策略選擇和計算中心點的并行滑動門方式，這些方法提高了GPU的利用率，降低了數(shù)據(jù)獲取延遲，充分發(fā)揮了CUDA架構的并行計算優(yōu)勢。通過與傳統(tǒng)規(guī)約并行方法相比，滑動門并行計算算法擁有更好的加速比。通過實驗證明了CUDA在并行計算方面的有效性，結果表明了在GeForce730上實現(xiàn)了算法2.8倍的加速比，有效地提高了K均值算法在計算機上的運行效率。

[1]LloydS.LeastsquaresquantizationinPCM[J].IEEETransactionsonInformationTheory,1982,28 (2): 129-137.

[2] 陶 冶,曾志勇,余建坤,等.并行k均值聚類算法的完備性證明與實現(xiàn)[J].計算機工程, 2010, 36 (22) :72-74.

[3] 喻金平,鄭 杰,梅宏標,等. 基于改進人工蜂群算法的K均值聚類算法[J]. 計算機應用,2014,34(4):1065-1069,1088.

[4] 霍迎秋，秦仁波，刑彩燕，等.基于CUDA的并行K-means聚類圖像分割算法優(yōu)化[J].農(nóng)業(yè)機械學報,2014,45(11):72-74.

[5]NguyenH.GPUGems[M].Boston:Addison-WesleyProfession,2008.L

[6]RezaR,DanielR,EllickC,etal.Aparallelimplementationofk-meansclusteringonGPUs[A].ProcofInternationalConferenceonParallelandDistributedProcessingTechniquesandApplications[C].Springer-Verlag, 2008:340-345.

[7]MarioZ,MichaelG.AcceleratingK-meansonthegraphicsprocessorviaCUDA[A].Procofthe1stInternationalConferenceonIntensiveApplicationsandServices[C]:IEEE,2009:7-15.

[8]BaiHT,HeLL,OuyangDT,etal.K-meansoncommodityGPUswithCUDA[A].ProcofWRIWorldCongressonComputerScienceandInformationEngineering[C].ACMPress, 2009:651-655.