高斯骨架差分进化算法与回声状态网络的结合应用

谢霖铨，曾孟麒，杨火根

(江西理工大学理学院，江西 赣州 341000)

时间序列普遍存在于人类社会的方方面面,比如在能源领域上能源每年的消耗量以及排放量等[1]；天气领域中城市的风速、降水量、平均气温等[2-3];工业领域中的铁块的凝结透点[4]、电力负荷序列[5]和原油物性预测[6]等。

随着人工智能研究的不断发展,人工神经网络、支持向量机[7]等机器学习方法逐渐成为非线性时间序列建模和预测[8]的主要工具。但是传统的人工神经网络如BP(Back Propagation)[9]存在训练复杂、网络结构难以确定、易陷入局部最优和收敛速度慢等问题。支持向量机(Support Vector Machine，SVM)若训练模型样本过大,存在训练时间长和训练效果不佳等问题。随着学者们的不断研究,发现回声状态网络[10](Echo State Network,ESN)在时间序列预测方面存在一定的优势,其训练算法简单,计算快速,且能保证解的全局最优性。虽然ESN具有以上优点,但是也存在一些问题,例如储备池的适应性问题、共线性问题、稳定性问题以及病态解的出现会影响网络的泛化能力,容易产生过拟合等。

近年来,学者们针对ESN存在的一些问题进行一些改进。文献[11]利用遗传算法(Genetic Algorithm,GA)直接优化回声状态网络,实现自适应非线性动力系统的控制。文献[12]采用人工蜂群算法(Artificial Bee Colony Algorithm,ABC)优化ESN,满足了聚氯乙烯聚合过程的实时控制要求。文献[13]利用粒子群(Particle Swarm Optimization,PSO)优化算法对ESN中未经训练的权重进行优化,提高网络的预测性能。文献[14]利用狼群优化算法(Grey Wolf Optimizer,GWO)直接优化ESN参数,避免了网络参数的设立,提高网络的精确性。

虽然已经有学者采用智能优化算法优化回声状态网络参数,但由于回声状态网络的储备池包含节点较多、神经元节点之间的连接状态较复杂,优化效果还有待提高。人工蜂群算法拥有较好的鲁棒性,但由于不能有效地平衡勘探和开采导致算法容易陷入局部极值。狼群优化算法因其参数少、易于实现、收敛速度快等受到广大学者喜爱,但其种群多样性差和后期算法收敛速度慢易陷入局部最优等,阻碍了算法寻优性能。粒子群算法由于参数的初始化会对算法性能产生较大影响,往往会使算法陷入局部最优。研究中发现差分进化算法[15](Differential Evolution,DE)结构简单、易于实现,受控参数少以及鲁棒性强,且具有较强的全局搜索能力,寻优速度快,更具实用性。该算法自提出以来,已被广泛用于约束优化、聚类优化、非线性优化、神经网络优化以及其他领域。综上,本文选用简单高效的差分进化算法。

为使差分进化算法可以更加高效地找到最优解以及适用性更强。本文提出一种改进的高斯骨架差分进化算法。结合三角高斯变异策略和标准差分进化算法中的DE/best/1变异策略,提出一种新的策略选择因子,同时交叉概率采用自适应生成。为了避免算法陷入局部最优,设立了算法停滞扰动策略来提升算法寻优性能。综上所述,本文采用改进的高斯骨架差分进化算法来优化回声状态网络中储备池的五个重要参数,通过策略选择机制来平衡差分进化算法的勘探和开采能力,从而加快该算法的收敛速度,更加准确地找到最优解,提高预测精度。

1 回声状态网络

回声状态网[10](Echo State Network,ESN)是Jaeger于2004年提出来的一种简化版的递归神经网络,由于其用稀疏连接的储备池代替全连接的隐含层,增强了对动态系统的建模能力,避免一般神经网络基于梯度下降法的收敛速度慢、易陷入局部最优问题。

1.1 ESN基本原理

基本ESN的拓扑结构如图1所示,回声状态网络是一种三层递归神经网络,由输入层、隐含层和输出层三部分组成,隐含层又称储备池,含有成百上千个稀疏递归连接的神经元,神经元之间的连接权值随机生成且固定不变。

图1 ESN拓扑结构示意图Figure 1 Schematic diagram of the ESN topology

ESN的状态方程和输出方程如下：

x(t)=φ{Winu(t)+Wxx(t-1)+bx}

(1)

y(t)=wTx(t)+b

(2)

上式中u(t)∈RM×1是M元输入向量,y(t)∈RM×1是M元输出向量,bx∈RN×1为输入偏置,b∈RM×1为输出偏置。由当前时刻t的输入u(t)和前一时刻t-1的储备池内部状态x(t-1)通过激活函数的映射得到当前时刻的状态x(t)∈RN×1。φ(?)为神经元激活函数,可以选取Sigmod函数或者tanh函数。输入连接矩阵Win∈RN×N的元素在区间[-1,1]取值。Wx∈RN×N为储备池内部连接矩阵,稀疏连接。通过伪逆求得输出连接矩阵w∈RN×M。输入连接矩阵Win和内部连接矩阵Wx随机生成,且在回声状态网络的训练阶段始终保持不变。网络只需求输出连接矩阵w,因此降低了计算复杂度。

1.2 网络关键参数

回声状态网络的核心是储备池,储备池性能的好坏取决于五个重要的参数：初始化储备池大小Nini、储备池大小N、储备池的普半径R、稀疏度D以及输入缩放因子S,如何选取以上五个参数至关重要。下面介绍储备池的参数选择对模型性能的影响。

(1)初始化储备池大小

初始化储备池大小的Nini是储备池初始激活神经元个数,让储备池中形成较好的回声,使ESN具有较好的动力特性。

(2)储备池大小

储备池的大小N即为储备池中神经元数目,是影响ESN预测性能最重要的参数。N值过大,容易造成过拟合；N值过小,又容易造成欠拟合。

(3)储备池的普半径

普半径R为内部连接矩阵Wx。的最大特征值的绝对值。R的取值范围一般为[0,1]之间,但对于不同的时间序列其取值将视情况而定。

(4)稀疏度

储备池内部神经元连接的稀疏程度称为稀疏度D,储备池的神经元之间不是全连接的,而是少部分连接。具体实现方法是使储备池的连接权重Wx中的大多数元素等于零。文献[16]中提到稀疏度D∈[0.01,0.04]之间即可保证储备池具有足够的动力特性。

(5)输入变换因子

输入变换因子S是指在信号输入储备池前缩放的比例因子,表征输入连接权值的取值范围。根据式(1)可知,其大小决定激活函数的工作区间,也决定了输入对储备池的状态变量作用的强度。

回声状态网络参数的选择往往会针对数据的不同特性而变化,因此如何选择适合不同数据的储备池参数是本文研究的重点。本文采用改进的差分进化算法来优化回声状态网络储备池的参数,得到适合当前数据的储备池参数,用于提高网络的预测性能。

2 改进差分进化算法

差分进化算法[15](Differential Evolution,DE)是基于群体智能的优化算法。由于该算法简单易于实现、鲁棒性好、搜索能力强等优点被广泛应用于工程问题的求解。

2.1 差分进化算法

DE采用实数编码方式,其算法原理与遗传算法十分相似,进化流程与遗传算法相同：变异、交叉和选择。DE算法中的选择策略通常为锦标赛选择,而交叉操作方式与遗传算法也大体相同,但在变异操作方面使用差分策略,即利用种群中个体间的差分向量对个体进行扰动,实现个体变异。DE的变异方式,有效利用群体分布特性,提高算法的搜索能力,避免遗传算法中变异方式的不足。

(3)

其中,xi(0)表示种群中第0代的第i个个体,xj,i(0)表示第0代的第i个个体的第j维。NP表示种群大小,rand(0,1)表示在(0,1)区间均匀分布的随机数。

2)变异操作。DE通过变异策略来产生新的个体,有利于后代进行搜索以便找到最优的解。常用的变异策略主要有以下6种:

DE/rand/1

vi(g)=xr1(g)+F×(xr2(g)-xr3(g))

(4)

DE/rand/2

vi(g)=xr1(g)+F(xr2(g)-xr3(g))+F(xr4(g)-xr5(g))

(5)

DE/best/1

vi(g)=xbest(g)+F·(xr1(g)-xr2(g))

(6)

DE/best/2

vi(g)=xbest(g)+F·(xr1(g)-xr2(g))+F·(xr3(g)-xr4(g))

(7)

DE/current-to-best/1

vi(g)=xi(g)+F·(xbest(g)-xi(g))

+F·(xr1(g)-xr2(g))

(8)

高斯三角差分进化策略

μ=(xr1(g)+xr2(g)+xr3(g))/3

σ=(|xr1(g)-xr2(g)|+|xr1(g)-xr3(g)|

+|xr2(g)-xr3(g)|)/3

vi(g)=N(μ,σ)

(9)

其中,F为缩放因子,xbest(g)表示第g代种群中最优个体,xi(g)表示第g代种群中第i个个体,i≠r1≠r2≠r3≠r4≠r5,μ、σ为高斯变异策略个体的均值和标准差,N(?)为服从均值为μ标准差为σ的正态分布。在进化过程中,为了保证解的有效性,必需判断个体中每一维的值是否满足边界条件,如果不满足边界条件,则按照式(3)重新生成个体。这样就得到变异的中间个体vi(g)。

(3)交叉操作。对第g代种群{xi(g)}及其变异的中间体{vi(g)}进行个体间的交叉操作:

(10)

其中,CR为交叉概率,jrand为[1,2,…,D]的随机整数。生成的试验个体为ui,j(g),{i=1,2,…NPj=1,2,…,D}

(4)选择操作,就是对新旧个体进行的一个淘汰制操作,其规则为比较新旧个体的适应度,适应度好的个体保留到下一代继续参与进化,适应度差的个体就被淘汰,以此来获得适应度较好的个体,选择公式如下:

(11)

式中xi(g)为第g代的第i个个体,如果个体xi(g)的适应度比实验个体ui,j(g)适应度更优,则个体xi(g)进入下一代参与进化,反之实验个体ui,j(g)进入下一代参与进化。

反复执行2)、3)、4)步操作达到设定好的最大代数G{g=0,1,2,…,G}得到一个最优种群NP,种群中适应度最好的个体就是最优个体xbest。

2.2 MGBDE算法

Wang等[17]于2013年提出了一种基于高斯变异的新型骨架差分进化算法(GBDE)及其改进算法MGBDE。其中,MGBDE算法将GBDE算法中的高斯变异策略和标准DE算法中的DE/best/1策略进行组合,种群中每一个个体随机等概率地从中选择一种变异策略进行进化,如式(15)所示。

(12)

其中,Mi表示个体xi执行的变异策略。

2.3 IGBDE算法

上述提到将高斯变异策略和标准DE算法中的DE/best/1变异策略进行组合得到一种新的变异策略Mi,但是这种变异策略Mi在进化过程中每个个体都是等概率的执行DE/best/1变异策略和高斯变异策略,完全没用考虑到DE算法在进化初期趋向于勘探,在进化后期趋向开采,故本文提出一种策略选择机制ZZ,使DE算法在进化初期大概率的选择利于勘探的高斯变异策略,在进化后期大概率的选择利于开采的DE/best/1变异策略。策略选择机制ZZ公式如下:

ZZi(g)=1-

(13)

其中ω是阈值,本文取值为0.001.ZZi表示为第i个个体的策略选择概率。两个个体xi、xbest在搜索空间中离得很远的情况下ZZi较大,则说明算法应该利于勘探,故选择勘探能力强的高斯变异策略,反之,两个个体xi、xbest在搜索空间离得很近的情况下ZZi较小,则说明算法应该利于开采,故选择开采能力强的DE/best/1变异策略。在进化过程中ZZi随着个体xi、xbest一起更新演变。

虽然高斯变异策略勘探能力强,但利用了xbest个体来指导进化,没有利用差分进化的核心思想,即随机选择的个体之间的差异进行高斯变异,使种群多样性降低,故本文采用Peng等[18]于2015年提出一种基于三角骨架的DE算法,该算法提出了一种新的变异策略——高斯三角变异策略。也就是本文的式(9)。将式(12)的高斯变异策略换成高斯三角变异策略,提高算法的前期的勘探能力。结合策略选择机制ZZ,形成新的变异策略如下：

(14)

其中μ,σ由式(9)得来的。

分析发现,在算法寻优初期,个体适应值与全局最优值相差较大,则ZZi较大,种群中绝大部分个体执行三角高斯变异策略,对解空间进行探索;随着进化的推进;个体不同程度的向最优个体靠近,此时根据个体适应值与最优适应值的比较,为每一个个体指定当前最合适的变异策略;在进化后期,大部分个体收敛于最优个体附近,此时ZZi较小,种群中绝大部分个体转为执行DE/best/1变异策略,算法聚焦于在最优个体附近开发。IGBDE算法继承了骨架算法家族近似无参的优点,新引进的控制参数ZZ随进化过程自适应调节,无需人工干预和任何先验知识。

此外,交叉概率在一定程度上会影响算法性能。为提高交叉操作的适应性,本文采用文献[19]的方法控制交叉概率的设置,即：

(15)

最后借鉴文献[20]提出的停滞扰动策略来防止算法陷入局部极值。停滞扰动策略公式如下：

(16)

根据文献[20]给出τ取值为0.2时对单峰问题和多峰问题均表现出良好的寻优性能,故本文τ取值为0.2。

3 改进高斯骨架差分进化算法优化储备池参数

ESN优点在于网络学习过程中仅仅需要更新储备池到输出层之间的连接权值,而输入连接权值和储备池中神经元内部的连接权值一般都是在随机赋值后保持不变。故储备池参数设置直接影响ESN的预测性能,人工调节参数既费时又不能选择出最佳的参数值,所以本文提出基于改进的高斯骨架差分进化算法优化回声状态网络(Improve Gaussian Bare-Bones Differential Evolution-Echo State Network,IGBDE-ESN)的参数。本文将实验数据分为训练集、测试集两部分,目标函数为训练集误差最小。算法伪代码如算法1所示,算法流程图如图2所示。

图2 IGBDE-ESN算法流程图Figure 2 IGBDE-ESN algorithm flowchart

算法1改进的高斯骨架差分进化算法优化回声状态网络

步骤1设置种群代数g=0,初始化交叉概率CR=0.9,最大迭代次数为G,设置储备池参数的取值范围。

步骤2设置个体参数,即每一维对应储备池的一个参数。初始化并评价种群NP,保存最优个体。

步骤3根据式(13)计算初代ZZi

步骤4迭代,判断g≤G时,设置监视器count=0,根据式(12)生成变异个体vi(g),对变异个体越界处理。

步骤5根据式(10)生成实验个体ui(g),评价ui(g)。

步骤6根据式(11)进行选择操作,如果生成新的个体,则count=count+1,反之count保持不变。根据式(15)更新CRi,根据式(13)更新ZZi。

步骤8输出最优个体,在测试集上验证。

4 仿真实验

为验证本文所提模型的有效性,本文选择Mackey-Glass时间序列、赣州月平均气温—降雨量数据集进行仿真实验,并提出2组分别适合大样本和小样本的ESN储备池参数取值范围。同时,在相同数据集上,与人工蜂群[11]优化ESN(ABC-ESN)、狼群算法[13]优化ESN(GWO-ESN)、粒子群[12]优化ESN(PSO-ESN)和极限学习机(ELM)预测模型的仿真结果进行比较。模型预测性能好坏的评价指标为均方根误差(Root Mean Square Error,RMSE)、平均对称绝对误差(Symmetric Mean Absolute Percentage Error,SMAPE)和标准化均方根误差(Normalized Root Mean Squared Error,NRMSE)。RMSE、SMAPE、NRMSE的计算公式定义如下:

(17)

(18)

(19)

4.1 Mackey-Glass时间序列

本文采用龙格库塔[21](Runge-Kutte)方法来生成Mackey-Glass混沌时间序列数据。Mackey-Glass(MG)时滞微分方程如下:

(20)

当θ>17时,整个序列是混沌的,无周期的,且不收敛不发散,本文取θ=17时产生的10 000组MG时间序列。

由于MG时间序列数据较大适合大样本的ESN储备池参数取值范围。改进高斯骨架差分进化算法优化ESN模型参数的设定:初始化储备池大小Nini范围设为[10,150]、储备池规模N范围设为[500,1500]、谱半径R范围设为[0.5,1.5]、稀疏度D范围设为[0.01,0.1]及输入变换因子S范围设为[0.1,1];种群大小NP设为25,最大迭代次数设为100,阈值ω设为0.001,τ值设为0.2.前9000条MG数据用于训练,后1 000条MG数据用于测试。人工蜂群优化ESN模型、灰狼优化ESN模型以及粒子群优化ESN模型的参数设置同改进的高斯骨架差分优化ESN模型一致。表1给出了对于MG时间序列IGBDE选出的最优参数值。

表1 Mackey-Glass时间序列:IGBDE-ESN模型参数Tab.1

表2给出不同模型对MG时间序列的预测结果,可以看出本文所提出的模型在RMSE、SMAPE、NRMSE方面均优于其他对比模型,且更具优势。其中狼群、人工蜂群优化ESN模型的预测精度不及粒子群优化ESN模型,由此可以看出在同时优化多个不同范围的参数时,狼群、人工蜂群优化ESN存在一定的局限性。

表2 Mackey-Glass时间序列:测试集仿真结果Table 2 Mackey-Glass time series:test set simulation results

图3给出了IGBDE-ESN对MG测试数据的预测曲线和误差曲线。从图中可以看出IGBDE-ESN能够很好的拟合Mackey-Glass时间序列,且误差值越来越趋近于零。

图3 IGBDE-ESN对MG测试数据的预测曲线和误差曲线Figure 3 Prediction curve and error curve of IGBDE-ESN on MG test data

图4给出对于MG时间序列,本文所提出的改进的高斯骨架差分进化算法优化ESN(IGBDE-ESN)、人工蜂群优化ESN(ABC-ESN)、灰狼优化ESN(GWO-ESN)以及粒子群优化ESN(PSO-ESN)在迭代过程中的适应度值(Fitness)的变化曲线图.为了更清楚地显示各个模型适应度曲线的差别,对适应度的值取以10为底的对数.从图中可以看出IGBDE-ESN模型在迭代过程中,误差越来越小,最终误差值在其他同类模型中最低。

图4 Mackey-Glass时间序列:不同模型的适应的曲线Figure 4 Mackey-Glass time series：Adaptation curves for different models

4.2 赣州月平均气温序列

在中国气象网上获取赣州市从1971年1月到2021年10月的月平均气温数据记录值,作为仿真实验数据,采样间隔为月.共计609组数据,其中前450条数据作为训练集,后159条数据作为测试集.该数据集采用小样本的ESN储备池参数取值范围.改进高斯骨架差分进化算法优化ESN模型参数的设定:初始化储备池大小Nini范围设为[10,50]储备池规模N范围设为[20,150]、谱半径R范围设为[0.1,1]、稀疏度D范围设为[0.01,0.04]及输入变换因子S范围设为[0.0001,0.1];种群大小NP设为25,最大迭代次数设为100,阈值ω设为0.001,τ值设为0.2.为证明本文所提模型的有效性,对于ABC-ESN模型、GWO-ESN模型以及PSO-ESN模型的参数设置同IGBDE-ESN模型一致.表3给出本文所提模型针对赣州月平均气温序列选出的最优储备池参数.

表3 赣州月平均气温：IGBDE-ESN模型参数Table 3 Monthly average temperature in Ganzhou:IGBDE-ESN model parameters

表4给出了不同模型对赣州月平均气温数据的预测结果,可以看出本文所提模型在RMSE、SMAPE、NRMSE方面均优于其他对比模型。其中狼群、人工蜂群优化ESN模型的预测精度不及粒子群优化ESN模型但也相差不大,由此可见粒子群、狼群、人工蜂群对赣州月平均气温的预测大体上差不多,但也不及本文所提模型。

表4 赣州月平均气温:测试集仿真结果Table 4 Monthly average temperature in Ganzhou:Test set simulation results

图5给出了IGBDE-ESN模型对赣州月平均气温序列测试数据的预测曲线和误差曲线。从图中可以看出,本文模型能够很好地拟合赣州月平均气温序列曲线,绝对误差较小。

图5 IGBDE-ESN对赣州月平均气温测试数据的预测曲线和误差曲线Figure 5 Prediction curve and error curve of IGBDE-ESN on monthly average temperature test data in Ganzhou

图6给出对赣州月平均气温数据集,IGBDE-ESN、ABC-ESN、GWO-ESN及PSO-ESN在迭代过程中的适应度值的变化曲线图,为了更清楚地显示各个模型适应度曲线的差别,对适应度的值取以10为底的对数,从图中可以看出IGBDE-ESN算法收敛速度快,且误差较小。

图6 赣州月平均气温序列:不同模型的适应的曲线Figure 6 Sequence of monthly average temperatures in Ganzhou：Adaptation curves of different models

5 结论

本文通过选用不同样本容量的时间序列,利用改进的高斯骨架差分进化算法来动态的选择回声状态网络参数,并找到了两组适合不同样本容量的储备池参数范围,以适应于不同的时间序列的动力学特性,从而提高ESN的预测精度和泛化性能。对差分进化算法改进主要有2点,一是提出了一种新的策略选择机制,使个体能够随着进化动态的选择合适的变异策略,以提高DE算法性能。二是将高斯三角变异策略替换掉策略选择机制中的高斯变异策略,进一步增强DE算法在进化初期的勘探能力。对两组时间序列进行预测分析,仿真实验结果表明,本文所提模型较其他预测模型更好,既具有较高的预测精度,又具有较快的收敛速度,在时间序列预测分析中具有实用性和有效性、普适性。