余 娜 張曉清 袁伏全 趙燕杰
(中國西寧 810000 青海省地震局)
近年來,基于統(tǒng)計物理學(xué)與地震學(xué)的圖像信息學(xué)(Pattern Informatics,簡稱PI)算法在地震較活躍地區(qū)的中長期預(yù)測研究中得到了廣泛應(yīng)用。Rundle 等(2002)將PI 算法應(yīng)用到美國南加州地區(qū)地震的預(yù)測研究中,認(rèn)為該方法對加州地區(qū)地震活動有較好的預(yù)測效果;國內(nèi)一些研究者也將PI 算法應(yīng)用到中國臺灣(Chen et al,2005)、中國大陸西部(蔣長勝等,2008;Zhang et al,2010,2013;Jiang et al,2010,2011,2013;張小濤等,2014;袁伏全等,2016)等地區(qū)地震的研究中,認(rèn)為PI 算法的預(yù)測效果顯著優(yōu)于隨機概率法。Zhang 等(2013)和張小濤等(2014)利用PI 算法分別對2008 年汶川8.0 級、2014 年于田7.3 級、2013 年蘆山7.0 級地震進(jìn)行地震危險性的回溯性研究,結(jié)果表明模型計算參數(shù)的選取對PI 算法的預(yù)測結(jié)果有一定影響,對于7 級以上大地震,選擇較大的網(wǎng)格尺度和較長的預(yù)測時間窗可能會取得更好的預(yù)測效果;袁伏全等(2016)利用PI 算法對青海地區(qū)中強地震進(jìn)行地震危險性的回溯性研究,通過采用不同的網(wǎng)格尺度進(jìn)行回溯性預(yù)測檢驗,認(rèn)為對于青海地區(qū)采用0.2°×0.2°網(wǎng)格尺度時預(yù)測效果更好。
印度板塊對歐亞板塊的擠壓作用使青藏塊體區(qū)域構(gòu)造運動強烈,地震頻發(fā),中國大陸7 級以上強震主要發(fā)生在青藏塊體及鄰近區(qū)域(吳哲,2018)。本文研究區(qū)(30.0°—41.0°N,88.0°—105.0°E)包括羌塘活動地塊、巴顏喀拉活動地塊、柴達(dá)木活動地塊和祁連活動地塊等4 個Ⅱ級活動地塊?;顒訅K體邊界帶是我國大陸強震的主體帶和集中區(qū),因此,研究活動地塊邊界帶的強震趨勢顯得尤為重要(張浪平等,2010)。本文在前人工作的基礎(chǔ)上,利用PI 算法對研究區(qū)MS≥6.0 地震進(jìn)行回溯性研究,對不同時間尺度下的預(yù)測結(jié)果進(jìn)行ROC 檢驗,并給出適合該地區(qū)的模型計算參數(shù)。
PI 算法的實質(zhì)是通過對研究區(qū)地震活動增強和平靜進(jìn)行分析,計算并預(yù)測在中長期時間尺度上可能發(fā)生中強地震的概率。主要步驟:首先,將研究區(qū)按一定尺度劃分空間網(wǎng)格,對落入網(wǎng)格的不小于截止震級的地震事件構(gòu)建時間序列;然后,對地震活動強度變化進(jìn)行歸一化處理,計算每個網(wǎng)格發(fā)震概率;最后,發(fā)震概率減去背景概率得到發(fā)震概率較高的區(qū)域,即“PI 熱點分布圖”(Tiampo et al,2002;Nanjo et al,2006;蔣長勝等,2008)。具體計算步驟如下。
(1)對研究區(qū)進(jìn)行空間網(wǎng)格劃分,每一網(wǎng)格為xi。
(2)對落入網(wǎng)格且不小于截止震級的地震事件的網(wǎng)格構(gòu)建1 個時間序列Ni(t),其中,Ni(t)為中心坐標(biāo)和相鄰的8 個格點(Moor 近鄰)(Moore,1962)在時間t 內(nèi)單位時間的地震數(shù)目;t0為研究區(qū)地震目錄的起始時刻;tb為滑動變化各時間序列的起始時刻;t1為地震活動中異常學(xué)習(xí)的起始時刻;t2和t3分別為預(yù)測時段的起始時刻和終止時刻(圖1)。
圖1 空間網(wǎng)格劃分及時間序列構(gòu)建Fig.1 Division of spatial grid and construction of time series
(3)地震活動強度函數(shù)Ii(tb,t)為從時刻tb到t 的單位時間內(nèi)發(fā)生在網(wǎng)格i 中的不小于截止震級的平均地震數(shù)
(4)對不同時間段的地震活動強度求平均,再除以標(biāo)準(zhǔn)偏差,進(jìn)行標(biāo)準(zhǔn)化處理
(5)計算地震活動強度函數(shù)平均變化量,以較少隨機擾動的影響
(6)未來強震發(fā)生概率Pi(t0,t1,t2)為地震強度函數(shù)平均變化量的平方。第i 網(wǎng)格的概率值減去所有網(wǎng)格概率的平均值得到強震發(fā)生在第i 網(wǎng)格的概率
將ΔPi(t0,t1,t2) >0的格點視為研究區(qū)的地震熱點。在PI 圖像中將地震活動性較強或地震頻次較高的30%格點畫出(陳建志等,2013),由其反映預(yù)測時間窗內(nèi)空間上發(fā)生目標(biāo)地震概率較高的區(qū)域,用lg(ΔP/ΔPmax)的值表示其發(fā)震概率。
ROC 檢驗方法是在國際合作項目(CSEP)中進(jìn)行科學(xué)預(yù)測試驗的方法之一,其目的是檢驗預(yù)測結(jié)果的優(yōu)劣程度。ROC 檢驗既要考慮“命中率”,同時也要考慮“虛報率”。有意義的預(yù)測必須是ROC 值大于0.5。ROC 值為ROC 曲線與隨機預(yù)測曲線所圍成的有效面積,ROC 值越大,預(yù)測效能越好。
選取青藏塊體為(30.0°—41.0°N,88.0°—105.0°E)研究對象(圖2)。歷史上該區(qū)域曾發(fā)生過多次7 級以上強震。本文使用中國地震臺網(wǎng)中心提供的1970—2019 年8 月全國ML≥2.0 地震目錄,共計20 896 次地震。其中,ML2.0—2.9 地震3 286 次;ML3.0—3.9 地震14 611 次;ML4.0—4.9 地震2 307 次;ML5.0—5.9 地震605 次;ML6.0—6.9地震74 次;ML≥7.0 地震13 次。Jiang 等(2011)基于川滇地區(qū)強震的回溯性研究討論了強余震對PI 算法的影響,認(rèn)為余震對PI 算法計算結(jié)果的影響時間不超過1 年;Tiampo 等(2002)認(rèn)為震后余震反映了區(qū)域高應(yīng)力的釋放。因此,本文沒有對地震目錄進(jìn)行剔除余震處理。
圖2 研究區(qū)空間位置矩形框為研究區(qū)域Fig.2 The spatial location of the study area
截止震級Mc是PI 算法中重要的參數(shù)之一,而截止震級Mc的選取與該地區(qū)的最小完備震級有關(guān)(Woessner et al,2005)。在計算中,若截止震級Mc選取過高,則參與計算的地震數(shù)目會減少;若截止震級Mc選取過低,則有些區(qū)域因監(jiān)測能力較弱而記錄不到地震。因此,本文將截止震級Mc設(shè)定為ML3.0,這既能保證充足的地震數(shù)據(jù)量,又能滿足截止震級應(yīng)至少小于目標(biāo)震級2 個震級單位的要求(Holliday et al,2005,2006)。
一 些 研 究 者(Rundle et al,2002;Chen et al,2005;Nanjo et al,2006;Jiang et al,2010;袁伏全等,2016)利用PI 算法對中國大陸西部地震進(jìn)行回溯性檢驗,認(rèn)為網(wǎng)格參數(shù)選取為0.2°×0.2°,預(yù)測時間窗選為8 a 時,PI 算法的預(yù)測效果較好。本文參照前人的研究結(jié)果(Rundle et al,2002;Chen et al,2005;Nanjo et al,2006;Jiang et al,2010;袁伏全等,2016),選取0.2°×0.2°空間網(wǎng)格尺度,以30 d 為步長連續(xù)向前滑動,目標(biāo)震級≥MS6.0。研究區(qū)域中強地震頻發(fā),2009 年8 月至2019 年8 月研究區(qū)共發(fā)生9 次MS≥6.0 地震(圖3,表1),為了考察研究區(qū)MS≥6.0 地震在不同模型參數(shù)下的回溯性預(yù)測效果,在研究區(qū)范圍、網(wǎng)格尺度、截止震級以及滑動步長不變的情況下,選取預(yù)測時間窗為10 a、5 a、3 a 進(jìn)行回溯性研究,模型參數(shù)設(shè)置結(jié)果見表2。
圖3 2009 年8 月至2019 年8 月研究區(qū)MS≥6.0 地震震中分布Fig.3 Distribution of MS≥6.0 earthquakes since August 2009
表1 2009 年8 月至2019 年8 月研究區(qū)MS≥6.0 地震參數(shù)Table 1 MS≥6.0 earthquakes parameters in the study area since 2009
表2 模型計算參數(shù)設(shè)置Table 2 Selection of calculation parameters
根據(jù)表2 的模型計算參數(shù),得到了3 個回溯性預(yù)測檢驗時段的PI 地震熱點分布圖(圖4)。
圖4 研究區(qū)回溯性PI 預(yù)測圖像及ROC 檢驗結(jié)果(a)、(b)2009 年8 月1 日—2019 年8 月1 日;(c)、(d)2014 年8 月1 日—2019 年8 月1 日;(e)、(f)2016 年8 月1 日—2019 年8 月1 日Fig.4 PI retrospective forecast and ROC test for the research area
選定研究區(qū)地震目錄的起始時刻 t0為1970 年1 月1 日,分別選取地震活動中異常學(xué)習(xí)的起始時刻t1為1999 年8 月1 日、2009 年8 月1 日和2013 年8 月1 日,預(yù)測時段的起始時刻t2分別為2009 年8 月1 日、2014 年8 月1 日 和2016 年8 月1 日,預(yù)測時段的終止時刻t3為2019 年8 月1 日,因此共3 個回溯性預(yù)測時間段。圖4(a)為10 a 預(yù)測時間窗(2009 年8 月1 日至2019 年8 月1 日)內(nèi)的地震熱點分布。由圖4(a)可見,異常出現(xiàn)在大柴旦—宗務(wù)隆山斷裂、鄂拉山斷裂、達(dá)布遜湖—霍布遜湖斷裂、祁連山北緣斷裂、岷江斷裂和汶川—北川斷裂附近。根據(jù)預(yù)測規(guī)則可知,這些區(qū)域是預(yù)測時間窗內(nèi)空間上相對危險的地區(qū),發(fā)生目標(biāo)地震的概率較高。實際上,該預(yù)測時間窗內(nèi)研究區(qū)共發(fā)生了9 次MS≥6.0地震,其中有4 次地震震中落在由PI 算法計算得出的地震熱點區(qū)域,分別為2009 年海西6.4級、2013 年蘆山7.0 級、2016 年門源6.4 級、2017 年九寨溝7.0 級地震。
從5 a 預(yù)測時間窗(2014 年8 月1 日至2019 年8 月1 日)內(nèi)的PI 預(yù)測圖像[圖4(c)]可見,2016 年門源6.4 級、2017 年九寨溝7.0 級地震震中均位于PI 預(yù)測圖像的地震熱點區(qū)域;從3 a 預(yù)測時間窗(2016 年8 月1 日至2019 年8 月1 日)內(nèi)的PI 預(yù)測圖像[圖4(e)]可見,2017 年九寨溝7.0 級地震震中位于PI 預(yù)測圖像的地震熱點區(qū)域。
采用ROC 方法回溯性檢驗不同預(yù)測時間窗的預(yù)測效果,發(fā)現(xiàn)10 a 時間尺度的ROC 值為0.610,5 a 時間尺度的為0.725,3 年時間尺度的為0.793??梢钥闯觯A(yù)測時間窗內(nèi)的強震震中基本都在PI 預(yù)測圖像的地震熱點區(qū)域,而且PI 算法預(yù)測結(jié)果明顯優(yōu)于隨機概率法[圖4(b),4(d),4(f)]。雖然可供參考的強震事件數(shù)較少,但在對PI 算法的回溯性檢驗中發(fā)現(xiàn),2009 年海西6.4 級、2016 年門源6.4 級、2013 年蘆山7.0 級、2017年九寨溝7.0 級地震震中都落在地震熱點的叢集區(qū),這表明PI 算法對青藏塊體的強震具有一定的預(yù)測效能,并且將預(yù)測時間窗設(shè)定為3 a 時間尺度時其預(yù)測效能較好。
2009 年8 月至2019 年8 月研究區(qū)共發(fā)生9 次MS≥6.0 地震,其中,6 次發(fā)生在活動塊體邊界帶上(圖3),分別為2009 年海西6.4 級、2010 年玉樹7.1 級、2013 年蘆山7.0級、2014 年康定6.3 級、2016 年門源6.4 級、2017 年九寨溝7.0 級地震??紤]到研究區(qū)強震大部分發(fā)生在活動塊體邊界帶,且活動塊體邊界帶也是我國大陸強震的主體帶和集中區(qū),因此定性分析了活動塊體邊界帶“目標(biāo)”地震與地震熱點叢集間的關(guān)系,發(fā)現(xiàn)青藏塊體MS≥6.0 地震與活動塊體邊界帶上地震熱點間對應(yīng)關(guān)系較好,活動塊體邊界帶的“目標(biāo)”地震都發(fā)生在地震熱點附近區(qū)域,而在活動塊體內(nèi)部的MS≥6.0 地震震前無地震熱點情況。僅就PI 算法而言,很難對該現(xiàn)象給出具體原因。但在地震中長期預(yù)測中應(yīng)用PI 算法時,將其計算結(jié)果與活動塊體邊界帶結(jié)合起來考慮強震相對危險區(qū)域,在地震預(yù)測研究中可能有一定參考價值。
本文利用中國地震臺網(wǎng)中心提供的1970 年以來全國地震目錄,在分析研究區(qū)最小完備震級的基礎(chǔ)上,應(yīng)用PI 算法,選取空間網(wǎng)格尺度0.2°×0.2°,截止震級為ML3.0,目標(biāo)震級≥ML6.0,將預(yù)測時間窗設(shè)定為10 a、5 a、3 a 不同時間尺度進(jìn)行回溯性研究,并對預(yù)測結(jié)果進(jìn)行ROC 檢驗,發(fā)現(xiàn)PI 算法對青藏塊體強震有一定預(yù)測效能,預(yù)測效果明顯優(yōu)于隨機概率法,并且將預(yù)測時間窗設(shè)定為3 a 時間尺度時預(yù)測效能較好。