0 引言

地震孕育实质是震源区介质变形、损伤进而导致失稳破裂的过程，可能与岩石圈非均匀脆性介质特定加载阶段有关。尹祥础等(1994)基于岩石应力-应变本构关系提出了加卸载响应比(Load/Unload Response Rate，LURR)理论。基于此理论，国内学者逐渐发展了一套具有明确物理意义的地震中短期预测方法(尹祥础，2004；余怀忠等，2010；余怀忠等，2026)，并在地震预测工作中取得较好的成效(宋治平等，1996；王海涛等，1998；张永仙等，1998；余怀忠等，2004；张浪平等，2006；刘月等，2023)。加卸载响应比方法通常选用一定时间窗和空间窗释放能量(即Benioff应变)作为响应量计算，近年来也尝试以地壳形变、地下流体、地磁、断层氢气等地球物理观测作为响应量进行检验(彭克银等，2000；张浪平等，2008；于晨等，2020；李泽平等，2024；贾东辉等，2024)，而以高精度超导重力作为响应量的研究目前仍是空白。

重力变化与地震活动存在密切关系(祝意青等，2012；陈石等，2015；胡敏章等，2021；郝洪涛等，2022)，高精度重力仪可观测地球内部物质特别是深部流体迁移活动，深部物质运移会导致介质内部应力集中，当达到临界状态时会触发地震，其孕震过程也可以用加卸载响应比理论解释。超导重力仪是目前性能最好的连续相对重力仪，具有高精度、高稳定性、低漂移率等特点，可用于监测各类地球物理现象特别是地震活动等(孙和平等，1997；张克亮等，2013；刘子维等，2018)。

依托“一带一路”地震监测台网项目，黑龙江省地震局在依兰—伊通断裂带北段依兰台布设超导重力仪iGrav-057，在嫩江断裂带附近碾子山台布设了超导重力仪iGrav-062，为区域孕震研究提供了高精度连续重力观测资料(图 1)。近期依兰—伊通断裂北段区域地震活动显著增强，接连发生2次M_L＞4.0地震，分别为2023年6月11日嘉荫M_L4.1和2025年1月28日汤旺M_L4.0地震。在嫩江断裂附近泰来地区，2025年11月2日和11月7日也相继发生2次M_L3.6地震。已有研究表明，加卸载响应比方法可以提取到M_L4.5地震前几个月至半年的前兆异常(张永仙等，1991)，在M0~2.0地震记录比较完整的条件下，也有可能预测M2.1以上、甚至更小地震(尹祥初等，2005)。本研究基于超导重力数据，对汤旺M_L4.0和泰来M_L3.6地震进行震前LURR异常回溯性检验，讨论超导重力数据结合LURR方法在本区域孕震分析中的适用性。

1 原理和方法

一般认为震源区介质孕震过程为非线性系统。孕震过程分三个阶段：弹性变形、损伤和亚失稳(失稳)阶段(图 2)。当系统处于弹性阶段，加载响应量和卸载响应量相当，孕震系统处于稳态。随着孕震区应力不断累积，岩石介质产生非弹性变形，内部微裂隙不断积累，进入损伤阶段。当系统处于损伤阶段，加载响应量大于卸载响应量，系统将失稳。

地壳应力包括构造应力和潮汐应力，潮汐应力量级小于构造应力，但是潮汐应力变化量远远大于构造应力变化量。地震成核失稳意味着地壳处于高应力阶段，此时若给予一个相对较小的潮汐应力扰动，会促使地壳介质发生较大改变(尹祥础，1987；刘晓霞，2010)。当断层接近临界状态时，可将潮汐应力引起的库仑破裂应力的增减作为加、卸载的判断条件(任建辉等，2015)。

LURR定义为

$ Y=\frac{X_{+}}{X_{-}} $

(1)

$ X=\lim\limits_{\Delta P \rightarrow 0}\left(\frac{\Delta R}{\Delta P}\right) $

(2)

其中，“+”表示加载，“-”表示卸载，P和R分别代表载荷和响应，ΔR为载荷P发生微小ΔP时响应R的变化，X为响应率。

以高精度超导重力、地壳形变等观测数据为响应量，计算公式可表示为

$ Y=\frac{\frac{\left(\sum\limits_{i=1}^{N_{+}} G_{\mathrm{i}}\right)}{\left(N_{+}\right)}}{\frac{\left(\sum\limits_{i=1}^{N_{-}} G_{\mathrm{i}}\right)}{\left(N_{-}\right)}} $

(3)

其中，G_i为第i个观测值，N₊和N_-分别为加载阶段和卸载阶段观测值总数。即一定时间窗内，响应量Y为加载阶段记录数据的平均值与卸载阶段记录数据的平均值之比。

2 超导重力数据及预处理

本文收集整理了2024年2月1日—2025年1月31日汤旺M_L4.0地震前的依兰超导重力数据，以及2025年2月1日—2025年9月31日泰来M_L3.6地震前的碾子山超导重力数据。以依兰台超导重力数据为例，进行数据预处理说明。

首先去除记录数据中存在尖峰、阶跃、间断等异常数据，数据预处理软件采用Tsoft软件(Vauterin，1998)和韦进等(2012)开发的DataPreProcess4 Singal软件。对于残差信号，预处理完成后将扣除的理论合成固体潮加回，得到“修复”的固体潮数据。随后进行各项改正，包括固体潮改正、海潮改正、气压改正、极移改正、漂移改正等后，得到包含可能与地震异常变化相关的重力残差信号。固体潮改正和气压改正需对“修复”的固体潮数据进行潮汐分析，采用Eterna调和分析软件，得到实际潮汐参数表和实测导纳值。

3 数据处理及计算分析 3.1 数据处理

在地震预测实践中，应依据潮汐力在地震破裂面上引起的库伦破裂应力变化判断加载和卸载过程，但地震破裂面在震前未知，通常可采用区域构造应力场作为震源参数。具体方法为参考历史上区域附近较大地震的震源机制解资料，以断层面走向、倾角和滑动方向3个参数间接获得介质构造应力状态，计算潮汐应力在破裂面上的有效剪应力(彭克银等，2000；于晨等，2020)。

参照Yin等(2008)和Yu等(2006、2020)构建的参数分区模型，设置本区域构造应力场震源参数为：走向30°，倾角85°，滑动角174°，震源深度10km，摩擦系数0.4，时间窗长30天，滑动步长7天，计算预处理超导重力值LURR时序变化。预处理重力值包含固体潮等潮汐成份，且存在正负值(图 3(a))。

计算结果(图 3(b))显示，在超导重力值同时存在正、负值情况下，会导致LURR为负值且变幅较大，因此需对相对重力值增加偏移量，使相对重力值平移后全部为正。偏移量取超导重力最小值的绝对值。对经平移后的重力数据(图 3(c))计算LURR时序曲线(图 3(d))，可以看出LURR在0.6左右变化，其原因可能为重力固体潮日变幅较大，约为200×10^-8~300×10^-8m/s²，导致LURR比值偏离1，微弱异常信号淹没在背景值中。由于固体潮理论值可根据公式预先计算得到，为非异常值，若要识别并突显LURR异常变化，需去除理论固体潮。

3.2 扣除理论固体潮重力残差LURR分析

设置区域应力场为震源参数，以去除理论固体潮的重力残差(图 4(a))作为响应量计算得到LURR时序曲线(图 4(b)中圆圈折线图)。结果显示，汤旺地震前1~3个月LURR值呈明显的峰值点变化特征。岩石力学实验也证明了LURR异常存在峰值点规律(尹祥础，2015)，即孕震初期LURR大多在1上下波动，随后逐渐上升至峰值，但地震往往不在峰值点发生，而是有所滞后。峰值点到发震时刻的这段时间，称为延迟时间，延迟时间与介质的均匀程度有很大关系(张磊等，2025)。

通过P波初动方法得到汤旺地震的震源机制解(图 1)。研究表明震源机制解能够比较准确地给出局部区域的构造剪切应力方向，而构造剪应力方向与地震断层面上的滑移矢量一致。地震的加卸载状态可通过潮汐有效剪应力在构造剪切应力方向上的投影，即潮汐应力在地震断层面上引起的库仑应力变化来判断(余怀忠等，2026)。以震源机制解(走向143°，倾角70°，滑动角20°，震源深度10km，摩擦系数0.4)为震源参数计算得到LURR曲线(图 4(b)中五角星折线图)，并与区域应力场参数计算结果(图 4(b)圆圈折线图)对比，结果显示以显著地震震源机制解为震源参数得到的异常更显著。

结合震源机制和地震地质背景，判断发震断层为乌伊岭断裂，其为区域主断裂依兰—伊通断裂西侧北西向次级断裂(图 1)。研究表明，许多地震并不在主断裂上发生，而往往发生在主断裂一侧的次级断裂上，此时地震将受这些小断裂与主断裂的共同影响(刘特培等，2003)。已有研究表明，依兰—伊通断裂北东向主断裂影响区域构造运动，北西向次级断裂控制地震活动(王智雄，2022)。超导重力计算结果显示，选用区域构造应力场(主断裂)和震源机制解(次级断裂)为震源参数均能观测到LURR异常变化，说明区域地震受主断裂构造运动与次级断裂地震活动共同影响，也与以往本区域形变场-块体相对运动研究结果一致。

3.3 精细重力残差LURR分析

通过调和分析计算各分潮波群潮汐参数，以此潮汐参数构建当地合成潮汐模型，剔除观测记录中的固体潮和海潮负荷影响。利用依兰超导重力仪同步气压观测数据，采用Eterna调和分析程序，计算大气导纳值为-0.2739 ×10^-8m/s²/hPa。利用该实测导纳值扣除气压影响，得到精细重力残差变化。残差时间序列中仍存在振幅2×10^-8~3×10^-8m/s²小台阶及诸多地震突跳变化，进一步扣除小台阶、突跳变化和扣除极移影响得到变化曲线(图 4(c)。从结果中可明显看出，依兰超导重力累积1年增大14×10^-8m/s²，随后1月28日发生汤旺M_L4.0地震。为观察震后重力是否继续累积，延长数据时段至震后1个半月(2025年3月17日)，处理得到精细处理结果(图 4(c))，可以看出震后重力时间序列基本稳定变化，1个半月时间仅增大约1×10^-8m/s²。

利用精细重力残差变化计算LURR时序，结果显示LURR在1.0附近变化，范围大多在0.01以内(图 4(d))，LURR异常消失，分析认为数据处理扣除的实际固体潮、海潮、气压改正、极移可能是介质达到临界状态后触发地震的影响因素。

通过依兰超导重力数据实例计算，认为以超导重力为响应量提取LURR异常，应采用预处理后的数据，首先将观测值进行平移为正数，然后仅扣除理论固体潮。基于以上处理策略，对距依兰超导台226km的伊春汤旺M_L4.0地震和距碾子山超导台140km的泰来M_L3.6地震计算LURR，并分别进行分析和讨论。

4 结果和讨论 4.1 汤旺M_L4.0地震分析 4.1.1 超导重力计算结果

依兰—伊通断裂北段分为北东向主断裂和北西向次级断裂群，主断裂和次级断裂为同一构造系统，为此次地震孕震区，依兰台位于主断裂带附近(图 1)。超导重力LURR异常反映孕震区局部出现高应力集中，进而导致孕震区次级断裂群岩石产生裂隙并进一步扩展，流体物质沿断裂带侵入，引起超导重力加卸载响应比值产生异常变化。参考本区域地震地质研究成果(余中元等，2017)，认为汤旺M_L4.0地震的发震机制为：在西太平洋板块俯冲作用下，依兰—伊通断裂北段两侧地块产生相对运动差异，与主断裂垂直分布的北西次级断裂断块间受剪切应力作用向断陷带中心滑动。在此过程中高应力区岩石产生裂隙，深部流体物质涌入，造成应力聚集，推动重力值持续增大14×10^-8m/s²，此过程被超导重力仪观测到(图 4(c))。在应力积累到临界状态时，微小的潮汐应力或气压影响会产生较大的系统响应，导致系统失稳，经历破裂成核阶段，最终引发地震。

地震发生后，更新超导重力数据至2025年9月30日，仅扣除理论固体潮进行处理(图 5(a))。可看出重力观测曲线存在明显下降趋势，通过对比同址观测的依兰量子绝对重力变化(图 5(b))，认为下降趋势为仪器漂移导致，计算依兰超导重力仪器漂移率为-8×10^-8m/s²/a。去除仪器漂移后，得到依兰台超导重力残差变化(图 5(c))，可以看出与同时间段量子绝对重力观测曲线(同样仅扣除理论固体潮)趋势较一致。

对比去仪器漂移与未去漂移的超导重力LURR时变曲线(图 5(d))，两者LURR曲线震后均已恢复至正常背景值附近，背景值约为0.998，震后平静时段接近一年，孕震系统暂处于稳定阶段。两者形态接近说明LURR异常变化基本不受仪器线性漂移的影响，这也是2024年超导重力精细处理后重力值呈类似线性增大形态(图 4(c))，而图 4(d)中LURR结果仍变化稳定，未探测到异常的原因。

4.1.2 基于地震目录的LURR时空演化特征

应用超导重力资料作为响应量探测震前加卸载响应比异常尚属首次，而以Benioff应变作为LURR响应量则为常用有效方法，已广泛应用于地震预测实践，本研究尝试结合Benioff应变作为响应量对汤旺震例进行对比验证分析。

从图 1可看出依兰—伊通断裂带北段及西侧区域地震丛集分布，且处于同一孕震环境。截取1973—2025年历史地震目录(图 1中蓝色圆圈)，计算整个区域加卸载响应比时序。为避免地震数据太少导致加卸载响应比时间序列强烈波动，计算时间窗内包含多个加、卸载循环过程(刘月，2014)。设置计算参数为：时间窗360天，滑动步长10天，最小地震数目为5，震级上限M_L4.0，下限为M_L0。设置震源参数为：走向143°，倾角70°，滑动角20°，震源深度10km，摩擦系数0.4，计算结果如图 6所示。

LURR时序扫描结果(图 6)显示，本区域在出现LURR高值后3~12个月时间内，往往对应发生M_L≥4.0地震，且异常值越大，持续时间越长，对应地震个数也越多，说明本区域LURR异常和地震孕育存在相关性。对于近期2次M_L≥4.0地震(2023年6月11日嘉荫M_L4.1和2025年1月28日汤旺M_L4.0地震)，震前均出现LURR异常变化，其变化特征符合“正常—高值—正常—发震”规律。但仅从时序曲线无法获得时空孕震演化过程，因此对近期2次地震进行空间扫描和孕震分析。

选取黑龙江及周边地区2009—2025年地震目录(图 1中浅色透明圆圈)，空间扫描参数设置经纬度范围：40°N~55°N，118°E~136°E，以0.25°步长为单位分别沿经度和纬度方向滑动，扫描半径200km，扫描区域内M_L0~4.0地震的Benioff应变为响应量，摩擦系数取0.4，计算时间窗18个月，最小地震个数20，以时间窗截止日期标记为扫描结果日期。

LURR空间演化特征显示，嘉荫M_L4.1地震前1年9个月(2020年4月1日—2021年9月3日)，在震中附近出现LURR高值异常，随后该区域异常幅值产生波动，先减弱后逐渐增强，至震前9个月(2021年4月1日—2022年9月30日)异常幅值达到最大，震前3个月消失，随后2023年6月11日发生嘉荫M_L4.1地震，震后异常区迁移至南部依兰超导重力台附近，该区域LURR值逐渐增大，至震前6个月达到最大，随后LURR值逐渐减弱，至2025年1月31日发生汤旺M_L4.0地震。目前整个依兰—伊通断裂北段LURR处于低值平静阶段(图 7)。

孕震区域陆续出现异常“集中—增强—减弱”过程，表明该地区在北东向挤压背景应力场控制下，应力应变不断累积增强导致震中及附近区域岩石损伤或弱化，位于嘉荫M_L4.1地震震中区应变能积累超过介质强度而引发地震。进而在依兰超导重力台附近区域再次出现应力集中区，时间段对应超导重力数据2024年重力值的持续增大，表明区域应力逐渐集中，当达到临界状态时触发汤旺M_L4.0地震。超导重力和Benioff应变加卸载响应比方法均监测到汤旺地震震前异常，且异常区一致。

4.2 泰来M_L3.6地震分析 4.2.1 超导重力计算结果

2025年10月2日和10月7日，黑龙江省齐齐哈尔市泰来县附近连续发生2次M_L3.6地震(分别位于46.38°N、123.70°E和46.40°N、123.69°E)。本区域历史地震并不十分活跃，附近主要断裂有嫩江断裂、滨州断裂、第二松花江断裂。

碾子山超导台距泰来M_L3.6地震震中140km。采用上述方法处理最新的超导重力资料(数据更新至2025年9月30日)，得到重力残差变化(图 8(a))。结果显示，2024年8月至2025年9月重力值增大了15×10^-8m/s²，可能为正向漂移，计算年漂移率为12×10^-8m/s²。去除漂移后，年重力值变化范围在10×10^-8m/s²以内(图 8(c))。

分别计算未去漂移和去除仪器漂移的重力残差LURR值平均值和中误差(图 8(b)、8(d))，设置2倍中误差为阈值，超过阈值范围为异常值。从图 8(b)可以看出，泰来M_L3.6地震前8个月，LURR值为1.0003，超出正常区间[0.9979，0.9996]，异常幅度为125%。说明震前碾子山台及附近区域存在应力扰动和临界地壳介质损伤，为本次地震的异常。去漂移重力残差LURR异常特征也与之相同。

4.2.2 GNSS面应变结果

为了与汤旺M_L4.0地震震前异常演化过程对比研究，尝试对泰来M_L3.6地震进行LURR空间扫描。但因泰来地震震中及周边区域地震活动性弱，LURR空间扫描难以满足得到合理结果的条件。而GNSS应变场时空变化反映的同样是应力应变积累和集中的时空演化过程，本文采用GNSS应变场云图结果进行替代。

收集震中周边省测绘局CORS站和“陆态网络”GNSS站资料，采用Gamit软件解算坐标时间序列，经预处理后设置18个月窗长、1个月步长，计算应变场变化云图(图 9)，具体计算方法见文献(常金龙等，2021)。从图 9可以看出，自2025年3月开始碾子山台东部、泰来地震震中东北部嫩江断裂、滨州断裂及第二松花江断裂逐渐形成挤压高值区，且逐月增强，同时东南部形成面膨胀高值区，随后在高值区边缘发生地震。孕震过程对应碾子山超导重力LURR异常“集中—增强—减弱”过程。综合分析结果，可解释本次地震的孕震过程为在地表(GNSS)和地壳应力(重力)共同作用下，震中区域介质受挤压应力集中失稳破裂而触发地震。

5 结论

超导重力资料精度高，连续性好且稳定，年漂移量小。与基于Benioff应变为响应量的LURR方法相比，以超导重力为响应量可以更好地认知少震地区地震孕育的物理过程。预处理需对相对重力增加偏移量，调整为正值，计算时仅需扣除理论固体潮，如果存在明确的线性漂移，应尽量扣除仪器漂移。

通过震例分析，依兰—伊通断裂北段和嫩江断裂区域震前可观测到较明显的LURR异常。综合分析超导重力、Benioff应变等不同响应量的LURR异常变化，辅以GNSS形变场变化特征，可揭示地震前震源区介质孕震时空演化过程。

汤旺M_L4.0地震前，依兰超导重力呈持续增大趋势，反映了依兰—伊通断裂及西侧孕震区岩石介质在北东向挤压应力作用下产生裂隙并扩展，流体在裂隙中运移引起较大范围质量迁移变化，流体迁移区位于震中与依兰超导台之间。随着深部流体物质沿区域断裂系统渗透，震中及附近同一构造区域应力累积增强，破裂带岩石出现损伤至应变能超过介质强度极限，触发地震。区域Benioff应变LURR时空扫描震前呈先增强后减弱特征，佐证了本次地震震前异常及孕震演化过程。

泰来M_L3.6地震前孕震区地表形变场呈现持续挤压变形，表明地壳存在应力应变积累，形变场异常时间段内超导重力LURR存在近同步异常变化形态。推断泰来地震孕震过程为地表和地下介质受挤压，挤压高值区边缘震中附近岩石介质达到临界状态后失稳破裂而发震。在地震较不活跃地区，可结合GNSS应变场和LURR时间扫描方法，共同研究地震孕育的演化过程，为震前异常识别提供特征参考。

基于超导重力的加卸载响应比方法在黑龙江省地震重点监视区震例回溯检验中，提取到中小震震前异常，刻画了时空孕震演化过程，应用效果显著。在众多的观测数据中，超导重力变化是一种对地震反映非常敏感的地球物理前兆异常，但受台站数、台站位置、区域主次断裂分布、震中矩等因素影响，各震例异常出现时间、持续过程和结束时间存在一定差异。目前两台超导重力仪器运行时间段尚短，还需要积累更多震例进行分析总结。

致谢: 中国地震台网中心余怀忠研究员团队为本研究提供了计算程序，审稿专家对论文提出了宝贵意见和建议，图件采用GMT绘制，在此一并表示衷心感谢。