2016年1月21日，青海门源发生M_S6.4地震，震中位于祁连山中东段的冷龙岭断裂附近，冷龙岭断裂是全新世的左旋走滑兼逆冲的断裂, 断裂的走向近北西向(袁道阳等，2004).这次地震发生在我国地震重点监视区的青藏高原东北缘地区, 该地区自2009年以来，每年开展2期流动重力观测，2013年7月甘肃岷县漳县M_S6.6地震后观测到祁连山中东段地区存在重力异常变化, 祝意青等(2014)曾指出，祁连山中东段的甘青交界地区重力场在时间演化上变化强烈，空间上差异变化剧烈，因此，存在中-长期强震的危险背景；在中国地震局第二监测中心2016年地震趋势研究报告“中国大陆及重点区域重力场动态变化综合分析”一节中进一步指出，门源、天祝、武威地区出现了重力剧烈变化及四象限分布特征，重力升、降差异达100×10^-8m·s^-2以上，祁连山中东段(震中位置在37.5°N, 102.2°E附近)有6级左右地震发生的可能.可以看出，基于流动重力观测分析做出的中期预测，与2016年1月21日青海门源发生的M_S6.4地震(37.68°N，101.62°E)对应较好，尤其是在地点的预测上，预测的危险区中心位置距离中国地震台网测定的门源M_S6.4地震震中相距只有55 km.在此之前的相关研究也表明，基于流动重力观测资料曾对2008年四川汶川M_S8.0、2013年四川芦山M_S7.0地震和2013年甘肃岷县漳县M_S6.6地震做出了一定程度的中期预测(祝意青等, 2008, 2013, 2014).这表明，依据区域重力场时-空变化的分析，可以开展强震中期预测的探索，尤其是强震可能发生地点的判定，这在地震三要素预测中尤为重要.因此，系统深入分析青藏高原东北缘重力观测资料，进而对2016年青海门源M_S6.4地震前重力场的时空动态演化特征与规律、地震孕育发生过程中的前兆表现有一个较为客观的认识，这无论是从地震预报角度来看，还是从地球动力学的角度来看，都是有实际意义的.

青藏高原东北缘是青藏高原向大陆内部扩展的前缘部位，该地区构造变形与地震活动显著, 一直受到国内外地学家高度重视.为了监测青藏高原东北缘地区地壳运动变化、发现可能的地震中短期前兆，中国地震局第二监测中心、甘肃地震局和宁夏地震局3个单位自20世纪80年代分别在该地区建立了地震重力监测网，通过不断的优化与整合，形成了目前的青藏高原东北缘整体重力监测网(图 1)，这有利于系统分析研究青藏高原东北缘重力场时空动态变化，探讨地面重力变化与地壳运动及构造活动和强震活动的关系(祝意青等, 2004, 2014).2016年门源M_S6.4地震发生在青藏高原东北缘重力监测网区域内.本文主要分析2011年以来青藏高原东北缘地区的重力场变化.

图 1

Fig. 1

图 1 青藏高原东北缘重力测量路线及活动断裂略图 Fig. 1 Map showing gravity survey routes and major active faults of the northeastern margin of Tibetan plateau

对资料的处理：①采用LGADJ程序进行拟稳平差，以获得统一起算基准下的重力变化(祝意青等，2004).②平差计算时，先对多期重力观测资料计算结果进行整体分析，初步了解各台仪器的观测精度后，合理确定各台仪器的先验方差，再重新平差计算，以得到最佳合理解算结果(祝意青等，2014).③对平差结果中极个别周围环境改变或不重合的测点进行了删除.④用Kriging对重力观测数据进行拟合推估，以便突出显示构造因素的重力效应.各期的观测资料及整体平差计算后重力点值精度统计情况见表 1.各测点重力测值的平均精度均优于12×10^-8m·s^-2，反映观测资料质量可靠.

表 1(Table 1)

表 1 青藏高原东北缘重力测量情况及精度统计 Table 1 Information of gravity survey in northeastern margin of Tibetan plateau 观测时间 仪器型号 点值最小精度 (10-8m·s-2) 点值最大精度 (10-8m·s-2) 点值平均精度 (10-8m·s-2) 2011-05—06 LCR-G793、Burris-054、CG5-509、CG5-524 4.6 12.2 8.3 2012-05—06 Burris-053、Burris-054、CG5-509、CG5-524 4.6 12.6 8.2 2013-05—06 Burris-053、Burris-054、CG5-511、CG5-207 5.3 15.5 10.0 2014-04—05 LCR-G853、LCR-G854、 LCR-G829、 LCR-G1103、Burris-054、Burris-057 5.2 11.9 8.2 2015-04—05 LCR-G853、 LCR-G829、 LCR-G854、 Burris-095、Burris-053、Burris-101 5.7 14.6 9.2 2015-09—10 LCR-G793、 LCR-G818、 LCR-G854、 Burris-095、Burris-050、Burris-054 4.9 19.7 11.7

表 1 青藏高原东北缘重力测量情况及精度统计 Table 1 Information of gravity survey in northeastern margin of Tibetan plateau

对2011年以来青藏高原东北缘地区获得的多期重力观测资料，进行年际及累积重力场动态变化分析(图 2中小圆点表示重力测点位置):

图 2

Fig. 2

图 2 青藏高原东北缘重力场动态变化图像(单位：10-8m·s-2) (a) 2011-05—2012-05；(b) 2012-05—2013-05；(c) 2013-05—2014-05；(d) 2011-05—2014-05. Fig. 2 Maps of gravity changes in northeastern margin of Tibetan plateau (contour unit: 10-8m·s-2)

(1) 2011年5月—2012年5月期间(图 2a)，重力变化总体趋势是自西南向东北出现由负向正逐渐增加，重力差异变化达80×10^-8m·s^-2以上，重力变化等值线走向总体与祁连—海原大断裂带的走向基本一致，并沿断裂构造线出现重力变化高梯度带，反映重力变化受区域应力场作用和深大断裂活动的控制.重力变化在自南向北由负向正的变化过程中于门源、天祝地区及其附近出现局部重力异常区及等值线的拐弯，异常区中心位于庄浪河断裂与毛毛山断裂、金强河断裂及冷龙岭断裂交汇区附近.

(2) 2012年5月—2013年5月期间(图 2b)，重力变化比较复杂，仔细分析具有以下3个特征：其一，甘东南地区重力出现显著差异变化，临夏重力正值变化，岷县、天水重力负值变化，重力差异变化达150×10^-8m·s^-2，并沿岷县—隆德形成重力变化高梯度带，2013年7月22日甘肃岷县M_S6.6地震发生在重力变化高梯度带上(祝意青等，2014)；其二，宁夏地区自南向北出现由负向正的趋势性变化，在海原、隆德出现与2011年5月—2012年5月(图 2a)反向变化的重力局部异常；其三，河西走廊地区在门源—天祝和祁连—民乐出现两个局部重力异常区，2013年9月20日门源M_S5.1地震发生在这两个局部重力异常区的过渡地带.

(3) 2013年5月—2014年5月期间(图 2c)，重力变化表现为一个新的态势，仔细分析具有以下3个特征：其一，甘东南地区重力出现差异变化较大，岷县重力负值变化、天水地区重力正值变化，重力差异变化达90×10^-8m·s^-2，并沿岷县—隆德形成重力变化高梯度带，但重力变化与2012年5月—2013年5月(图 2b)反向，岷县M_S6.6地震发生在重力变化反向恢复过程中；其二，宁夏地区自南向北出现由负向正的趋势性变化，但重力变化平缓；其三，河西走廊地区在武威—张掖及西宁地区出现两个局部重力正变化异常区，祁连及门源—天祝出现两个局部重力负变化异常区，门源附近形成重力变化四象限分布特征，2013年9月20日青海门源5.1级地震发生在重力变化四象限中心附近.

(4) 2011年5月—2014年5月期间(图 2d)，3年尺度的累积重力变化总体表现为自南向北由负向正的趋势性变化，重力差异变化达100×10^-8m·s^-2以上，重力变化与布格重力变化背景场一致(祝意青等，2012).重力变化具有以下3个特征：其一，2013年岷县M_S6.6地震发生在重力变化异常区及其伴生的重力变化高梯度带上；其二，河西走廊地区沿祁连山断裂带出现重力变化高梯度带，并在门源、天祝发生转折弯曲，2016年门源M_S6.4地震发生在与北西向冷龙岭断裂走向基本一致的重力变化高梯度带零值线及梯度带的拐弯部位.其三，宁夏北部重力变化相对平缓.

上述分析中提到2012年5月—2013年5月期间重力差异变化达150×10^-8m·s^-2，2011年5月—2014年5月期间重力差异变化达100×10^-8m·s^-2，后者时间长于前者，但量级更小.这是由于2012年5月—2013年5月期间出现显著变化的重力测点位于甘东南地区，主要反映2013年7月岷县M_S6.6地震前震中附近的重力变化前兆信息(祝意青等, 2014)，2013年5月—2014年5月期间岷县M_S6.6震中附近重力变化仍剧烈，但重力变化与2012年5月—2013年5月(图 2b)反向，岷县M_S6.6地震发生在重力变化反向恢复过程中.由于2011年5月—2014年5月期间的重力变化包含了岷县M_S6.6震后的反向变化，因此，甘东南地区的重力变化量级小于岷县M_S6.6地震前的短期变化, 但祁连山中东段的累积重力变化明显高于2012年5月—2013年5月期间, 较好地反映了门源M_S6.4地震前的重力累积变化特征.这说明研究区域重力场变化应动态跟踪，不同时段的重力场变化速率是不一致的.

为了研究门源M_S6.4地震前区域重力场的短期前兆变化特征，我们绘制了地震前1年尺度和半年尺度的区域重力场动态变化图像(图 3).分析图 3表明：

图 3

Fig. 3

图 3 门源MS6.4地震前重力场动态变化图像(单位：10-8m·s-2) (a) 2014-05—2015-05；(b) 2015-05—2015-09. Fig. 3 Maps of gravity changes before the Menyuan MS6.4 earthquake (contour unit:10-8m·s-2)

(1) 2014年5月—2015年5月期间(图 3a)，重力变化具有以下3个特征：其一，甘东南地区重力变化平缓，重力变化在-10×10^-8～＋30×10^-8m·s^-2之间，主要表现为2013年7月22日甘肃岷县M_S6.6地震后新的准均匀态重力变化；其二，宁夏地区自南向北出现由正向负的趋势性变化，表现为与2013年5月—2014年5月期间反向的重力变化；其三，河西走廊地区在门源—天祝一带出现新的局部重力异常，门源、天祝、武威地区出现了重力剧烈变化及四象限分布特征，重力差异变化达100×10^-8m·s^-2以上；2016年1月21日门源M_S6.4地震发生在重力差异变化剧烈的四象限中心附近，与断裂走向基本一致的重力变化高梯度带零值线上.

(2) 2015年5月—2015年9月期间(图 3b)，区域重力场变化比较剧烈，河西走廊地区在门源—天祝一带重力变化表现为与上期(图 3a)反向的变化特征，沿北西向的冷龙岭断裂及北北西向的武威—天祝断裂出现重力变化高梯度带，门源M_S6.4地震发生在与北西向的冷龙岭断裂走向基本一致的重力变化高梯度的拐弯部位.

分析1年尺度的区域重力场动态图像可以发现，2011年5月—2012年5月(图 2a)自西南向东北出现由负向正逐渐增加的急剧重力变化及沿祁连山主构造断裂出现重力变化高梯度，而且区域重力场异常变化形态与其布格重力异常的空间分布具有很大程度的相关性，这可能是区域应力增强引起的大空间尺度的趋势性变化(祝意青等，2012；Chen et al., 2016)；2012年5月—2013年5月(图 2b)祁连山中东段在门源—天祝和祁连—民乐出现两个局部重力异常区，2013年9月门源M_S5.1地震发生在这两个局部重力异常区的过渡地带，甘东南出现的剧烈重力异常变化，较好地突出了岷县M_S6.6地震前的重力变化前兆信息(祝意青等, 2014)；2013年5月—2014年5月(图 2c)重力变化在门源附近出现四象限的分布特征，它虽较好地对应了2013年9月门源M_S5.1地震，但该地区重力异常持续时间长、重力变化幅度大，重力差异变化达100×10^-8m·s^-2以上，根据以往的震例研究(顾功叙等，1997；申重阳等，2009；Zhu et al., 2010)，该地区应具有更强地震发生的可能；2014年5月—2015年5月(图 3a)门源、天祝、武威地区出现了与2013年5月—2014年5月反向的重力剧烈变化及四象限分布特征，2016年1月门源M_S6.4地震发生在重力差异变化剧烈的四象限中心附近，与断裂走向基本一致的重力变化高梯度带零值线上.年际尺度的重力场动态图像较好地反映了门源M_S6.4地震前震中附近的重力变化，是一个由“区域性重力异常→局部重力异常→四象限分布特征性异常→反向变化发震”的系统演化过程(祝意青等, 2004, 2013).另外，区域重力场的变化对测区南部2013年发生的甘肃岷县漳县M_S6.6地震也有较好的反映(祝意青等，2014).

分析区域重力场累积动态图像(图 2d)可以看出，2011年5月—2014年5月期间的累积重力场的异常变化可分为三级.一级变化为自西南向东北出现由负向正的趋势性变化，主要反映门源M_S6.4地震前区域应力场增强引起的大空间尺度重力场的有序性变化；二级变化表现为区域重力场趋势变化中的大型突变，即研究区内沿祁连—海源大断裂出现的延伸长、变幅大的重力变化梯度带.三级异常变化为沿祁连—海源大断裂出现重力变化梯度带(二级变化)外，还在断裂带附近出现门源和天祝2个局部的、不同值的重力负异常区(三级变化)，门源震中位于门源负重力异常区及与祁连—海源断裂带走向基本一致的重力变化高梯度带上、重力变化等值线拐弯的地区，较好地反映了强震中期危险地点与区域重力场的局部异常、高梯度带及其拐弯、交汇部位有关(Zhu et al., 2015).

分析地表重力变化对深部物质运动信息的反映，应分析地表变形运动对地表重力变化的影响.一般来讲，地表重力变化直接受地表垂直运动的影响，每抬升(或下降)1 cm，将引起测点1.9×10^-8～2.0×10^-8m·s^-2的下降(或上升)重力变化.目前高精度垂直运动观测主要通过水准观测来获取，但因其观测周期长，与重力观测不同步，故只能利用有关水准测量成果粗略估算地表垂直运动对重力变化的影响.1970—2011年期间垂直形变速度场图(图 4)表明，相对于稳定的华南地台，青藏高原东北缘现今总体上呈现差异性的隆升运动(郝明，2012；王双绪等，2013).其中：①西秦岭—六盘山地区是该区域上升速率较快的地带，西秦岭北缘、六盘山断裂带附近的隆升速率达5～6 mm·a^-1；②祁连山东段的天祝隆起区上升速率约为3～4 mm·a^-1；③2016年门源M_S6.4地震震中附近隆起速率为1 mm·a^-1.总的来说，青藏高原东北缘垂直形变年速率不超过6 mm·a^-1，地表垂直运动对年际重力变化的贡献低于2×10^-8m·s^-2，对1至3年不同时间尺度的重力变化来说，垂直运动对其影响在观测精度范围内.进一步对比2011—2014年重力变化与垂直形变速率(图 5)，可以看出：①垂直形变上升最剧烈的地区重力负值变化最为显著，这可用膨胀扩容和质量迁移模式来解释(Li and Fu, 1983; Chen et al., 1979; Kuo et al., 1999), 其主要观点认为区域应力场的增强变化引起了地壳介质的裂隙增大和贯通，并引起了深部地壳或上地幔热物质的上涌侵入，进而引起地表重力场的异常变化；②震中区域的东西两侧垂直差异运动较为显著，震中区域也是重力场变化最为剧烈的地方，这可用断层位错和蠕动模式(Okubo，1991；Sun，2004；Fu and Sun, 2008)来解释，其主要观点认为，断层的震前蠕动和同震错动将引起地表变形与重力场发生相应的变化.

图 4

Fig. 4

图 4 青藏高原东北缘垂直形变速率(1970—2011)(单位：mm·a-1) Fig. 4 Map of vertical deformation rates in northeastern margin of Tibetan plateau for the period 1970 to 2011 (contour unit: mm·a-1)

图 5

Fig. 5

图 5 青藏高原东北缘重力变化(等值线表示)与垂直形变速率(色标表示)对照图 Fig. 5 Comparison between gravity changes and vertical deformation rates in northeastern margin of Tibetan plateau

GNSS测量显示，在青藏地块NE向运动的环境动力作用下，青藏高原东北缘GNSS水平运动偏向北东(江在森等，2009；王双绪等，2013).重力变化也表现为：在青藏高原东北缘呈现自西南向北东由负向正的重力变化(图 2d).即沿GNSS水平运动的方向重力增加，说明重力变化受地下致密作用(密度增加)比地表隆升作用更占优势地位.由2011—2014年GNSS水平运动获得的面膨胀率分布(图 6)进一步反映出青藏高原东北缘地区呈现挤压收缩的特征，在民乐、门源及武威一带面收缩率达到峰值(-25×10^-9a^-1)，2016年门源地震震中位于面收缩率峰值附近，2013年7月岷县M_S6.6地震发生在武都及岷县北的两个面收缩率峰值之间的过渡地带.

图 6

Fig. 6

图 6 青藏高原东北缘面膨胀率分布(2011—2014)(单位：10-9a-1) Fig. 6 Distribution of surface expansion rates in northeastern margin of Tibetan plateau for the period 2011 to 2014 (unit: 10-9a-1)

为便于对比分析青藏高原东北缘重力和GNSS观测反映的面应变情况，我们把2011—2014年期间GNSS水平运动获得的面膨胀率与重力变化绘制在同一张图上，将重力变化等值线标示在面膨胀率彩色图上，得到较为直观的重力变化和面膨胀率变化图像(图 7).分析对比图 7所示的青藏高原东北缘2011—2014年期间的重力变化与面膨胀率图像可以看出，祁连山中东部重力上升变化强烈的地带面收缩率也最为显著，门源M_S6.4地震发生在重力上升变化的高梯度带与面收缩率峰值附近，门源震中附近的重力变化高梯度走向与面收缩率走向基本一致.

图 7

Fig. 7

图 7 青藏高原东北缘重力变化(等值线表示)与面膨胀率(色标表示)对照图 Fig. 7 Comparisons between gravity changes and surface expansion rates in northeastern margin of Tibetan plateau

综合上述分析认为，青藏高原东北缘近期显著的重力变化是该地区深部壳、幔物质运移作用(滕吉文等，2008；陈运泰等，2013；李德威等，2013)下引起的地表重力变化效应.已有研究表明，强震易发生在与构造活动有关联的重力变化四象限中心地带或正、负异常区过渡的高梯度带上(Zhu et al., 2015)，也易发生在沿活动性断裂的断块垂直差异运动强烈或兼有强走滑运动的地方(江在森等，2009；郝明，2012).2016年门源M_S6.4地震发生在祁连山地震带的冷龙岭断裂附近，2011—2014年重力变化和GNSS及较长时期的水准观测成果表明，门源震中附近是重力变化高梯度带拐弯的地区，也是面膨胀率变化过渡带和垂直形变高梯度带的地区；临震前震中附近区域重力场出现四象限变化特征(图 3).门源震中地区重力异常变化的几何形态与面膨胀率空间分布如此密切相关，进一步证实了青藏高原东北缘存在深部壳、幔物质运移深层次的物质与能量的交换和动力作用，即深部壳、幔边界与上地幔物质和能量尚在进行强烈交换，引起活动断层物质变迁和构造变形，在地表产生相应的重力变化.

本文系统分析研究2016年青海门源M_S6.4地震前青藏高原东北缘区域重力场动态变化特征，获得的主要认识有：

(1) 青藏高原东北缘重力网的流动重力观测资料显示2016年门源M_S6.4地震前出现了较好的中期前兆性变化图像，即区域重力场变化先呈现大尺度空间范围的有序性及与祁连—海原断裂构造带走向基本一致的重力变化高梯度带(图 2d)、后呈现震中附近特征性异常及与地震孕育发生有关的局部重力异常区(图 3a)，2008汶川地震、2013芦山地震之前也出现了这种类似现象(祝意青等, 2008, 2013).

(2) 年际尺度的重力场动态图像较好地反映了门源M_S6.4地震前震中附近的重力变化，是一个由“区域性重力异常→局部重力异常→四象限分布特征性异常→反向变化发震”的系统演化过程.另外，区域重力场的变化对测区南部2013年7月发生的甘肃岷县漳县M_S6.6地震也有较好的反映.

(3) 流动重力资料对2016年门源M_S6.4震中地点的判定，进一步证实区域重力场观测对未来强震震中位置的判定具有独到的优势.强震易发生在重力变化四象限分布中心附近或正、负异常区过渡的高梯度带上，并须结合地震构造活动情况.门源M_S6.4地震发生在重力差异变化剧烈的四象限中心附近，与断裂走向基本一致的重力变化高梯度带零值线上.

(4) 重力场的空间分布及其随时间变化与地壳垂直与水平运动及地质构造活动等观测结果有一定的对应关系，强震易发生在沿活动性断裂的断块垂直差异运动强烈或兼有强走滑运动的地方.2016年门源M_S6.4地震发生在祁连山地震带的冷龙岭断裂附近，门源震中附近是重力变化高梯度带拐弯的地区，也是面收缩率最为显著和垂直形变高梯度带的地区.

感谢参与青藏高原东北缘重力与GNSS测量的全体工作人员，他们精心观测，获得了高质量资料.感谢两位匿名审稿专家的宝贵意见.