王秀英 张聪聪 王成亮
1 中国地震局地壳应力研究所,北京市安宁庄路1号,100085
2 邯郸市地震局,邯郸市丛台路496号,056008
地震前兆观测数据是进行地震预报和各种地球科学研究的基础。随着数据应用研究的深入,前兆观测系统不断改善,前兆观测数据的观测精度和采样率不断提高,导致观测数据量激增[1]。由于前兆观测过程中观测数据会受到各种干扰和影响,工作人员需对由干扰导致的“问题数据”进行预处理。另外,在进行数据的常规方法或特定方法应用之前,也需对问题数据进行先期处理,以满足不同方法对数据的要求。由于前兆观测数据量巨大,问题数据的定位靠人工检查,不仅缺少判别的客观性,效率也很低。由于观测数据的连续变化,基于这些数据提取的常用统计特征量均值、方差等也在不断变化,变化的特征量不能反映不同时段数据总体特征,也无法利用它们进行比较、分析和判断,常规的均值、方差等特征量数据判别方法在用于地震前兆观测数据分析时也不适用。因此,需要研究更有效的方法来进行此项工作。
本文基于大数据挖掘的思路,设计一种利用信息熵值作特征量检测地震前兆观测数据的方法,可以快速从大量数据中检测出问题数据,大大缩短数据预处理中非正常数据的定位时间。
数据特征量提取方法设计中,借鉴信息领域信息熵的概念,通过前兆观测数据信息熵的定义来提取前兆观测数据的特征值。
Shannon[2]将信息熵定义为:离散随机事件出现概率的不确定性。一个系统越有序,信息熵就越低;越混乱,信息熵就越高。所以,信息熵是系统有序化程度的一个度量。
假设X是一个离散随机变量,它的取值范围R={x1,x2,…,xn}是有限可数的。设pi=P{X=xi}为事件xi的发生概率,则X的信息熵定义为[2]:
由信息熵定义知[3],对于地震前兆观测时间序列,如果数据变化完全随机无规律,则每个事件出现的概率大致相等,其信息熵值最大;相反,如果数据变化简单规律,则某类事件出现的概率较大,其他类事件出现的概率较小,表现为信息熵越小;如果观测数据为一个常数,则某类事件的出现概率为1,其他类事件的出现概率为0,最终的信息熵值为0。因此,信息熵可以反映观测数据的变化特性,通过检验信息熵数据的变化可以检测数据的变化。这样,通过一个特征量使数据整体的变化特性得以反映,从而简化时间序列数据的检测。
地震前兆观测数据大多序列较长,重复数据较少。因此,在进行信息熵特征量提取之前必须对观测时间序列进行有效降维[4]。
对时间序列X用符号化方法降维及计算信息熵的具体操作步骤如下:
1)计算序列X的均值。
2)将序列中每一个观测值xi与序列均值比较,xi≥时,取值1;xi<时,取值0。
3)得到一个与原序列等长的0、1序列。
4)将序列符号化。根据符号化字符数长度的要求,将01序列按不同长度截取。如符号化字符数取16时,0、1序列截取长度为4,每个长度为4的0、1序列可以确定一个字符,这样原来长度为1 440的时间序列,降维后成为长度为360的字符序列。
5)统计各字符出现的概率,按式(1)计算信息熵。
通过符号化方法可降低数据序列的维度。至于序列长度降低的程度,需要依据数据序列的特征及要提取的特征选择合适的取值。
为说明信息熵对数据整体情况的反映,以地震前兆台网产出的观测数据为例进行验证。
图1为云南省云龙地震台水平摆倾斜观测北南分量2008~2012得到的信息熵曲线,信息熵计算中字符长度取值4。可以看到,该测项的信息熵值大致在0.7~0.8,存在一些信息熵值很大和很小的点;2008、2009和2011年信息熵数值变化较平稳,2010年变化则较为剧烈,2012和2013年更甚,尤其是2013年,信息熵数据极不稳定,变化非常剧烈,偏离正常值更小的情况非常多。
图1 云南省云龙地震台水平摆倾斜观测北南分量2008~2012年观测数据信息熵曲线Fig.1 Informational Entropy Curve of the north-south tilt observation of horizontal pendulum from 2008to 2012at Yunlong station,Yunnan province
对这些信息熵值明显偏离正常变化范围的点,利用其日变时序曲线进行检查。
筛选2008~2013年得到的信息熵数据,信息熵值≤0.6的数据共有44条,其中2008年2条、2010年2条、2011年4条、2012年6条、2013年30条。这44条信息熵值对应的观测数据时序曲线全部对应单点大幅突跳或大幅短时台阶的情况(图2)。
其中图2(a)的短时段数据跳跃从时间上与汶川地震发震时间对应,应是汶川地震的地震波造成的观测数据异常,而且幅度非常大。图2(b)~(d)存在明显的单点突跳或短时段数据台阶,由于幅度较大,导致原始曲线形态无法反映。作为对比,图2(e)是该测项正常形态的日变曲线。图2(e)与(a)~(d)中各条曲线对比可以确定,图2(a)~(d)中的时序曲线确实存在偏离正常情况的变化,而这种变化通过信息熵数据也得以反映。图2(e)正常曲线形态的观测数据计算所得信息熵为0.73,属于图1中大多数信息熵值的变化范围;另外几个信息熵值都远远小于正常信息熵值的变化范围,图2(a)~(d)曲线反映了这些数据的确存在异常变化。
经逐一对比,在2008~2013年所得信息熵数据中筛选44个明显偏小的信息熵值,对应日观测时序数据全部存在明显的非正常变化。
筛选2008~2013年计算信息熵中数值≥0.9的数据,共382个,其中2008年14个、2009年10个、2010年39 个、2011年11 个、2012年27 个、2013年90个。382个信息熵对应的观测数据时序曲线全部存在偏离正常的情况。偏离情况大致可以分为:1)正常形态上存在短时较大幅度干扰;2)短时大幅度干扰导致原曲线形态压缩改变;3)有异于常规的形态呈现;4)观测数据有连续高频干扰存在;5)连续随机干扰及较大幅度的突跳,导致曲线形态改变;6)基本形态存在连续小幅度干扰和相对较大幅度的突跳。图3为几类比较典型的观测数据异常曲线形态。
由图3可知,信息熵值较大时,观测数据异于正常形态的情况多种多样。这说明信息熵值增加时,数据中无规律变化的成分增加。信息熵值的突然大幅变化可以反映其代表的观测时序数据发生较大变化或受到较大干扰。
图2 信息熵异常小值及正常值对应观测数据时序曲线Fig.2 Daily observation data curves corresponding to the very small and normal informational entropy values
图3 信息熵异常大值对应观测数据时序曲线Fig.3 Daily observation data curvescorresponding to the very large informational entropy values
由计算实例可知,对于信息熵值明显偏离正常变化范围的情况,对应的数据都存在比较严重的问题,这说明信息熵的确有反映原观测时间序列某些特性的能力。对于信息熵值介于明显偏大和明显偏小区间的数据,即信息熵值介于[0.6,0.9]的数据,大部分表现为正常形态,少量处于正常范围内的数据有局部小的干扰。虽然这种干扰从曲线形态上可以明显看到,但由于幅度较小,持续时间较短,无法在信息熵中得以反映。但总的表现规律为,信息熵较大时(大于0.9),数据曲线形态表现较复杂;信息熵较小时,曲线形态更趋简单,符合信息熵对数据特性的描述。
分析表明,信息熵具有反映数据总体变化特征的能力,利用信息熵可以快速发现存在较大异常变化的数据。对2008~2013年6a每日观测数据的信息熵计算提取,在Intel Core 2 Quad CPU,2.66GHz主频的计算机上,用时约2′29″。对于如此巨大的数据量,如果以人工逐日查看曲线的方式进行检测,在短时间内根本无法完成。
信息熵之所以能够反映观测数据的某些形态变化,与该方法中数据的降维符号化方法相关。当数据中有幅度特别巨大的突跳时,即使持续时间很短,也可能造成均值的改变,进而影响符号化过程,最终在信息熵数值中得以反映。另外,如果数据较原有形态发生了较大改变,也会导致符号化后各符号出现的概率发生变化,最终使得信息熵值改变。所以,为了降维所选取的符号化方法在信息熵计算过程中非常关键。计算实例所用方法是一种最简单的方法。具体应用时,可以结合应用目的和所用数据的特征设计不同的符号化方法[5-6],以反映不同的曲线变化特征,从而达到快速检测某些数据特性变化的目的。
信息熵值可以作为一个数据序列的特征量使用。地震前兆观测测项众多,数据在量纲和量级上都存在很大差异,观测数据无法直接比较,采用统计中的均值或方差等特征量也不能解决这个问题。另外,由于观测数据量级的差异,即使同一观测测项不同,观测点的数据也不能直接比较。而信息熵值无量纲,且能反映数据样本的某些变化特征,可以利用信息熵特征量,对不同前兆观测测项或者同一观测测项不同观测地点观测数据进行直接比较,解决前兆观测数据无法直接对比的问题。
文中给出的信息熵计算是基于原始观测数据序列的,为了突出某些变化特征,可以先行对原始序列进行转换,如差分等操作,再进行信息熵计算提取,以突出反映更多的数据内部特性。
应用示例中仅给出了信息熵值明显偏离正常变化范围的情况,对应的观测数据都存在比较明显的非正常变化。对于信息熵值介于其间的情况,其中也存在一些数据异常的案例,但由于异常数据持续时间较短、幅度较小,通过信息熵数据无法反映。所以,利用信息熵值法处理一些特别严重的问题数据,效果较好。
信息熵方法对于前兆观测时序数据具有比较好的检测效果,可以快速定位问题数据,而且该方法对于所有前兆观测数据都适用,不存在学科或观测测项处理方法上的差别,尤其适用于目前前兆数据中心数据量巨大的情况。
[1]周克昌,赵刚,王晨,等.中国地震前兆台网观测技术系统整合[J].中国地震,2013,29(2):270-275(Zhou Kechang,Zhao Gang,Wang Chen,et al.Upgrade and Integration of the Precursor Observation Network of China[J].Earthquake Research in China,2013,29(2):270-275)
[2]Shannon C E.A Mathematical Theory of Communication[J].The Bell System Technical Journal,1948,27(7):379-423
[3]王栋,朱远甡.信息熵在水系统中的应用研究综述[J].水文,2001,21(2):9-14(Wang Dong,Zhu Yuanshen.Informational Entropy and the State-of-the-Art of Its Application in Hydrology,Water Resources and Water Environment[J].Hydrology,2001,21(2):9-14)
[4]李海林,杨丽彬.时间序列数据降维和特征表示方法[J].控制与决策,2013,28(11):1 718-1 722(Li Hailin,Yang Libin.Method of Dimensionality Reduction and Feature Representation for Time Series[J].Control and Decision,2013,28(11):1 718-1 722)
[5]任江涛,何武,印鉴,等.一种时间序列快速分段及符号化方法[J].计算机科学,2005,32(9):166-169(Ren Jiangtao,He Wu,Yin Jian,et al.A Fast Time Series Segmentation and Symbolization Method[J].Computer Science,2005,32(9):166-169)
[6]钟清流,蔡自兴.基于统计特征的时序数据符号化算法[J].计算机学报,2008,31(10):1 857-1 864(Zhong Qingliu,Cai Zixing.The Symbolic Algorithm for Time Series Data Based on Statistic Feature[J].Chinese Journal of Computers,2008,31(10):1 857-1 864)