目前,在地震活动分析中存在许多预报指标,如b值、C值、D值、Y值、η值等,这些指标从不同侧面反映了地震活动时间、空间和强度特征。目前这类参数较多,同时它们之间还可能存在一定的相关性。在一些中强以上地震前有些参数会出现较明显的异常,而另一些参数并不出现异常,这些都给实际预报带来困难。王炜等(2005)、李滔(2007)用主成分分析法对地震活动性参数进行了约简,同时将该方法应用到江苏溧阳6.0级地震中。本文利用主成分分析法,对1970年以来东北14次5级以上浅源地震进行了研究,实现了多项预报指标的简化。

主成分分析方法是将描述样本特征的多个可能有一定相关性的指标(陆远忠等,1999; 韩渭宾,2003)化为少数几个综合指标的一种统计分析方法。主成分分析法在最大程度保留原始数据信息的基础上,还能对高维变量空间进行综合与简化,并能够客观地确定各个指标的权数,避免了主观随意性(方开泰,1989; 唐启义,冯光明,2002; 吕悦军等,1997)。在原始数据的基础上应用主成分分析法,可以找出由若干个指标线性组合而得到的综合指标,即若干个主成分。这些主成分可以最大化地反映原来指标的信息,同时彼此之间相互独立。

本文选用1970年以来,东北地区14次5级以上浅源地震前后震中附近地区的地震资料,选择与地震活动强度有关的6个参量,即3级以上地震频次N、b值、η值、M_f值、A(b)值和A_C值进行主成分分析。这6个参量之间具有一定的相关性,在不同时段各参量的变化各有所异,预报效果并不理想。但是根据主成分分析方法可以得到反应地震强度异常特征的综合指标W,发现在地震前W均出现明显的异常变化,多数震例震后异常变化恢复,表明综合指标W可以较好地反映地震活动的异常特征。主成分分析方法是简化预报参量的有效工具,在地震预报中具有一定的应用前景。

主成分分析法在力保原始数据信息丢失最小的情况下,对高维变量空间进行降维处理,经过线性组合,以少数的综合变量取代原有的多维变量。

设原始变量为x₁,x₂,…,x_p,进行主成分分析后得到的主成分(综合变量)为z₁,z₂,…,z_m(m<p),它们是x₁,x₂,…,x_p的线性组合。综合变量z₁,z₂,…,z_m构成的坐标系是在原坐标系经平移和正交旋转后得到的,故称z₁,z₂,…,z_m空间为m维主超平面。在主超平面上,第一主成分z₁对应于数据变异(贡献率e₁)最大的方向。对于z₂,z₃,…,z_m,依次有e₂≥e₃≥…≥e_m。因此,z₁是携带原始数据信息最多的一维变量,而m维主超平面是保留原始数据信息最大的m维子空间。

主成分分析法的步骤如下:

(1)为了排除数量级和量纲不同带来的影响,首先对原始数据进行标准化处理:

x'_ij=(x_ij-x^-_i)/(σ_i), i=1,2,…,n.(1)

式中,x_ij为第个指标第个样本的原始数据,x^-_i和σ_i分别为第i个指标的样本均值和标准差。

(2)根据标准化数据表(x'_ij)_p×n,计算相关系数矩阵R=(r_ij)_p×n。其中:

r_ij=1/n∑((x_ij-x^-_i)(x_ij-x^-_j))/((σ_iσ_j)).(2)

(3)计算R的特征值和特征向量。根据特征方程│R-λI│=0,计算特征根λ_i,并使其从大到小排列:λ₁≥λ₂≥…≥λ_p,同时可得对应的特征向量u₁,u₂,…,u_p,它们标准正交。u₁,u₂,…,u_p称为主轴。这里,I为单位矩阵。

(4)计算贡献率e_i 和累计贡献率

e_i=λ_i〖JB<1*/〗∑^pλ_k,(3)

E_m=∑^mλ_k〖JB<1*/〗∑^pλ_k.(4)

(5)计算主成分

z_m=∑^pu_mjx_j.(5)

(6)综合分析。一个m维主超平面究竟以多大的精度来代替原始变量系统,才能尽可能多的保留原始数据信息?这可以通过求累计贡献率E_m来判断。一般取E_m>85%的最小变量m(m<p),就可确定主超平面的维数m,从而可对m个主成分进行综合分析。

(7)根据主成分z_i分析得到的主成分和相应的权值(贡献率)e_i,计算本文定义的反映地震活动时间、空间和强度异常特征的综合指标:

W=^me_iz_i.(6)

因为m个主成分已基本保留了这些预报参数的信息,所以综合指标W包含了这些参数从不同侧面反映地震活动时间、空间和强度异常的基本特征。本文选择了从不同侧面反映地震活动时间、空间和强度特征的6个参量:地震频次N(M_L≥3.0)、b值、η值、M_f值、A(b)值、A_C值进行主成分分析。

为说明主成分分析方法如何用于实际地震预报中,下面以1975年2月4日海城M_S7.3地震为例。图1为7.3级地震发生前后震中周围300 km范围内的一些地震活动性参数随时间的变化曲线。除地震频次N选取的起始震级为M_L3.0外,其它参数计算所选取的起始震级均为M_L2.0。参数计算累计时间为1年,滑动步长为1个月。从图1可以看到,虽然地震前这些参数变化各有所异,但是大多参数变化并不显著。

表1为通过主成分分析方法得到的上述参量在各主成分中的系数(特征向量)、特征值与贡献率。当取3个主成分时,累计贡献率达到93%; 当取4个主成分时,累计贡献率达到98%。这表明取前3个主成分已包含了样本中的绝大部分信息量。由表1可知,主成分1中的A(b)、M_f值系数较大,是构成主成分1的主要参数; 主成分2中的地震频次b值、η值系数较大,是构成主成分2的主要参数; 主成分3中的A_C值、N值系数较大,是构成主成分3的主要参数。

本文取前3个主成分根据式(6)计算地震综合指标W。图2为海城7.3级地震前后震中附近地区地震活动的综合指标W随时间的变化。可以看到,在1975年海城7.3级地震前6个月时间,震中附近地区W值出现明显抬升的形态,异常幅度大于1.0。这表明地震活动的综合指标W可以很好地反映地震活动增强特征(王炜等,1994,1999)。

图1 海城地震参数随时间变化曲线
Fig.1 Curve of the seismic activity parameter with time

表1 各参量在各主成分中的系数(特征向量)、特征与贡献率
Tab.1 Parameters in primary component of the coefficient(eigenvector)、features and contribution rate

图2 地震强度综合指标W值随时间的变化
Fig.2 Curve of earthquake intensity synthetic index W with time

本文选取1970年以来东北地区(38°～55°N,117°～135°E)小震目录,1970年～2008年东北地区共发生M_S≥5.0浅源地震14次(剔除了海城、岫岩地震强余震,不含唐山老震区地震),分别取上述6个参数,进行主成分分析,进而研究这些中强震前W值的异常变化。考虑到地震孕育过程中可能的影响范围及东北地地震活动实际情况,对5～6级地震一般选取200～250 km左右的圆域; 对7级以上地震选取300 km左右的圆域。除地震频次N选取的起始震级为M_L3.0外,其它参数的计算选取起始震级M_L2.0,并剔除余震,资料计算时间窗为12月,滑动步长为1个月。

采用主成分分析方法,取累计贡献率>85%的最小的m(m<p)确定主成分个数m。在此基础上,对m个主成分使用(b)式进行综合分析,得到反映地震活动时间、空间和强度异常特征的综合指标W,并计算其在中强以上地震前后随时间的变化曲线。图3为表2中各次地震前后的变化曲线,由图3可知,综合指标表现为W值大体在0以下波动,在中强震发生前0.5～3 a的中期阶段,综合指标W都出现明显的增高异常,异常高值通常超过1.0。有的震例是在高值发震,如海城、岫岩地震; 有的震例是高值恢复后发震,如前郭、鄂伦春地震。这表明,综合指标W可以较好地反映地震活动异常的综合特征,它在震中附近地区震前的异常变化明显。

表2 1970年以来东北MS≥5.0地震
Tab.2 MS≥5.0 earthquakes in the Northeast China since 1970

图3 东北14次MS≥5.0浅源地震前W值异常
Fig.3 W value anomaly before 14 shallow earthquakes(MS≥5.0)in the Northeast China

图4为对与地震活动强度有关的6个参量进行因子分析时,得到的4个公共因子得分随时间变化曲线,可以看到图2中综合指标W值的基本形态主要由主因子1确定,这是由于主因子1的贡献率最大,达到58%。而因子分析的结果显示,因子1主要反映了异常相对较好的M_f值、A(b)值2个参数的信息。其它2个因子得分各有所异,且异常变化特征不甚显著,但由于它们的贡献率较小,因此综合指标W值可以较好地反映出震前的异常变化。由此可以看出,在因子分析中可以通过研究各参量的相关矩阵内部结构,得到各参量的因子载荷和因子贡献率,将一些包含异常信息较大的参量和因子赋予较大的权值,从而使得综合指标W值的异常变化明显。

图4 主成分1(a)、2(b)、3(c)随时间的变化
Fig.4 Curve of the primary component 1,2 and 3 with time

目前在地震预测中,反映地震活动性的参数较多,这些参数从不同侧面反映了地震在时间、空间和强度方面的特征(李献智,侯建盛,2000)。由于各个指标之间有一定的相关性,因而所得的统计数据反映的信息在一定程度上重叠,而且各个参数的形态变化各不相同,这不但给分析预测工作带来不必要的麻烦,还可能使预测工作存在主观随意性。主成分分析法可在保证原始数据信息丢失最少情况下,对高维变量空间进行降维处理,经过线性组合,以少数的综合变量取代原有的多维变量。本文所使用综合反映地震活动时间、空间和强度特征的W值,包含了上述6个参数中的绝大部分信息,可以较好地反映地震活动的时间、空间和强度异常变化特征。通过主成分分析方法得到的反映地震异常特征的综合指标W在东北14次5级以上中强浅源地震前出现一定的异常变化,表明综合指标W可以较好地反映地震活动异常特征,在地震预报中具有一定的应用前景。