混沌参数优化RBF算法的震前ENPEMF信号强度趋势预测

2020-12-18郝国成谭淞元曾佐勋

东北大学学报（自然科学版） 2020年12期

郝国成，锅娟，谭淞元，曾佐勋

(1.中国地质大学(武汉) 机械与电子信息学院，湖北武汉 430074； 2.中国科学院测量与地球物理研究所大地测量与地球动力学国家重点实验室，湖北武汉 430077； 3.中国地质大学(武汉)复杂系统先进控制与智能自动化湖北省重点实验室，湖北武汉 430074；4.中国地质大学(武汉) 智能地学信息处理湖北省重点实验室，湖北武汉 430074； 5.中国地质大学(武汉) 地球科学学院，湖北武汉 430074)

地球天然脉冲电磁场(the Earth’s natural pulse electromagnetic field, ENPEMF)是指可在地球表面接收的由天然场源产生的综合电磁总场[1-2].地震、滑坡等地质灾害现象可在地表产生甚低频(very low frequency, VLF)信号脉冲波动，“微破裂机-电转换”机制和“地壳波导”是上述电磁现象的机理之一[3].ENPEMF信号中的电磁异常信息具有潜在变化趋势及典型的非平稳特征，可反映地质活动的孕育发展趋势，可用于震前电磁异常监测分析[4].

随着信息技术、人工智能及机器学习理论的不断发展，电磁信息预测模型逐渐成为前沿热点.由于影响ENPEMF信号的场源较多，具有非平稳信号的特点，很难根据脉冲波形归纳其规律.非线性模型中的神经网络方法是基于经验风险最大化原则的机器学习算法，非线性拟合能力较强，可以对采集到的震前ENPEMF信号进行建模，并拟合其强度趋势的变化特点.

径向基函数(radial basis function, RBF)神经网络可逼近任意的非线性函数，其优点是学习速度快、非线性逼近能力强，具有良好泛化能力，可为ENPEMF信号构建非线性预测模型.本文提出基于混沌参数优化RBF算法的预测模型对其强度趋势进行预测，为数据分析和灾害监测提供支持.首先采用混沌理论对实测ENPEMF数据进行分析，其中假邻近法(false nearest neighbor, FNN)及自相关函数法分别求得嵌入维数和延迟时间等混沌特征参数，并优化RBF神经网络；采用训练完成的混沌参数优化RBF神经网络模型对ENPEMF数据进行预测，并与传统RBF神经网络进行比较.结果表明，混沌参数RBF算法可基于混沌理论确定动态系统的混沌特性，从采集到的震前ENPEMF信号强度数据中找到其变化态势，并预测14 d(4月7～20日)的ENPEMF数据强度趋势，且预测效果及精度均优于传统RBF预测模型.

1 混沌理论

由于震前ENPEMF信号的产生机理复杂、孕育过程非线性，其信号强度数据具有非平稳特点和混沌特性.因此，本文引入混沌理论对其数据内部特征进行挖掘，找到隐藏的混沌特点.假设一段时间内采集的震前ENPEMF信号数据为{x(tj),j=1,2,…,n}，其中n表示采集的数据点数，通过混沌理论中的假邻近法及自相关函数法对数据进行处理，得到数据变化形式：

X(t)=[x(t),x(t+τ),…,x(t+(m-1)τ)] .

(1)

式中，τ和m分别表示震前ENPEMF信号数据的延迟时间和嵌入维数，用于描述该信号隐藏的混沌特征并为RBF神经网络输入节点提供判断依据.

1.1 假邻近法

时间序列的本质是将系统高维空间坐标的运动轨迹投影到低维空间.当嵌入维数较小时，系统空间轨道中本来相距很远的相点相互挤压折叠，未能充分展开，这些点为假邻近点[5].随着嵌入空间维数的增加，轨道逐渐展开，投影到低维空间的假邻近点随之分离.当所有的假邻近点消失时所对应的最小嵌入空间维数即为最佳嵌入维数.给定正整数m，可构造m维重构向量：

ym(n)=(x(n),x(n+τ),…,
x(n+(m-1)τ))T.

(2)

(3)

将维数从1维增加到m+1维.m+1维空间重构向量如式(4)所示：

ym+1(n)=(x(n),x(n+τ),…,
x(n+(m-1)τ),x(n+mτ))T.

(4)

(5)

(6)

(7)

(8)

对所有的重构向量，利用判据找出邻近点中的假邻近点，并记录下所有假邻近点的数目FN(m).继续增加维数，当找到一个整数mε使得FN(mε)=0时，mε即为所求嵌入维数.当假邻近点所占比率即假邻近率随着嵌入维数的增加趋于平稳不再降低时，所对应的嵌入维数m为最佳嵌入维数.本文通过统计假邻近点数的比率随嵌入维数升高逐渐减小，最后维持不变的情况，确定最优嵌入维数.

1.2 自相关函数法

由于实际时间序列长度有限且存在噪声，选取合适的延迟时间至关重要.延迟时间τ过小，将使重构的系统由于相关性较强造成相空间的挤压，不能充分展示系统的动力特征；延迟时间τ太大，会造成相邻两时刻的动力学形态剧烈变化，使构造的相空间比实际空间复杂.

自相关函数法可在降低相关性的同时保证原动力学的系统信息不丢失，使重构相空间能充分展现系统拓扑性质和几何性质.首先写出时间序列的自相关函数，然后做出自相关函数随时间变化的函数图，找到自相关函数首次达到零点时对应的时间即为延迟时间τ.自相关函数定义为

(9)

自相关函数值随时间变化逐步下降，当其下降到初始值的(1-1/e)时对应的时间为所求延迟时间.

2 RBF神经网络算法

径向基函数(RBF)神经网络是一种包含输入层、隐含层和输出层的三层神经网络[6].其中输入层与隐含层之间为非线性变换，从隐含层到输出层为线性变换.

在RBF神经网络中，输入层仅作为通道传输信号.隐含层中神经元的变换函数为径向基函数，通过非线性变换可将信号从输入层传递到隐含层.输出层是对输入信号的响应.RBF神经网络结构可根据具体问题在训练阶段进行自适应调整.

在网络中，X=(x1,x2,…,xn)T为输入样本，Y=(y1,y2,…,yn)T为输出响应.RBF神经网络算法需要求解基函数的中心Ci，基函数的宽度Di以及隐含层到输出层的权值Wi三个参数.

RBF神经网络的训练过程分为两步：首先进行无监督学习，计算输入层与隐含层之间的和，得到隐含层输出：

(10)

其中：i=1,2,…,N;‖X-Ci‖为欧氏范数.

然后采用最小二乘法求隐含层与输出层之间的权值ωi.最终得到RBF神经网络的输出y：

(11)

将经过参数训练的RBF神经网络用于预测混沌时间序列.可利用混沌时间序列数据校正上述权值参数，提高神经网络的非线性泛化能力.其中，当输入层节点个数为混沌时间序列的嵌入维数m时，预测结果较好；隐含层节点数目根据实验实时调整确定：预先设定RBF神经网络的精度值，隐含层节点个数递增，当神经网络达到预设精度时，该节点个数即为神经网络的隐含层节点数.

本文提出的基于混沌参数优化RBF神经网络算法的预测模型工作流程如图1所示.

图1 本文算法预测模型工作流程

首先对实测ENPEMF数据进行混沌分析，然后分别采用假邻近法和自相关函数法求得最优嵌入维数和延迟时间等混沌特征参数，将得到的参数作为确定输入层节点个数的依据并优化RBF神经网络；对参数优化的RBF神经网络进行训练，学习其内部混沌特征，最后用训练完成的混沌参数优化RBF算法预测ENPEMF信号的强度趋势.

3 实验仿真

隐藏在混沌时间序列内部的某些特征信息可以通过延迟时间及嵌入维数表现出来.本文用经典的Rossler混沌时间序列验证假邻近法及自相关函数法获得嵌入维数和延迟时间等参数的可行性.

Rossler系统可用式(12)微分方程组进行描述：

(12)

选取参数a=b=0.2，c=5，初值x(0)=-1，y(0)=0，z(0)=1，积分时间步长h=0.05，生成长度为3 000的连续混沌时间序列如图2所示.

图2 Rossler混沌时间序列

Rossler混沌时间序列在有限区域内运动时趋向于一个稳定的点，完全展开系统内部的混沌特性.嵌入维数的取值范围为[1,8]，阈值的判别门限范围为[2,15]，Rossler时间序列长度为3 000，采用假邻近法计算时间序列x分量的嵌入维数m，结果如图3所示.

当嵌入维数从1增加到4时假邻近率急速下降；当嵌入维数达到5时，假邻近率趋于平缓，此时的嵌入维数达到理想值，为所求最佳嵌入维数，即Rossler混沌时间序列的嵌入维数为5.

采用自相关函数法求延迟时间τ，选择x分量序列计算延迟时间τ，仿真结果如图4所示.

图4中，直线为初始值的(1-1/e)，曲线为自相关函数曲线，选取值为自相关函数下降到初始值的(1-1/e)时所对应的时间，即Rossler混沌时间序列的延迟时间τ为16 s.

图3 假邻近法求嵌入维数

图4 自相关函数法求延迟时间

取迭代后长度为2 000的时间序列，其中前500个数据进行训练，后1 500个数据进行预测.以嵌入维数m=5作为RBF神经网络的输入节点判断依据，采用训练完成的混沌参数优化RBF算法对数据进行预测.图5为预测模型的仿真结果.“×”表示真实值，“-”表示预测值.图6为预测结果的绝对误差.

图5 Rossler混沌时间序列的预测结果

图6 预测结果的绝对误差

图5和图6中，混沌参数优化RBF算法可以较好拟合经典的Rossler混沌时间序列，且预测误差较小.

4 混沌参数优化RBF算法的预测研究

本文使用俄罗斯科学院托木斯克分院GR-01型设备接收ENPEMF信号.在武汉九峰地震台放置了3台设备，接收方向为W-E和N-S，3个通道为CN1，CN2，CN3.设备记录了ENPEMF信号的AH数据(超过设定阈值的脉冲幅度)及NH数据(超过设定阈值的脉冲个数)，可表征地表天然磁场的强弱[7].设备的工作频率为甚低频段：5～25 kHz，在武汉设置的接收频率为14.5 kHz.

ENPEMF信号为非周期、非平稳信号，具有明显混沌特性，设备输出为数字量化后的信号，数据存储格式为：时间-幅度-脉冲数(t-AH-NH)[8].幅度单位仅为信号包络大小变化的参照量，已不具有原来量纲的直接意义.

对采集的ENPEMF数据进行混沌特性分析，找到其内部隐藏的混沌特征及趋势变化特点，结合RBF神经网络算法对其信号强度进行预测，识别孕震信息[9].本文采用Grassberger和Procaccia提出的关联维数算法(G-P算法)求解ENPEMF信号关联维数，可判断信号是否具有混沌特性[10].

对于时间序列x(1),x(2),…,x(t)，其长度为M，对其进行相空间重构，得到向量：X(t)=[x(t),x(t+τ),…,x(t+(m-1)τ)],其中t=1,2,…,N,N=M-(m-1)τ.给定正数ε足够小，当空间向量间的距离小于ε时，向量关联.关联向量的关联积分表达式为

(13)

其中θ(·)为Heaviside阶跃函数，满足式(14)：

(14)

当时间序列的长度N→∞，半径ε→0时，关联积分与半径的关系为

(15)

其中，D为所求关联维数，变形后得到式(16)：

(16)

给定一系列半径ε和嵌入维数m，作半径随嵌入维数变化的关联积分图组，用最小二乘法对图中lbC(N,ε)～lbε最接近直线的一段拟合最佳直线，该直线斜率即所求关联维数D.

本文数据于2013年4月20日在中国四川省芦山地震期间收集，地震的位置和地震台的位置如图7所示[11].

图7 芦山地震的位置和地震台的位置

图8为4月10～20日通道2的AH数据，箭头指向中国四川省7.0级芦山地震发生时间.

图8 4月10～20日通道2的AH数据

地震发生前的11 d内，ENPEMF数据从14日到15日有较大的峰值变化，在16日回落至正常水平.在17日和18日观测到显著的峰值变化，在19日信号又跌落至正常.因此，在地震发生前ENPEMF信号脉冲强度会发生剧烈变化.

本文选择经过平滑及归一化的20 d ENPEMF数据作为实验数据(4月1～20日).其中前6 d数据(4月1～6日)为模型的训练样本，后14 d数据(4月7～20日)作为模型的预测样本.利用G-P算法计算4月1～6日ENPEMF信号的关联维数.信号的关联积分组在一定范围内呈近似直线分布；随着嵌入维数的增加，直线斜率增大，且最后关联维数趋于稳定，说明ENPEMF信号具有混沌特性.

对选取的ENPEMF信号中前6 d数据(4月1～6日)进行预处理，采用假邻近法计算嵌入维数m，自相关函数法计算延迟时间τ，得到图9和10，从中了解震前ENPEMF信号数据的混沌特性.

图9 假邻近法求嵌入维数

图10 自相关函数法计算延迟时间

图9中，随着嵌入维数的增加，当嵌入维数为4时，假邻近率趋于平稳，此时的嵌入维数即ENPEMF信号数据的最佳嵌入维数.图10中，信号自相关函数达到初始值的(1-1/e)时τ=5 s，即ENPEMF信号数据的延迟时间为5 s.用得到的参数确定RBF神经网络的输入节点个数为4，进而将训练完成的混沌参数优化RBF算法用于数据预测.

选择4月1～6日数据作为训练样本，训练混沌参数优化RBF神经网络，其中输入层有4个节点，输出层1个节点，隐含层有6个节点.选择径向基高斯函数作为隐含层神经元传递函数，输出为线性函数.最后，利用训练完成的混沌参数优化RBF神经网络预测模型和传统的RBF神经网络预测模型分别实现对4月7～20日ENPEMF数据的单步预测.图11为所提混沌参数优化RBF神经网络预测模型的结果，图12为传统RBF神经网络预测模型结果.

图11 混沌参数优化RBF预测结果

图12 传统RBF预测结果

图11和12中，分别采用混沌参数优化RBF预测模型与传统RBF预测模型对4月7～20日ENPEMF信号强度数据进行预测，两种模型均可以模拟出采集到的地震前14 d(7～20日)实际ENPEMF信号强度的波动.对于整体数据范围，传统 RBF 神经网络模型不能较好地跟踪实际值的变化，而本文所提预测模型对ENPEMF信号具有较好的跟踪拟合性能；对于17日的数据剧烈波动时刻，混沌参数优化RBF预测模型相较于传统的RBF预测模型拟合效果更好，具有较好的预测结果，误差较小，预测优势明显.为更精确评估所提预测模型的预测效果，选取绝对误差作为ENPEMF数据预测精度评价指标，结果如图13所示.