8 基于主成分分析和BP神经网络的因素与地震及其预测方法探讨

第28章基于主成分分析和BP神经网络的因素与地震及其预测方法探讨——因素与地震分析28.1案例背景地震是地壳快速释放能量过程中造成的振动。虽然预测地震是世界性难题，但当今科学界普遍认为，有可能反映地震前兆特征的指标不少于10个。已经有专业仪器在多个定点实时按秒记录这些指标的数据，期望通过对记录数据的分析研究找到地震的前兆特征。现已采集到某地2005年1月1日至2010年6月30日按小时观测的10多个指标的数据，和该地区该时期内已发生地震的时刻、经纬度、震级及震源深度的数据。这些数据中隐藏着地震发生的前兆特征。科学地截取这些数据的有用片段，对数据进行合理地预处理，用数学方法揭示地震前兆的数据特征，是一项很有意义的研究工作。案例所给数据中的这10多个指标，究竟哪些与地震的发生有关，有何种关系，是单一关系还是复合关系；除这10多个指标外还有哪些因素，以及题给指标中的哪些指标的哪种数学模型更能反映地震的前兆特征等等，人们迄今仍不很清楚，需要进行深入研究。地震数据的观测是持续进行的，随着时间的推移数据的规模会不断扩大。从中挖掘地震的前兆特征，必须有合理的数学模型，也必须有科学高效的算法分析平台。因此，本案例要求结合附件中给出的实际记录数据，尝试完成以下任务。任务一：分析数据特征，建立数学模型以度量各指标对地震发生的敏感程度。任务二：构造由某些或全部指标构成的综合指标，使其尽可能地集中反映地震发生前的数据特征的统计规律。任务三：结合题给数据，广泛查阅与地震相关的其它指标的数据和分析方法，建立数学模型来研究地震发生前的数量特征。28.2.模型建立本案例的研究背景是地震预报，这是一个目前的世界范围内尚未解决的科学难题。当前人们的科学技术水平尚未达到完全掌握地震的自身规律的程度。目前的地震预报主要是在对地震活动及相关测震指标的历史资料的基础上进行模型建立与分析预报的。28.2.1任务一指标对地震敏感程度分析将指标的敏感度定量为有震报准率，是许绍燮提出的R值评分多指标的其中之一，但又有所改进，具体认为，有震报准率c值如下：c预测有效的数据总数异常数据总数(28-1)本案例中，数据的异常与平均值x相比较，将各数据标准化如下：·2·数学建模30个案例分析min()max()min()xxxxxvvv(28-2)式中：x——指标标准化的数值；x——原始数据；min()xv——一列数据的最小值；max()xv——一列数据的最大值。本章利用VB程序（程序详见28.3.1）统计在平均值上下10%~100%之外数据个数1n，除以有效数据总数N，若满足如下公式，则认为异常：120%nN(28-3)统计满足此异常标准的幅度，即偏离平均的百分比，由于各个指标的分布情况不同，故不同，选取满足20%的百分比。在异常指标中找寻预测有效的数据，与异常数据之比即为最终敏感度，由于选取不同，c的结果不同，根据实际情况，地震预报的不准确性，选择c值较小的，比较各个数据的报准率，得出结论。最终统计结果如表28.1。表28.1：12个指标对地震敏感度的统计表指标标准化的平均值有效数据总数异常数据总数预测有效的数据总数异常度报准率幅度电磁波幅度EW0.7610594816841310.0085740.00242110%电压0.659286463567513710.1620720.0094510%气压0.46155248168495120.0102770.02424260%水温0.033444813923160.0047990.02597410%水位0.6739694816867862810.1408820.04140940%电磁波幅度NS0.64390448168573340.0118960.05933710%雨量0.103117481681600.00033200%气温0.4990124816848813780.1013330.07744370%气氡0.064624816864368320.0172730.12927380%地温0.84935748165175240.0036330.13714340%倾斜仪NS0.4516541947530465840.1564060.191727100%倾斜仪EW0.370851071617984870.1677860.270857100%根据以上分析得出则各指标对地震发生的敏感度比较排序如下：倾斜仪EW倾斜仪NS地温气氡气温电磁波幅度NS水位水温气压电压电磁波幅度EW雨量。28.2.2任务二基于主成分分析法的综合指标模型对题中所给的有可能反映地震活动前兆特征的12个指标（电压、电磁波幅度EW、电磁波幅度NS、地温、水位、气温、气压、水温、气氡、雨量、倾斜仪NS、倾斜仪EW）的分析发现，这些指标可以从不同的侧面反映了地震活动的相关特征。这12个指标中，两两之间还可能存在一定的相关性。此外，在实际的地震预报计算分析过程中发现，常常有些指标参数在某些地震前会出现比较明显的异常，而另一些指标参数并并没有表现出较大的异常；而在正常情况下（无震），也常常会发生某些指标参数出现明显异常的现象，而其他一些指标参数并不出现异常变化。正是指标测量数据的这些特点，给地震预报分析带来的很大的困难。为解决指标参数数据异常变化不一致而带来预报意见不同的影响，构造由全部指标构成一个综合预报指标，并尽可能反映地震发生前指标数据数量特征的统计规律，本章拟采用主成分分析法简化多项预报指标参数，构建地震预报研究的综合指标。第2章·3·28.2.2.1主成分分析法原理主成分分析法可以把多个描述样本特征的、并且可能具有一定相关性的多个指标参数化简为少数几个综合指标的一种统计分析方法。因为每个变量都在不同程度上反映了所研究问题的某些信息，并且指标之间彼此有一定的相关性，因而所得的统计数据反映的信息在一定程度上有重叠。故主成分分析法利用降维的思想，能够在力保原始数据信息，最大限度地减少信息丢失的基础上，对高维变量空间进行降维处理，能够对原始变量系统进行最佳的综合与简化，它可以客观地确定各个指标参数的权重，避免主观判别带来的随意性。在原始数据的基础上，利用主成分分析法经过线性变换和舍弃部分信息，可以找出由若干指标组合而成的综合指标，即若干个主成分。而这些主成分就可以尽可能地反映原来指标的特征与信息，同时彼此间相互独立。设原始变量为12,,,pxxx…，进行主成分分析后得到的主成分（综合变量）为12,,,mzzz…（mp），它们是12,,,pxxx…的线性组合。新变量12,,,mzzz…构成的坐标系是在原坐标系经平移和正交旋转后得到的，称12,,,mzzz…空间为m维主超平面。在主超平面上，第一主成分1z对应于数据变异（贡献率1e）最大的方向。对于23,,,mzzz…，依次有2mee…。因此，1z是携带原始数据信息最多的一维变量，而m维主超平面是保留原始数据信息量最大的m维子空间。主成分分析法的步骤如下：（1）为了排除数量级和量纲不同带来的影响，首先对原始数据进行标准化处理，标准化处理公式如下：'()ijiijixxx（1,2,,;1,2,,ipjn……）（28-4）式中ijx位第i个指标第j个样本的原始数据；ix和i分别为第i个指标的样本均值和标准差。（2）根据标准化数据表'ijpnx（），计算相关系数矩阵ijpnRr=（），其中1()()1nkiikjjijkijxxxxrn（28-5）（3）计算R的特征值和特征向量。根据特征方程0RI，计算特征根i，并使其从大到小排列：12p…，同时可得对应的特征向量12,,,puuu…，它们标准正交。12,,,puuu…称为主轴。这里的I为单位矩阵。（4）计算贡献率ie和累计贡献率mE1iipiie（28-6）11mjjmpiiE（28-7）（5）计算主成分mZ1pmmjjjZux（28-8）这时，各主成分相互独立。·4·数学建模30个案例分析（6）综合分析。一个m维主超平面究竟以多大的精度来近似代替原始变量系统，才能确保尽可能多的原始数据信息？这可以通过求累计贡献率mE来判断。一般取85%mE的最小变量m（mp），就可以确定主超平面的维数m，从而可对m个主成分进行综合分析。（7）根据主成分iz分析得到的主成分和相应的权值（贡献率）ie，计算本文定义的反映地震活动表现异常特征的综合指标：1miiiWeZ（28-9）由于m个主成分已基本保留了这些预报指标参数的数据信息，所以综合指标W包含了这些指标参数从不同方面反映地震活动的基本特征。本章选择题给全部指标参数（12个）进行主成分分析计算。28.2.2.2基于主成分分析法的综合指标模型由于题中所给测量数据的不完整性，本章选取2005年10月29日1:00—2008年12月31日23：00的12个指标的27839组数据进行分析计算，此段之间内发生地震6次。做出12个指标数据随时间的变化过程曲线，结合6次地震发生时间，可以看到虽然地震发生前某些指标参数出现异常，但并不是每次地震前均会出现明显异常，而且各指标参数的异常变化各有所异，大部分指标参数的变化也并不显著。运用主成分分析法的理论，本章编写Matlab程序（程序详见28.3.2）对附件2（12个指标的27839组数据）进行主成分分析计算，得到12个指标参量在个各主成分中的系数（特征向量）及特征值与贡献率，见表28.2。图28.1为6次地震发生前后的一些参量随时间变化的曲线，参数计算的累积时间为30个月，滑动步长为1个月，从图中和一些文献可以看到，一些参量出现较为明显的异常，另一些参量的异常变化则不明显。图28.16次地震气温、倾斜仪EW、电压、地温随时间的变化表28.2为通过主成分分析法得到的上述参量在各主成分中的系数（特征参量）即特征值与贡献率。表5-2各参量在各主成分的系数（特征向量）、特征与贡献率主成分1主成分2主成分3主成分4主成分5主成分6电压-0.01703-0.275080.0556690.6301050.155404-0.38662电磁波幅度EW0.360851-0.44905-0.28043-0.020840.0448260.108391电磁波幅度NS0.364195-0.4391-0.27993-0.004590.0472040.101636第2章·5·地温0.2437290.294381-0.06040.1817930.2460450.36811水位0.3564470.379839-0.099650.2242820.033099-0.07456气温0.158494-0.071990.68592-0.0846-0.07488-0.08253气压-0.323860.379344-0.374770.174520.0836820.052405水温-0.10236-0.118640.3127290.2649110.3548080.711377气氡0.3109220.1841550.1869920.4402020.042456-0.24689雨量0.0206790.0367740.019224-0.365780.868816-0.32486倾斜仪NS-0.42272-0.196920.1908570.0828210.037748-0.05093倾斜仪EW0.3662270.2484580.214058-0.27738-0.11580.051783特征值2.9322.20191.74011.2010.97850.8853贡献率24.43%18.35%14.50%10.01%8.15%7.38%累积贡献率24.43%42.78%57.29%67.29%75.45%82.83%续表5-2主成分7主成分8主成分9主成分10主成分11主成分12电压0.103456-0.55580.156938-0.041360.064989-0.009917721电磁波幅度EW0.003580.0919190.0