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