利用概率密度分布提取流体观测资料中的高频异常信息——以2008年汶川8.0级地震为例

孙小龙， 王广才， 晏锐

1 中国地质大学(北京)， 北京　100083 2 中国地震局地壳应力研究所， 北京　100085 3 中国地震台网中心，北京　100045

利用概率密度分布提取流体观测资料中的高频异常信息

——以2008年汶川8.0级地震为例

孙小龙1,2， 王广才1， 晏锐3*

1 中国地质大学(北京)， 北京100083 2 中国地震局地壳应力研究所， 北京100085 3 中国地震台网中心，北京100045

摘要随着地下流体观测技术的提高，尤其是地下流体观测数字化改造以后，观测资料的采样频率明显提高，这些高频采样的观测资料中蕴含着丰富的构造信息，如何从分钟值甚至更高频率采样的观测资料中提取有用的异常信息，是目前从事地震地下流体资料分析人员最为关注的科学问题之一.本文引入概率密度分布法，分析了2008年汶川8.0级地震前南北地震带及其附近区域72个测点的数字化水位和水温分钟值采样高频观测资料，结果显示：汶川8.0级地震前有16个测点水位和14个测点水温出现高频异常信息，出现高频信息异常的观测点多集中在滇西南构造带，异常出现的时间呈现出由南向北推移的特征.据此认为概率密度分布法在流体资料的高频异常信息提取方面具有一定的可靠性与适用性，可为数字化高频观测资料异常提取提供借鉴.

关键词流体资料； 概率密度分布； 高频异常信息； 汶川8.0级地震

1引言

21世纪初，针对观测井孔和观测仪器老化、观测技术水平跟不上技术现代化发展要求的现状，中国地震局在“九五”、“十五”期间实施了前兆台站观测技术改造工程，对全国前兆观测井进行了数字化、网络化改造.经历了近10年的实践与资料积累，数字化观测的优势逐渐被行业内同行所认可，与之前的模拟观测相比，数字化观测具有数据传输速度快、数据信息量丰富、人为误差少等优点，从而使地下水位的高频、短周期信息大量增加，为捕捉地震孕育及发生过程中的前兆异常信息提供了有利条件.流体观测技术经历了数字化改造后，观测设备采样率明显提高，除了之前常采用的月均值、旬均值和日均值外，还出现了整点值、分钟值、秒值以及更高频的数据.

但是，在日常的数据处理与震情跟踪工作中，业内普遍沿用先前针对模拟观测所提出的各类分析方法(张素欣等，2002；杨明波等，2009)，数据类型的选取多采用月均值、日均值和整点值，数字化观测资料中的分钟值、秒值以及更高频的数据信息，仅限于应用在流体对远场地震的同震响应分析(刘耀炜等，2005；孙小龙和刘耀炜，2008；曹玲玲和高安泰，2010； Sun and Liu,2012)，更普遍的现象是在某些时候，工作人员为了异常分析需要，还须将数字化观测资料中的分钟值、整点值等转化成日均值、月均值进行异常提取，变相增加了分析人员的工作量.出现以上现象的原因，并非是数字化观测技术或观测资料本身存在某种缺陷，而是由于在数字化改造实施后，尚缺乏与之相适应的分析处理方法.因此，在数字化观测技术改造完成后，摆在地震监测与预报人员面前的首要问题是尽快提出与数字化观测资料相匹配的数据处理方法，即数字化观测资料中高频信息的异常提取方法研究.

近几年来，科研人员针对数字化观测资料的特征与其所面临的科学问题，探索性地提出了一些适合数字化资料的数据处理与异常分析方法.例如，杨丛杰等(2005)运用小波分析方法对1990 年常熟、1995 年苍山和1996 年南黄海三个中强地震前江苏地区井水位固体潮的变化特征进行了研究，发现井水位M2波潮汐因子在地震前几个月几乎都出现一个幅度较大、周期为半月或一个月左右的异常信号，表明小波分析方法在处理和分析井水位潮汐资料方面可能是一种有效的方法.邱泽华等(2010，2012)应用小波分解-超限率法处理并分析了汶川8.0级地震前姑咱台和宁陕台钻孔应变高频数字化资料在震前出现短周期“毛刺”的异常现象，分析认为这种高频的异常是震前构造运动的表现.刘琦和张晶(2011)引入时频分析中的S变换法，对姑咱台四分量钻孔应变观测曲线中出现的“压脉冲”和“潮汐畸变”异常信号进行了分析，结果表明在背景信号的基础上，周期约为10～60 min的信号在汶川地震前开始大量出现，震后逐渐衰减.晏锐等(2011)基于临界慢化概念，运用表征临界慢化现象的自相关系数和方差法，分析了汶川8.0级地震前氡浓度观测资料，结果表明不同台站的水氡浓度震前都存在明显的临界慢化现象.Manshour等(2009，2010)曾用概率密度分布法成功提取了全球12次强震发生前高频垂直速度观测数据中的异常现象，结果表明该方法在提取高频异常信息时效果明显.

由此可见，数字化观测资料的高频成份中隐含着大量的前兆异常信息，运用适当的数据处理方法，是可以识别出这些前兆异常的.

2概率密度分布法识别异常信息

随着地震观测数据采样率的不断提高，其数据量也急剧增大.如何从大量的观测数据中提取出分析所需的有用信息，是当前数字化观测资料分析研究的重点.目前，在提取数字化观测资料中的高频信息时，普遍采用高通滤波法(邱泽华等，2010，2012)、小波分析法(Daubechies,1988；晏锐等，2011；邱泽华等，2012)、多项式拟合法(Manshour et al.,2009，2010)、经验模态分解法(Huang et al., 1998；孙小龙等，2011)等，这些方法既能有效地提取出观测数据的低频趋势变化项，也能提取出高频的短周期信息.

本文在地震流体数字化观测资料的高频信息成份提取过程中，采用了三阶多项式法，该方法在日常的数据处理中多用来提取流体观测资料中的长趋势变化，而本文主要利用原始数据在消除趋势变化后的剩余信息，即短周期的高频信息成份.三阶多项式法主要是利用多项式拟合法对原始数据进行趋势拟合，得到一条最接近原始曲线趋势变化的三阶多项式变化曲线，利用原始曲线数据减去趋势变化曲线数据，得到原始数据中的高频信息成份.

对一组时序数据D(t)，以s为窗长进行三阶多项式拟合，得到趋势项拟合曲线为

T(t)=at3+bt2+ct+d，

(1)

其中，a、b、c和d均为拟合系数，t为时间.高频信息成份Z(t)=D(t)-T(t)，给定窗长s后依次滑动进行拟合求参数，并最后得到原始数据中的高频信息Z(t).

图1为四川南溪台水位分钟值曲线及其利用三阶多项式(窗长600 min)提取的高频信息，其中图1a和图1b分别为一个窗长的原始观测曲线和剔除趋势后的高频信息曲线(a中虚线为依据原始观测值利用三阶多项式拟合得到的趋势变化线)，图1c和图1d分别为2007年9月1日至2008年5月11日的原始观测曲线和剔除趋势后的高频信息曲线.从图可以看出：

(1) 在一定窗长范围内，三阶多项式方法基本上能拟合得到观测值的趋势变化线(图1a)，进而可去除趋势得到观测曲线的高频信息.如果观测曲线中有明显的潮汐或气压变化形态，如部分水位观测曲线，该方法一定程度上能消除掉观测曲线中固体潮或气压的半日潮变化(图1b)，若想达到较好的拟合效果，窗长应小于半天.

(2) 从年尺度和月尺度的曲线对比来看，三阶多项式方法能很好地剔除掉原始观测曲线(图1c)中的长期趋势变化，得到原始观测曲线中的高频信息(图1d).

(3) 图1d中剔除掉趋势变化后的高频信息曲线，可看出类似于“水位固体潮半月潮”的波动变化.从图1a可看出，虽然三阶多项式法可很好地拟合得到观测值的趋势变化，但从局部来看，也有不吻合的部分(图1b中3∶00—5∶00时段)，导致其去趋势残差(高频信息)值较大.另外，水位半月潮是由于水位固体潮响应幅度和相位的规律性变化引起，在利用三阶多项式进行拟合去趋势时，这种规律性的变化也会直接影响到拟合残差的波动变化.因此，结合图1a和图1b，这种类似于“水位固体潮半月潮”的信息，并非真正的水位半月潮信息，而是受其影响的拟合残差的波动信息.

图1　四川南溪台水位分钟值观测曲线及其高频信息Fig.1　Curves of water level recorded (a—b) and high-frequency information (c—d) at Nanxi station, Sichuan

地震流体观测资料，其长周期的趋势变化受构造环境和水文地质条件的影响较大(郭增建等，1974；杨明波等，2009)，而其短周期的高频信息在无其他干扰因素的条件下，多表现为一种服从正态分布的随机信号.从统计的角度分析，在去除趋势变化后，且无其他影响因素存在的条件下，观测资料变化值服从正态分布X(t)～N(0,σ).但是，如果存在某种干扰(如构造活动或应力场变化)，σ值就会发生波动，并且σ值的波动会依干扰强度的变化而发生变化.从能量传递的角度考虑，其波动值服从对数正态分布lnσ～N(0,λ)(Castaing et al., 1990)，Z(t)值的概率密度分布可表示为

(2)

式中，σ和λ分别为Z(t)和lnσ的标准差，σ0为σ的数学期望值.

如图2所示，分别服从正态分布和对数正态分布的两组信号X(图2a)和Y(图2b)，其概率密度分布分别为图2d和图2e，由二者乘积合成的信号Z=XY(图2c)的概率密度分布为图2f.由图可以看出，单纯服从正态分布的原始信号X，在加入干扰信号Y后，其合成信号Z的概率密度分布P(Z)与X的概率密度分布P(X)形态有了明显的不同，并且这种形态上的不同会依据干扰信号的强度而发生变化.如图3所示，同一个原始信号X，在不同的干扰强度(干扰强度决定了λ值的大小)下，其概率密度分布表现出了明显的不同，干扰强度越大，其合成信号的概率密度分布P(Z)与正态分布P(X)在形态上的偏离程度越大.

因此，在理想状态下观测资料中含有的高频信息其概率密度分布应符合正态分布N(0,σ)，高频信息的标准差σ相对稳定，干扰信息的λ值接近于零.

图2　随机信号及其概率密度分布(a) 服从正态分布的随机信号； (b) 服从对数正态分布的随机信号； (c) 合成信号； (d) 正态分布的概率密度； (e) 对数正态分布的概率密度； (f) 合成信号的概率密度.Fig.2　Random signals and their probability density distribution(a) Normal distribution of random signals； (b) Logarithm normal distribution of signals； (c) Synthetic signals； (d) Normal distribution of probability density； (e) Logarithm normal distribution of probability density； (f) Probability density of synthetic signals.

图3　随机信号的叠加及其概率密度分布Fig.3　Stack of random signals and their probability density distributions

如果存在干扰因素，如地震孕育期间构造应力状态的变化，σ值就会有发生变化，那么观测资料中含有的高频信息其概率密度分布就会偏离正态分布，干扰信息的λ值就会大于零.对于一组类似于图1中所提取的高频信号Z(t)，由式(2)可知， Castaing等人(1990)将其概率密度分布的数学表达式表示为

(3)

式中，Z表示观测数据，P为该数据的概率密度分布，F为符合正态分布的随机信号，其标准差为σ，而G为未知的干扰信号，其标准差为λ.F和G的数学表达式分别为

(4)

(5)

实际应用时，可先统计得出原始数据的概率密度分布，然后依式(3)求得干扰G的标准差λ，一般用最大似然法或最小二乘法拟合求得.为得到最佳拟合值，用表达式为

(6)

式中，P(z)表示观测数据的概率密度分布，其标准差为σ，PCastaing为利用式(3)计算得到的概率密度分布，σCastaing为其相应的标准差.

原始数据的概率密度分布和式(3)计算得到的概率密度分布，依据式(6)可求得二者之间的相似程度，当式(6)中的2为最小时，表示二者相似程度最高，进而可得到干扰项G的标准差λ.图4为利用上述方法统计并计算的四川南溪台水位分钟值不同时段的高频信息概率密度分布，可以看出水位高频信息的概率密度分布值，在几乎不存在干扰信号的时段(2007-09-12—2007-10-12)能很好地吻全正态分布(即λ值接近于0，λ=0.03)，而在临近地震发生时段(2008-04-12—2008-05-12)，明显偏离正态分布(即λ=0.32)，说明此时段存在某种干扰信号，这种干扰有可能与周边构造活动或地震孕育有关.

图5为利用以上方法处理的2007年9月12日—2008年5月11日四川南溪台水位分钟值高频信息干扰信号强度时序值曲线(30天窗长，2日滑动，短横线表示拟合误差)，处理结果显示在2008年1月上旬前，该台水位受干扰强度很小(λ值接近于0)，而在2008年1月上旬后干扰强度逐渐增强，直到2008年5月12日汶川8.0级地震发生，截止地震发生前λ2值超过0.6.这种干扰信号的增强，预示着观测点周边构造活动或应力积累水平的增加.由此可见，流体观测资料高频信息中隐含着与构造活动相关的前兆信息，利用概率密度分析方法在一定程度上可以识别出这种前兆信息.

图4　四川南溪台水位分钟值高频信息概率密度分布(不同时段)Fig.4　Probability density distributions of high-frequency information for water-level minute values at Nanxi station of Sichuan (in different time intervals)

图5　四川南溪台水位分钟值高频信息及其概率密度Fig.5　High-frequency information from water-level minute values at Nanxi station of Sichuan and its probability density

图6　汶川8.0级地震前流体高频信息异常曲线Fig.6　Curves of high-frequency fluid anomaly information before Wenchuan MS8.0 earthquake

3汶川8.0级地震前流体异常提取

本文利用上述概率密度分析方法，处理了中国大陆地区南北地震带上72个流体观测点的水位、水温观测数据(共计122条数据)，数据时间段为2008年汶川8.0级震前一年(2007年9月12日—2008年5月11日).处理结果显示，共有30条数据在地震前出现了异常变化，即干扰信号的λ值明显增强.异常测项中水位出现异常的有16项，水温有14项，个别观测点水位和水温均出现了概率密度分布异常(如云南开远观测点).图6为其中部分异常变化的λ2值曲线，可以看出水位、水温的高频信息异常多出现在地震前3个月内，并且这种异常多以震前上升和高值为主.另外，各观测点异常出现的时间和λ2最高值不尽相同，在一定程度上反映了观测点周边应力水平受强震孕育影响的程度也各不相同.部分测点的λ2值在地震前出现了起伏变化(如图6中的云南姚安水位)，这也反映了局部构造活动或区域应力场在受强震孕育影响后的反复调整过程.

图7a所示为本文所选南北地震带流体观测点中出现高频信息异常变化和未出现异常变化测点的空间分布图，由图可以看出，出现异常变化的测点多分布在滇西南构造带，震源区以北出现异常的测点较少.图7b所示为出现异常测点的λ2高值的大小及其出现时间的空间分布图，虽然这些异常点的λ2高值在空间分布上无明显规律，但其出现时间有明显的由南向北推移的特征：如图7b所示，较早出现的异常多集中在南北地震带南段的云南南部区域，其异常出现时间多在地震前90天以上(蓝色虚线以下，截止2008年5月12日12时)，之后出现的异常多集中在云南北部地区，其异常出现时间多为60—90天(蓝色虚线与红色虚线之间)，最后出现的异常点多集中在靠近震中的四川地区，其异常出现时间在震前30天左右(红色虚线以上).

利用概率密度分布方法提取的汶川8.0级地震前流体高频信息异常变化特征显示： (1) 地震前出现高频信息异常的观测点多集中在滇西南构造带； (2) 异常出现的时间有由南向北推移的特征； (3) 异常的λ2值在空间分布上无明显的变化规律.由概率密度分布法的分析原理可知，流体观测资料高频信息中的λ值的变化在一定程度上反映了观测点受周边构造应力调整的作用，结合以上异常变化特征，分析认为在汶川8.0级孕育过程中，滇西南构造带的构造应力作用或构造活动受孕震过程的影响程度较高，并且这种影响有由南向北推进的过程，这种影响是否于大区域的构造动力学背景有关，值得进一步深入分析.另外，异常的λ2值大小无明显空间分布的特征表明，λ2值的大小可能与观测点本身对构造应力作用的敏感程度有关，如不同的构造位置或不同的水文地质条件等.

4结论与讨论

作为一种统计学方法，概率密度分布法曾多次用在地震前兆异常识别(Manshour et al., 2009，2010)、太阳风扰动分析(Carbone et al., 1996； Burlaga and Lazarus, 2000；Sorriso et al., 2001)、外汇市场分析(Ghashghaie et al., 1996)以及人体心率异常检测(Kiyono et al., 2005)等领域，本文将其引入地震地下流体数字化资料分析中，期望能为资料分析和异常识别提供一种新的手段.

当观测资料中包含的地震异常比较显著时，可从原始观测曲线中直接识别(何案华等，2012)，但也有一些异常是隐含在观测资料中的，需要用一些数学处理方法来识别(邱泽华等，2010,2012；刘琦和张晶，2011；晏锐等，2011).概率密度分布法识别地震异常前，需先进行观测数据的高频信息提取，提取方法也有多种，如三阶多项式拟合去趋势法、小波分析法等、经验模态法等，提取效果均可满足实际工作需要.图8为分别利用三阶多项式和小波分析方法提取的四川南溪台水位分钟值观测数据中的高频信息及其概率密度分布曲线，从图可以看出，尽管两种方法的处理结果稍有不同，但利用其处理出的高频信息均可有效地识别出地震前的异常信号，即地震前概率密度分布曲线的λ2值出现了明显的高值变化.

本文利用概率密度分布方法提取的汶川8.0级地震前流体高频信息异常变化结果表明，地震前存在一批流体高频信息异常，这些异常多集中在滇西南构造带，并且异常出现的时间有由南向北推移的特征，这些异常有可能与大区域的构造动力学背景有关.

图7　汶川8.0级地震前流体高频信息异常空间分布Fig.7　Spatial distributions of high-frequency fluid anomaly information before Wenchuan MS8.0 earthquake

图8　三阶多项式与小波分解法对比Fig.8　Comparison of three-order polynomial and wavelet decomposition methods

从板块运动的角度分析，汶川大地震是印度板块推挤欧亚板块的结果，其最根本动力来源是青藏高原和华南地块之间的相对运动在断裂带上产生的能量积累和释放(Wang et al., 2001；Zhang et al., 2004；张培震等，2008).这种大尺度范围的构造动力学作用，也会影响到地壳内部岩体孔隙介质的变形，进而引起地下水流状态(如流速、流向、水温等)的变化.结合概率密度分布法自身的原理，当流体敏感观测点所在区域存在明显的构造动力作用时，原本服从正态分布规律的高频信息，就会偏离正态分布，指示干扰强度的λ值就会出现高值异常变化.

从异常的空间分布来看，虽然震源区以北的西北地区和以南的西南地区背景观测点数量相当，但绝大多数异常集中分布于震源区以南的滇西南构造带，而在震源区以北的区域分布很少.张培震等的研究表明，汶川8.0级特大地震孕育和发生是不同力学性质、不同变形方式的多个地质单元相互作用的结果，涉及3个地质单元：川西高原由于地壳结构的软弱而发生强烈震前变形，是孕震的变形单元；龙门山断裂带的组成物质强度很高、断裂产状不利于滑动，震前的变形很缓慢，是孕震的闭锁单元；四川盆地由于刚性好不易变形而对青藏高原的向东扩展起着阻挡作用，是孕震的支撑单元.震前变形主要发生在变形单元(川西高原)，闭锁单元(龙门山断裂带)只发生缓慢的变形和运动，支撑单元由于起着阻挡作用而与闭锁单元之间几乎不发生相对运动或变形(张培震等，2009).同时，川滇地区GPS现代地壳运动速度场和应变率场分析结果也表明，汶川8.0级地震震源区以南的滇西南构造带的地壳运动速度和应变率明显高于震源以北区域(张培震等，2009； Wu et al., 2011).可见，流体高频信息异常空间分布与周边地壳运动场具有很好的一致性，因此认为二者之间有一定的关联性.

从异常出现的时间来看，异常最先出现在距离震源区较远的滇南地区，之后是川滇交界及其附近区域，最后向震源区推进.参考邱泽华等(2010)对异常形成机理的分析与Lockner等(1991)、Lockner(1996)的实验结果，在断层成核之前，声发射源在一个相当广大的区域上发展，伴随成核过程开始，声发射源向断层部位集中，周围远处的声发射源趋于消失.这种观点可以解释高频信息的异常演化特征：汶川8.0级地震发生前，也即断层破裂前，距离震中较远的区域异常开始出现，之后逐步向震中位置集中.与本文研究结果类似，晏锐等的结果表明汶川震前水氡浓度都存在明显的临界慢化现象(晏锐等，2011).

图9为云南高大井水温高频信息概率密度自观测以来至今λ2值时序曲线与周边地震的对应图，从图中可以看出在周边发生的多次中强地震前，该井水温概率密度λ2值均出现了不同幅度的升高(λ2值高于0.02).对比异常信息表(表1)，认为异常幅度值与地震震级存在一定的正相关，即震级大越大异常幅度越高.同时，异常幅度值与震中距也有一定的关系，除了在周边中强地震前λ2值出现异常升高外，在周边发生较近区域个别中小地震发生前也出现了明显的异常升高现象.毕竟概率密度分布法是一种统计学手段，虽然能有效地识别出部分地震前的异常信号，但在某些λ2值异常升高时段周边并没有发生地震，如2008年9—11月λ2值的异常升高并无相应地震对应，这可能与不同地震的孕震机理不同相关.

表1　云南高大井周边地震信息表

识别地震前兆信息是一项较困难的工作，如何界定某个信息是正常还是异常更加不易.在积累探索过程中，多以未发生地震前的信息为背景值，当信号值与该背景值之间出现明显差异可认为是异常信号.但是，由于每个台站自身的井孔条件、信息捕捉能力以及所处的构造位置不同，其地震前出现异常的形态和幅值也不尽相同，如图6中各台点高频信息的概率密度分布λ2值，在汶川8.0级地震前其异常幅度也各不相同.如何在地震前有效地识别出异常信息，这需要一个不断积累的过程，须对每个台站、每个测项的历史资料震例进行总结，进而归纳提炼出其异常指标，方可在日后的地震预报应用中加以应用.从图9可知，本文中云南高大井水温的高频信息的λ2峰值高于0.02即可认为是异常信号，且其异常峰值多出现在震前45天内.进一步的统计结果显示，该井水温震前的高频信息异常除了与震级、震中距有关外，还与地震的震源机制密切相关，如图10所示.图10b中方块大小表示有地震对应的高频信息λ2峰值时间与发震时间的时间差(dT，单位为天，地震序号与表1相对应)，方块颜色表示时间差的大小，实线为时间差约为10天的分界；图10a为各地震的震源机制解及其空间分布，蓝色震源球表示异常峰值时间与发震时间的差值dT小于10天，红色震源球表示dT大于10天.从图10可以看出，表1所列21个地震中，异常峰值时间与发震时间的差值大于10天的地震几乎全为走滑型地震(2008年汶川8.0级地震EQ1例外)，小于10天的地震多为非走滑型地震.这种震源机制解与异常持续时间之间的相关性，一定程度上反映了地震在孕震、发震过程中的区域构造应力作用或大范围的构造动力学背景.

图9　云南高大井水温高频信息概率密度时间序列Fig.9　Temporal sequence of probability density of high-frequency information from water temperature data of Gaoda well in Yunnan

图10　云南高大井水温高频信息异常对应的地震信息(a)地震震源机制解及其空间分布；(b)地震震级、震中距及异常持续时间.Fig.10　High-frequency anomalies in water temperature of Gaoda well in Yunnan and corresponding earthquake data(a) Focal mechanism solutions of earthquakes and its spatial distribution； (b) Magnitude, epicenter distance and duration of abnormities.

因此，与晏锐等(2011)的研究类似，本文的研究是回溯性的个案研究，虽然概率密度分布法可有效地识别出数字化流体观测资料中的高频信息异常，但仅是一种数学处理方法，对高频信息异常的形成机理以及与地震孕育过程相互关系的分析仅限于初步讨论，有待于积累更多的震例分析和进一步的深入研究.

致谢本项工作得到中国地震局监测预报司预报管理处的大力支持，中国地震局预测研究所邵志刚博士在研究过程中给予了悉心指导和帮助，审稿专家的修改建议使作者对该项研究的内容有了更深层次认识，在此表示衷心的感谢.

References

Burlaga L F, Lazarus A J. 2000. Lognormal distributions and spectra of solar wind plasma fluctuations: Wind 1995—1998.JournalofGeophysicalResearch, 105(A2): 2357-2364.

Cao L L, Gao A T. 2010. Coseismic response characteristics of digital records of water level and water temperature in Gansu caused by WenchuanMS8.0 earthquake.ActaSeismologicaSinica(in Chinese), 32(3): 290-299.

Carbone V, Veltri P, Bruno R. 1996. Solar wind low-frequency magnetohydrodynamic turbulence: extended self-similarity and scaling laws.NonlinearProcessesinGeophysics, 3: 247-261.

Castaing B, Gagne Y, Hopfinger E J. 1990. Velocity probability density functions of high Reynolds number turbulence.PhysicsD, 46(2): 177-200.

Daubechies I. 1988. Orthonormal bases of compactly supported wavelets.CommunicationsonPureandAppliedMathematics, 41(7): 909-996.

Ghashghaie S, Breymann W, Peinke J, et al. 1996. Turbulent cascades in foreign exchange markets.Nature, 381(6585): 767-770.

Guo Z J, Qin B Y, Feng X C. 1974. Discussion on the change of ground-water level preceding a large earthquake from an earthquake source model.ActaGeophysicaSinica(in Chinese), 17(2): 99-105.

He A H, Zhao G, Liu C L, et al. 2012. The anomaly characteristics before Wenchuan earthquake and Yushu earthquake in Qinghai Yushu and Delingha geothermal observation wells.ChineseJ.Geophys. (in Chinese), 55(4): 1261-1268, doi: 10.6038/j.issn.0001-5733.2012.04.021.

Huang N E, Shen Z, Long S R, et al. 1998. The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis.Proc.R.Soc.LandA, 454(1971): 903-995.Kiyono K, Struzik Z R, Aoyagi N, et al. 2005. Phase transition in a healthy human heart rate.PhysicalReviewLetters, 95(5): 058101.

Liu Q, Zhang J. 2011. Application of S transform in analysis of strain changes before and after Wenchuan earthquake.JournalofGeodesyandGeodynamics(in Chinese), 31(4): 6-9.

Liu Y W, Yang X H, Liu Y M, et al. 2005. The characters of groundwater responses induced by Sumatra Ms 8.7 earthquake. In: Survey and Forecast Department of CEA ed. Effect in Chinese Mainland Induced by Sumatra 8.7 Earthquake in Indonesia, 2004 (in Chinese). Beijing: Earthquake Press, 131-258.

Lockner D A. 1996. Brittle fracture as an analog to earthquakes: Can acoustic emission be used to develop a viable prediction strategy.J.Acous.Emission., 14(3-4): S88-S101.

Lockner D A, Byerlee J D, Kuksenko V, et al. 1991. Quasi-static fault growth and shear fracture energy in granite.Nature, 350(6313): 39-42.

Manshour P, Ghasemi F, Matsumoto T, et al. 2010. Anomalous fluctuations of vertical velocity of Earth and their possible implications for earthquakes.PhysicalReviewE, 82(3): 036105.Manshour P, Saberi S, Sahimi M, et al. 2009. Turbulencelike behavior of seismic time series.PhysicalReviewLetters, 102(1): 014101.

Qiu Z H, Tang L, Zhang B H, et al. 2012. Extracting anomaly of the Wenchuan earthquake from the dilatometer recording at NSH by means of wavelet-Overrun Rate Analysis.ChineseJ.Geophys. (in Chinese), 55(2): 538-546, doi: 10.6038/j.issn.0001-5733.2012.02.016.

Qiu Z H, Zhang B H, Chi S L, et al. 2011. Abnormal strain changes observed at Guza before the Wenchuan earthquake.ScienceChinaEarthSciences, 54(2): 233-240.

Sorriso-Valvo L, Carbone V, Giuliani P, et al. 2001. Intermittency in plasma turbulence.PlanetaryandSpaceScience, 49(12): 1193-1200.

Sun X L, Liu Y W. 2008. Characteristics and mechanism of co-seismic responses in the Tayuan well.EarthquakeResearchinChina(in Chinese), 24(2): 105-115.

Sun X L, Liu Y W. 2012. Changes in groundwater level and temperature induced by distant earthquakes.GeosciencesJournal, 16(3): 327-337.Sun X L, Liu Y W, Yan R. 2011. Application of empirical mode decomposition method to groundwater data.JournalofGeodesyandGeodynamics(in Chinese), 31(2): 80-83.

Wang Q, Zhang P Z, Freymueller J T, et al. 2001. Present-day crustal deformation in China constrained by Global Positioning System measurements.Science, 294(5542): 574-577.

Wu Y Q, Jiang Z S, Yang G H, et al. 2011. Comparison of GPS strain rate computing methods and their reliability.GeophysicalJournalInternational, 185(2): 703-717, doi: 10.1111/j.1365-246X.2011.04976.x.

Yan R, Jiang C S, Zhang L P. 2011. Study on critical slowing down phenomenon of radon concentrations in water before the WenchuanMS8.0 earthquake.ChineseJ.Geophys. (in Chinese), 54(7): 1817-1826.Yang C J, Feng Z S, Song D W, et al. 2005. Preliminary application of wavelet analysis method to extract the characters of tidal factor in well water level before earthquakes.NorthwesternSeismologicalJournal(in Chinese), 27(2): 163-167. Yang M B, Kang Y H, Zhang Q, et al. 2009. Tendencious fall of groundwater table in Beijing region and recognition of earthquake precursor information.ActaSeismologicaSinica(in Chinese), 31(3): 282-289.

Zhang P Z, Shen Z K, Wang M, et al. 2004. Continuous deformation of the Tibetan Plateau from global positioning system data.Geology, 32(9): 809-812.

Zhang P Z, Wen X Z, Xu X W, et al. 2009. Tectonic model of the great Wenchuan earthquake of May 12, 2008, Sichuan, China.ChineseSci.Bull. (in Chinese), 54(7): 944-953.

Zhang P Z, Xu X W, Wen X Z, et al. 2008. Slip rates and recurrence intervals of the Longmen Shan active fault zone, and tectonic implications for the mechanism of the May 12 Wenchuan earthquake, 2008, Sichuan, China.ChineseJ.Geophys. (in

Chinese), 51(4): 1066-1073.Zhang S X, Zhang Z G, Liu J M, et al. 2002. Application and analysis of digitized observation data of underground water level in Hebei Province.Earthquake(in Chinese), 22(4): 89-93.

附中文参考文献

曹玲玲, 高安泰. 2010. 汶川MS8.0地震引起的甘肃数字化水位、水温同震响应特征分析. 地震学报, 32(3): 290-299.

郭增建, 秦保燕, 冯学才. 1974. 从震源孕育模式讨论大震前地下水的变化. 地球物理学报, 17(2): 99-105.

何案华, 赵刚, 刘成龙等. 2012. 青海玉树与德令哈地热观测井在汶川与玉树地震前的异常特征. 地球物理学报, 55(4): 1261-1268, doi: 10.6038/j.issn.0001-5733.2012.04.021.

刘琦, 张晶. 2011. S变换在汶川地震前后应变变化分析中的应用. 大地测量与地球动力学, 31(4): 6-9.

刘耀炜, 杨选辉, 刘永铭等. 2005. 地下流体对苏门答腊8.7 级地震的响应特征.∥中国地震局监测预报司编. 2004年印度尼西亚苏门答腊8. 7 级大地震及其对中国大陆地区的影响. 北京: 地震出版社, 131-258.

邱泽华, 张宝红, 池顺良等. 2010. 汶川地震前姑咱台观测的异常应变变化. 中国科学: 地球科学, 40(8): 1031-1039.

邱泽华, 唐磊, 张宝红等. 2012. 用小波-超限率分析提取宁陕台汶川地震体应变异常. 地球物理学报, 55(2): 538-546, doi: 10.6038/j.issn.0001-5733.2012.02.016.

孙小龙, 刘耀炜. 2008. 塔院井水位和水温的同震响应特征及其机理探讨. 中国地震, 24(2): 105-115.

孙小龙, 刘耀炜, 晏锐. 2011. 经验模态分解法在地下水资料处理中的应用. 大地测量与地球动力学, 31(2): 80-83.

晏锐, 蒋长胜, 张浪平. 2011. 汶川8.0级地震前水氡浓度的临界慢化现象研究. 地球物理学报, 54(7): 1817-1826.

杨丛杰, 冯志生, 宋德伟等. 2005. 小波分析方法在提取井水位潮汐因子震前变化特征的初步应用. 西北地震学报, 27(2): 163-167.

杨明波, 康跃虎, 张庆等. 2009. 北京地下水位趋势下降动态及地震前兆信息识别. 地震学报, 31(3): 282-289.

张培震, 闻学泽, 徐锡伟等. 2009. 2008年汶川8.0级特大地震孕育和发生的多单元组合模式. 科学通报, 54(7): 944-953.

张培震, 徐锡伟, 闻学泽等. 2008. 2008年汶川8.0级地震发震断裂的滑动速率、复发周期和构造成因. 地球物理学报, 51(4): 1066-1073.

张素欣, 张子广, 刘俊明等. 2002. 数字化水位观测资料的应用研究. 地震, 22(4): 89-93.

(本文编辑张正峰)

Extracting high-frequency anomaly information from fluid observational data:a case study of the WenchuanMS8.0 earthquake of 2008

SUN Xiao-Long1,2, WANG Guang-Cai1, YAN Rui3*

1ChinaUniversityofGeosciences(Beijing),Beijing100083，China2InstituteofCrustalDynamics，Beijing100085,China3ChinaEarthquakeNetworksCenter,Beijing100045,China

AbstractWith improvement of underground fluid observation technology, the sampling frequency of data has been obviously improved, especially after the digital transformation. There is a great deal of tectonic information in high-frequency observational data. It is a key problem how to extract the useful high-frequency information from theses data. We introduced the PDF (probability density distribution) method, and analyzed the observed fluid data of water level and water temperature at 72 stations in the area of the north-south seismic zone before the Wenchuan MS8.0 earthquake on 12 May, 2008. The analysis results show that 16 water level and 14 water temperature data had the high-frequency anomalies before the MS8.0 earthquake, and the points with abnormal information were concentrated in the tectonic belts of southwest Yunnan. The high-frequency anomalies appeared from south to north over time. On the basis of these results we considered that the PDF method has certain reliability and applicability to extract the high-frequency anomaly information from observational data of fluid. It can provide a reference for information extraction of high-frequency digital observation.

KeywordsFluid observational data； Probability density distribution； High-frequency abnormal information； Wenchuan MS8.0 earthquake

基金项目地震行业科研专项经费项目(201408019-03)和国家自然科学基金资助项目(41502239)共同资助.

作者简介孙小龙,男,1981年生,副研究员，主要从事地震地下流体研究及地震预测研究.E-mail:xlsun04@163.com *通讯作者晏锐，男，1978年生，副研究员，主要从事地震地下流体研究与地震预测研究.E-mail:ray_2005@163.com

doi:10.6038/cjg20160512 中图分类号P315

收稿日期2015-06-20，2016-01-12收修定稿

孙小龙， 王广才， 晏锐. 2016. 利用概率密度分布提取流体观测资料中的高频异常信息——以2008年汶川8.0级地震为例.地球物理学报，59(5):1673-1684,doi:10.6038/cjg20160512.

Sun X L, Wang G C, Yan R. 2016. Extracting high-frequency anomaly information from fluid observational data: a case study of the WenchuanMS8.0 earthquake of 2008.ChineseJ.Geophys. (in Chinese),59(5):1673-1684,doi:10.6038/cjg20160512.