地球物理学报  2012, Vol. 55 Issue (8): 2551-2560   PDF    
利用Yabuki & Matsu'ura反演方法计算 2011年日本东北地区太平洋海域Mw9.0级地震同震滑动分布
王阅兵 , 金红林 , 付广裕 , 孟国杰     
中国地震局地震预测研究所,北京 100036
摘要: Yabuki & Matsu'ura反演方法是利用ABIC最佳模型参数选取方法和平滑的滑动分布作为约束条件,由形变观测数据计算发震断层滑动分布.本文基于日本列岛同震GPS观测数据和发震断层曲面构造模型,利用Yabuki & Matsu'ura反演方法计算2011年日本东北地区太平洋海域Mw9.0级地震的发震断层同震滑动分布.反演结果表明,断层面上的最大滑动量为35 m,较大滑动分布在浅于30 km的震源中心上部,最大破裂集中在20 km深度的地方.其地震矩约为3.63×1022N·m,对应的矩震级为Mw9.0.模拟结果显示Yabuki & Matsu'ura反演方法更适用于倾角低于40°的断层模型反演.最后,本文基于上述方法获得的发震断层滑动模型,利用地球体位错理论正演计算该地震在中国及其邻区产生的远场形变,正演计算结果基本可以解释由中国GPS陆态网络观测到的同震形变.
关键词: Yabuki &      Matsu'ura反演方法      Akaike's Bayesian Information Criterion (ABIC)      断层滑动分布模型      日本Mw9.0级地震     
Estimation of co-seismic slip distribution of the 2011 Tohoku-Oki Mw9.0 earthquake using Yabuki & Matsu'ura's inverse method
WANG Yue-Bing, JIN Hong-Lin, FU Guang-Yu, MENG Guo-Jie     
Institute of Earthquake Science,China Earthquake Administration,Beijing 100036,China
Abstract: Yabuki&Matsu'ura's inverse method is to find optimum model parameters from geodetic data, combining prior information of both the smoothness of fault slip and maximum likelihood method with Akaike's Bayesian Information Criterion (ABIC). We apply Yabuki & Matsu'ura's inverse method to estimate the co-seismic slip distribution of the 2011 Tohoku-Oki Mw9.0 earthquake on a curve fault model from GPS observation data. The estimated maximum slip is 35 m. The strong slips are delineated in an area with depth less than 30 km, shallower than the seismic focus. The biggest ruptures concentrate in an area with depth about 20 km. The total moment is about 3.63?1022 N穖. The corresponding moment magnitude is Mw9.0. Our simulation results confirm that the Yabuki & Matsu'ura's inversion method is suitable for faults with dip angle less than 40?. Finally, we compute the far-field deformations in Northeast China of the 2011 Tohoku-Oki Mw9.0 earthquake using spherical dislocation theory of Sun et al.. Our theoretical predictions agree well with the observed co-seismic displacement from China GPS Observation Network..
Key words: Yabuki &      Matsu'ura's inverse method      Akaike's Bayesian Information Criterion (ABIC)      Fault slip distribution model      The 2011 Tohoku-Oki Mw9.0 earthquake     
1 引言

发生在板块俯冲区域的大地震(Mw9.0 以上)引起的震动和海啸一般都会造成巨大的损失,如1960年智利Mw9.5级地震,1964年阿拉斯加Mw9.2级地震以及2004年苏门答腊Mw9.2级地震.2011年3月11日,世界时05∶46,日本本州东部地区太平洋海域发生了Mw9.0级地震,此次地震引发了巨大海啸,造成了很大的经济损失和人员伤亡.这次地震发生在太平洋板块和北美板块间的俯冲带内,属于板块边界带的逆冲断裂型地震.

日本岛的构造背景非常复杂,处于三大板块的三联点上,其北半部分属于北美板块,东侧承受太平洋板块俯冲,而西侧则以接力方式,向欧亚板块亦即中国的东北地区俯冲,同时,其南半部分与中国大陆地区相连,构成欧亚板块的东缘部分,承受来自于太平洋与菲律宾板块的接力俯冲.在日本东北地区,太平洋板块以83mm/a速率向日本列岛俯冲[1],由此引发该地区的形变运动.至1995年日本关西地区发生了7.3级的阪神大地震后,为了监测日本岛上的地壳形变运动,日本国土地理院(The Japanese Geographical Survey Institute (GSI))建立了一个由一千多个GPS连续观测站组成的观测网GEONET[2],以20km 左右的间隔覆盖了整个日本半岛,进行实时形变监测.这次日本东北地区太平洋海域Mw9.0级地震发生后,GEONET 观测网成功地监测到地震的同震形变信号,为计算地震滑动分布以及地震发生机制等研究提供了宝贵的观测数据.

GPS测量到的同震形变数据,揭示了日本东北地区太平洋海域Mw9.0 地震对日本岛的地壳形变的影响,运用GPS形变数据可以研究和估计地震活动引起的断层滑动分布模型.地震发生后,国内外研究人员利用GPS观测到的地震同震形变数据,反演了该地震的断层滑动分布[3-9].他们的反演计算都用到了同震GPS形变数据,但因给出的断层构造和反演方法的不同,估算出大同小异的滑动分布结果.

本文将利用Yabuki & Matsu′ura[10]反演方法和曲面断层构造模型数据,计算断层面上的滑动量分布情况.该反演方法不同于目前普遍使用的Okada[11]的矩形断层模型计算,它是一种点源计算方法,以断层滑动分布为平滑的先验条件,结合观测数据,利用ABIC 和最小二乘法,找出最佳参数(parameter),计算断层滑动分布.文中首先通过模拟计算方法探讨该方法的可靠性和局限性,然后利用GPS同震形变数据进行反演计算.最后通过计算得到的滑动量分布模型,利用球体位错理论[12]计算这次日本东北地震对中国及其邻区形变场的影响.

2 Yabuki & Matsu′ura反演算法基本原理

在均匀弹性半无限空间中,利用Maruyama[13]提出的弹性动力学表示理论得到地表的形变和地下断层滑动量的关系,从而得到需要求解的观测方程:

(1)

其中,H为一个N×M维系数矩阵,d为观测数据,a为模型参数,e为误差,该误差由测量误差和模型误差两部分组成.Yabuki & Matsu′ura[10]假设误差e服从高斯分布,即:

(2)

其中,σ2 是一个未知的协方差尺度因子,E为单位协方差矩阵.当模型误差比观测误差小时,Yabuki & Matsu′ura[10]认为这个假设是成立的.

Yabuki & Matsu′ura[10]通过引入先验约束来光滑断层面上的滑动量,并假设断层面上的扭曲形变并不强烈,从而得到先验约束的向量形式:

(3)

其中,GM×M维对称矩阵,a为模型参数,r为常量.

由观测方程(1)和服从正态分布的随机误差(2)式可得到使观测数据和模型参数相关联的随机模型,从而得到观测数据d关于模型参数和尺度因子σ2 的似然函数.由先验约束(3)可得到断层面上滑动分布光滑约束的概率密度函数,通过未知参量ρ2 来调节断层面上滑动分布的光滑性.利用贝叶斯理论,将似然函数和概率密度函数合并可得到一个贝叶斯模型,从而给出观测数据的后验概率密度方程:

(4)

其中,

(5)

为了得到参数的最佳估计,Yabuki & Matsu′ura[10]利用Akaike[14]提出的基于熵最大原理的贝叶斯信息量(ABIC),定义如下:

(6)

Yabuki & Matsu'ura[10]给出ABIC的最终形式:

(7)

(4) ~(7)式中:H为格林函数组成的系数矩阵,NM分别为矩阵H的行维和列维,E为单位协方差矩阵,G为约束断层模型M×M的稀疏矩阵,P为矩阵G的秩,ΛP为矩阵G的特征值矩阵,C为与α2 有关的常量;α2=σ2/σ2σ2 为协方差尺度因子,σ2为控制断层面上滑动分布光滑的参量;a* =HTE-1H+α2G)-1HTE-1dd为观测数据组成的列向量.

利用数值的方法搜索α2 使得ABIC 的值达到最小,从而可得到参数α2 的最佳估计值珔α2 和模型参数的最佳评估.更详细的推导过程见Yabuki & Matsu′ura[10]一文.

3 Yabuki & Matsu′ura反演方法适用性分析

Yabuki & Matsu′ura[10]反演方法的基本思想是通过构造一系列的带有参数的贝叶斯模型(这些参数可以通过约束的权重来调节),利用ABIC 方法从概率论的角度确定一个最佳模型,并在线性约束的条件下找到模型参数的客观解.Yabuki & Matsu′ura[10]反演方法是在计算格林函数的时候利用了Maruyama[13]方法,该方法与Okada[11]方法区别在于前者在计算格林函数的时候不要求计算区域为矩形,是点源的计算,更具有灵活性.计算复杂的曲面断层构造反演时,只要给出等值线的数据和破裂范围就可以直接进行反演,获取平滑的滑动分布.而Okada[11]方法对处理简单的平面矩形断层更为方便,计算复杂的断层构造面破裂时,由不同倾角的子断层拼接为复杂曲面,在子断层的拼接处容易出现空白或重叠现象,对计算结果带来影响.

为了探讨Yabuki & Matsu′ura[10]反演方法的可靠性和局限性,首先在研究区域给出简单的平面断层模型和理论断层滑动分布数据,并利用Okada[11]给出的断层面上的滑动分布与地表形变的关系的正演计算方法,计算给定断层理论滑动分布在地表实际GPS台站位置所产生的理论形变量.基于以此形变量,利用Yabuki & Matsu′ura[10]反演方法计算断层模型的滑动量分布,并与给出的理论断层滑动分布进行比较,分析该反演方法可靠性和局限性.

图 1a为理论断层滑动分布模型,假定滑动方向一致,rake角为45°,即断层的上盘相对于下盘逆冲方向与断层走向的逆时针夹角为45°.断层模型的走向为201°,倾角为8°,断层长度为600km, 断层宽度为280km, 断层深度是从5km 开始计算.

图 1 断层面上的滑动量分布特征 (a)给定的断层面上的滑动量分布特征;(b)通过Yabuki & Matsu′ura[10]方法反演得到的断层滑动量分布;颜色代表滑动分布量值的大小. Fig. 1 The sliding distribution on the fault surface (a)Preset the sliding distribution on the fault surface; (b)Estimate the fault sliding distribution from random noise displacement data; Color represents the value of sliding distribution.

通过给定的滑动分布模型,利用Okada[11]半无限平面位错计算方法,计算地表形变特征,在计算的过程中,选取日本列岛上具有GPS观测站的点作为计算点,得到在这些站点上地表形变量.如图 2(a, b)表示水平和垂直形变分布.

图 2 利用Okada[11]方法计算得到的地表形变特征 (a)水平形变特征;(b)垂直形变特征. Fig. 2 The surface displacement using Okada method. (a) The horizontal displacement; (b) The vertical displacement.

利用上面得到的地表形变特征作为观测数据,考虑到观测数据具有随机的误差,将得到的观测数据加入随机误差,加入的随机误差为ee服从均值为0,标准差为0.003的高斯分布.

利用带有随机误差的数据来反演断层滑动分布.根据Yabuki & Matsu′ura[10]反演方法,在计算的过程中,我们将给出的断层模型划分为26×15个小块,每个小块用其中心点来表示,如图 2 所示,蓝色的点为断层模型通过划分后的小块的中心点.为了得到最佳模型参数,即要找到一个合适的α2 值使得ABIC 的值达到最小值,我们通过搜索法,给出一系列α2 值,计算ABIC 的值,从而可以得到α2 与ABIC 之间的关系(图 3).

图 3 α2 值与ABIC 之间的关系 Fig. 3 The relationship between α2 and AIBC

图 3 中可以看出,当α2 =1.1×10-5,ABIC达到最小值-10830,这时α2 为最佳值,通过取α2最佳值,计算断层滑动模型参数的最佳估计值,如图 1b所示.

对比图 1a图 1b的滑动分布可以看出,尽管局部区域反演获得的滑动量要比给定的值偏小,但估算的滑动分布结果与给出的理论滑动分布在整体上非常相似,能够很好的反应出给出的滑动分布特征.局部地区的微小偏差是由地表数据分布的不均匀和不充分以及给出的理论滑动分布数据不光滑所致,即图 1b中的反演计算结果,仅仅依靠于断层上盘部分的GPS观测台站点的形变观测数据;给出的理论滑动分布数据在计算过程中平滑掉.基于覆盖整个断层面的均匀分布地表形变数据进行反演计算,发现所获的估算值和理论值相似度明显提高,可见地表观测数据不均匀分布对反演结果的影响比较明显.模拟计算结果表明我们可以利用Yabuki & Matsu′ura[10]反演方法,通过地表形变数据来反演得到实际地震断层滑动量的分布特征.当然,利用GPS形变观测数据进行反演计算时,考虑到观测数据分布特征,反演结果中会存在一定的误差.

假设Xi为给出的理论断层模型的滑动量分布,Yi为通过Yabuki & Matsu′ura[10]反演方法计算出来的滑动分布量,其中i为模型划分成子块体的编号.给出下列公式来计算反演得到的滑动分式中,N布量与给出的滑动分布量的不吻合度:

式中,N表示划分的子块体的总数.可利用Misfit值的大小判断反演结果的误差.

为进一步分析Yabuki & Matsu′ura[10]反演方法在不同倾角地震断层构造反演计算中的影响,我们计算了每个倾角所对应的Misfit值.Misfit值与倾角之间的关系如图 4所示.

图 4 倾角与Misfit之间的关系 Fig. 4 The relationship of dip and Misfit

图 4中可以发现,Yabuki & Matsu′ura[10]方法对断层倾角具有敏感性,当倾角增大时,该反演方法得到的结果与真实的滑动量分布之间的Misfit是逐渐增加的.而且在计算过程中发现,当数据带有误差的时候,用该方法对高角度断层的反演结果很不稳定,与真实的滑动分布结果有较大的差异;而对于低角度的断层影响比较小,通过该方法能够得到可靠的滑动分布的特征.因此,Yabuki & Matsu′ura[10]方法适用于倾角为40°以下的低角度断层的反演.对于2011年日本东北地区太平洋海域Mw9.0地震来说,其断层的模型属于低角度的,可以用Yabuki & Matsu′ura[10]方法来反演此次地震的滑动分布量特征.

在模拟反演的过程中已知断层的模型为平面模型,而在实际的情况中,断层面往往是比较复杂的曲面.用平面模型代替曲面断层必然会产生模型误差,对于这个模型误差,在Yabuki & Matsu′ura[10]方法中,并没有单独的对其进行处理,而是将模型误差和随机误差放在一起计算,假设这两个误差总体成高斯分布,这个假设的前提是模型误差要小于随机误差,因此,在实际的计算过程中,断层模型参数的选取非常重要,断层模型参数的好坏直接影响反演结果的质量.

4 日本列岛GPS 同震形变观测数据和发震断层构造

2011年3月11日,日本大地震发生以后,日本国土地理院公开提供了日本列岛上GEONET GPS连续点的同震形变数据,美国JPL 数据处理中心则利用GIPSY 数据处理软件进行处理后,给出了由地震引起的三维同震形变位移场,结果公布在JPL 网站(ftp://sideshow.jpl.nasa.gov/pub/usrs/ARIA/ARIA_coseismic_offsets.v0.3.table).本文采用了地震发生前后05∶40-05∶55 UTC 之间5 min 方程解的处理结果,因此这个结果主要包含同震形变信息,其震后形变信息可以忽略.考虑到震中距离、观测形变量大小以及观测误差等原因,选取了282个站点作为同震形变观测数据进行了计算.选取站点位置和其同震形变位移如图 5所示.

图 5 2011年日本东北地区太平洋海域Mw9.0级地震的地表同震位移 (a)中橘黄色的箭头表示GPS观测到的水平同震位移;(b)中橘黄色的箭头表示GPS观测到的垂直同震位移. Fig. 5 Surface displacement of the 2011 off Pacific Coast of Tohoku Mw9.0 earthquake (a) The orange arrows is the horizontal displacement; (b) The orange arrows is the vertical displacement.

通过数据可以看出,日本东北部海岸向东移动的最大值大约有4.3m, 海岸线平均下沉的最大值大约0.5m.从图 5中还可以看出,同震位移的垂直变化趋势较水平变化要迅速.

基于Nakajima等[15-16]和Kita 等[17]通过对地震速度结构和俯冲型地震分布的研究给出的在日本东部地区太平洋板块俯冲深度分布数据,我们构造了曲面的地震断层模型.如图 5所示,等值线表示俯冲板块的深度,以10km 为单位的等间隔线划分.通过对太平洋板块日本东部地区的俯冲深度数据进行二次B 样条插值,对断层的划分进行一个平滑约束,得到更加真实的曲面断层模型.根据断层破裂范围和海沟深度,给定破裂面的长度范围为600km, 破裂的宽度范围为280km, 最浅深度为6km.

5 发震断层同震滑动反演结果

在计算过程中,将整个曲面断层通过二次B 样条插值划分成24×12个子块体,实际计算的模型参数个数为22×10×2共440个参数.

此次地震是发生在板块边界,被认为是由板块相对运动引起的构造应力的长期积累,达到临界后集中释放.因此,假设板块间的滑动量与板块间的相对运动在同一个方向上是合理可信的.利用这个先验信息,可以对断层面上的滑动分布方向进行约束.2011年日本东部地区太平洋海域Mw9.0级地震发生是由太平洋板块向日本岛俯冲运动引起,其俯冲速率大约8~8.5cm/a[17],根据这个信息,可以把断层面上的同震滑动方向约束在N85°E±45°范围内.

为了找到模型参数的最大相似评估,利用了图 5(a, b)中给出的同震水平位移和垂直位移观测数据,通过给定α2 值来计算ABIC 的值,当某个α2 能够使ABIC 值达到最小时,得到的模型参数为最大相似评估.通过搜索法,给出一组α2 的值来计算ABIC的值,得到α2 与ABIC之间的关系,如图 6所示.

图 6 α2 与ABIC 之间的关系 Fig. 6 The relationship of a2 and ABIC

图 6α2 与ABIC 的关系可知,当α2 =0.0066时,ABIC 有最小值-2207,计算此时的模型参数值,可得到模型参数的最佳评估,即得到断层面上的滑动分布的最佳评估,如图 7所示.

图 7 反演得到断层面上的滑动分布特征 Fig. 7 The slip distribution on the fault surface

图 7中,黑白花球是USGS公布的震中的位置和震源机制解,箭头指示断层面上滑动矢量方向,颜色代表滑动量的值,最大滑动量达到了35m.利用反演得到的断层面上的滑动分布,取μ=4×1010Pa, 可以计算地震发生时释放的地震矩为3.63×1022N·m.

本文反演得到的断层面上最大滑动量35 m 和其分布,基本与Iinuma等[6]和Miyazaki等[8]反演结果相似,两者都使用了相同的断层构造数据和类似的约束方法.但比邵志刚等[5]利用GPS数据反演得到的最大滑动量25.8m 大,小于Ito[7]等联合利用日本岛上和海底GPS 观测数据反演得到的最大滑动量65m.前者是由于反演过程中选取的约束和断层构造不同所造成的,后者是由于增加了断层近场的海底GPS形变观测数据造成.得到的地震矩比Yagi等[18]的4.5×1022N·m 要小,该差异可能的原因也是选取的断层模型的不同,反演的过程中我们利用了光滑的约束条件.计算得到的地震矩与Ozawa[3]得到的3.43×1022N·m 一致,但Ozawa[3]得到的最大滑动量为27m.Hashimoto 等[19]利用GPS陆地形变观测数据和海底形变观测数据[20]联合反演得到的最大滑动量为32m、地震矩为3.8×1022N·m;Koketsu 等[21]用强震数据、远震数据和GPS形变数据反演得到的最大滑动量为40m、地震矩为3.8×1022N·m;Yokota等[22]利用强震数据、远震数据、大地测量数据和海啸数据联合反演得到的最大滑动量为35m 和地震矩为4.2×1022N·m.各种数据联合反演结果与本文结果也有一些相似.

通过反演得到的滑动量分布来计算观测点上的GPS站点形变值,图 8给出了观测点上形变计算值和测量结果的比较,可以看出反演得到的断层面上的滑动分布在地表观测点上引起的形变计算值与观测值,垂直和水平分量都符合得比较好.

图 8 GPS观测点上同震形变位移计算值与测量值对比 (a)蓝色箭头观测得到的地表水平位移,红色箭头为计算值;(b)蓝色箭头为观测得到的地表垂直位移,红色箭头为计算值. Fig. 8 The comparison between the computed and observed c〇-seismic displacement. (a) Blue rows is the observed surface horizontal displacement, red rows is the computed displacement; (b) Blue rows is the observed surface vertical displacement, red rows is the computed displacement.
6 日本Mw9.0 级地震引起的远场同震变形

基于反演得到的滑动分布结果,利用Sun等[12]提出的球体位错形变计算方法计算了2011 年日本东北地震在中国大陆东北地区的影响,并与中国大陆陆态网GPS连续站观测到的同震的位移[23]相比较.在图 9中,红色的箭头是通过Sun等[12]方法计算得到的理论形变场,蓝色箭头是GPS 观测结果,两者在量级和方向上是比较一致的,但也存在细微差别,在JLYJ站点达到了8 mm, 引起这些误差的原因大致有以下几个方面:(1)Sun 等[12]提出分层球体位错模型只在径向(垂直)上考虑了介质的分层特性,而忽略了断层、初始应力、横向介质不均匀等因素所造成的影响,导致远场理论形变与实际观测之间的差异;(2)GPS 观测本身存在2 mm 左右的观测误差;(3)计算中利用的地表同震形变数据,只局限在较远于破裂面的日本岛上数据,因此反演结果中反映的断层浅部滑动量可能小于实际滑动量.将两者之间的差值大于2mm 的台站点挑出做出表格,如表 1所示.远场理论位移与观测位移的一致性从一个独立的角度证明了本文反演模型的可靠性.另外,可以通过大地震的滑动反演结果,研究地震的远场效应.

表 1 计算得到的形变场与观测值相减的差值 Table 1 The Difference of the computed and observed displacement
图 9 远场同震形变红色的箭头表示利用反演得到的滑动分布量来计算的远场同震形变,蓝色箭头是由GPS连续观测得到的远场同震形变. Fig. 9 The inversed model predict the Far field co-seismic deformation, red arrows is the computed surface displacements, blue arrows is the observation using GPS.

同时,利用计算获得的远场同震形变数据估算了相应的应变场影响.如图 10 所示的GPS 观测点上主压应变分布中可看到,日本东北地震在中国东北地区吉林和辽宁地区的产生的影响比较大,将改变该区域原有的应变分布状况.是否加强或削弱该地区地震危险性,这将在以后有待开展研究.

图 10 日本东北地区太平洋海域Mw9.0地震在中国东北地区产生的主压应变 Fig. 10 2011 off Pacific Coast of Tohoku Mw9.0 earthquake effects strain on the northeast China
7 总结与讨论

在运用Yabuki & Matsu′ura[10]反演方法进行计算的过程中,本文以先验条件方式给定滑动分布矢量的变化范围,并利用平滑的滑动分布作为约束条件和Akaike′s BayesianInformation Criterion(ABIC)最佳模型参数选取方法,反演2011 年日本东部地区太平洋海域Mw9.0 级地震的同震滑动分布.反演过程中,日本列岛观测到的GPS 水平位移与垂直位移资料同时参与了计算.本文依据一个曲面断层模型进行反演计算,该模型的建立依据是Nakajima[15-16]和Kita[17]等通过对地震速度结构和俯冲型地震分布研究给出的日本东部地区太平洋板块的俯冲深度数据.反演计算结果显示,断层面上的最大滑动量为35m, 较大滑动量主要分布在浅于30km 的震源中心上部,最大破裂集中在20km 深度的地方,可见震前该地方耦合率比较高,随着应力积累的不断增加,最终突破临界值而发生地震,发生较大的同震滑动.

本文的模拟计算证明Yabuki & Matsu′ura[10]反演方法比较适用于断层倾角小于40°的低角断层构造模型的反演计算.2011年日本东北地区太平洋海域Mw9.0地震是发生在低倾角断层上的巨型地震,适合运用Yabuki & Matsu′ura[10]反演方法研究该地震的发震断层滑动分布量特征.本文反演得到的断层滑动模型的地震矩为3.63×1022N·m, 比Yagi等[18]的4.5×1022N·m 要小,该差异可能与断层模型的选取有关.另外,由于本文在计算过程中运用了B样条插值方法对滑动分布量进行拟合,而由于插值的光滑和边缘效应,消减了断层面上的滑动分布量,从而导致地震矩的减小.反演得到的断层面上最大滑动量要比邵志刚等[5]利用GPS 数据反演得到的最大滑动量25.8m 要大,这是由于反演过程中选取的约束不同所造成的.我们计算得到的最大滑动量和地震矩,与利用类似的断层构造和反演方法的Iinuma等[6]、Miyazaki等[8]、Ozawa 等[3]得到的结果基本一致,同时与Hashimoto 等[19]、Koketsu 等[21]、Yokota 等[22]等利用多种数据联合反演得到的结果也是一致的.

依据日本列岛水平形变和垂直形变数据分布(图 8)可知,地震断层的滑动分布对地表水平形变响应范围较大,而对垂直形变响应范围较小,该特征定性说明断层面上的较大滑动量应分布在较浅的位置,这个结果与Yagi等[18]结论一致,同时也证明了垂直形变资料参与反演计算的重要性:单纯利用水平位移资料进行反演分析必然会导致比较大的误差.本文在考虑发震构造曲面分布的同时,利用GPS水平位移与垂直位移资料联合反演2011年日本东北地区太平洋海域Mw9.0 地震的发震断层滑动分布模型,从而具有较高的可信度.本文研究结果与其它部分国内外研究结果的比较分析提示,形变反演中观测数据分布、断层构造以及反演计算方法差异都可导致最后滑动分布大同小异的结果,具体分析导致这些区别的原因是我们将在下一步研究中展开的问题.

可利用同震形变数据反演断层滑动模型,然后利用Sun等[12]的球体位错理论研究远场形变.本文依照上述思路计算了2011 年日本东北地区太平洋海域Mw9.0地震在华北地区产生的远场同震形变,该结果基本可以解释GPS观测结果.

致谢

感谢东京大学地震研究所青木阳介博士提供了日本东北部地区太平洋板块边界电子数据.感谢日本国土地理院和美国JPL 公开提供了日本列岛GEONET GPS连续站观测数据和同震形变观测数据.

参考文献
[1] DeMets C, Gordon R G, Argus D F. Geologically current plate motions. Geophysical Journal International , 2010, 181(1): 1-80. DOI:10.1111/gji.2010.181.issue-1
[2] Sagiya T. A decade of GEONET: 1994-2003-The continuous GPS observation in Japan and its impact on earthquake studies. Earth Planets and Space , 2004, 56: xxix-xli. DOI:10.1186/BF03353077
[3] Ozawa S, Nishimura T, Suito H, et al. Coseismic and postseismic slip of the 2011 magnitude-9 Tohoku-Oki earthquake. Nature , 2011, 475(7356): 373-376. DOI:10.1038/nature10227
[4] Simons M, Minson S E, Sladen A, et al. The 2011 magnitude 9.0 Tohoku-Oki earthquake: mosaicking the megathrust from seconds to centuries. Science , 2011, 332(6036): 1421-1425. DOI:10.1126/science.1206731
[5] 邵志刚, 武艳强, 江在森, 等. 基于GPS观测分析日本9.0级地震同震位错与近场形变特征. 地球物理学报 , 2011, 54(9): 2243–2249. Shao Z G, Wu Y Q, Jiang Z S, et al. The analysis of coseismic slip and near-field deformation about Japanese 9. 0 earthquake based on the GPS observation. Chinese J. Geophys. (in Chinese) , 2011, 54(9): 2243-2249.
[6] Iinuma T, Mako O, Yusaku O, et al. Coseismic slip distribution of the 2011 off the Pacific coast of Tohoku Earthquake (M9.0) estimated based on GPS data—Was the asperity in Miyagi-oki ruptured?. Earth Planets Space , 2011, 63: 643-648. DOI:10.5047/eps.2011.06.013
[7] Ito T, Kazuhiro O, Tsuyoshi W, et al. Slip distribution of the 2011 off the Pacific coast of Tohoku Earthquake inferred from geodetic data. Earth Planets Space , 2011, 63(7): 627-630. DOI:10.5047/eps.2011.06.023
[8] Miyazaki S, McGuire J J, Segall P. Seismic and aseismic fault slip before and during the 2011 off the Pacific coast of Tohoku Earthquake. Earth Planets Space , 2011, 63(7): 637-642. DOI:10.5047/eps.2011.07.001
[9] Pollitz F. http://supersites.earthobservations.org/Politz_slip-results.jpg/
[10] Yabuki T, Matsu'ura M. Geodetic data inversion using a Bayesian information criterion for spatial distribution of fault slip. Geophysical Journal International , 1992, 109(2): 363-375. DOI:10.1111/gji.1992.109.issue-2
[11] Okada Y. Surface deformation due to shear and tensile faults in a half-space. Bulletin of the Seismological Society of America , 1985, 75(4): 1135-1154.
[12] Sun W, Okubo S, Fu G, et al. General formulations of global co-seismic deformations caused by an arbitrary dislocation in a spherically symmetric earth model-applicable to deformed earth surface and space-fixed point. Geophysical Journal International , 2009, 177(3): 817-833. DOI:10.1111/gji.2009.177.issue-3
[13] Maruyama T. Statical elastic dislocation in an infinite and semi-infinite medium. Bulletin of the Earthquakes Research Institute , 1964, 42(1): 289-368.
[14] Akaike H. Likelihood and the Bayes procedure. Trabajos de Estadisticay de Investigacion Operative , 1980, 31(1): 143-166. DOI:10.1007/BF02888350
[15] Nakajima J, Hasegawa A. Anomalous low-velocity zone and linear alignment of seismicity along it in the subducted Pacific slab beneath Kanto, Japan: Reactivation of subducted fracture zone?. Geophysical Research Letters , 2006, 33(16): L16309. DOI:10.1029/2006GL026773
[16] Nakajima J, Hirose F, Hasegawa A. Seismotectonics beneath the Tokyo metropolitan area, Japan: Effect of slab-slab contact and overlap on seismicity. Journal of Geophysical Research , 2009, 114(B8): B08309.
[17] Kita S, Okada T, Hasegawa A, et al. Anomalous deepening of a seismic belt in the upper-plane of the double seismic zone in the Pacific slab beneath the Hokkaido corner: Possible evidence for thermal shielding caused by subducted forearc crust materials. Earth and Planetary Science Letters , 2010, 290(3-4): 415-426. DOI:10.1016/j.epsl.2009.12.038
[18] Yagi Y, Fukahata Y. Introduction of uncertainty of Green's function into waveform inversion for seismic source processes. Geophys. J. Int. , 2011, 186(2): 711-720. DOI:10.1111/j.1365-246X.2011.05043.x
[19] Hashimoto C, Noda A, Matsu'ura M. The Mw9.0 northeast Japan earthquake: total rupture of a basement asperity. Geophysical Journal International , 2012, 189(1): 1-5. DOI:10.1111/gji.2012.189.issue-1
[20] Sato M, Ishikawa T, Ujihara N, et al. Displacement above the hypocenter of the 2011 Tohoku-Oki earthquake. Science , 2011, 332(6036): 1395. DOI:10.1126/science.1207401
[21] Koketsu K, Yokota Y, Nishimura N, et al. A unified source model for the 2011 Tohoku earthquake. Earth and Planetary Science Letters , 2011, 310(3-4): 480-487. DOI:10.1016/j.epsl.2011.09.009
[22] Yokota Y, Koketsu K, Fujii Y, et al. Joint inversion of strong motion, teleseismic, geodetic, and tsunami datasets for the rupture process of the 2011 Tohoku earthquake. Geophysical Research Letters , 2011, 38: L00G21. DOI:10.1029/2011GL050098
[23] 王敏, 李强, 王凡, 等. 全球定位系统测定的 2011 年日本宫城 MW9.0 级地震远场同震位移. 科学通报 , 2011, 56(20): 1593–1596. Wang M, Li Q, Wang F, et al. GPS measurements of co-seismic displacements at far field caused by the 2011 magnitude 9.0 Tohoku-Oki earthquake. Chinese Science Bulletin (in Chinese) , 2011, 56(20): 1593-1596.