基于半監(jiān)督核均值漂移聚類的地震相識(shí)別研究

郝茜茜 ，周亞同 ，任婷婷

0 引言

地震相識(shí)別是地震層序劃分的前提，是油氣藏勘探和儲(chǔ)層預(yù)測(cè)的基礎(chǔ)[1].在地震剖面上包含著豐富的地質(zhì)信息，已有很多聚類方法被用于地震相識(shí)別.例如K均值[2]，模糊C均值聚類[3]，DBSCAN聚類[4]，但上述都是一些常規(guī)的聚類方法.

近年又涌現(xiàn)了一些性能優(yōu)良的聚類方法.例如SOM聚類[5]，次勝者懲罰競(jìng)爭(zhēng)學(xué)習(xí)聚類[6]和均值漂移聚類.均值漂移[7]是一種基于非參數(shù)估計(jì)的密度聚類算法，通過(guò)迭代搜索特征空間中的樣本點(diǎn)，使聚類中心始終向密度最大的方向移動(dòng).該算法迭代速度快，無(wú)需人為規(guī)定聚類個(gè)數(shù)且可以對(duì)任何的集群結(jié)構(gòu)聚類，Subbarao等[8]和Vedaldi等[9]實(shí)現(xiàn)了均值漂移在核空間的聚類.此方法被廣泛應(yīng)用于目標(biāo)跟蹤[10]、圖像分割[11-12]、圖像去霧[13]和廣播音頻[14]等方面，但目前為止還未被用于地震相識(shí)別.

半監(jiān)督學(xué)習(xí)是近年來(lái)較受關(guān)注的方法之一，通過(guò)將半監(jiān)督學(xué)習(xí)與一些常規(guī)聚類方法結(jié)合，應(yīng)用少量的先驗(yàn)信息去指導(dǎo)聚類過(guò)程以使聚類結(jié)果更準(zhǔn)確.例如林超[15]通過(guò)對(duì)半監(jiān)督學(xué)習(xí)方法與k均值聚類結(jié)合，解決了算法的約束違反問(wèn)題.Kulis等[16]將半監(jiān)督學(xué)習(xí)與圖譜聚類結(jié)合優(yōu)化了圖形數(shù)據(jù)聚類.尹學(xué)松等[17]先根據(jù)成對(duì)約束得到投影空間后在投影空間進(jìn)行k均值聚類，然后在由線性判別法得到的子空間中再次聚類.Tuzel等[18]和Anand等[19]實(shí)現(xiàn)了半監(jiān)督學(xué)習(xí)與核聚類法的結(jié)合.

核均值漂移聚類不局限于集群結(jié)構(gòu)的類型，適用于復(fù)雜多變的地震數(shù)據(jù)結(jié)構(gòu)，半監(jiān)督學(xué)習(xí)又可以根據(jù)已知的先驗(yàn)信息指導(dǎo)聚類過(guò)程.基于以上考量，本文研究基于半監(jiān)督核均值漂移聚類（SKMS）的地震相自動(dòng)識(shí)別算法.理論數(shù)據(jù)模型和實(shí)際地震數(shù)據(jù)聚類均表明SKMS是一種有效的地震相劃分方法.

1 核均值漂移（KMS）聚類原理

在核均值漂移（KMS）聚類中，歐式空間擴(kuò)展為一般的內(nèi)核空間.令χ為輸入空間，則有n個(gè)樣本xi∈χ，i=1，…，n.假設(shè)空間 χ表示為 Rd，x 通過(guò)映射函數(shù) φl(shuí)，l=1，…，dφ映射到 dφ維特征空間 H，即

定義對(duì)角帶寬矩陣為hiIdφ×dφ，i=1，…，n，y∈H，則在特征空間H中的核密度估計(jì)為

對(duì)公式（2）求梯度可得空間H中的均值漂移向量為

定義ei為特征空間H中的第i維標(biāo)準(zhǔn)基向量，有ei∈Rn，則φ（xi）=Φei，將此式與式（4）代入式（3）中有

從而可以得到均值漂移向量，該方法同樣收斂于局部分布模式.這樣通過(guò)賦予合適的映射函數(shù)就可以實(shí)現(xiàn)核均值漂移聚類.

2 半監(jiān)督均值漂移（SKMS）聚類原理

在核均值漂移算法的基礎(chǔ)上，該算法采用成對(duì)約束來(lái)指導(dǎo)聚類過(guò)程.在聚類之前，需要根據(jù)先驗(yàn)經(jīng)驗(yàn)采集must-link和cannot-link約束組成成對(duì)約束.聚類過(guò)程為：首先將所有點(diǎn)通過(guò)核函數(shù)映射到高維空間中，然后在高維核空間中對(duì)所有成對(duì)約束進(jìn)行線性變換，最后將特征點(diǎn)全部投影到約束向量零空間，使兩點(diǎn)之間距離符合距離目標(biāo)參數(shù)的聚為一類.

2.1 通過(guò)更新核函數(shù)實(shí)現(xiàn)線性變換

對(duì)樣本進(jìn)行變換，就是將樣本從特征空間向約束向量的零空間投影的過(guò)程，此過(guò)程可以通過(guò)更新式（5） 中的核矩陣隱形實(shí)現(xiàn).定義（j1，j）2為成對(duì)約束，表示j1，j2被強(qiáng)制成為一對(duì)，既可以是must-link成對(duì)約束也可以是cannot-link成對(duì)約束，有.若給定nc個(gè)成對(duì)約束NC，dφ維的約束向量可以表示為其中n維向量zj表示為第j個(gè)成對(duì)約束的指標(biāo)向量，則含nc個(gè)約束向量的約束矩陣A=ΦZ，其中Z= [z1，z2，…，znc]為n×nc階指標(biāo)矩陣.定義線性變換矩陣為

其中：s為a的縮放因子.當(dāng)s=1/aTa時(shí)，該變換就變?yōu)閺奶卣骺臻g向約束向量a的零空間投影；當(dāng)0〈s〈2/aTa時(shí)，該變換減小成對(duì)點(diǎn)之間的距離；當(dāng)s〈0或s>2/aTa時(shí)，該變換增大成對(duì)點(diǎn)之間的距離.

令成對(duì)點(diǎn)的距離為d>0，則有

將約束向量a=Φz(mì)代入到公式（8）中得

2.2 通過(guò)logdet布雷格曼散度實(shí)現(xiàn)核函數(shù)更新

布雷格曼散度是一種類似距離度量的方式.若有n×n維的矩陣X和Y，則有布雷格曼散度公式如下

其中：X和Y為半正定矩陣且X和Y的秩不大于n.對(duì)X和Y進(jìn)行奇異值分解，有X=V∧VT，Y=UΘUT，進(jìn)而可以求得logdet布雷格曼散度

現(xiàn)通過(guò)logdet布雷格曼散度實(shí)現(xiàn)核函數(shù)更新.給定m個(gè)must-link成對(duì)約束集M和n個(gè)cannot-link成對(duì)約束集C，有m+c=nc.must-link約束的目標(biāo)距離為dm，cannot-link約束的目標(biāo)距離為dc，則最終的更新核矩陣問(wèn)題轉(zhuǎn)化為求logdet布雷格曼散度最小化問(wèn)題，即目標(biāo)函數(shù)為

由于在logdet布雷格曼散度中的第一個(gè)參數(shù)X要求必須是凸的，用logdet散度最小化更新核矩陣可以保證算法收斂到全局最優(yōu)解.

3 應(yīng)用SKMS進(jìn)行地震相識(shí)別的步驟

用SKMS進(jìn)行地震相識(shí)別的步驟如圖1所示.對(duì)于給定地震數(shù)據(jù)，首先進(jìn)行地震屬性提取和歸一化，然后對(duì)優(yōu)選后的地震屬性采用SKMS聚類法得到地震相識(shí)別結(jié)果.

在圖1中需要用到SKMS聚類，它的具體步驟為：

圖1 用SKMS進(jìn)行地震相識(shí)別的步驟Fig.1 The steps of seismic facies identification using SKMS

步驟1：各參數(shù)初始化.初始化約束距離參數(shù)dm和dc，dm為must-link約束初始距離，dc為cannot-link約束初始距離；初始化成對(duì)約束集M和C；

步驟2：計(jì)算初始核矩陣k；

步驟3：采用logdet散度更新核矩陣k?；

步驟4：對(duì)于地震屬性的n個(gè)樣本點(diǎn)i=1，2，…，n

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

為了驗(yàn)證SKMS聚類效果，分別對(duì)理論模型和實(shí)際地震數(shù)據(jù)進(jìn)行處理，將處理結(jié)果與k均值聚類、核k均值聚類（KK聚類）、譜聚類[20]、均值漂移聚類（MS聚類）、全局核k均值聚類（GKK聚類）[21]和自組織神經(jīng)網(wǎng)絡(luò)聚類（SOM聚類）等算法對(duì)比.

4.1 理論數(shù)據(jù)模型實(shí)驗(yàn)

在地質(zhì)結(jié)構(gòu)中通常含有褶皺、尖滅和套疊等結(jié)構(gòu)單元.首先理論模擬出這3種結(jié)構(gòu)，并采用SKMS聚類，然后將這3種理論結(jié)構(gòu)模型采用上述5種算法對(duì)之聚類.聚類結(jié)果如圖2～圖4所示.

實(shí)驗(yàn)1：褶皺結(jié)構(gòu)如圖2所示.

圖2 褶皺結(jié)構(gòu)的各種算法聚類結(jié)果Fig.2 The clustering results of various algorithms for folds

實(shí)驗(yàn)2：尖滅結(jié)構(gòu)如圖3所示.

圖3 尖滅結(jié)構(gòu)各種算法聚類結(jié)果Fig.3 The clustering results of various algorithms for pinch-out

實(shí)驗(yàn)3：套疊結(jié)構(gòu)如圖4所示.

圖4 套疊結(jié)構(gòu)的各種算法聚類結(jié)果Fig.4 The clustering results of various algorithms for telescope

表1 理論數(shù)據(jù)模型的詳細(xì)信息Tab.1 Detailed information on theoretical data models

表2 各算法運(yùn)行時(shí)間對(duì)比Tab.2 Comparison of running time of each algorithm

表3 各聚類算法準(zhǔn)確率對(duì)比（CA）Tab.3 Accuracy ratio of each algorithm （CA）

從表2看出，SKMS總體比KK、MS、K的運(yùn)行時(shí)間長(zhǎng)，但與GKK相比，總體運(yùn)行時(shí)間短.在樣本個(gè)數(shù)為202時(shí)，譜聚類的時(shí)間比SKMS長(zhǎng)，隨著樣本個(gè)數(shù)的增長(zhǎng)，SKMS比譜聚類的運(yùn)行時(shí)間的增長(zhǎng)速度快.對(duì)聚類結(jié)果的評(píng)價(jià)標(biāo)準(zhǔn)，本文采用準(zhǔn)確性（cluster accuracy，簡(jiǎn)稱CA）度量.從表3可以觀察到：SKMS的聚類準(zhǔn)確性明顯優(yōu)于其它的幾種算法，SKMS在不同的理論模型上的準(zhǔn)確性均達(dá)到了90%以上.而且SKMS在計(jì)算時(shí)均能正確估計(jì)聚類個(gè)數(shù).

4.2 實(shí)際地震數(shù)據(jù)實(shí)驗(yàn)

本實(shí)驗(yàn)采用荷蘭北海F3地震數(shù)據(jù).在opendtect6.0和Matlab2013a軟件平臺(tái)下采用傾角控制中值濾波后的地震數(shù)據(jù)進(jìn)行運(yùn)算.選取聯(lián)絡(luò)測(cè)線1 000縱剖面，剖面范圍為主測(cè)線的450-550道，時(shí)間線time的1724-1820部分，該縱剖面存在明顯的波形反射構(gòu)型.

從地震屬性中選取相關(guān)系數(shù)較小的瞬時(shí)振幅、瞬時(shí)頻率、瞬時(shí)相位屬性，將這3種地震屬性作歸一化處理.應(yīng)用這3種屬性，采用SKMS對(duì)屬性聚類，選取了3類共30個(gè)點(diǎn)如下所示，將這些點(diǎn)組成成對(duì)約束并采用SKMS聚類，效果如圖5所示.

圖5 SKMS聚類效果Fig.5 SKMS clustering result

在聯(lián)絡(luò)剖面1 000中，共標(biāo)記了3類數(shù)據(jù)，如圖5a）所示，因截取的剖面數(shù)據(jù)為25×101，即共有2 525個(gè)地震數(shù)據(jù)，所以標(biāo)簽數(shù)據(jù)占總地震數(shù)據(jù)數(shù)的1.19%，在SKMS聚類中，可以生成個(gè)mustlink成對(duì)約束，然后再構(gòu)造同樣個(gè)數(shù)的cannot-link成對(duì)約束，選取高斯核函數(shù)，其中σ為0.5，SKMS聚類將地震相劃分為6類，如圖5b）所示.

同樣應(yīng)用上述3種屬性，分別采用不同的方法聚類，結(jié)果如下所示.

在KK聚類、譜聚類和GKK聚類中，均采用同一個(gè)高斯核函數(shù)，σ取0.5，聚類個(gè)數(shù)設(shè)為6類，如圖6 a） ～圖6 c） 所示.在MS聚類中，帶寬值為2.5，聚類結(jié)果為6類，如圖6 e） 所示.k均值聚類結(jié)果如圖6 f）所示.SOM聚類為商業(yè)軟件的聚類效果.SKMS與其他算法的聚類效果相比，層次更分明，層與層之間的邊界清晰，在綠圈范圍內(nèi)，能夠?qū)⒁恍┪⑿〉貙訁^(qū)分出來(lái).

圖6 各算法聚類效果Fig.6 clustering results of each algorithm

5 結(jié)論

本文采用SKMS對(duì)地震屬性聚類，利用已知的少量先驗(yàn)信息對(duì)聚類過(guò)程約束，達(dá)到提高地震相劃分結(jié)果的精確性的目的.均值漂移屬于密度聚類，能夠自動(dòng)優(yōu)化判斷聚類個(gè)數(shù)，且可以對(duì)任意數(shù)據(jù)結(jié)構(gòu)都有效.SKMS聚類法融合了MS聚類和半監(jiān)督聚類的優(yōu)勢(shì)，將SKMS聚類與其他聚類算法相比，比無(wú)監(jiān)督聚類的聚類結(jié)果準(zhǔn)確度有了很大提高，劃分地震相的層次更分明.

[1] Robert E.Sheriff.Structural interpretation of seismic data[M].American：American Association of Petroleum Geologists，1982：14.

[2] 龐銳，魏嘉.利用K均值聚類方法進(jìn)行地震相識(shí)別[C]//臧紹先.中國(guó)地球物理學(xué)會(huì)第二十四屆年會(huì)論文集.北京：中國(guó)地球物理學(xué)會(huì).2008：132.

[3] 張陽(yáng)，邱隆偉，李際，等.基于模糊C均值地震屬性聚類的沉積相分析[J].中國(guó)石油大學(xué)學(xué)報(bào)自然科學(xué)版，2015，39（4）：53-61.

[4] 楊瑞超.DBSCAN算法在地震相劃分中的應(yīng)用[D].西安：西安科技大學(xué)，2011：1-43.

[5] 張龔，鄭曉東，李勁松，等.基于SOM和PSO的非監(jiān)督地震相分析技術(shù)[J].地球物理學(xué)報(bào)，2015，58（9）：3412-3423.

[6] Zhan Shifan，Lei Li，Wei Xiong，et al.Automatic geological body identification using the modified rival penalized competitive learning clustering algorithm[C]//Seg Technical Program Expanded.USA：Society of Exploration Geophysicists.2011：4424

[7] Cheng Yizong.Mean Shift，Mode seeking，and clustering[J].Pattern Analysis&Machine Intelligence IEEE Transactions on，1995，17（8）：790-799.

[8] Subbarao R，Meer P.Nonlinear mean shift for clustering over analytic manifolds[C]//Jean-Philippe Tardif.IEEE Computer Society Conference on Computer Vision and Pattern Recognition.USA：IEEE Computer Society，2006：1168-1175.

[9] Vedaldi Andrea，Soatto Stefano.Quick shift and kernel methods for mode seeking[M].France：Computer Vision-ECCV 2008，2008：705-718.

[10]馬麗，常發(fā)亮，喬誼正，等.基于改進(jìn)的均值漂移算法的目標(biāo)跟蹤[J].計(jì)算機(jī)工程，2006，32（24）：175-177.

[11]伍艷蓮，趙力，姜海燕，等.基于改進(jìn)均值漂移算法的綠色作物圖像分割方法[J].農(nóng)業(yè)工程學(xué)報(bào)，2014，30（24）：161-167.

[12]白培瑞，李良，趙奇，等.基于均值漂移的醫(yī)學(xué)超聲圖像分割改進(jìn)算法[C]//中國(guó)智能自動(dòng)化會(huì)議.南京：中國(guó)自動(dòng)化學(xué)會(huì)，2009：1426-1431.

[13]陸海俊，汪榮貴，楊娟，等.基于均值漂移的暗原色先驗(yàn)圖像去霧算法[J].合肥工業(yè)大學(xué)學(xué)報(bào)自然科學(xué)版，2016，39（9）：1205-1210.

[14] 鄭繼明，俞佳.基于 Mean-Shift的廣播音頻聚類算法[J].計(jì)算機(jī)應(yīng)用，2009，29（10）：2741-2743，2750.

[15]林超.基于成對(duì)約束的半監(jiān)督聚類算法研究及其并行化實(shí)現(xiàn)[D].西安：西南交通大學(xué)，2013：1-51.

[16]Kulis Brian，Basu Sugato，Dhillon Inderjit，et al.Semi-supervised graph clustering：a kernel approach[J].Machine Learning，2009，74（1）：1-22.

[17]尹學(xué)松，胡恩良，陳松燦.基于成對(duì)約束的判別型半監(jiān)督聚類分析[J].軟件學(xué)報(bào)，2008，19（11）：2791-2802.

[18]Tuzel O，Porikli F，Meer P.Kernel methods for weakly supervised mean shift clustering[C]//IEEE，International Conference on Computer Vision，ICCV 2009，Kyoto，Japan，September 27-October.DBLP，2009：48-55.

[19]Anand Saket，Mittal Sushil，Tuzel Oncel，et al.Semi-supervised kernel mean shift clustering[J].IEEE Transactions on Pattern Analysis&Machine Intelligence，2014，36（6）：1201-15.

[20]Choromanska A，Jebara T，Kim H，et al.Fast spectral clustering via the Nystr?m method[M].Germany：Algorithmic Learning Theory.Springer Berlin Heidelberg，2014：367-381.

[21]Chen W Y，Song Y，Bai H，et al.Parallel spectral clustering in distributed systems[J].IEEE Transactions on Pattern Analysis&Machine Intelligence，2011，33（3）：568-586.