近年来，中国铁路在六次大提速的背景下，客货运均呈现出高速化、重载化的趋势.由此导致的轨道质量日益恶化也影响着行车安全和乘车体验^[1-2].传统的轨道养护方法主要分为“故障修”和“周期修”^[3]，已经无法满足铁路高速运行的需要.所以，利用轨道动态检测数据拟合轨道质量历史发展趋势^[4-5]，并对其未来的发展趋势进行有效预测^[6]，对于合理安排养护周期、提高养护效率具有重大意义.

轨道质量的发展具有趋势性和随机性.国内外的学者对轨道质量的预测进行了大量研究.日本学者三和雅史^[7-8]等根据实际养护作业的经验提出基于Logistic分布的特征值推移法来指导大型养路机械(Multiple Tie Tamper，MTT)的作业；日本铁道技术研究所(RTRI)研究人员KAWAGUCHI A^[9]等以轨向不平顺的标准偏差来建立非线性退化模型，并以此为依据确定工务部门的养护计划；曼彻斯特大学(University of Manchester)的HYSLIP J^[10]提出以轨检车得到的轨道几何数据为基础，使用分形法对检测数据进行分析并做出预测.国内学者则主要围绕TQI开展研究.曲建军^[11]将灰色系统不确定性理论引入轨道质量的预测领域，建立了轨道不平顺TITCGM(1，1)-PC灰色非线性预测模型；许玉德^[12]提出了利用动态轨检数据建立线性预测模型，然后利用TQI预测值来编制养护计划；韩晋^[13]提出了基于非等间距灰色理论和BP神经网络残差校正的预测方法.

但是，现有的预测方法均存在缺陷.例如，线性预测模型对异常点无法进行较好拟合，鲁棒性差；而以灰色理论为基础的BP神经网络残差校正模型，因BP网络初始权值和阈值均为随机赋值，导致预测结果缺乏稳定性，预测效果欠佳.

基于此，本文提出一种将非等间距灰色模型^[14-16]和遗传算法^[17-18]优化的Elman神经网络^[19-20]相结合的预测方法(GM-GA-Elman).首先通过非等间距灰色模型GM(1，1)得到TQI序列的大致发展趋势和残差序列，然后利用遗传算法优化Elman神经网络的权值和阈值，最后使用优化后的Elman神经网络进行残差校正，最终得到一个更加准确的TQI预测序列.

1 非等间距灰色模型 1.1 非等间距序列的生成

设原始TQI序列值为

$ {X^{\left( 0 \right)}} = \left\{ {{x^{\left( 0 \right)}}\left( {{t_1}} \right),{x^{\left( 0 \right)}}\left( {{t_2}} \right), \cdots ,{x^{\left( 0 \right)}}\left( {{t_n}} \right)} \right\}. $

(1)

序列x⁽⁰⁾(t_i)所对应的检测时间为t_i，故原始TQI序列值的时间差为

$ \Delta {t_i} = {t_i} - {t_{i - 1}},i = 2,3, \cdots ,n. $

(2)

由于TQI序列的检测时间并非等差分布，即

$ \Delta {t_i} \ne {\rm{const,}}i = 2,3, \cdots ,n. $

(3)

式中const表示常量，X⁽⁰⁾为非等间距序列.

对原始TQI序列做一阶累加生成，并考虑到检测时间差Δt_i，则有

$ {x^{\left( 1 \right)}}\left( {{t_i}} \right) = \sum\limits_{k = 1}^i {{x^{\left( 0 \right)}}\left( {{t_k}} \right)\Delta {t_k}} ,i = 2,3, \cdots ,n. $

(4)

x⁽¹⁾的初值一般取x⁽¹⁾(t₁)=x⁽⁰⁾(t₁)，但是预测曲线不一定经过原始序列初值点，且由于原始序列中所有的点均与预测模型相关，所以本文将初值选取为原始序列均值，即

$ {x^{\left( 1 \right)}}\left( {{t_1}} \right) = \frac{1}{n}\sum\limits_{k = 1}^n {{x^{\left( 0 \right)}}\left( {{t_k}} \right)} . $

(5)

1.2 非等间距灰色模型的建模

非等间距序列模型得到的一阶累加生成序列X⁽¹⁾为

$ {X^{\left( 1 \right)}} = \left\{ {{x^{\left( 1 \right)}}\left( {{t_1}} \right),{x^{\left( 1 \right)}}\left( {{t_2}} \right), \cdots ,{x^{\left( 1 \right)}}\left( {{t_n}} \right)} \right\}. $

(6)

由一阶累加生成模块建立GM(1, 1)模型，白化形式微分方程为

$ \frac{{{\rm{d}}{x^{\left( 1 \right)}}\left( t \right)}}{{{\rm{d}}t}} + a{x^{\left( 1 \right)}}\left( t \right) = u,t \in \left[ {0,\infty } \right). $

(7)

将上式离散化，即对上式左右两边在[t_i-1, t_i]上进行积分：

$ {x^{\left( 0 \right)}}\left( {{t_i}} \right)\Delta {t_i} + a\int\limits_{{t_{i- 1}} }^{{t_i}} {{x^{\left( 1 \right)}}\left( t \right){\rm{d}}t} = u\Delta {t_i},i = 2,3, \cdots ,n. $

(8)

称

$ {z^{\left( 1 \right)}}\left( {{t_i}} \right) = \int\limits_{{t_{i- 1}} }^{{t_i}} {{x^{\left( 1 \right)}}\left( t \right){\rm{d}}t} ,i = 2,3, \cdots ,n. $

(9)

为x⁽¹⁾(t_i)在区间[t_i-1, t_i]上的白化背景值，一般将其简化^[15]为

$ {z^{\left( 1 \right)}}\left( {{t_i}} \right) = \frac{{\left( {{x^{\left( 1 \right)}}\left( {{t_i}} \right) + {x^{\left( 1 \right)}}\left( {{t_{i - 1}}} \right)} \right)}}{2}\Delta {t_i}. $

(10)

离散化差分方程为

$ {x^{\left( 0 \right)}}\left( {{t_i}} \right)\Delta {t_i} + a{z^{\left( 1 \right)}}\left( {{t_i}} \right) = u\Delta {t_i}. $

(11)

式中：a为发展系数，用来控制系统发展态势的大小，u为灰色作用量，用来反映数据变化的不确切关系.

利用最小二乘法求得待辨识参数a，u为

$ {\left( {a,u} \right)^{\rm{T}}} = {\left( {{\mathit{\boldsymbol{B}}^{\rm{T}}}\mathit{\boldsymbol{B}}} \right)^{ - 1}}{\mathit{\boldsymbol{B}}^{\rm{T}}}\mathit{\boldsymbol{Y}}. $

(12)

其中：

$ \mathit{\boldsymbol{B}} = \left[ {\begin{array}{*{20}{c}} { - {z^{\left( 1 \right)}}\left( {{t_2}} \right)}&1\\ { - {z^{\left( 1 \right)}}\left( {{t_3}} \right)}&1\\ \vdots&\vdots \\ { - {z^{\left( 1 \right)}}\left( {{t_n}} \right)}&1 \end{array}} \right], $

(13)

$ \mathit{\boldsymbol{Y}} = \left[ {\begin{array}{*{20}{c}} {{x^{\left( 0 \right)}}\left( {{t_2}} \right)}\\ {{x^{\left( 0 \right)}}\left( {{t_3}} \right)}\\ \vdots \\ {{x^{\left( 0 \right)}}\left( {{t_n}} \right)} \end{array}} \right]. $

(14)

将参数a，u代入微分方程解得时间响应函数为

$ \begin{array}{l} {{\hat x}^{\left( 1 \right)}}\left( {{t_i}} \right) = \left( {{x^{\left( 0 \right)}}\left( {{t_1}} \right) - \frac{u}{a}} \right){{\rm{e}}^{ - a\left( {{t_i} - {t_1}} \right)}} + \frac{u}{a},\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;i = 2,3, \cdots ,n. \end{array} $

(15)

从而得到GM(1, 1)初步预测值为

$ {{\hat x}^{\left( 0 \right)}}\left( {{t_i}} \right) = {{\hat x}^{\left( 1 \right)}}\left( {{t_i}} \right) - {{\hat x}^{\left( 1 \right)}}\left( {{t_{i - 1}}} \right),i = 2,3, \cdots ,n. $

(16)

1.3 背景值的优化

由式(12)~(15)可知背景值z⁽¹⁾(t_i)的构造形式决定灰色模型的预测精度.由式(10)可知，在区间[t_i-1, t_i]上使用梯形公式来代替曲线x⁽¹⁾(t_i)与时间横轴围成的面积造成精度误差.为了减少此误差，利用文献[21]中的方法，即把背景值表达式优化为

$ {z^{\left( 1 \right)}}\left( {{t_i}} \right) = \frac{{{x^{\left( 1 \right)}}\left( {{t_i}} \right) - {x^{\left( 1 \right)}}\left( {{t_{i - 1}}} \right)}}{{\ln {x^{\left( 1 \right)}}\left( {{t_i}} \right) - \ln {x^{\left( 1 \right)}}\left( {{t_{i - 1}}} \right)}} \cdot \Delta {t_i}. $

(17)

2 GA-Elman残差校正模型 2.1 Elman神经网络

Elman神经网络是J.E.Elman于1990年首先提出来的一种典型的局部回归神经网络. Elman网络具有与多层前向网络(如BP神经网络)相似的多层结构.但是不同于传统的BP神经网络，Elman网络中多出的承接层，增加了网络处理动态信息的能力和适应时变特性的能力^[19-20].其结构主要包括输入层、隐含层、输出层、承接层，网络模型见图 1.

2.2 GA-Elman残差校正模型

遗传算法(Genetic Algorithm, GA)是一种经典的启发式搜索算法，通用性强，能对多个领域的问题进行全局最优解搜索^[17-18].

轨道质量状态具有很强的趋势性和随机性，本文将遗传算法和Elman神经网络结合起来形成GA-Elman模型，即利用遗传算法对Elman神经网络的初始权值和阈值进行优化，优化后的网络在训练过程中能对突变数据进行充分学习并消除预测过程中随机因素的干扰，从而提高了预测精度. GA-Elman残差校正模型的原理如下：

1) 通过改进的非等间距灰色模型得到初步预测值并求得残差序列

$ \delta \left( {{t_i}} \right) = {x^{\left( 0 \right)}}\left( {{t_i}} \right) - {{\hat x}^{\left( 0 \right)}}\left( {{t_i}} \right),i = 1,2, \cdots ,n. $

(18)

2) 为了避免不同量级的数据相互影响，对残差序列δ(t_i)进行归一化处理如下：

$ y = \frac{{\left( {{y_{\max }} - {y_{\min }}} \right) * \left( {\delta - {\delta _{\min }}} \right)}}{{{\delta _{\max }} - {\delta _{\min }}}} + {y_{\min }}. $

(19)

式中：y_min=-1，y_max=1.

3) Elman神经网络参数的确定：

神经网络的训练效果对参数极为敏感，具体参数的选取对实际情况的依赖度较高，采用经验值的场合较多^[22-23].其中，隐含层节点数一般遵循公式${n_{\rm{h}}} = \sqrt {{n_{\rm{i}}} + {n_{\rm{o}}}} + m$，其中n_h为隐含层节点数，n_i为输入层节点数，n_o为输出层节点数，m为1~10之间的整数.通过多次实验发现选取隐含层节点数为3时，网络结构较为稳定.同时，由于本文待处理数据量较小，如果算法在短时间内快速收敛则容易陷入局部最优，故选取学习速率为0.1，附加动量因子取0.9以保证算法能高效运行又不至于陷入局部最优.另外，当隐含层节点数较小时，选取网络激活函数为双曲正切S型函数${\rm{tansig}}\left( \mathit{n} \right) = \frac{2}{{1 + {{\rm{e}}^{ - 2n}}}} $-1以及输出层传递函数为线性传递函数purelin(n)=n可以保证系统稳定运行；选取的迭代次数要足够大保证网络能充分训练以接近全局最优，故本文选取迭代次数为4 000次.

4) 遗传算法参数初始化：

待处理数据量较小，经多次实验发现，当Elman网络的隐含层节点数为3时，种群在30代以内基本接近收敛，故设置进化代数为30代，种群规模为10；同时，经多次测试发现在GA-Elman模型运算过程中，染色体的平均适应度容易快速下降接近至最佳适应度导致模型快速收敛，故选用较低的交叉概率p_cross=0.2以降低收敛速度，选用较高的变异概率p_mutation=0.1来提高种群产生更优个体的能力，保证了最终结果平稳地收敛至最优解；

5) 染色体编码：

实数编码方式可以扩大寻优范围，同时无需进行频繁地编解码，相较于二进制编码，运算效率显著提高^[24]，故本文选用实数编码方式.

神经网络在初始化权值和阈值系数时一般是在(-1, 1)的范围内随机产生.适当扩大编码范围可适当提高染色体多样性，产生新的优秀个体的概率也随之变大.此范围过大则会使得权值和阈值过于分散，不利于种群进化.故本文设置权值和阈值系数的范围为(-3, 3)，在此范围内随机产生10组权值和阈值系数构成一条染色体，该染色体包含了构成一个完整神经网络的全部权值和阈值系数；

6) 确定适应度函数：

适应度函数的选取直接决定了遗传进化的效果.将每个个体携带的权值和阈值赋入神经网络进行训练并预测输出，把预测输出${{\hat y}_{\rm{f}}}$和期望输出y_f之间的误差绝对值和E作为个体的适应度F_fitness，由此可知，适应度值越小，个体越优秀.

$ \begin{array}{*{20}{c}} {E = {k_1}\sum\limits_{i = 1}^N {\left( {\left| {{y_{\rm{f}}}\left( i \right) - {{\hat y}_{\rm{f}}}\left( i \right)} \right|} \right),} }\\ {{F_{{\rm{fitness}}}} = E.} \end{array} $

式中:k₁为系数，N为网络输出节点数.

7) 遗传算法选择操作：

选择操作一般都有轮盘赌法、随机遍历抽样法、锦标赛法等^[25].选用轮盘赌法.对于个体i，其被选择的概率是p_i，则：

$ {f_i} = {k_2}/{F_i} $

$ {p_i} = \frac{{{f_i}}}{{\sum\limits_{j = 1}^N {{f_i}} }} $

式中：F_i为个体i的适应度值，k₂为取倒系数，N为种群规模.

8) 交叉操作：

采用实数交叉，第k个染色体a_k和第l个染色体a_l在第j位的交叉^[26]操作如下

$ \left\{ \begin{array}{l} {a_{kj}} = {a_{kj}}\left( {1 - b} \right) + {a_{lj}}b,\\ {a_{lj}} = {a_{lj}}\left( {1 - b} \right) + {a_{kj}}b. \end{array} \right. $

9) 变异操作：

选取第i个个体的第j个基因a_ij进行变异^[26]，方法如下：

$ {a_{ij}} = \left\{ \begin{array}{l} {a_{ij}} + \left( {{a_{ij}} - {\rm{HB}}} \right) * f\left( g \right),r > 0.5;\\ {a_{ij}} + \left( {{\rm{LB}} - {a_{ij}}} \right) * f\left( g \right),r \le 0.5; \end{array} \right. $

$ f\left( g \right) = {r_2}{\left( {1 - \frac{{{G_{{\rm{cur}}}}}}{{{G_{\max }}}}} \right)^2}. $

式中：HB和LB分别为基因编码值的上下限，G_cur为当前迭代次数，G_max为总迭代次数，r和r₂为[0, 1]间的随机数.

根据以上参数设置建立GA-Elman残差校正模型.当种群进化到每一代染色体满足以下条件之一时，视为满足终止条件：

a) 每一代平均适应度值都收敛于最佳适应度值时；

b) 当连续多代染色体的平均适应度与最佳适应度差值较小，但又不收敛于最佳适应度时，为了提高算法运算效率，终止进化.

最后，将满足终止条件得到的最优种群解码为网络权值和阈值系数，并赋值给Elman网络进行训练得到GA-Elman残差校正模型.

10)残差校正：

利用训练好的Elman网络对残差序列δ(t_i)进行预测得到最终残差输出$\hat \delta \left( {{t_i}} \right)$，将其与初步预测值$ {{\hat x}^0}\left( {{t_i}} \right)$相加得到网络最终输出

$ x\left( {{t_i}} \right) = {{\hat x}^{\left( 0 \right)}}\left( {{t_i}} \right) + \hat \delta \left( {{t_i}} \right),i = 1,2, \cdots ,n $

3 GM-GA-Elman预测模型

根据以上理论分析，GM-GA-Elman预测模型原理如下利用非等间距灰色模型对TQI序列进行预测，把得到的初步预测值与真实值进行比较得到残差序列.最后将残差序列输入到GA-Elman残差校正模型中，得到本文的最终预测值.

具体预测步骤如下：

1) 在进行灰色建模之前，要对输入数据进行级比检验以判断其是否适用于灰色模型建模.方法如下

$ \lambda \left( {{t_i}} \right) = \frac{{{x^{\left( 0 \right)}}\left( {{t_{i - 1}}} \right)}}{{{x^{\left( 0 \right)}}\left( {{t_i}} \right)}},i = 2,3, \cdots ,n $

若对∀λ(t_i), i=2, 3, …, n都有

$ \lambda \left( {{t_i}} \right) \in \left( {{{\rm{e}}^{ - \frac{2}{{n + 2}}}},{{\rm{e}}^{\frac{2}{{n + 2}}}}} \right),i = 2,3, \cdots ,n $

则输入数据可以直接使用灰色模型建模.否则需要选取适当常数c对原始数据序列X⁽⁰⁾进行适当平移

$ {X^{\left( 0 \right)}} = {X^{\left( 0 \right)}} + c. $

2) 将步骤1)中处理完以后的数据做累加处理，得到一阶累加序列X⁽¹⁾；

3) 将一阶累加序列X⁽¹⁾输入至非等间距GM(1，1)模型，同时利用优化后背景值表达式(17)代替传统背景值表达式(10)进行灰色建模，求得辨识参数a，u；

4) 由a，u得到最终的时间响应函数并求得初步预测值${{\hat X}^{\left( 0 \right)}}$和残差序列δ(t_i)；

5) 将δ(t_i)输入至GA-Elman残差校正模型中得到残差校正值$\hat \delta \left( {{t_i}} \right)$；

6) 由残差校正值$\hat \delta \left( {{t_i}} \right)$和初步预测值${{\hat X}^{\left( 0 \right)}}\left( {{t_i}} \right)$得到最终预测值

$ x\left( {{t_i}} \right) = {{\hat x}^{\left( 0 \right)}}\left( {{t_i}} \right) + \hat \delta \left( {{t_i}} \right),i = 1,2, \cdots ,n $

7) 若在步骤1)中进行了平移，则反平移得到最终的TQI预测序列.流程图见图 2.

4 TQI预测效果分析

选取提速干线沪昆线上行K226.0~K226.2以及K226.4~K226.6的轨检车TQI数据对GM-GA-Elman预测模型进行验证.为了验证本文模型的可靠性，选取参考文献[11~13]的方法进行对比.

决定系数(R-square)、相关系数以及均方根误差(Root Mean Square Error, RMSE)是统计学中常用的拟合优劣判定标准^[26].

4.1 K226.0~K226.2区间段预测结果分析

选取K226.0~K226.2线路段2007年9月~2008年8月期间共19组TQI序列作为原始训练序列，预测2008年9月~2008年12月共7组TQI值.

图 3为经过非等间距灰色模型预测后的结果.

由图 3可知，单一的非等间距灰色模型只能反映TQI序列的大致变化趋势，难以体现局部变化.然后，利用GA-Elman残差校正模型进行残差校正，得到的结果见图 4.

种群进化曲线见图 5，当进化至第17代时，种群趋于收敛，说明此时模型已接近最优解.

本文预测结果与参考文献[11~13]的方法预测结果的对比见表 1和表 2.

结果表明，文献[11]中建立的TITCGM(1，1)-PC模型使用周期性函数对原始数据列中隐含周期成分和一定的随机成分进行拟合，平均相对误差为5.46%，RMSE较高，达到了0.430，模型预测效果一般；由于线性模型结构过于单一，所以文献[12]建立的线性预测模型预测效果较差，均方根误差过高；文献[13]中利用BP神经网络进行残差校正，虽然可以较好地拟合训练区间，但由于模型未对初始权值和阈值进行优化，导致预测结果不够理想，RMSE达到0.224；本文预测方法充分发挥了神经网络的优势，对原始TQI序列拟合效果较好，相对平均误差低至1.86%，同时，R-square、相关系数、均方根误差等三项统计指标均优于文献[11]、^[13]的方法.虽然从相关系数来看，文献[12-13]的方法与本文方法差别不大，但是从R-square来看，本文方法的优势十分明显.

4.2 K226.4~K226.6区间段预测结果分析

另选取K226.4~K226.6线路段2007年9月~2008年8月期间共19组TQI序列作为训练数据，预测2008年9月~2008年12月共7组TQI值.

经过非等间距灰色模型预测后的结果见图 6.

通过GA-Elman残差校正模型进行残差校正后得到的结果见图 7.

遗传算法种群进化曲线见图 8，由图可知当种群进化至30代时，种群趋于收敛，说明此时模型已接近最优解.

结果表明，文献[12]中的方法由于未考虑突变数据点的影响，使得线性模型中的斜率过大，导致预测结果几乎不具备研究价值；预测序列平均相对误差低至1.89%，低于文献[11]的5.74%和文献[13]的2.59%；均方根误差显著下降，低至0.165；确定系数也明显优于文献[11]、[13].虽然在个别点的精度不如文献[11]、[13]，但是整体来说，本文方法的多项统计性能均优于文献[11]、[13]，且在不同区段上均保持稳定，鲁棒性更好.

本文预测结果与其他方法的对比见表 3和表 4.

5 结论

1) 对非等间距灰色模型进行背景值优化，利用积分的几何意义，将累加生成序列视为非齐次指数函数形式，重构非等间距灰色模型的背景值，有效提高了模型的预测精度值；

2) 利用遗传算法进行全局寻优得到权值和阈值的最优解，相对于其他方法中对初始权值和阈值随机赋值，本文方法稳定性更高；

3) 由于Elman神经网络具有动态反馈特性，对于发展趋势具有典型随机性的TQI序列来说，Elman神经网络能更好地拟合原始数据，预测效果更佳；

4) 使用多种统计指标评价本文方法，提高了评价的准确性和可靠性；

5) 本文验证了Elman神经网络在小样本、贫信息的预测问题中，可以对短期内的数据进行较好地预测，但是如何将其引入中长期预测问题中，还有待进一步研究.