綜合分析方法在中強地震前后電離層TEC異常分析中的應用＊

王 剛 廖洪月 唐好叢

1)西安市地震監(jiān)測中心，西安 710007

2)西安市地震局，西安 710007

引言

中強地震具有較強的破壞性，給人們生命財產帶來巨大的損失。然而，地震產生的機理相當復雜，世界各國對地震的預測預報研究一直處于探索階段。如何提升地震監(jiān)測預報的有效性一直是各國研究人員面對的難題，眾多的地震工作者不斷努力的利用各種分析方法在不同領域進行探索。

自1964年美國阿拉斯加大地震時，Leonard等[1]發(fā)現(xiàn)地震期間電離層有異常擾動現(xiàn)象后，各國學者便開始從事電離層的變化與地震之間關系的研究工作，并逐漸成為當前研究的熱點之一。Antsilevich[2]發(fā)現(xiàn)1966年塔什干地震期間電子含量有明顯增加的現(xiàn)象；Weaver等[3]發(fā)現(xiàn)1969年千島群島地震期間電離層也出現(xiàn)了異常擾動。大量研究發(fā)現(xiàn)，大于5級的地震發(fā)生前幾天到幾個小時，震中附近上空電離層都會出現(xiàn)異?，F(xiàn)象[4-6]。

近年來，我國有眾多相關領域的研究人員在關于電離層TEC變化與地震的相關性方面做了大量的研究工作。馬一方等[7]研究發(fā)現(xiàn)蘆山地震前在震中附近上空東向和北向梯度VTEC均出現(xiàn)了負異常。姚宜斌等[8]利用震中附近GNSS觀測數(shù)據(jù)，采用滑動窗口法發(fā)現(xiàn)2011年3月11日本地震前電離層有擾動現(xiàn)象。張小紅等[9]利用IGS提供的震中附近4個格網點TEC數(shù)據(jù)，采用時間序列法得到2012年1月10日蘇門答臘島7.2級地震前13天電離層TEC時間序列，并分析了時間序列法、傳統(tǒng)滑動時窗法和四分位法在預測電離層TEC參考背景值和精度方面的差異。姚璐等[10]和張明敏等[11]先后研究了玉樹地震和九寨溝地震前后電離層TEC的變化異常。

本文采用迭代法和滑動四分位距法對同一震例進行電離層TEC數(shù)據(jù)分析。說明不同方法在具體震例中的應用及各自的特點，為通過分析電離層TEC異常預測地震發(fā)生積累經驗。

1 數(shù)據(jù)震例及分析方法

1.1 數(shù)據(jù)及震例

本文采用經NASA處理后的IGS全球電離層格網數(shù)據(jù)進行電離層異常分析，其空間分辨率為5°（經度）×2.5°（緯度），時間分辨率為 2 小時。

選取震例為北京時間2019年10月12日22時55分發(fā)生在廣西壯族自治區(qū)玉林市北流市（22.18°N，110.51°E）的5.2級地震（以下簡稱北流地震），震源深度10 km，圖1為震中及格網點示意圖。

圖1 北流地震震中及格網點示意圖Fig.1 Schematic diagram of epicentral and grid points of Beiliu earthquake

選取該震例主要是為了便于在迭代法研究預測發(fā)震區(qū)域的基礎上，采用不同方法對于其預測區(qū)域發(fā)生地震震中周邊電離層TEC數(shù)據(jù)進行數(shù)據(jù)分析，在較短時間序列內發(fā)現(xiàn)電離層TEC變化異常，從而可以組成中長期預測指標和短臨預測指標相結合的前兆異常分析方法。

1.2 分析方法

1.2.1 迭代法

迭代法是一種不斷用變量的舊值遞推新值的過程。使用迭代法進行電離層TEC數(shù)據(jù)分析過程中，空間分辨率為 5°（經度）×2.5°（緯度），時間分辨率為 2 h，連續(xù)對 20 年來全球經度（?180°—180°），緯度（?87.5°—87.5°）范圍內的 5183（73×71）個格點 TEC 數(shù)據(jù)進行系統(tǒng)處理分析。作者研究中為了更好的在較大區(qū)域范圍內找到顯著的電離層TEC異常，經過多次嘗試后最終選擇7×7組成的大格網（經度跨度35°，緯度跨度17.5°）區(qū)域進行分析應用。在數(shù)據(jù)分析中，將各格點TEC值數(shù)據(jù)轉化成數(shù)據(jù)變化率，在大格網區(qū)域內，以本格網為中心，將本格網變化率數(shù)據(jù)與周邊48個網格單元對比，若本格網數(shù)據(jù)全部大于或小于周邊網格，則判斷當前格網出現(xiàn)了一次異常，然后以365日為滑動窗口計算格網的異常頻次年均線，最后將年均值與歷史最低年均值比較，當年均線值比歷史最低年均線值增大1倍時視為均線異常。

由文獻[12]可知，在TEC數(shù)據(jù)分析中，以年均值線的突變現(xiàn)象作為判斷標準，突升是TEC活躍性增強異常，反之則是TEC活躍性減弱異常，突變指數(shù)越大異常越明顯。經對比研究發(fā)現(xiàn)：日本3?11大地震、汶川地震、玉樹地震、雅安地震、九寨溝地震等地震震前震中區(qū)或周鄰地區(qū)均出現(xiàn)了明顯的TEC活躍性變化特征。以2008年5月12日發(fā)生的汶川8.0地震（31°N，103°E）為例，地震前后震中所在區(qū)域TEC活躍性發(fā)生明顯變化，從2005年開始緩慢上升，2007年上升加劇，短時間內就達到分析方法設定的異常判斷標準，即增大1倍，2007年底達到峰值，隨后呈下降趨勢。之所以出現(xiàn)這種情況，根據(jù)地震發(fā)生的機理判斷可能是因為2007年地下板塊運動達到頂峰，可運動空間變小，導致活動性減弱所致。TEC活躍性至地震發(fā)生時仍處于異常判斷標準之上，地震發(fā)生后其活躍性快速下降至異常判斷標準之下（圖2a）。分析以汶川為中心的大區(qū)域TEC活躍性，發(fā)現(xiàn)汶川周圍地區(qū)震前出現(xiàn)高值異常，且這種現(xiàn)象在地震前持續(xù)時間超過1年（圖2b）。

圖2 汶川地震 TEC 活躍性指數(shù)變化圖（a）與 TEC 活躍指數(shù)時段平面圖（2008 年 5 月） （b）Fig.2 Changes of TEC activity index in Wenchuan earthquake （a） and plan diagram of TEC activity index period（May 2008） （b）

1.2.2 滑動四分位距法

作者在進行迭代法研究電離層變化異常工作的同時，提出采用滑動四分位距法對預測區(qū)域發(fā)生地震前后電離層數(shù)據(jù)進行分析的思路，力求尋找電離層TEC變化短臨異常與中強地震的相關性。

四分位距是一種穩(wěn)健統(tǒng)計技術中用于表示數(shù)據(jù)離散度的一個量，常用來檢查數(shù)據(jù)的異常情況。四分位距法在數(shù)據(jù)分析中具有能較為準確地獲得預測參考背景值，并可以比較準確的計算出背景參考值上下限的特點。所謂的四分位數(shù)就是將數(shù)列分為4部分，一個數(shù)列可以有3個四分位數(shù)，即下四分位數(shù)、中位數(shù)和上四分位數(shù)。以一個16個數(shù)的數(shù)列為例，將其從小到大排列為x1，x2，……，x16，則

其中，Q1表示在該數(shù)值以下的數(shù)據(jù)占總數(shù)的25%；Q2表示在該數(shù)值以下的數(shù)據(jù)占總數(shù)的50%；Q3表示在該數(shù)值以下的數(shù)據(jù)占總數(shù)的75%。

在統(tǒng)計學上IQR=1.34σ，即四分位距的期望值是標準差的1.34倍。在本文數(shù)據(jù)分析中將Q2±1.5IQR作為判定TEC是否為異常值的標準控制值，即上限=Q2+1.5IQR，下限=Q2?1.5IQR。其閾值約為標準差的2倍，異常檢驗置信度約為95%。

滑動四分位距法可以看做滑動窗口法與四分位距法的結合，就是對一段時間序列內的數(shù)據(jù)采用四分位距法進行數(shù)據(jù)處理，從而得到較為合理的控制限值用來判讀變化異常。由于電離層變化受季節(jié)影響較大，選擇時間窗口不宜過長。因此，在具體數(shù)據(jù)分析中以15天作為時間窗口，即分析每天數(shù)據(jù)時取當天及前14天數(shù)據(jù)組成數(shù)列，進行數(shù)據(jù)分析，實際觀測值超出上下限值即為異常值。

2 TEC 數(shù)據(jù)分析處理

北京時間2019年10月12日22時55分在廣西玉林市北流市（22.18°N，110.51°E）發(fā)生 5.2 級地震，震源深度 10 km。

迭代法分析結果顯示，自2018年中廣西北部灣出現(xiàn)TEC異常且持續(xù)，直至2019年7月異常達到近10年來最大值，隨后異常出現(xiàn)降低，但仍處于分析方法設定值以上，直至地震發(fā)生，地震發(fā)生后TEC異常呈快速降低形態(tài)，2020年初基本回歸原活動范圍（圖3a）。分析以震中區(qū)域為中心的大區(qū)域TEC活躍性，發(fā)現(xiàn)震中周圍地區(qū)震前出現(xiàn)高值異常，且這種現(xiàn)象在地震前持續(xù)時間超過1年（圖3b）。

圖3 北流地震 TEC 活躍性指數(shù)變化圖（a）與 TEC 活躍指數(shù)時段平面圖（2019 年 7 月） （b）Fig.3 Variation chart of TEC activity index of Beiliu earthquake （a） and time interval plan of TEC activity index（July 2019） （b）

選取2019年9月1日—10月16日期間，地震發(fā)生時間前15天，發(fā)生后3天組成的時間序列，以15日為時間窗口應用滑動四分位法進行數(shù)據(jù)分析處理。下面圖中依次顯示了震中周圍各格網點的TEC異常變化情況。通過數(shù)據(jù)分析可以得出，各格網點TEC異常變化具有較好的一致性，震前10—14天、7天、3天、2天、當日均出現(xiàn)正值異常，震前1天出現(xiàn)負值異常（圖4）。

圖4 北流地震震中周圍格網點TEC時間序列（黑色柱線為發(fā)震時間）Fig.4 TEC time series of grid points around the epicenter of the Beiliu earthquake （the black column is the time of seismogenesis）

根據(jù)相關研究文獻[13]可知，電離層異常變化易受到如太陽活動變化、地磁活動異常、天氣變化因素等多種因素的影響。因此，研究地震前后電離層異常變化需要考慮時間窗口內影響因素的變化情況。在下面的異常分析中將結合F10.7指數(shù)、Ap指數(shù)、Kp指數(shù)、Dst指數(shù)等主要影響因素進行分析。相關太陽活動和地磁活動數(shù)據(jù)均來源于國家科學院空間環(huán)境預報中心（圖5）。

圖5 9月27日—10月15日北流地震前后太陽活動和地磁變化情況Fig.5 Changes of solar activity and geomagnetism before and after the Beiliu earthquake from September 27 to October 15

根據(jù)國家科學院空間環(huán)境預報中心對2019年9—10月份的太陽活動和地磁觀測數(shù)據(jù)報告可知，9月份太陽活動水平較低，無C級以上級別耀斑產生；地磁活動9月27—30日出現(xiàn)小磁暴現(xiàn)象，對TEC變化有所影響。10月份太陽活動水平依然較低，無C級以上級別耀斑產生；地磁活動10月9日Kp指數(shù)、Ap指數(shù)處于高值，10月10—11日，太陽風速上升至510 m/s左右，地磁活動水平較高，會對TEC變化有所影響。

因此，可以判斷震前11—14天，震前1—3天TEC變化異常受到地磁影響因素較大，不宜作為地震前兆異常。

通過剔除對電離層TEC異常影響因素后，獲得電離層TEC異常后還有震前10天、7天及地震當日3組異常，即10月2日、10月5日、10月12日。

3 TEC 空間分布異常分析

為了進一步分析震前3組TEC異常是否與此次地震有關，本文分別根據(jù)數(shù)據(jù)分析結果給出了3組異常高值前后（5°—40°N，75°—125°E）范圍內電離層TEC異常分布圖。

根據(jù)數(shù)據(jù)分析結果給出的10月5日異常高值出現(xiàn)在14:00UT，據(jù)此給出10月5日10:00UT—20:00UT之間的TEC空間分布圖，發(fā)現(xiàn)震中附近區(qū)域僅在18:00UT出現(xiàn)較短時間異常（圖6）。因此，可以認為該數(shù)據(jù)分析異常不宜作為地震前兆異常。

圖6 10月5日數(shù)據(jù)異常高值時間TEC空間分布圖（黃色標記為震中位置）Fig.6 Spatial distribution of TEC with abnormally high data on October 5 （epicenter location is marked in yellow）

根據(jù)數(shù)據(jù)分析10月2日異常高值時間給出10月2日2:00UT—12:00UT間隔兩小時震中周邊區(qū)域TEC空間分布圖?？梢悦黠@看出，TEC異常區(qū)域自2:00UT出現(xiàn)在震中東南方向，隨時間推移逐漸沿22°N由東向西移動，8:00UT異常區(qū)域出現(xiàn)在震中附近，直至12:00UT遠離震中區(qū)域（圖7a）。

圖7 10月2日、12日TEC空間分布圖（黃色標記為震中位置）Fig.7 The spatial distribution of TEC on October 2 and 12 （the epicenter is marked in yellow）

根據(jù)數(shù)據(jù)分析10月12日異常高值給出10月12日2:00UT—12:00UT間隔兩小時震中周邊區(qū)域TEC空間分布圖?？梢悦黠@看出，TEC異常區(qū)域變化同10月2日情況基本一致，除異常區(qū)域自東向西移動外，6:00UT異常區(qū)域出現(xiàn)在震中附近，直至8:00UT異常區(qū)域遠離震中區(qū)域（圖7b）。

綜合考慮可以認為，10月2日與10月12日的電離層TEC異常與廣西北流地震有關。

因此，對中強地震發(fā)生前后數(shù)天內電離層TEC數(shù)據(jù)采用滑動四分位距法進行分析，盡管易受到太陽活動和地磁變化的影響，但將去除干擾因素后的異常作為地震前兆短臨判斷指標，仍具有一定的有效性。這一結論與其他研究人員分析地震前電離層TEC變化異常的研究結果基本一致[14-17]。

4 結語

本文采用迭代法和滑動四分位距法對同一震例中電離層TEC異常進行了分析。迭代法具有分析較長時間內大區(qū)域電離層TEC異常的特點，滑動四分位距法可以對單個格網較短時間序列內電離層TEC異常進行分析。在地震監(jiān)測預報日常工作中可以嘗試將迭代法和滑動四分位距法分析電離層TEC異常的方法綜合使用來進行地震前兆判斷研究，以達到中長期預報和短臨預報相結合的目的。但由于地震發(fā)生的機理復雜多變，利用電離層擾動對地震進行預測存在明顯的局限性，還有待于參考多種不同分析方法進行研究和震例印證，來提高地震預測的有效性。

致謝

感謝陜西省地震局方偉高工給予的指導和支持，以及西安市地震監(jiān)測中心同事在研究過程中提供的幫助。