王小娟 趙 亮 蔡 潤 周 坤 尹欣欣
1 甘肅省地震局,蘭州市東崗西路450號,730000
2 中冶成都勘察研究總院有限公司,成都市三色路199號,610023
地震活動目錄是地震臺網(wǎng)運行產(chǎn)生的綜合數(shù)據(jù)庫,主要包括地震位置、震源時間和震級信息。這些數(shù)據(jù)是分析地震活動的空間、時間和規(guī)模分布、地震活動變化以及危險評估的基礎(chǔ)。隨著地方、國家和全球地震監(jiān)測網(wǎng)絡(luò)逐步改善,全世界在發(fā)震地區(qū)記錄到的地震事件越來越多,使得地震目錄得以完善,基于地震目錄的研究在地震學(xué)研究中具有重要意義。研究地震目錄的目的是調(diào)查特定研究區(qū)域內(nèi)地震的時間、空間和規(guī)模分布特征,進(jìn)而估算一些基本地震目錄相關(guān)參數(shù),如古登堡-里克特(G-R)b值和完整性震級(MC)等,是地震學(xué)常用的分析方法[1-2]。在研究局部地震活動時,通常需要將地震進(jìn)行分區(qū)。聚類分析是一種多變量方法,通過將觀測值分組,在數(shù)據(jù)集中搜索相關(guān)模式。該方法的目標(biāo)是找到一個最佳分組,每個組內(nèi)的觀測值或?qū)ο笙嗨?同質(zhì))。然而,這些集群彼此不同(異構(gòu)),數(shù)據(jù)之間的距離決定了數(shù)據(jù)的相似程度。數(shù)據(jù)之間的距離小表示數(shù)據(jù)的高相似性水平,相反,數(shù)據(jù)之間的距離大表示數(shù)據(jù)的低相似性水平[3]。常規(guī)的地震區(qū)劃方法是通過對構(gòu)造特征、地殼特征等進(jìn)行主觀分析,來了解地震特征,這種分析容易出現(xiàn)主觀判斷錯誤[4]。Zamani等[5]提出一種利用層次聚類分析進(jìn)行地震區(qū)劃的替代方法。這種方法是地震區(qū)劃的起點,有可能改進(jìn)和重新定義新的數(shù)據(jù)集。該方法也可用于研究新構(gòu)造、地震構(gòu)造、地震分區(qū)和發(fā)震區(qū)的危險性評估。地震資料是識別地震構(gòu)造區(qū)最重要的信息來源之一,歷史和儀器地震數(shù)據(jù)的模式識別為從大量數(shù)據(jù)中提取有用信息提供了一種更穩(wěn)健、更合適的工具。Ansari等[6]提出基于最大似然估計模糊修正的聚類方法,聚類分析結(jié)果與地震構(gòu)造模型之間的比較表明,如果將地震事件的空間分布(震中)劃分為不同的區(qū)域,主要事件的聚類將取得最佳結(jié)果。
為探索聚類分析方法在地震空間上快速掃描地震群的實用性,本文使用K-means和高斯混合模型(Gaussian mixture model,GMM)兩種方法分別對甘東南及鄰近地區(qū)的原始目錄以及精定位目錄進(jìn)行聚類分析,結(jié)果發(fā)現(xiàn),結(jié)合精定位分析能讓地震聚類更好地將地震劃分成地震叢集,為研究地震群識別、區(qū)域地震危險性判斷提供參考。
甘肅東南及鄰近地區(qū)(32°~36°N,102°~106°E)處于青藏高原東北緣與東緣,構(gòu)造活動強烈,主要有祁連山地震帶和南北地震帶北段。研究區(qū)構(gòu)造上主要受印度板塊北推碰撞歐亞大陸的主動力控制,在往東北向擠壓過程中受到相對穩(wěn)定的阿拉善地塊及鄂爾多斯地塊的阻擋,使得該區(qū)域成為構(gòu)造活動與應(yīng)力場變化的敏感地區(qū)。研究區(qū)發(fā)育多條活動斷裂,地震活動率較高,歷史上曾發(fā)生過多次中強以上破壞性地震,近10 a來主要發(fā)生2013年岷縣MS6.6地震與2017年九寨溝MS7.0地震,造成重大的人員傷亡與財產(chǎn)損失,因此長期以來該地區(qū)一直是地學(xué)工作者研究的熱點區(qū)域之一。本文從中國地震臺網(wǎng)中心下載了研究區(qū)范圍內(nèi)2013-01-01~2021-11-01的地震臺網(wǎng)觀測報告,為保證地震震中位置的可靠性,選取最小臺站記錄數(shù)為6,按照如上挑選條件共下載地震事件15 590次。另外,由于研究區(qū)屬于甘肅、四川、陜西交界處,地震目錄中摻雜重復(fù)地震事件,因此首先去除不同臺網(wǎng)間記錄的邊界重復(fù)地震事件,在處理上以地震位置所在省份為準(zhǔn),只保留各臺網(wǎng)給出的自己省份地震事件,最終得到一共11 659次地震事件,地震震中分布以及時間分布如圖1(a)所示。地震臺網(wǎng)給出的地震目錄震源位置存在一定誤差,通常用地震精定位方法來減小該定位誤差。為提高地震定位精度,使用Waldhauser等[7]提出的雙差定位方法(HypoDD)對原始地震目錄進(jìn)行精定位處理,進(jìn)而對比精定位前后地震聚類分析結(jié)果的變化。在地震對分配過程中,以20 km作為兩個地震對之間的最大距離,超過20 km的不予考慮。由于事件震中分布比較集中,因此設(shè)定單一地震最多可以和20個地震組成地震對,速度模型選取甘青一維速度走時模型[2,8]。震中距范圍選取500 km作為閾值,震相走時經(jīng)過篩選后如圖1(b)所示。數(shù)據(jù)預(yù)處理完成后,使用HypoDD精定位方法進(jìn)行重定位,得到地震目錄6 654條,精定位結(jié)果地震分布見圖1(c)所示。
聚類技術(shù)是檢查數(shù)據(jù)、基于檢查進(jìn)行預(yù)測以及消除其中差異的重要方法。聚類用于識別和分組每天生成的不斷增長的數(shù)據(jù),并生成可以進(jìn)一步利用的模式和知識。K-means算法是數(shù)據(jù)挖掘中最常用的聚類方法之一。該算法以一個聚類的平均值作為聚類的質(zhì)心,通常根據(jù)數(shù)據(jù)之間的歐氏距離對數(shù)據(jù)進(jìn)行分組。K-means聚類分析是一個迭代過程,因為它是一種硬劃分算法,試圖將一組多元數(shù)據(jù)中的n個個體劃分為k個聚類,其中數(shù)據(jù)集中的每個個體分配給一個特定的聚類[9]。
高斯混合模型(GMM)是一種使用廣泛的聚類算法,該方法將高斯分布作為參數(shù)模型,并利用期望最大(expectation maximization,簡稱EM)算法進(jìn)行訓(xùn)練。K-means算法在特定約束條件下可以被看作是GMM的一種特殊形式。GMM是一種軟分類方法,只給出項目是某個類的概率值。利用概率密度函數(shù)表示類別的信息量比一個直接的分類結(jié)果要多。每個GMM由K個高斯分布組成,每個高斯分布稱為一個“component”,這些component線性加成在一起就組成了GMM的概率密度函數(shù):
(1)
根據(jù)式(1),若要從GMM的分布中隨機抽取一個點,可以通過2步來進(jìn)行:首先隨機在這K個component之中選取一個,每個component被選中的概率實際上就是其系數(shù)π(k);選中component之后,再單獨考慮從這個component的分布中選取一個點就可以了——這里已經(jīng)回到了普通的Gaussian分布,轉(zhuǎn)化成已知的問題。接下來就是聚類步驟,根據(jù)樣本數(shù)據(jù),假定這些數(shù)據(jù)是由GMM產(chǎn)生的,那么可以根據(jù)數(shù)據(jù)推算出GMM的概率分布,而GMM的K個component實際上就對應(yīng)了K個cluster。根據(jù)數(shù)據(jù)來推算概率密度通常被稱作概率密度估計。特別地,當(dāng)我們在已知(或假定)了概率密度函數(shù)的形式,而要估計其中的參數(shù)的過程被稱作參數(shù)估計。
根據(jù)對先驗信息的使用情況,分類技術(shù)可以分為監(jiān)督、半監(jiān)督和無監(jiān)督[10],高斯混合模型GMM和K-means屬于無監(jiān)督技術(shù),經(jīng)典的線性和二次辨別分析屬于監(jiān)督技術(shù)。無監(jiān)督學(xué)習(xí)最大的優(yōu)點是不需要先驗信息即可完成數(shù)據(jù)的分類。為選擇最佳聚類個數(shù),學(xué)者們提出許多信息準(zhǔn)則,通常通過加入模型復(fù)雜度的懲罰項來避免出現(xiàn)過擬合問題,本文選擇了常用的兩種模型選擇方法——赤池信息準(zhǔn)則(Akaike information criterion,AIC)和貝葉斯信息準(zhǔn)則(Bayesian information criterion,BIC)。在節(jié)點選擇上,使用了拐點法,如圖2所示,原始地震目錄和精定位目錄數(shù)據(jù)的最佳聚類個數(shù)分別為6和14。AIC和BIC的計算公式如下:
AIC=-2ln(L)+2k
(2)
BIC=-2ln(L)+kln(n)
(3)
式中,L為該模型下的最大似然,n為數(shù)據(jù)數(shù)量,k為模型變量個數(shù)。
根據(jù)AIC和BIC確定的最佳聚類數(shù)結(jié)果,分別對本文研究數(shù)據(jù)原始地震目錄和精定位目錄進(jìn)行聚類分析,在聚類個數(shù)選擇上分別選取6個和14個,具體分布見圖3和4。從結(jié)果上看,K-means聚類的地震群主要呈圓形或塊狀分布,而GMM聚類的地震群主要呈條帶狀分布,與地震主要沿斷層分布的特征相對應(yīng)。相對原始地震目錄,地震精定位結(jié)果對地震的分組更精確,更符合區(qū)域斷層幾何分布特征。研究范圍內(nèi)2次6級以上強震分別是2013年岷縣MS6.6地震與2017年九寨溝MS7.0地震,根據(jù)GMM方法聚類結(jié)果,只有精定位目錄將這兩次地震序列進(jìn)行了區(qū)分,原始目錄的聚類中岷縣MS6.6地震序列因為與周圍地震連接緊密以及最佳聚類數(shù)限制,沒有聚類成單獨的地震群。精定位后九寨溝地震序列空間分布更為收縮,與周邊零散地震鄰近距離拉長,聚類精度更高(圖3(b)、圖4(b)綠色框處)。為了確定本文聚類方法對震群識別的準(zhǔn)確性,對比前人的地震精定位結(jié)果[11-13],本文得到的兩次地震序列震群空間展布形狀均與前人的結(jié)果相一致,這也驗證了本文方法的可靠性。
從兩種聚類方法的原理上看,K-means聚類主要是計算聚類中心到各點的歐氏距離,適合于圓形分布的數(shù)據(jù)集群。地震主要是斷層活動產(chǎn)生,而從已有的斷層數(shù)據(jù)來看,多數(shù)斷層呈條帶狀分布,因此K-means聚類在地震聚類中未能取得理想的結(jié)果。而GMM聚類方法基于概率密度判斷相似性,假設(shè)數(shù)據(jù)點是呈高斯分布(橢圓形)的,這與現(xiàn)實地震條帶狀分布更接近,因此該方法相對于K-means方法來說更適用于地震聚類分析。地震精定位處理通常會將地震的空間位置更好地聚集起來,從而更容易發(fā)現(xiàn)斷層的幾何展布形態(tài)。另外,對相對孤立的地震事件,按照設(shè)定的地震配對標(biāo)準(zhǔn)不符合條件者會予以去除,這些孤立事件散布在各地震群之間,增加了地震空間數(shù)據(jù)模型的復(fù)雜度,使得對最佳聚類數(shù)的選擇存在精定位前后不一致的結(jié)果,進(jìn)而對地震聚類分析產(chǎn)生一定干擾。因此結(jié)合地震精定位和GMM聚類方法能夠提高地震震群識別,這也解釋了形狀更接近圓形的岷縣地震序列在沒有進(jìn)行精定位處理前GMM聚類也無法將其識別區(qū)分的原因。
從本文結(jié)果來看,地震空間數(shù)據(jù)的聚類結(jié)果與現(xiàn)實相吻合,并且可以通過使用聚類算法對地震數(shù)據(jù)進(jìn)行建模來確定區(qū)域內(nèi)地震活動特征。隨著計算機性能的逐步提高,大數(shù)據(jù)高效分析也得到全方位的普及,而聚類分析是一個歷史悠久且易于理解的過程,適用于從大型不相干數(shù)據(jù)集中提取有意義的特征。為了將地震目錄的空間分布劃分地震群,本文使用了K-means和GMM兩種聚類方法,分別對2013~2021年研究區(qū)的11 659條地震目錄及6 654次精定位結(jié)果進(jìn)行聚類分析,獲得如下認(rèn)識:
1)將地震分組為不同的集群可用于改進(jìn)區(qū)域內(nèi)地震活動的機制識別和模式識別。本文展示了地震空間聚類流程,可作為分析地震數(shù)據(jù)集的工具以及用于解釋結(jié)果的可視化。
2)在地震群個數(shù)的選擇中,可以通過計算AIC和BIC來確定最佳聚類數(shù),其中,本文所用數(shù)據(jù)中原始地震目錄和精定位目錄的最佳聚類數(shù)分別為6和14,結(jié)合區(qū)域斷層分布和數(shù)量來看,地震精定位對聚類分析有較好的促進(jìn)作用。
3)GMM聚類的地震群主要呈條帶狀分布,更符合實際地震主要沿斷層分布的普遍情況。根據(jù)本文結(jié)果,精定位地震目錄和GMM聚類方法可以從眾多地震分布中突出強震(如九寨溝MS7.0和岷縣MS6.6地震)地震序列的震群分布。
此外,與其他無監(jiān)督算法(如自組織映射)相比,K-means和GMM兩種無監(jiān)督技術(shù)較為直接和高效。另外,這兩種聚類方法并不是分類不同地震事件的最佳方法,而是可能方法,也適用于其他具有更復(fù)雜數(shù)據(jù)集的地震活動地區(qū)。額外的判別器會導(dǎo)致更多維度的問題,而本研究只處理二維問題,計算簡單、效率高。根據(jù)本文研究結(jié)果,可以在沒有先驗信息的地區(qū)使用無監(jiān)督聚類算法來對地震群進(jìn)行快速識別分析。
隨著地震臺網(wǎng)的進(jìn)一步密集化以及地震預(yù)警臺網(wǎng)的建立,長期連續(xù)觀測的地震目錄與日俱增,本文使用的聚類分析方法在大數(shù)據(jù)處理時具有明顯優(yōu)勢,可以第一時間提取地震的空間分布特征,為實時地震趨勢判斷提供參考。另外,本文提出的地震最佳聚類方法給地震集群以及活動性分析提供了新思路。
致謝:甘肅、四川測震臺網(wǎng)提供地震目錄數(shù)據(jù),本文作圖使用了Python和GMT軟件,在此一并表示感謝。