一種基于K-Means均值聚類算法的大地測(cè)量數(shù)據(jù)粗差探測(cè)方法

楊玉，吳云龍，單維鋒，石江浩

(1.防災(zāi)科技學(xué)院，河北省地震動(dòng)力學(xué)重點(diǎn)實(shí)驗(yàn)室，河北 三河 065201；2.中國(guó)地震局地震研究所，中國(guó)地震局地震大地測(cè)量重點(diǎn)實(shí)驗(yàn)室，湖北 武漢 430071； 3.北京航空航天大學(xué)，北京 100191)

1 引 言

歐洲空間局(ESA)于2009年3月17日發(fā)射了地球重力場(chǎng)和海洋環(huán)流探測(cè)衛(wèi)星[1](GOCE，Gravity field and steady-state Ocean Circulation Explorer)，該衛(wèi)星以100km空間分辨率、1mGal的精度測(cè)定全球地球重力場(chǎng)，并以1-2cm的精度測(cè)定全球大地水準(zhǔn)面[2]。GOCE衛(wèi)星搭載的高靈敏度重力梯度儀可對(duì)地球重力場(chǎng)每分鐘的變化進(jìn)行三維測(cè)量[3]，獲得的數(shù)據(jù)可更深入了解地球內(nèi)部結(jié)構(gòu)，觀察海洋和氣候變化。由于衛(wèi)星重力梯度在測(cè)量時(shí)受到衛(wèi)星系統(tǒng)周圍環(huán)境、檢校、定姿等問(wèn)題的強(qiáng)烈干擾，造成衛(wèi)星重力梯度數(shù)據(jù)中不可避免地出現(xiàn)粗差[4]。因此在進(jìn)行重力梯度數(shù)據(jù)的預(yù)處理時(shí)，必須進(jìn)行相應(yīng)的粗差探測(cè)與剔除來(lái)獲得數(shù)據(jù)保障。

粗差探測(cè)一直是國(guó)內(nèi)外大地測(cè)量數(shù)據(jù)處理的熱點(diǎn)問(wèn)題之一，許多專家學(xué)者提出不同角度的粗差探測(cè)方法。常見(jiàn)的一類方法是研究粗差數(shù)據(jù)套入模型中，如暴景陽(yáng)等[5]基于逆?zhèn)鞑ド窠?jīng)網(wǎng)絡(luò)對(duì)多波束測(cè)深數(shù)據(jù)中的粗差進(jìn)行剔除，該方法用實(shí)驗(yàn)?zāi)Ｐ蛯?duì)比分析擬合圖，驗(yàn)證了逆?zhèn)鞑ド窠?jīng)網(wǎng)絡(luò)方法的有效性；基于數(shù)學(xué)統(tǒng)計(jì)原理，閆廣峰等[6]人利用L1范數(shù)良好的抗差特性，在矩陣初等變換的理論上，求證了L1范數(shù)具有粗差定位的能力；為了豐富粗差探測(cè)理論，崔太岷等[7]基于半?yún)?shù)的方差膨脹模型，推導(dǎo)了Score檢驗(yàn)公式，從而提高參數(shù)估計(jì)值接近真實(shí)值的精度。另一類方法是分析觀測(cè)值之間的相關(guān)性，如張建等[8]人利用基于卡方檢驗(yàn)的粗差抗差算法減少觀測(cè)值由于相關(guān)性引起的粗差誤判，保證粗差定位的準(zhǔn)確性；劉根友[9]模擬觀測(cè)值，采用最小二乘原理確定粗差；李玉芝等[10]人在不等權(quán)的條件下，以模擬觀測(cè)值和含粗差觀測(cè)值之間出現(xiàn)分群現(xiàn)象為目標(biāo)，對(duì)粗差進(jìn)行定位。

作為機(jī)器學(xué)習(xí)算法中的常用算法，聚類算法通過(guò)對(duì)衛(wèi)星數(shù)據(jù)進(jìn)行數(shù)據(jù)挖掘和可視化分析，探究海量數(shù)據(jù)中有價(jià)值的特征信息。作為當(dāng)前最流行、使用簡(jiǎn)單廣泛、易于實(shí)現(xiàn)且計(jì)算效率較高的聚類算法，K-Means均值聚類算法對(duì)于大數(shù)據(jù)集，具備較小的時(shí)間和空間復(fù)雜度，該方法在地球科學(xué)領(lǐng)域中也得到驗(yàn)證并使用。余德清[11]應(yīng)用K-Means均值聚類算法對(duì)水面信息進(jìn)行提取，通過(guò)算法本身的收縮運(yùn)算，結(jié)合觀測(cè)數(shù)據(jù)進(jìn)行對(duì)比，研究湖水資源量的變化。

本文基于K-Means均值聚類算法對(duì)數(shù)據(jù)進(jìn)行聚類，分析了重力梯度數(shù)據(jù)之間的相似度，提出了基于K-Means均值聚類算法對(duì)衛(wèi)星重力梯度變量這樣的大數(shù)據(jù)集有高效的聚類效果，并設(shè)計(jì)了實(shí)驗(yàn)，對(duì)比原始數(shù)據(jù)與聚類結(jié)果，驗(yàn)證算法性能。K-Means均值聚類算法的研究未來(lái)還有更廣闊的空間，更多優(yōu)化算法會(huì)在粗差探測(cè)領(lǐng)域中發(fā)揮作用。

2 K-Means均值聚類算法

K-Means均值聚類算法屬于機(jī)器學(xué)習(xí)中的數(shù)據(jù)挖掘技術(shù)，它等價(jià)于數(shù)據(jù)集的分組，給定一個(gè)數(shù)據(jù)集，聚類算法可將數(shù)據(jù)集中的每個(gè)數(shù)據(jù)劃分成不同類型的組。理論上，同一組中的數(shù)據(jù)具有相似的屬性或特征，而不同組中的數(shù)據(jù)具有高度不同的屬性或特征。

典型的聚類算法內(nèi)部計(jì)算過(guò)程，如圖1所示：

圖1 典型的聚類算法內(nèi)部計(jì)算過(guò)程

K-Means均值聚類算法是一種基于迭代求解、距離的無(wú)監(jiān)督聚類分析方法，其運(yùn)算速度快，操作簡(jiǎn)單，常被應(yīng)用在連續(xù)型的序列數(shù)據(jù)領(lǐng)域。K-Means均值聚類算法最關(guān)鍵的技術(shù)就是選取參數(shù)k，一般有兩種方法選擇最優(yōu)的k值，分別為：手肘法，輪廓系數(shù)法。

方法一：手肘法

手肘法的核心指標(biāo)是誤差平方和(SSE，Sum of the Square errors)越小，即聚合程度越高。

(1)

當(dāng)k小于最佳聚類中心個(gè)數(shù)時(shí)，由于k的增大會(huì)提高每個(gè)簇的聚合程度，故SSE的值會(huì)迅速下降；而當(dāng)k接近最佳聚類中心個(gè)數(shù)時(shí)，隨著增加k值，聚合程度的效果不大，SSE下降的趨勢(shì)減小，k越增大而SSE越趨于平緩，即SSE和k的關(guān)系圖是一個(gè)手肘的形狀，隨著分類的類別數(shù)增加，這個(gè)肘部對(duì)應(yīng)的k值就是數(shù)據(jù)的估計(jì)聚類數(shù)。實(shí)驗(yàn)開(kāi)始假設(shè)k>0，通過(guò)畫(huà)出K-SSE曲線，找出下降途中的拐點(diǎn)，即可較好地確定K值，曲線的第一個(gè)(或最顯著的)拐點(diǎn)表示估計(jì)的簇?cái)?shù)。

方法二：輪廓系數(shù)法

a(i)是某樣本xi與同簇的其他樣本的平均距離，稱為凝聚度，b(i)是某樣本xi與最近簇中所有樣本的平均距離，稱為分離度。

首先計(jì)算樣本i到同簇其他樣本的平均距離a(i)，a(i)越小，說(shuō)明樣本i越應(yīng)該被聚類到該簇。

然后計(jì)算樣本i到其他某簇的所有樣本的平均距離b(i)，b(i)越大，說(shuō)明樣本i越不屬于其他簇。

最后根據(jù)樣本i的凝聚度a(i)和分離度b(i)，定義樣本i的輪廓系數(shù)，輪廓系數(shù)S(i)越大，聚類效果越好[12]。表達(dá)式(2)定義如下：

(2)

按照最鄰近原則把待分類樣本點(diǎn)分到各個(gè)簇。然后按平均法重新計(jì)算各個(gè)簇的質(zhì)心，從而確定新的簇心。一直迭代，直到簇心的移動(dòng)距離小于某個(gè)給定的值。

3 數(shù)據(jù)和評(píng)估指標(biāo)

3.1 模型重力梯度數(shù)據(jù)

采用300階的EMG96模型，基于正常橢球GRS80模擬了徑向重力梯度觀測(cè)值Vzz，時(shí)長(zhǎng)為1天，共 17 280個(gè)數(shù)據(jù)，將其視為“純凈”數(shù)據(jù)。表1是模擬衛(wèi)星重力梯度數(shù)據(jù)采用的相關(guān)參數(shù)，表2是數(shù)據(jù)中常見(jiàn)的統(tǒng)計(jì)值[13]。

數(shù)據(jù)模擬采用的相關(guān)參數(shù) 表1

模擬重力梯度數(shù)據(jù)的統(tǒng)計(jì)值(單位：10-9s-2) 表2

3.2 模擬生成衛(wèi)星重力梯度數(shù)據(jù)

計(jì)算整個(gè)數(shù)據(jù)序列的標(biāo)準(zhǔn)差σ(如表2)，以0為期望、0.01σ為標(biāo)準(zhǔn)差模擬生成含白噪聲的序列，再向白噪聲序列中加入169個(gè)的粗差，其中149個(gè)粗差以隨機(jī)分布的形式加入整個(gè)數(shù)據(jù)序列中；20個(gè)粗差以連續(xù)的形式加入數(shù)據(jù)序列，最終得到含有白噪聲和粗差的衛(wèi)星重力梯度數(shù)據(jù)。

3.3 粗差探測(cè)評(píng)定指標(biāo)[14]

在衛(wèi)星重力梯度數(shù)據(jù)粗差探測(cè)的模擬研究中，通常采用以下兩個(gè)指標(biāo)來(lái)評(píng)價(jià)粗差探測(cè)的效果：

(1)成功率：

ORS=ns/no

(3)

式中：ns為探測(cè)到的粗差的個(gè)數(shù)，no為數(shù)據(jù)集中粗差的總個(gè)數(shù)。

(2)失敗率：

ORF=nf/n

(4)

式中：nf為有效數(shù)據(jù)被錯(cuò)誤地探測(cè)為粗差的個(gè)數(shù)，n為總的數(shù)據(jù)個(gè)數(shù)。

需要指出的是，GOCE實(shí)測(cè)數(shù)據(jù)中的粗差個(gè)數(shù)是未知的。本文通過(guò)模擬計(jì)算和分析，利用上述兩個(gè)指標(biāo)對(duì)粗差探測(cè)方法進(jìn)行評(píng)價(jià)，可為GOCE實(shí)測(cè)數(shù)據(jù)的粗差探測(cè)提供參考依據(jù)。

4 訓(xùn)練策略

4.1 用訓(xùn)練集對(duì)K-Means均值聚類算法進(jìn)行訓(xùn)練

K-Means均值聚類算法是將數(shù)據(jù)集按照樣本之間的距離大小，自動(dòng)劃分為k組，選定k個(gè)初始的簇中心，按照以下步驟迭代細(xì)化[15](如圖2所示)：

圖2 K-Means均值聚類算法實(shí)驗(yàn)流程

訓(xùn)練策略中，數(shù)據(jù)樣本處理步驟為：

(1)輸入數(shù)據(jù)集，選定好k個(gè)簇中心。隨機(jī)選擇數(shù)據(jù)集中的k個(gè)點(diǎn)分別作為k個(gè)簇的初始中心，然后計(jì)算其他點(diǎn)到這k個(gè)簇初始中心點(diǎn)的距離，從而決定這些其他點(diǎn)都屬于那個(gè)簇中心。比如某個(gè)點(diǎn)A到這k個(gè)點(diǎn)的距離產(chǎn)生k個(gè)數(shù)值，最小的那個(gè)值就是點(diǎn)A對(duì)應(yīng)的簇；

(2)從數(shù)據(jù)集中隨機(jī)選擇k個(gè)分類對(duì)象作為初始化的k個(gè)聚類質(zhì)心。簇中所有數(shù)據(jù)的均值通常被稱為這個(gè)簇的“質(zhì)心”。在一個(gè)二維平面中，一簇?cái)?shù)據(jù)點(diǎn)的質(zhì)心的橫坐標(biāo)就是這一簇?cái)?shù)據(jù)點(diǎn)的橫坐標(biāo)的均值，質(zhì)心的縱坐標(biāo)就是這一簇?cái)?shù)據(jù)點(diǎn)的縱坐標(biāo)的均值；

(3)將帶有聚類的數(shù)據(jù)放到一個(gè)聚類集合中，利用歐式距離計(jì)算相似度；

歐式距離表達(dá)式：

(5)

(4)根據(jù)聚類結(jié)果，不斷迭代，更新聚類質(zhì)心，使得類內(nèi)的相似度最大，類間的相似度最小。

(5)計(jì)算每個(gè)簇包含全部點(diǎn)的坐標(biāo)平均值作為新的質(zhì)心，進(jìn)行2次迭代，直到誤差平方和最小，即結(jié)果趨于平穩(wěn)收斂。

K-Means均值聚類算法對(duì)于大數(shù)據(jù)集樣本處理具有較高的效率且是可伸縮性的，該算法的迭代優(yōu)化功能改善了初始監(jiān)督學(xué)習(xí)樣本識(shí)別不合理的地方。當(dāng)樣本結(jié)果密集時(shí)，該算法的效果也較好，是其他聚類算法(如：譜聚法，高斯混合模型算法等)的基礎(chǔ)算法。

4.2 訓(xùn)練所得模型對(duì)數(shù)據(jù)集進(jìn)行處理

該方法的聚類階段如圖3，圖4所示(圖中x和y分別代表數(shù)據(jù)點(diǎn)的橫縱坐標(biāo)值)。該聚類算法包含兩個(gè)部分：其一是k，k代表著類的數(shù)目，實(shí)驗(yàn)中將數(shù)據(jù)聚為5類。其二是means，表示每次計(jì)算聚類中心的時(shí)候采取的是計(jì)算平均值。

圖3 原始數(shù)據(jù)

圖4 聚類結(jié)果(k=5)

從圖中可以看出，離群點(diǎn)是數(shù)據(jù)集中部分區(qū)域的粗差，在訓(xùn)練時(shí)無(wú)法在聚類過(guò)程中統(tǒng)一到集合中。

4.3 K-Means均值聚類模型下的粗差探測(cè)效果

如圖5是數(shù)據(jù)集在K-Means均值聚類模型下的粗差探測(cè)效果，圖中清晰地標(biāo)注了粗差的位置，可以看出，聚類算法能夠準(zhǔn)確地識(shí)別粗差。實(shí)驗(yàn)中，將訓(xùn)練所得數(shù)據(jù)與輸入的原始數(shù)據(jù)(含粗差數(shù)據(jù))進(jìn)行對(duì)比處理后，以此來(lái)初步判斷網(wǎng)絡(luò)模型的有效性。

圖5 在K-Means均值聚類模型下的粗差探測(cè)效果

5 計(jì)算結(jié)果與分析

K-Means均值聚類算法在重力梯度異常點(diǎn)檢測(cè)數(shù)據(jù)處理中具有顯著的優(yōu)勢(shì)。表3為應(yīng)用K-Means均值聚類算法對(duì)重力梯度數(shù)據(jù)粗差探測(cè)的結(jié)果，圖6展示了應(yīng)用K-Means均值聚類算法前后的粗差探測(cè)效果對(duì)比。結(jié)果顯示，經(jīng)過(guò)模型訓(xùn)練后的數(shù)據(jù)具有顯著的粗差識(shí)別能力。

重力梯度數(shù)據(jù)粗差探測(cè)結(jié)果 表3

圖6 重力梯度數(shù)據(jù)粗差探測(cè)前后的對(duì)比(片段)

實(shí)驗(yàn)證明，基于K-Means均值聚類算法，在選定聚類中心個(gè)數(shù)(簇的個(gè)數(shù))、質(zhì)心、迭代次數(shù)后，數(shù)據(jù)與數(shù)據(jù)之間區(qū)別較明顯，且質(zhì)心大小相近時(shí)，其聚類結(jié)果較理想。對(duì)于處理大數(shù)據(jù)集合，該算法非常高效，適應(yīng)性較好。

6 結(jié) 論

本文基于K-Means均值聚類算法進(jìn)行詳細(xì)的闡述，應(yīng)用該方法對(duì)重力梯度變量的模擬數(shù)據(jù)進(jìn)行粗差探測(cè)，結(jié)果顯示，K-Means均值聚類算法方法能夠高效地應(yīng)用在衛(wèi)星重力梯度數(shù)據(jù)的預(yù)處理中。未來(lái)的研究工作可將改進(jìn)的聚類算法與深度學(xué)習(xí)結(jié)合進(jìn)行探測(cè)，以探索得到更適用于粗差探測(cè)的方法。