引言

地震是地球上频繁发生的一种自然灾害，是不可避免的，给人类社会带来了巨大的损失（杨桐等，2013）。近年来的大量研究表明，电磁异常对地震十分敏感（赵国泽等，2017）。国内外学者研究发现电离层观测到的多种数据，例如电场、磁场、等离子体密度、空间粒子通量变化等，在地震前后均出现异常变化，这些电磁异常现象与地震活动都具有相关性（安张辉等，2011；万剑华等，2012；泽仁志玛等，2012；Zhang et al，2012；闫相相等，2014）。

自2004年第一颗专门用于探测地震电离层扰动的电磁卫星（Detection of Electro-Magnetic Emissions Transmitted from Earthquake Regions，缩写为DEMETER）发射以来，科学工作者们在地震电离层方面作了大量研究。由于DEMETER是太阳同步轨道卫星，卫星飞过地面不同区域的地方时是固定的，昼侧和夜侧的轨道数据对应的地方时分别为10:30和22:30左右（杨牧萍等，2018），因此对DEMETER卫星电磁数据的背景场研究主要分昼侧和夜侧。杨牧萍等（2018）利用五年的DEMETER电场功率谱数据，分昼侧和夜侧构建了东北亚地区电离层电场背景场并分析其动态变化。Bertello等（2018）分析了八年的DEMETER卫星电场和磁场数据，构建了拉奎拉地震区昼侧和夜侧的电磁背景场，并针对2009年拉奎拉地震提取异常轨道。

蜂群（Swarm）卫星一天约飞行16个圆轨，每个圆轨的飞行时间大约为1.5个小时。其飞过同一区域的地方时是不断变化的，而不同地方时下磁场数据的背景值相差几十甚至上百倍，因此对Swarm卫星建立昼侧、夜侧背景场并不适用，还会影响提取的异常轨道的准确性。目前针对Swarm卫星磁场数据的研究主要是利用滑动时窗等方法进行单轨道分析，例如：Marchetti和Akhoondzadeh （2018）利用Swarm卫星磁场数据对墨西哥地震进行异常提取，在震前130天出现了异常累计加速；Akhoondzadeh等（2018）利用Swarm卫星磁场数据对2017年伊朗地震进行分析，发现在震前8—11天出现了地震前兆异常，然而单轨道分析提取出的异常也有可能由区域内其它短期局部性非震因素引起。此外，卫星观测的磁场数据受地磁活动的影响，当地磁活动很强时，电离层的扰动可能会掩盖因地震引起的磁场变化（李凯艳，2020）。目前排除地磁活动影响的常规作法是选择低地磁活动日的卫星轨道、舍弃高地磁活动日的卫星轨道（泽仁志玛等，2012；Zhang et al，2012；杨牧萍等，2018；Marchetti et al，2019；Ouyang et al，2020），但这样会造成卫星观测数据应用不充分。

考虑到上述不足，本文拟利用Swarm卫星磁场数据建立长期的高时间分辨率的时变背景场，去除地方时对磁场数据的影响，以及区域性的短期非震因素对结果的影响。此外，通过变分模态分解（variational mode decomposition，缩写为VMD）去除地磁活动对数据的影响，使卫星全部的观测数据可以应用于地震分析中。利用建立的时变背景场对2020年M_W7.7牙买加地震影响区域内的磁场Y分量（东向分量）数据进行异常提取，并与传统的昼、夜侧背景场的异常提取结果进行对比分析。此外，对岩石层、大气层以及电离层的多参数进行震前异常分析，并对三个圈层的异常出现时间进行解释。

1.   地震与数据

1.1   牙买加地震

2020年1月28日19:10:24，在古巴南部和牙买加西北部加勒比海地区发生了M_W7.7强震（USGS，2020a），震中位于（19.419°N，78.756°W）（图1），震源深度为14.9 km。根据Dobrovolsky公式R＝10^0.43M计算地震影响区域的半径R （Dobrovolsky et al，1979），其中，M为震级，R为孕震区半径。针对2020年牙买加M_W7.7地震计算得到R＝2046.4 km，对应于图1所示区域（0.58°S—39.42°N，58.76°W—98.76°W）。

图  1  2020年牙买加 M_W7.7地震的地理位置和地震影响区（USGS，2020a）

Figure  1.  Location of the 2020 Jamaica M_W7.7 earthquake and related earthquake affected areas （USGS，2020a）

下载: 全尺寸图片 幻灯片

1.2   电磁卫星数据

2013年11月22日，欧洲空间局发射了Swarm卫星群（Friis-Christensen et al，2006），星群包括Alpha，Bravo和Charlie三颗卫星。Swarm卫星群能实现对地球磁场高精度、高分辨率的测量，目前已经广泛应用到地震前兆异常提取中。由于中低纬地区地震磁场Y分量受外源磁场扰动小（在高纬场向电流等也会影响Y分量），更有可能受到岩石圈活动的影响（Pulinets，Ouzounov，2011），因此，本文选择Swarm卫星群Alpha卫星的磁场Y分量数据进行研究，数据采样率为1 Hz。考虑到Alpha卫星完成一个轨道大概需要90 min，地方时变化一个小时需要大概11天，大约133天可以完成地方时的重访（李凯艳，2020），本文选用2017年1月1日至2020年6月30日的磁场Y分量数据构建牙买加地震影响区域内的时变背景场，将该区域内震前90天至地震当天的轨道数据作为异常提取对象。

1.3   太阳与地磁活动指数

太阳是地磁场变化的根本来源（丁鉴海等，2011），太阳射电流量F_10.7能表征太阳活动性，一般F_10.7低于90 sfu的年份是太阳活动低年（李凯艳，2020）。

卫星观测数据受地磁活动的影响，而地磁指数是对空间电流体系引起的地磁扰动的定量描述，可以表征地磁活动水平（王楚钦等，2015；金巍等，2017）。本文选用的地磁活动指数为a_p指数，该指数又被称为“等效的行星性3小时幅度”，是描述全球3小时时段内地磁活动的指标（丁鉴海等，2011）。

2.   时变背景场建立方法

时变背景场的建立分为两个步骤：首先利用太阳活动指数F_10.7排除太阳活动性高的数据，并基于VMD将地磁活动对磁场数据的影响去除；然后利用区域网格化建立时间分辨率为1 h的时变背景场，去除地方时对磁场数据的影响。

利用VMD去除数据的部分分量，降低地磁活动对卫星磁场数据的影响。VMD是Dragomiretskiy和Zosso （2014）提出的一种自适应算法，该方法可以自适应地将地磁活动信号从卫星磁场数据中有效分离，这样既可以降低地磁活动这一非震因素的干扰，又可以利用更多的卫星磁场数据。

VMD是通过迭代搜寻约束变分模型的最优解，自适应地将信号分解成多个有限带宽的固有模态函数，每个模态函数的中心频率及带宽在迭代求解变分模型的过程中不断更新（Dragomiretskiy，Zosso，2014；池成全，2020；赵亚军等，2020），其约束变分模型表达式为

式中，{u_k}＝{u₁， u₂， ···， u_k}为分解后的模态，{ω_k}＝{ω₁， ω₂， ···， ω_k}为各模态对应的中心频率，∂_t为随时间的导数，δ（t）为狄拉克方程，$* $为卷积运算，|| ||₂为向量的二范数，f为分解出的k个模态之和。

选用2017年1月1日至2020年6月30日期间Swarm星群Alpha卫星的磁场Y分量数据，对数据作预处理，减去CHAOS-7模型（Finlay et al，2020）以去除主磁场、岩石圈磁场以及外源磁场。对预处理后的每一条轨道数据进行VMD，分解成两个模态。通过计算分解出模态的平均能量与a_p指数的相关系数，定量分析两个模态与地磁指数的关系。两个模态的平均能量φ_IMFi的计算公式为

式中，i代表第i个模态，j代表第j个数据点，N代表该轨道的数据点数。

计算每一条轨道分解出的两个模态的平均能量φ_IMFi与a_p指数的相关系数，结果如图2所示。图2a中蓝线和红线分别是每条轨道分解出的第一个模态的平均能量φ_IMF1和对应的滑动平均值$\overline \varphi $_IMF1；图2b中蓝线和红线分别是每条轨道分解出的第二个模态的平均能量φ_IMF2和对应的滑动平均值$\overline \varphi $_IMF2；图2c中蓝线和红线分别是每条轨道对应的a_p指数和对应的滑动平均值$\overline a $_p。

图  2  两个模态的平均能量与对应的$a_p $指数

（a） 轨道第一个模态的能量；（b） 轨道第二个模态的能量；（c） 轨道对应的$a_p $指数

Figure  2.  Results of average energy of the two modes and their corresponding $a_p $ index

（a） The average energy of the first mode of all the tracks；（b） The average energy of the second mode of all the tracks； （c） The $a_p $ index at the moment corresponding to the time of the track

下载: 全尺寸图片 幻灯片

由图2可知，在$\overline a $_p取得高值时，$\overline \varphi $_IMF1也取得高值，在$\overline a $_p取得低值时，$\overline \varphi $_IMF1也取得低值，二者的变化趋势很接近，而$\overline \varphi $_IMF2与$\overline a $_p的变化趋势没有明显的相似之处。为了进一步判断$\overline \varphi $_IMF1，$\overline \varphi $_IMF2以及$\overline a $_p之间的关系，利用滑动窗分别计算$\overline \varphi $_IMF1与$\overline a $_p的相关系数λ₁和$\overline \varphi $_IMF2与$\overline a $_p的相关系数λ₂，分别计算相关系数的分布，结果如图3所示，图中红点是$\overline \varphi $_IMF1与$\overline a $_p的相关系数，占比大的值集中在0.7—0.9之间；图中蓝点是$\overline \varphi $_IMF2与$\overline a $_p的相关系数，较均匀地分布在0.1—0.6之间。其中$\overline \varphi $_IMF1与$\overline a $_p的相关系数中值为0.800 3，而$\overline \varphi $_IMF2与$\overline a $_p的相关系数中值为0.305 1。

图  3  模态的平均能量与对应的$\overline a_p $的相关系数

Figure  3.  The correlation coefficients between $\overline \varphi_{\rm{IMF1}} $ and $\overline a_p $ and those between $\overline \varphi_{\rm{IMF2}} $ and $\overline a_p $

下载: 全尺寸图片 幻灯片

$\overline \varphi $_IMF1与$\overline a $_p的相关性较高，而$\overline \varphi $_IMF2与$\overline a $_p相关性不高，故第一个模态受到地磁活动的影响较大，第二个模态受到地磁活动的影响较小。因此，剔除第一个模态以降低地磁活动对磁场数据的影响，接下来只对去除第一个模态后的轨道数据进行分析。

对2020年牙买加地震影响区域进行网格化，划分成5°×5° （地理经纬度）的64个小网格，去除该区域内M6.0以上地震（USGS，2020b）分别在其Dobrovolsky范围内（Dobrovolsky et al，1979）震前90天的数据，以排除地震对背景场的影响；接着把剩余的轨道数据按照地方时（local time，缩写为LT）划分成24个地方时范围，以1 h为间隔，分别是LT0 （00:00—00:59），LT1 （01:00—01:59），···，LT22 （22:00—22:59），LT23 （23:00—23:59）。根据地磁季节划分的方法，把3—4月份划分为春季，5—8月份划分为夏季，9—10月份划分为秋季，11—12月份以及来年1—2月划分为冬季（杨牧萍等，2018）。由于2020年牙买加地震震前90天正处于冬季，因此，本文建立的是冬季时变背景场。

计算每个小网格内24个地方时的数据点的能量，计算方法为

式中，f （t）是预处理后的数据，IMF1（t）是f （t）分解出的第一个模态。

分别以网格内所有数据点能量的中值作为对应地方时的背景值来建立冬季时变背景场，如图4所示（图中白色区域代表该区域无数据）。从图中可以看出，不同地方时的背景值相差很大，其中，地方时为白天（0—5时，18时—23时）的背景值远大于地方时为夜晚（6时—17时）的背景值。

图  4  牙买加地震影响区域在不同地方时（LT）的时变背景场

Figure  4.  The time-varying background field at different local time inside the earthquake affected areas of the Jamaica earthquake

下载: 全尺寸图片 幻灯片

选取整个牙买加地震影响区域，分别计算时变背景场里整个研究区域每个地方时范围内的中值并将其作为区域背景值，进一步分析该区域时变背景场在不同地方时的背景值差异。将该区域24个地方时的区域背景值进行拟合，结果如图5所示。当地方时为10时，背景值最大，为1.059 nT²；当地方时为4时，背景值最小，为0.005 9 nT²，二者相差180倍。由此可见，地方时对Swarm卫星磁场数据的影响非常大，建立时变背景场去除地方时的影响是非常有必要的，这会影响后续异常提取的准确性。

图  5  牙买加地震影响区域内24个地方时的平均背景值

Figure  5.  The average background values of 24 local times inside the earthquake-affected area of the Jamaica earthquake

下载: 全尺寸图片 幻灯片

3.   震前卫星磁场异常提取结果与讨论

3.1   基于时变背景场与昼侧、夜侧背景场的异常提取结果

本文通过相同的方法建立昼侧、夜侧背景场，选取2017年1月1日至2020年6月30日期间牙买加地震影响区域内的磁场数据，首先将数据划分到四季，然后将冬季数据按照地方时分成昼侧和夜侧的数据集，即地方时为0—5时，18时—23时期间的数据划分到夜侧数据集；地方时为6时—17时期间的数据划分到昼侧数据集，利用区域网格化和中值建立冬季昼侧、夜侧背景场。

分别使用时变背景场与昼侧、夜侧背景场提取牙买加地震震前异常轨道。其中基于时变背景场提取异常轨道是以每个网格内对应地方时的时变背景场为背景值，利用滑动四分位（Liu et al，2000）计算网格内对应地方时的所有数据点能量的上四分位Q₁、下四分位Q₃和四分距Q₁—Q₃，以2.5倍的四分距作为阈值偏差。滑动四分位的计算方法如式（4）—（6）所示。

将一个网格内一个地方时的所有数据点的能量进行排序，结果为

当一条轨道经过某一网格的数据能量的中值与对应地方时的背景值偏差超过2.5倍四分距时，即认为这条轨道为异常轨道。基于昼侧、夜侧背景场提取异常轨道与基于时变背景场提取异常轨道的方法相同。

分别将基于时变背景场与昼侧、夜侧背景场提取出的异常轨道进行累计，当第k天提取出n条异常轨道时，异常轨道累计数目N（k）在N（k－1）的基础上加n。异常轨道数目累计结果如图6所示。由图6a可知：震前90天至震前50天，基于时变背景场提取出的异常轨道累计数目呈线性增长；震前50天至震前43天，异常轨道累计数目加速增长；震前43天至地震当天，缓慢恢复到线性增长。由图6b可见：震前90天至震前54天，基于昼侧、夜侧背景场提取出的异常轨道累计数目呈现线性增长；震前54天至震前43天，异常轨道累计数目加速增长；震前43天至地震当天，缓慢恢复至基本不变化。

图  6  基于时变背景场（a）与昼侧、夜侧背景场（b）分别提取的异常轨道累计结果

Figure  6.  The cumulative number of anomalous tracks over time based on the time-varying background field （a） and the day-side and the night-side background field （b）

下载: 全尺寸图片 幻灯片

3.2   异常提取结果对比与分析

为了说明时变背景场相较昼侧、夜侧背景场的优势，本文将从背景场和震前异常提取两方面进行比较分析。考虑到地方时为白天的昼侧背景场与地方时为夜晚的夜侧背景场差异很大，将研究区域内24个地方时的时变背景场分成昼侧和夜侧并分别与昼侧、夜侧背景场进行对比，如图7所示。可见：对于地方时为白天的时变背景场，地方时为10时的背景值最大，为1.059 0 nT²，地方时为17时的背景值最小，为0.030 5 nT²，二者相差35倍；对于地方时为夜晚的时变背景场，地方时为20时的背景值最大，为0.018 5 nT²，地方时为0时的背景值最小，为0.005 9 nT²，二者相差3倍；而该区域内昼侧背景值为0.273 7 nT²，夜侧背景值为0.010 3 nT²。由此可见，不管是白天还是夜晚，在不同的地方时，时变背景场的背景值相差均很大，有的会高达几十倍。因此，建立时变背景场可以有效地把不同地方时的背景值差异凸显出来，这对于后续准确提取异常轨道是非常重要的。

图  7  地方时为白天（a）和夜晚（b）时昼侧、夜侧的两个背景场对比

Figure  7.  Comparison between dayside and night side background field in the day time （a） and the night time （b）

下载: 全尺寸图片 幻灯片

为了分析基于时变背景场与昼侧、夜侧背景场提取出异常轨道差异的原因，对两种背景场提取异常的阈值作对比分析，如图8a所示，以一个小网格为例。图8b是两个背景场阈值差异对比，若一条轨道在某一网格的能量值落在区域A，则该轨道在时变背景场与昼侧、夜侧背景场都被判定为非异常的轨道；若落在区域B，则该轨道是昼侧、夜侧背景场提取异常少于时变背景场的轨道，集中在白天和夜晚低背景值的地方时范围内；若落在区域C，则该轨道是昼侧、夜侧背景场提取异常多于时变背景场的轨道，集中在白天和夜晚高背景值的地方时范围内。

图  8  时变背景场与昼侧、夜侧背景场提取异常的阈值（a）及其对比图（b）

Figure  8.  The threshold （a） of the time-varying background field and the day-side and the night-side background field used to extract anomalous tracks and its comparison （b）

下载: 全尺寸图片 幻灯片

昼侧、夜侧背景场比时变背景场多提取出的异常轨道都是集中在背景值高的地方时，这是由于昼侧、夜侧背景场的时间分辨率比时变背景场低，不管是白天还是夜晚，高值背景值会被低值背景值拉低，导致部分不是异常的轨道被错误地识别为异常轨道。同理，昼侧、夜侧背景场比时变背景场少提取出的异常轨道均集中在背景值低的地方时，低值背景值会被高值背景值拉高，导致部分异常轨道不能被识别。因此，建立高时间分辨率的时变背景场相较昼侧、夜侧背景场能更加准确地识别并提取异常轨道。

3.3   震前多圈层异常提取结果与分析

本文分析了震前岩石层、大气层和电离层的参量，对地震影响范围内M＞4.6浅源地震（深度＜50 km）进行震前能量累计（Kanamori，Anderson，1975；Mignan et al，2007；Marchetti et al，2020）。图9a是岩石层的震前地震能量累计结果，对岩石层应变采用了贝尼奥夫应变分析方法累计地震前能量，从图中可以看出，在震前35天出现加速现象；图9b是大气层的总水气柱异常累计结果，大气层的总水气柱（total column water vapour，缩写为TCWV）异常提取对象为震中3°×3°范围内震前90天的TCWV数据，根据地震前兆识别的气候分析算法（climatological analysis for seismic precursor identification，缩写为CAPRI）提取大气参数地震前兆异常（Piscini et al，2019），图中显示震前45天出现异常累计加速；图9c是电离层的磁场异常轨道累计结果，震前50天出现异常轨道累计加速

图  9  岩石层能量（a）、大气层总水气柱（b）和电离层磁场异常轨道（c）累计结果

Figure  9.  Accumulative number of anomalies of the lithospheric energy （a）， the total column water vapourof atmosphere （b） and the ionosphere magnetic field （c）

下载: 全尺寸图片 幻灯片

自震前78天，岩石层能量开始缓慢累计，震前50天，在与牙买加地震相同的断层上发生M_W5.0地震（2019年12月9日，19.08°N，80.44°W），这可能是主震前的微破裂，这种微破裂可能产生超低频（ultra-low frequency，缩写为ULF）电磁波（Molchanov，Hayakawa，1995）。由于电磁波在电离层的扰动，我们推测震前50天观测到电离层磁场数据的异常累计加速现象可能与电磁波的扰动有关。Ventura和Di Giovambattista （2013）发现地震前可能存在流体（或气体）的地下运动，我们推测流体的地下运动可能与微破裂同时发生，可能需要几天的时间从断层向上传播到达海洋底部，并在水中释放一些气体（Ouzounov et al，2018），这些气体可能会导致大气参数的异常，因此，震前45天观测到总水气柱异常累计加速现象可能是由气体的释放导致的。马瑾和郭彦双（2014）提出，断层带上存在相对弱和相对强的部位，前者往往表现为断层预滑、慢地震或弱震，后者是应力锁定部分和快速失稳区域（Noda et al，2013）；在断层的黏滑过程中实际存在两次失稳，前者与弱部位的释放有关，后者与强部位的快速释放有关，表现为强震。因此本文推测岩石层能量累积在震前35至震前34天的加速可能与弱部位的释放有关，在震前22天至震前20天的加速可能与强部位的快速释放有关。岩石层、大气层和电离层中不同参数的异常累积结果表明，这些异常可能与牙买加地震的孕育有关。

4.   讨论与结论

本文针对蜂群卫星地方时不断变化的特点，建立了时间分辨率高的卫星磁场数据时变背景场，有效去除了地方时对数据的影响。基于建立的时变背景场和传统的昼、夜侧背景场，分别对2020年牙买加M_W7.7地震进行了震前异常提取，两者的对比分析结果显示：利用时变背景场可以准确地提取出异常轨道，而昼、夜侧背景场由于时间分辨率低、背景值不准确，会导致异常轨道的提取结果不准确。因此，本文提出的时变背景场相较传统的昼、夜侧背景场能更准确地提取出异常轨道，提高了地震震前异常提取的准确度，为地震前兆异常的研究提供了新思路。

通过时变背景场提取的异常轨道累计结果在震前50至43天出现加速增长，联合岩石层和大气层的研究，牙买加地震先在电离层出现异常，然后在大气层出现异常，最后在岩石层出现异常，结合Molchanov和Hayakawa （1995）、Ventura和Di Giovambattista （2013）、Ouzounov等（2018）以及马瑾和郭彦双（2014）的研究成果，我们对牙买加地震震前岩石层、大气层、电离层参数的异常出现时间的机理提出了合理的猜测：电离层的异常可能是由主震前微破裂产生的ULF电磁波导致的，而大气层的异常可能是由于流体的地下运动产生的，岩石层的异常则可能是由断层的黏滑过程中存在的两次失稳导致的。

本文仅对牙买加地震进行了分析，为了获取更加令人信服的结果，需要进一步对多个震例进行统计分析，特别是M7.0以上的大地震。此外，虽然本文联合了岩石层、大气层和电离层的异常提取结果进行了简单的分析，但是如何验证提取的异常与震前活动是否相关尚需进一步的探索和研究。