多策略融合学习萤火虫算法在年径流预测中的应用

谢智峰，吴润秀，吕 莉

(1.江西科技学院 信息工程学院，江西 南昌 330098；2.南昌工程学院 信息工程学院，江西 南昌 330099)

径流预测是根据径流形成规律，借助已知信息对未来某一段时间内的径流量做出定性或定量的预测，其在防汛抗旱、水资源配置和水库调度等方面发挥着重要作用[1]。径流形成受太阳活动、大气环流和人类活动等诸多因素影响，因此径流预测属于受多因素影响的非线性复杂问题。目前，径流预测的方法主要有小波分析(Wavelet Analysis，WA)[2]、时间序列分析(Time Series Analysis，TSA)[3]、人工神经网络(Artificial Neural Networks，ANN)[4]和支持向量回归(Support Vector for Regression，SVR)[5-7]等。其中，支持向量回归因在解决非线性复杂问题方面表现突出，已成为水文预测领域的常用方法之一。

群智能算法(Swarm Intelligence Algorithm，SIA)是一种新兴的仿生进化算法，具有模型简单、实现容易、鲁棒性强等特点。典型的群智能算法包括粒子群算法(Particle Swarm Optimization，PSO)[8-9]、人工蜂群算法(Artificial Bee Colony，ABC)[10-11]和萤火虫算法(Firefly Algorithm，FA)[12-14]等。由于这类算法不受优化问题的连续或可微的限制，已广泛应用于无线传感器优化[15-17]、图像压缩[18-19]和径流预测[20-22]等领域。

为获取最佳的SVR核参数组合，很多学者将群智能算法应用于SVR的参数寻优，不仅提高了参数寻优效率，而且获得了较好的泛化性能。王宏伟[23]等将遗传算法(Genetic Algorithm，GA)嵌入SVR，构建了GA-SVR预报模型，取得了较好的预报效果；畅明琦[24]等设计了基于不同预报因子的两组SVR模型，实验结果表明，多预报因子的SVR模型预测效果更佳；蒋林利[25]等提出基于粒子群和模拟退火的混合算法，实现SVR参数优化选取，实验验证了混合算法模型具备较好的泛化性能，预测精度高且稳定。

Yang[26]于2008年提出萤火虫算法。该算法思想来源于自然界中发光萤火虫向比自身更亮的萤火虫飞行的生物习性。萤火虫离散分布于一定的活动范围内，每只萤火虫因体内不同量的荧光素而发出强弱不一的光，随着时间推移，该范围内的萤火虫会逐渐聚集于最亮萤火虫附近，形成若干个相似亮度的聚集中心，类似于寻优过程。

为兼顾种群全局探索和局部开采的平衡，提高算法求解精度，提出多策略融合学习萤火虫算法(Multi-Strategy Fusion Learning Firefly Algorithm，MSFLFA)。首先，算法为最优萤火虫随机选择两只不同萤火虫作为学习对象，制定新更新方式，进行定量次数的进化，提高算法探索全局最优解的能力；然后，在非最优萤火虫的既有更新策略上添加新的扰动项，动态识别搜索域，增强算法的局部寻优能力；最后，不同萤火虫交替迭代进化，完成寻优。为验证MSFLFA性能，设计基准测试函数的仿真实验，结果表明，MSFLFA寻优精度高且综合性能优于文中涉及的比较算法。

1 支持向量回归和萤火虫算法

1.1 支持向量回归

支持向量回归常用的核函数[27]主要包括线性核函数、sigmoid 核函数、多项式核函数和径向基核函数。线性核函数适合线性可分问题；sigmoid 核函数类似于多层神经网络；多项式核函数和径向基核函数均可将低维非线性问题转化为高维线性问题。对于受多因子影响的非线性径流预报问题，线性核函数明显不适用；基于sigmoid 核函数的SVR回归预报结果不理想；径向基核函数相较于多项式核函数涉及的参数少，计算复杂度低。基于复杂度和易操作性原则，本文选取径向基核函数作为SVR核函数。

SVR模型可抽象为如下形式：

y=f(x|C，ε，σ)，

(1)

其中C为惩罚系数；ε为不敏感损失系数；σ为径向基核参数。

1.2 萤火虫算法

根据萤火虫算法寻优原理，建立算法优化模型，其相关数学表达式如下：

定义1萤火虫的光亮强度I：

I=I0e-γrij.

(2)

定义2萤火虫的吸引度β：

β=β0e-γrij，

(3)

其中I0和β0分别表示萤火虫的固有光亮强度和吸引度；γ表示光吸收系数；rij表示萤火虫i和萤火虫j之间的欧氏距离，见定义3。

定义3萤火虫i和萤火虫j的欧氏距离rij：

(4)

其中xi和xj分别表示萤火虫i和萤火虫j的搜索空间位置；xid(t)和xjd(t)分别为萤火虫i和j的第t代第d维位置；D表示搜索空间维度。

定义4萤火虫i向萤火虫j的学习后的位置xid(t+1)：

xid(t+1)=xid(t)+β(xjd(t)-xid(t))+αi(t)ε，

(5)

其中αi(t)表示萤火虫i的第t代步长因子；ε服从均匀分布，取值范围为[-0.5，0.5]。

2 多策略融合学习萤火虫算法

在标准萤火虫算法中，每只萤火虫均向其余萤火虫依次学习，采用这种全吸引模型更新方式，既损耗更新资源，又易使种群出现振荡现象，导致种群过早收敛，陷入局部最优[28]。为有效利用更新资源，避免种群“早熟”，本文为不同萤火虫制定了2种学习策略。文献[29]为广义中心粒子设计了深度学习策略，在优胜劣汰的评估机制下，不断挖掘更为优异的广义中心粒子，引导种群进化。为此，本文将深度学习策略应用于最优萤火虫，提高最优萤火虫的质量，加快完成全局挖掘最优解任务；非最优萤火虫采用新添加的领域搜索项更新模式，在随机吸引模型[28]的指导下，完成局部开采种群优势信息任务。最优萤火虫和非最优萤火虫分别采用不同的学习策略，构成一个进化单元，迭代寻优。

2.1 最优萤火虫学习策略

最优萤火虫随机选择两只不同的萤火虫作为学习对象，根据式(6)进行Deep_count次深度学习，使其最大限度挖掘全局最优信息。

(6)

种群进化前期，种群处于发散状态，萤火虫之间位置相距较远，较大的搜索步长，可有效提高种群全局探索最优解能力；种群进化后期，种群处于收敛状态，萤火虫位置相距较近，较小的搜索步长，可有效提高种群局部开采最优解能力。为此，通过柯西函数动态调整步长，并采用贪婪更新策略评估学习状态。评估结果为优，则更新状态；反之，保留原有状态。最优萤火虫学习策略定义如下：

(1)记录最优萤火虫当前位置XBest和适应值f(XBest)；

(3)按照式(6)更新最优萤火虫当前位置；

(4)评估最优萤火虫，结果为优则记录新位置new_XBest和新适应值new_f(XBest)，反之保留更新之前的位置及适应值，同时记录评估次数FEs=FEs+1；

(5)判断最优萤火虫是否完成Deep_count次深度学习，完成跳转至(6)，反之跳转至步骤(1)；

(6)记录最优萤火虫最终位置XBest和适应值f(XBest)，最优萤火虫更新完成。

2.2 非最优萤火虫学习策略

位置更新公式(5)由3项组成，第1项表示萤火虫Xi当前位置，第2项表示萤火虫Xi向萤火虫Xj移动的距离，第3项表示萤火虫Xi完成学习任务后，以自身为中心，搜索附近最优信息移动的步长。在整个种群进化过程中，前两项组合保证了萤火虫算法特有的学习方式，而第3项直至种群寻优任务结束，萤火虫始终做轻微扰动，无法根据当前搜索状态，有效地适应不同搜索域，易陷入局部最优，导致种群出现停滞。

为有效控制萤火虫在更新过程中陷入局部最优，同时适应不同的搜索空间，本文将式(7)作为新位置更新方式，以提高萤火虫逃离局部最优能力。

xid(t+1)=xid(t)+β(xjd(t)-xid(t)) +αi(t)ε(xjd(t)-xid(t)).

(7)

式(7)第3项在原有基础之上乘以两只萤火虫的实时位置差，未额外添加参数。在进化前期，萤火虫i和萤火虫j之间的距离较大，萤火虫j对萤火虫i吸引力很弱，第2项使得萤火虫i向萤火虫j移动距离较小，学习程度不够，为此，第3项结合萤火虫i和萤火虫j的实际距离，增大学习步长，增强学习程度，探索能力显著加强；在进化后期，两萤火虫之间相距很近，第2项和第3项值均相对较小，有效保证萤火虫进行局部开采；因此，式(7)动态地兼顾了种群在整个搜索过程的探索和开采能力。同时，该更新策略适用于不同搜索域，有效地提高了算法的泛化能力。

最优萤火虫完成任务之后，非最优萤火虫启动更新进程。为了充分利用更新资源，非最优萤火虫使用随机吸引模型替代全吸引模型，且每只萤火虫使用式(7)完成位置的搜索。至此，由最优和非最优萤火虫学习策略组合构成种群进化单元，不同萤火虫按照进化单元既定规则协作迭代，完成种群寻优任务。非最优萤火虫学习策略定义如下：

(1)记录除最优萤火虫外N-1个萤火虫当前位置Xi和适应值f(Xi)；

(2)设置随机吸引模型作为N-1个萤火虫位置更新方式；

(3)按照式(7)更新N-1个萤火虫当前位置；

(4)评估最优萤火虫，记录新位置new_Xi和新适应值new_f(Xi)，同时记录评估次数FEs=FEs+1；

(5)N-1个萤火虫完成更新任务并记录最终位置Xi和适应值f(Xi)，非最优萤火虫更新完成。

2.3 算法流程

本文将评估次数作为算法终止条件。其中，种群大小为N，总评估次数为MAX_FEs，当前评估次数为FEs。MSFLFA流程如下：

Step1：初始化种群，记录FEs=N；

Step2：执行最优萤火虫学习策略设定的进化流程；

Step3：执行非最优萤火虫学习策略设定的进化流程；

Step4：若FEs>MAX_FEs，则跳转至Step5，反之，跳转至Step2；

Step5：输出种群最优值和最优解。

3 实验仿真及分析

3.1 基准测试函数及参数设置

在仿真实验中，采用了经典基准测试函数[30]对MSFLFA进行测试。测试集函数包含12个，f1～f7是单峰函数，f8～f12是多峰函数，所有测试函数维度D设置为30。

分别选取同类算法FA，MFA，CFA，WSSFA，VSSFA，RaFA，ApFA算法与本文提出的MSFLFA进行实验比较。所有算法的最大评估次数均设置MAX_FEs为5e5。MSFLFA的萤火虫种群大小N为20，步长因子α为0.2，初始吸引度β0为1.0，光吸收系数γ为1.0，深度学习次数Deep_count为2 400(该参数设置依据参考文献[29])，其它比较算法的参数设置及仿真实验数据均来自参考文献[30]。

3.2 实验结果分析

表1给出了8种算法在基准测试函数上独立运行20次的平均适应值。分析结果可知，MSFLFA在12个测试函数上的精度均优于FA，VSSFA和WSSFA算法；对于单峰函数，MSFLFA在5个测试函数上的精度优于所有算法；对于多峰函数，MSFLFA在4个测试函数上的精度优于所有算法；对于单峰函数f6，MSFLFA，MFA，CFA，RaFA和ApFA均搜索到了理论最优值；MSFLFA在多峰函数f12上的精度显著优于所有算法。

表2分类给出了8种算法在各类型函数集上的Friedman值，结果表明MSFLFA在所有函数、单峰和多峰函数的综合性能均拥有最小Friedman值，算法综合优化性能最佳。

表1 8种算法在基准试函数上的平均适应值

表2 8种算法在基准测试函数上的Friedman值

为了更加直观观察算法的收敛速度和求解精度，图1显示了基准测试函数的所有算法收敛曲线图。从图1(a)，(b)，(c)，(d)，(f)和(i)中可以看出，MSFLFA收敛曲线在整个寻优进程中一直呈线性递减的趋势，表明该算法不断地搜索全局最优值，未出现陷入局部最优的情况。虽然图1(e)，(h)，(j)和(l)中的MSFLFA收敛曲线在搜索进程后期趋于平缓，但其都处在所有算法的收敛曲线下方，表明MSFLFA搜索精度均优于其它算法。图1(f)和(i)中的MSFLFA收敛曲线在寻优进程前期已收敛完成，表明该算法提前完成寻优任务并已搜索到全局最优值。

图1 8种算法在基准测试函数的收敛曲线图

4 MSFLFA在年径流预测中的应用

4.1 MSFLFA与SVR预测模型融合

选择以实测值和预测值的均方误差(Mean Squared Error，MSE)作为MSFLFA的目标函数，其定义如下：

(8)

确定MSFLFA目标函数之后，对{C，ε，σ}参数进行寻优的描述如下：

(1)确定模型的训练和测试样本，设置{C，ε，σ}参数的搜索范围；

(2)初始化种群数量N、步长因子α、初始吸引度β0、光吸收系数γ；

(3)设定搜索范围，初始化位置，将{C，ε，σ}参数与萤火虫的位置关联，即xi=(xi1，xi2，xi3)=(Ci1，εi2，σi3)；

(4)执行最优和非最优萤火虫学习策略设置的进化流程。若MSFLFA达到终止条件，则跳转至步骤(5)，否则，重复执行步骤(4)；

(5)输出最优萤火虫位置信息，即SVR模型的{C，ε，σ}参数值，SVR利用参数值进行模型预测。

4.2 MSFLFA-SVR预测模型应用

实验数据引自参数文献[23]，即以陕西府谷县黄甫川水文站1979—2008年的实测径流值为预测对象，同时以年降水量和蒸发量为年预测因子。

实测径流样本数量共30 a，将前25 a实测径流值作为MSFLFA-SVR模型的训练样本，后5 a实测径流值作为MSFLFA-SVR模型的预测样本，同时以实测径流值和预测径流值的MSE作为MSFLFA-SVR的评价函数。设萤火虫种群大小N为20，迭代次数ITERs为5000，步长因子α为0.2，初始吸引度β0为1.0，光吸收系数γ为1.0，深度学习次数Deep_count为24，惩罚系数C、不敏感损失系数ε和核参数σ的取值区间分别为[0.01，500]、[0，1]和[0.1，100]。

4.3 MSFLFA-SVR预测模型性能评价

通过MSFLFA对SVR模型进行参数优化，得到一组最优的参数组合{C，ε，σ}={37.87，0.02，4.14}，使得MSFLFA-SVR模型预测效果最佳。获得最佳模型之后，将30 a的预报因子输入到模型中，得到最佳的年预测径流值。同时，将年预测径流值与FA-SVR、WSSFA-SVR、MFA-SVR、RaFA-SVR及文献[24]中BP-ANN-SVR和PPR-SVR模型结果进行比较分析，预测结果见图2。

图2 1979—2008年各算法模型年径流值训练和预测结果

为对MSFLFA-SVR模型预测性能进行评价，选取平均绝对误差R1、平均相对误差R2%和确定性系数DC作为模型评价指标[31]。经计算，获得上述7种模型性能评价指标值，计算结果见表3。

表3 7种预报模型的拟合和预报值

分析表3，通过比较3个指标，对于训练样本，FA-SVR模型的3个指标均优于其它6种模型，说明该模型对数据拟合结果较好，MFA-SVR、WSSFA-SVR、RaFA-SVR和MSFLFA-SVR模型的3个指标值近似相等，与FA-SVR模型相比较，拟合效果差于FA-SVR模型，PPR-SVR和BP-ANN-SVR模型的3个指标基本一致，相较于其它5种训练结果，该2种模型拟合效果较差。对于测试样本，FA-SVR模型的3个指标发生了较大变动，R1和R2%分别增大了0.44和0.10，DC减小了0.80，表明该模型预测性能极差，其它6种模型的3个评价指标波动范围相似，相较于训练拟合结果，预测性能均有所下降。

为综合评估7种模型的性能，对上述评价指标进行秩均值检验(训练和测试样本各3个)，秩均值越小，表明模型综合预报性能好，测试结果见表4。通过秩均值排序。由表4数据可知，MSFLFA-SVR模型的秩均值最小，其综合预报性能均优于其它6种模型。

表4 7种模型评价指标的秩均值

5 结束语

MSFLFA通过为最优萤火虫和非最优萤火虫分配不同的学习策略，实现了种群的全局探索和局部开采，且在测试函数上寻优性能表现优异。为验证MSFLFA-SVR模型的性能，将其应用于陕西府谷县黄甫川年径流量的预测，与其它6种预测模型对比，MSFLFA-SVR模型取得了更好的预测效果，为径流预报提供了一种有效的方法。MSFLFA仅从更新策略改善算法性能，未对算法固有参数做优化，同时，算法应用方面仅进行了年径流预报，预报粒度较大，无法满足不同预报需求。因此，算法固有参数设置、月径流和日径流预报是下一步研究工作重点。