張麗娟 萬永革 王福昌 靳志同 崔華偉

1)防災(zāi)科技學(xué)院，三河 065201 2)河北省地震動力學(xué)重點實驗室，三河 065201 3)山東省地震局，濟南 250102

0 引言

活動斷層的幾何形狀在地球動力學(xué)以及地質(zhì)構(gòu)造研究中具有非常重要的作用(王鳴等，1992； 萬永革等，2008)。針對發(fā)震斷層面的研究發(fā)現(xiàn)斷層面的空間幾何分布特征往往非常復(fù)雜，且介質(zhì)的性質(zhì)在大范圍內(nèi)具有非均勻性和地質(zhì)構(gòu)造多變性的特點，一般而言并非一個確定的簡單平面。對于出露在地表之上的斷層，可以測量出露斷層面的幾何形態(tài)以對地震深部的幾何形態(tài)進行推斷(Xuetal.，2002，2009)。然而，來自地表的外營力，如地面的剝蝕、 雨水的沖刷作用、 地面生物生長等會給斷層地質(zhì)調(diào)查帶來很大的不確定性。對于隱伏斷層，可通過開挖探槽等地質(zhì)調(diào)查的方法厘定其空間位置(冉勇康等，1997，2018； 張培震等，2003)。然而，地質(zhì)調(diào)查僅能得到淺部的斷層形態(tài)，而斷層深部和淺部的形態(tài)有時差別較大，如張先康等(2002)通過深部地震反射發(fā)現(xiàn)1679年的三河-平谷大地震的深部與淺部構(gòu)造存在巨大差別。人工地震測深是一種確定深部活動斷層幾何形態(tài)的方法。然而，該手段需要使用大量人力物力對深部地震剖面進行探測，且地球物理反演的多解性也會給資料解釋帶來諸多不確定性?？紤]到地質(zhì)調(diào)查只能看到地表的破裂，對于深部破裂信息可通過余震分布來確定(萬永革等，2008)。

通常地震常發(fā)生在預(yù)先存在或新激活的斷層面上(Waldhauseretal.，2002； Waldhauser，2009； Rossetal.，2020)，因而最直接通常也最為有效的方法是通過對震源位置進行精準定位從而刻畫斷層面的精細結(jié)構(gòu)。王鳴等(1992)對1989年大同地震中部分余震目錄的空間分布特征和斷層面的空間結(jié)構(gòu)進行了研究和討論。 其采用最小二乘法基于余震分布求解斷層面，并用高斯-牛頓法進行了迭代，然而這種方法對初始參數(shù)的選擇比較敏感，且研究結(jié)果也并未給出斷層面的空間位置等參數(shù)。王福昌等(2010)為避免高斯-牛頓迭代算法過程中的缺陷并提高運算速度，采用主成分分析的數(shù)值分析方法對斷層面的傾角和走向進行了分析和討論，但依然存在未解決的問題，如地震斷層面的各參數(shù)估計值可能會受離群震源點的嚴重影響，從而導(dǎo)致估計結(jié)果失真。有學(xué)者通過建立模型并考慮定位過程中的誤差，基于應(yīng)力場及反演理論等地震學(xué)有關(guān)知識，采用模擬退火算法對1976年唐山地震進行了討論(萬永革等，2008)，給出了斷層面的滑動角、 傾角及走向等參數(shù)，并進行了誤差分析。以上方法一般只能基于存在一個主斷層面的前提下進行擬合，而不能對含有多個子斷層面的復(fù)雜地震破裂結(jié)構(gòu)進行處理。有學(xué)者對華北地區(qū)2001—2009年的地震進行精確定位(張廣偉等，2011)，給出了小地震分布的一些分布規(guī)律，并采用小地震空間分布估算得出斷層面的各個參數(shù)。該方法未明確給出劃分小地震群的數(shù)學(xué)模型和具體實現(xiàn)方法，從側(cè)面對萬永革等(2008)所給出方法的實用性及有效性進行了驗證。Ouillon等(2008)雖然給出了一種比較繁瑣的劃分小地震的方法，并將其應(yīng)用于Landers地震的研究中，但其計算速度仍然較慢。王福昌等(2013)在主成分分析的基礎(chǔ)上利用幾何模糊聚類算法給出了同時推測多個子斷層面的統(tǒng)計聚類算法，極大提高了計算速度和地震斷層面精度。

目前關(guān)于地震斷層面走向等空間精細的定量描述大多給出的都是主要斷層的走向，不能同時估計多個斷層面的空間分布，對于地震斷層面空間位置的刻畫依然有很大的改進空間。隨著數(shù)字地震臺站在中國大陸的廣泛布設(shè)和雙差定位方法(Waldhauseretal.，2002； Waldhauser，2009)及主事件定位方法(周仕勇等，1999)的提出，測定的震源數(shù)據(jù)精確度越來越高，能夠測定的地震震級下限也越來越低。在地震數(shù)量急劇增加的實際背景下，地震震源目錄資料中蘊含的大量關(guān)于地震所處斷層的幾何特征信息，給精準的空間定位的數(shù)學(xué)算法提供了很好的使用機會。萬永革(2022)基于密度的空間聚類方法對2021年云南漾濞地震序列的震源機制數(shù)據(jù)進行了處理，得到2簇節(jié)面中心與地震序列的空間密度線性分布較為一致的分支斷層，提供了獨立于地震波、 大地測量、 地質(zhì)資料等之外的確定斷裂帶形狀分布的方法，但并沒有對聚類結(jié)果來自于哪個地震進行統(tǒng)計分析，也未給出斷層面的空間定位。本文采用王福昌等(2013)提出的模糊聚類算法，對云南漾濞地震的震源數(shù)據(jù)進行聚類分析，并對斷層面區(qū)域進行了最優(yōu)化處理，在符合假設(shè)條件的前提下，給出了一系列最優(yōu)解，并通過一系列研究對地震序列的幾何形態(tài)進行精準推測和評估。

1 資料與算法

1.1 資料與數(shù)據(jù)來源

中國地震臺網(wǎng)正式測定，2021年5月21日21時48分云南大理白族自治州漾濞彝族自治縣發(fā)生了M6.4地震，震中位置位于(25.67°N，99.87°E)(圖 1)。震后，云南省地震局立即成立指揮部并啟動三級響應(yīng)，分配工作隊奔赴現(xiàn)場，開展震情監(jiān)視與災(zāi)害損失評估等工作。災(zāi)區(qū)附近的電力、 交通、 水利等基礎(chǔ)設(shè)施均受到了一定程度的影響。經(jīng)初步統(tǒng)計，地震有感范圍約 1.2 萬平方千米，大理市地區(qū)及漾濞周邊各縣普遍有震感。截至2021年5月27日18時，漾濞地震序列已致全縣共15786人受災(zāi)，共造成經(jīng)濟損失約2.021億元。

圖 1 云南漾濞地震序列的精定位數(shù)據(jù)分布圖Fig. 1 The relocation results of the Yunnan Yangbi earthquake sequence.云南漾濞的精地震序列定位用圓圈表示； 5級地震用紅色五角星表示； 主震用黃圓圈表示。圖中用不同顏色的圓圈表示地震深度。 黑色實線F1為維西-喬后斷裂。虛線多邊形為研究區(qū)所在的主要區(qū)域，橢圓形區(qū)域內(nèi)為第2類數(shù)據(jù)區(qū)域

地震發(fā)生后，許多國內(nèi)外研究者對此次地震的震源機制進行了研究(龍鋒等，2021； 余海琳等，2021)，對余震進行了精準定位(Yangetal.，2021)，并開展了GNSS與InSAR破裂滑動分布特點及同震變形場特征(于書媛等，2021； 張克亮等，2021)、 震源區(qū)電性結(jié)構(gòu)特征(葉濤等，2021)、 震源區(qū)地殼速度結(jié)構(gòu)(杜廣寶等，2021； 李大虎等，2021)等研究，為后續(xù)的防震減災(zāi)及地震研究等領(lǐng)域提供了基礎(chǔ)數(shù)據(jù)。位于維西-喬后及巍山附近的一條隱伏的斷裂帶(雷興林等，2021)是本次地震的發(fā)震斷層。此次云南漾濞地震主要在四川與云南塊體偏西的邊緣區(qū)域孕震，為該區(qū)有記載以來發(fā)生的最大地震。本次地震以走滑運動為主，主震之后余震活動頻繁，其地震序列發(fā)生在地震監(jiān)測臺網(wǎng)較為密集的區(qū)域，通過密集的地震觀測臺站可以較好地確定震級下限較低的地震目錄。梁姍姍等(2021)采用雙差定位方法對此次地震開展了研究，結(jié)果顯示地震的前震和余震序列存在顯著的時空差異性。其中，前震呈條帶狀，主要分布在NW-SE向，且前震震中的空間位置還存在遷移和往返的特點。而余震呈現(xiàn)出不對稱和共軛形態(tài)，主要分布在NW和NE向(圖 1)。地震帶南端NW向還有多組余震沿不同走向分布，揭示了震源區(qū)的介質(zhì)特性和斷層幾何結(jié)構(gòu)方面的復(fù)雜性。這種復(fù)雜特征不能僅模擬為單獨一條斷層，這為本文采用地震精定位結(jié)果同時求解多條斷層的方法應(yīng)用提供了絕好的機會。

本文基于中國地震臺網(wǎng)中心給出的(99°～101°E，24.5°～26.5°N)空間范圍內(nèi)的觀測報告結(jié)果，使用Waldhauser等(2002)提出的雙差定位法，對2021年5月18日—7月5日云南漾濞共7092個地震事件的101806條Pg和Sg到時數(shù)據(jù)進行了處理，對這些地震事件進行重定位。定位時，設(shè)置Pg的權(quán)重為1.0，Sg的權(quán)重為0.8。用10km作為相鄰地震事件的最大搜索閾值，相鄰地震事件的數(shù)量最多為20個，嚴格限制地震事件，并選擇至少6個震相記錄實現(xiàn)定位，最終得到了5827個地震事件的精確定位結(jié)果(圖 1)。其中，90.23%的事件被重新定位，其結(jié)果見表1。

表 1 雙差重定位使用的一維速度模型Table1 Double-difference relocation results by 1-D velocity model

定位結(jié)果顯示，漾濞地震序列整體呈NW-SE向分布，長約25km，位于維西-喬后斷裂西南側(cè)3～10km處，整體走向約為135°；M6.4主震震中位于(25.688°N，99.877°E)，震源深度約9.6km，M>5.0地震3次，分布在余震區(qū)的不同部位，余震區(qū)中段的最大前震為M5.6，北西端的主震為M6.4，余震區(qū)南東段的最大余震為M5.2。綜合地震序列深度剖面和震源機制解結(jié)果可知，發(fā)震斷裂應(yīng)為NW走向，并且可見整個序列地震的密集區(qū)域主要限定在主震的SE側(cè)，僅有少部分的地震分布在NW側(cè)，整體上空間展布呈現(xiàn)北西較窄、 南東較寬的特點。因此，可認為單側(cè)破裂特征較為明顯，主震方向為SE(雷興林等，2021)。此外，地震序列在空間結(jié)構(gòu)整體散落的疏密度上也有顯著差異性： 南東段相對較為分散，而北西段則更加密集，且在離開序列主干部分的北西段存在2個小地震叢集，似乎有與整體序列的NW-SE走向相偏離的趨勢(龍鋒等，2021)。這些復(fù)雜的地震序列分布為同時求解多斷層的參數(shù)提供了應(yīng)用契機。

1.2 模糊聚類算法

本文采用GK模糊聚類算法(王福昌等，2013)，使用自適應(yīng)距離對模糊算法進行了修改，每個子類對應(yīng)距離矩陣Ai用于描述數(shù)據(jù)的局部結(jié)構(gòu)。在FCM算法中將目標函數(shù)修改為GK-FCM算法的目標函數(shù)：

(1)

(2)

在給定的聚類數(shù)目下，滿足10標記終止誤差值，符號A表示距離矩陣取值，采用隨機或指定方式對初始矩陣進行劃分U(0)，按照以下步驟進行處理：

第1步，計算類中心：

(3)

第2步，計算模糊協(xié)方差矩陣：

(4)

并加入調(diào)整單位陣：

Fi=(1-γ)Fi+γdet(F0)1/pI

(5)

Fi=[φi1，φi2，…，φip]diag(λi1，λi2，…，λip)(φi1，φi2，…，φip)-1

(6)

第3步，計算距離：

(7)

第4步，更新劃分矩陣：

(8)

直至‖U(l)-U(l-1)‖<ε使目標函數(shù)最優(yōu)，從而得到每個類的中心位置以及每個樣本點對于每個類的隸屬程度，把所有滿足隸屬度大于閾值T的震源點全部歸并到同一類中。

利用不同的隸屬度取值和數(shù)據(jù)集的集散程度構(gòu)造有效函數(shù)，用有效函數(shù)最優(yōu)的方法確定初始聚類數(shù)，再結(jié)合主成分分析給出95%置信水平的子斷層面長度和寬度、 斷層面走向角和傾角等參數(shù)的估計值。

2 確定2021年云南漾濞地震序列的斷層面

從2021年云南漾濞地震序列的雙差定位結(jié)果在地表的投影圖(圖 1)可以看出，虛線多邊形區(qū)域和實線橢圓形區(qū)域地震序列明顯屬于不同的地震斷層，將2部分數(shù)據(jù)分開進行處理。處理實線區(qū)域內(nèi)的數(shù)據(jù)時，誤差值取ε=10-6。 處理虛線多邊形區(qū)域內(nèi)的數(shù)據(jù)時，取聚類個數(shù)為3，將實線橢圓形區(qū)域單獨處理，聚為一類。

圖 2 確定的斷層面和地震位置的平面圖(a)和立體圖(b)Fig. 2 Determined fault planes and hypocenters in map view(a)and in 3-D(b).方框為確定的斷層面； 圓圈為聚類所用的屬于斷層面的數(shù)據(jù)； 藍色叉號為聚類過程中判定的離群事件

圖 2 中用叉號表示離群的震源點，其余的各子類震源點以及未被指定給任何子類的震源點均采用圓圈標記。圖中還給出了漾濞地震序列各震源點的空間分布及斷層面在地面上的投影情況。用主成分方法對這些屬于某條斷層的事件數(shù)據(jù)進行分析并確定矩形特點的斷層面，將各子斷層投影至地表上，均為四邊形，并用A、 B、 C對其進行標記。沿NW-SE走向斷裂的3個斷裂面的走向分別為134.22°、 132.65°、 129.45°，斷裂面的傾角分別為87.14°、 81.96°、 74.77°，斷層面的范圍與余震分布能夠較好吻合。

本研究采用模糊聚類算法識別出3條斷裂面，與Yang等(2021)推斷的“漾濞地震破裂的東南端有至少3條斷裂”的結(jié)果相吻合。不同之處是Yang等(2021)識別的第3條斷裂是由西南側(cè)的一簇地震導(dǎo)致的，但那一簇地震并未形成條狀地震帶，究竟其是不是一條斷層，還需要通過其他地球物理探測手段證實，而本文識別出的3條斷裂帶中的2條位于主要斷裂的方向。萬永革(2022)采用震源機制節(jié)面聚類方法對2021年云南漾濞地震序列中的震源機制節(jié)面進行聚類，識別出2條斷裂面。斷裂面的分布與維西-喬后斷裂帶的方向幾乎平行，且向SE按照地震分布的密集程度逐漸分為2支斷裂。西側(cè)分支對應(yīng)本文的斷裂面A、 B，應(yīng)為主斷裂； 東側(cè)分支對應(yīng)本文的斷層面C，與本文的研究結(jié)果基本一致。

第2類數(shù)據(jù)(圖 1 中的實線橢圓形區(qū)域)位于維西-喬后斷裂及巍山斷裂附近，與主斷裂面斜交，小地震的叢集深度為7～10km。該區(qū)域共有276條數(shù)據(jù)，分支看起來較為明顯，但其面狀結(jié)構(gòu)不明顯，且震源點較少，故對其進行單獨處理。采用模糊聚類算法識別為“斷層面”。該“斷層面”的走向為235.66°，傾角為66.30°，傾向NW(圖 2)。從斷層面位置的分布來看，這簇地震帶位于維西-喬后斷裂帶附近，需要指出的是仍需通過其他地球物理探測手段確定該分布究竟是否為一條地震斷層面，本文稱其為偽斷層面。由位于維西-巍山-喬后斷裂附近的小地震集中活動，推測該地區(qū)后期存在持續(xù)地震活動的可能。

表2和圖3總結(jié)了前人得到的斷層面形狀結(jié)果，包括基于InSAR、 GPS及地震波數(shù)據(jù)等不同資料利用不同方法得到的結(jié)果。前人所給出的斷層面傾向既有SWW的結(jié)論，也有NEE的結(jié)果。從本文得到的3個主要斷裂面的結(jié)果來看，可得到以下結(jié)論： 3個斷裂面總體幾乎垂直于水平面，與已有結(jié)果較為吻合，這驗證了本研究結(jié)論的正確性。

表 2 不同研究組或機構(gòu)給出的漾濞地震斷層面的幾何參數(shù)Table2 The seismic fault planes with its geometry parameters determined by the previous researches

圖 3 其他作者得到的斷層幾何形狀與本研究結(jié)果的比較Fig. 3 Comparison of the research results of other researches with that of this paper.其他作者給出的斷層形狀采用點弧線表示，其中心形狀采用粗弧線表示，實弧線為本研究通過聚類分辨的3個斷層面的形狀。黑點表示其他作者給出的斷層面極點位置

為更好地分析研究區(qū)所在的主要區(qū)域的斷層面走向，取數(shù)據(jù)區(qū)域(圖 1 中的虛線多邊形所圍成的區(qū)域)的聚類個數(shù)為4，得到如下結(jié)果(圖 4)： A、 B、 C、 D 4個斷層面的走向分別為146.01°、 135.55°、 112.09°、 134.994°，傾角分別為84.03°、 88.03°、 82.77°、 83.82°。從結(jié)果來看，B、 D斷層面為主斷層面，與本文前部分主斷層面的方向較為接近，從一定意義上驗證了漾濞地震的主斷層面方向具有一致性，且從另一側(cè)面驗證了算法對主斷層識別的可靠性。在東側(cè)數(shù)據(jù)密集的區(qū)域識別出的2個斷層面的走向與前文虛線多邊形所圍研究區(qū)所在的主要區(qū)域的斷層面走向略有差異?？紤]到大地震破裂一般具有優(yōu)勢斷層斷裂面，這會導(dǎo)致余震發(fā)生的微小裂紋至少有一部分與優(yōu)勢斷層的斷裂方向具有一致性。本研究中的模糊聚類算法的主要目的之一是精確刻畫地震目錄所呈現(xiàn)的優(yōu)勢斷層面的空間分布，并探討其他斷層面的空間分布。即使某一部分地震數(shù)據(jù)識別出“斷層面”也不一定表明其為真正的地震斷裂面。這是由于主震破裂導(dǎo)致鄰域內(nèi)的應(yīng)力場發(fā)生較大改變，致使相同幾何形狀裂紋的錯動方向有較大的變化范圍，從而導(dǎo)致地震斷裂面的形狀出現(xiàn)較大的變化范圍。

圖 4 按4個斷層面處理的主破裂地震數(shù)據(jù)確定的斷層面和地震位置的平面圖(a)和立體圖(b)Fig. 4 Determined fault planes and hypocenters in map view(a)and in 3-D(b) in main rupture area and 4 fault planes.圖例同圖3

3 結(jié)論與討論

本研究結(jié)合主成分分析，利用GK模糊聚類方法給出了2021年云南漾濞地震活動斷層的網(wǎng)絡(luò)三維空間結(jié)構(gòu)。以簇來標記密度足夠高的區(qū)域，依據(jù)成叢的小地震在斷層面附近分布且在附近密集的原則，得到了與該地震序列分布走向大體一致的3個主要分支斷層面和1個隱伏斷裂面。3個主斷層面的走向分別為134.22°、 132.65°、 129.45°，傾角分別為87.14°、 81.96°、 74.77°，所得結(jié)果與地震震源數(shù)據(jù)的一致性較高。由研究結(jié)果推測2021年云南漾濞地震主要斷層面的分布與維西-喬后斷裂近平行，并逐漸分為3個分支，西側(cè)分支的傾向大體為SWW，對應(yīng)2個主斷裂面，并向SE逐漸分為3個斷裂分支。在維西-喬后斷裂帶附近還有一條斷層，走向為235.66°，傾角為66.30°，稱之為偽斷層。

本研究中使用的模糊聚類法適用于分支斷裂面比較復(fù)雜的斷裂帶，可同時識別多條斷層面，是基于大量地震震源斷裂面數(shù)據(jù)進行處理的一種聚類分析算法，從理論上可能得到統(tǒng)計意義的空間斷裂帶形狀。使用該方法分析地震斷層面時，如能配合地震精定位的地震形態(tài)展布、 地質(zhì)等其他手段進行綜合分析，將顯著增強研究結(jié)果的應(yīng)用價值。本文使用模糊聚類法在漾濞地震序列中得到了4個地震斷層面，位于維西-喬后斷裂帶附近的斷層面與主要斷裂面的差異較大，考慮可能是地震發(fā)生時與地下斷層和地質(zhì)結(jié)構(gòu)存在某種關(guān)系。同樣，必須指出的是，本文所述方法在云南漾濞地震序列研究中得出的主斷層面的空間位置、 2個斷層面的傾向分布、 2條斷裂的分支斷裂之間的連接性及其之間如何連接等相關(guān)科學(xué)問題仍須借助如地球物理探測技術(shù)等其他手段作進一步研究和考證。

另外，本文所得的結(jié)果難以評估地震對結(jié)果的貢獻度，也不能分辨哪個地震將對結(jié)果造成主要影響，而是根據(jù)每個震源點的數(shù)據(jù)隸屬度按照模糊聚類原則整體確定的斷裂面。算法程序在普通計算機上的運行速度已經(jīng)較快，且可以實現(xiàn)在短時間內(nèi)進行多次重復(fù)計算，如果再考慮隨機給出初始劃分矩陣，得到的分類結(jié)果可能會略有不同，可根據(jù)其他探測技術(shù)、 地質(zhì)條件以及地質(zhì)結(jié)構(gòu)進行分析討論，以得到更加準確的結(jié)果，這一點在王福昌等(2013)關(guān)于Landers地震的處理過程中也曾有描述。這種方法可以很好地依據(jù)余震的空間分布特點確定大地震子斷層空間的定位及其參數(shù)，進而可成功解讀斷層的幾何學(xué)特征、 地震序列活動的活動范圍，并可在無需人工指定的前提下自動完成對地震斷層面的中心位置以及整個斷層面空間方位的刻畫。若根據(jù)小地震分布預(yù)先通過人工給定一些斷層面的中心并構(gòu)造初始劃分矩陣，收斂速度將會更快，并得到更好、 更加精準的定位結(jié)果。

文中所提模糊聚類求解斷裂帶形狀的方法，快速計算分析得到與地震分布大體一致的分支斷裂帶的形狀及其空間分布，進一步驗證了該方法在地震斷層面空間定位的研究領(lǐng)域具有一定應(yīng)用價值和參考意義，可將這種基于數(shù)據(jù)分析的算法應(yīng)用于對地震斷裂帶的形狀分析和預(yù)判。在地震發(fā)生后進行災(zāi)害評估及開展地震災(zāi)前預(yù)測時，對斷層結(jié)構(gòu)的網(wǎng)絡(luò)性研究結(jié)果將對地震附近區(qū)域的應(yīng)力傳遞檢驗及相關(guān)假設(shè)提供幫助(Wanetal.，2010)。

致謝本文的繪圖使用GMT及MATLAB軟件繪制完成； 審稿專家為本文提出了建設(shè)性修改意見，對稿件的質(zhì)量提升幫助很大。在此一并表示感謝！