CN101592673B - 铁路沿线风速预测的方法 - Google Patents

Abstract

一种铁路沿线风速预测的方法，包括如下步骤：(1)选择小波分析法的分解与重构算法对拟预测风速序列进行多层分解与重构计算，实现将原始非平稳风速序列分解为多层平稳的风速序列；(2)利用改进时间序列分析法对分解重构后的各层风速序列分别建立相应的时序预测模型；(3)利用各层所建的对应预测模型对分解层风速序列进行超前多步预测计算，然后对各层风速预测值加权计算，获取原始风速序列的相关预测值。本发明可用于不同预报场合、可获得高精度短期预测的优化算法。

Description

铁路沿线风速预测的方法

技术领域

本发明属于大风恶劣环境下铁路安全行车的监测和控制技术领域，具体的涉及一种根据铁路沿线的测风站所测定的原始风速时间序列，精确的进行短时间内风速预测的方法。

背景技术

大风是影响铁路运输安全的主要气象灾害之一，在其作用下列车绕流流场明显改变，气动性能恶化，严重影响到列车的侧向稳定性。在风口区域的特大桥梁、高路堤、丘陵以及曲线等一些特殊路段，环境风导致列车气动力显著增大，与列车动作用力相互叠加，使列车脱轨、倾覆的可能性大大增加。

由大风引发的行车事故在世界各国时有发生。大风对铁路行车安全构成了严重威胁。为防止大风给列车运行安全带来的危害，日本、德国、法国等一些铁路运输发达国家主要采取的方法是在强风常发的个别区段设置观测点，当区段最高风速超过某一限值时，发出列车限速或停运指令。我国兰新铁路亦已初步建立了大风监测***，但该***仅能提供简单的气象监测功能，无法满足在大风环境下实时地根据线路和当前运行列车信息为安全行车提供限速指令或指挥调度决策的迫切要求。

沿线风速受高原特殊气流和地形地貌等诸多因素影响，具有突发多变的不确定性特征。因此，要真正实现对沿线列车的准确预警，***不仅要将列车车型载重、路况参数、测风站实测风速等信息作为发出限速预警指令的决定性因素，同时应该考虑沿线风速的未来变化趋势。只有这样才能提高预警***所发出列车限速预警指令的准确性，同时给被限速列车以一定的应急处理时间。因此，当前迫切需要对已经研建的青藏铁路大风监测预警***升级风速数值预报功能，实现***能够根据沿线测风站k时刻以前的实测历史风速，自动计算该测风站k+1及以后时刻的风速预报值。

铁路沿线风速时间序列与风电场风速、机械故障振动、铁路货运量等其他 时间序列均为同一类非平稳随机时间序列。它们在建模预测方面具有较大的相似性，均涉及序列非平稳判断处理、模型辨识、模型定阶、参数估计、预测计算等重要建模环节。

国外发生多起由于大风所引发的铁路重大安全事故，引起了国外相关研究机构和学者的高度重视。日本等多国学者对大风影响列车行驶安全的机理进行深入研究，并提出在沿线设置挡风墙等措施以减小风速危害。

对风速对象的预测研究也发展非常迅速。许多学者提出了各种风速预测的有效方法。德国学者Uwe Hopp mann提出一种线性外推法实现铁路沿线风速预测。该方法通过在铁路沿线设置气象站采集实时风速，对所采集的风速序列进行线性外推计算以获取未来时刻的风速预测值。美国学者Xie Lian提出一种基于非线性控制理论的飓风表面风预测方法。美国学者Milli-gan Michael利用时间序列法建立风速ARMA模型，实现某风电场的风速预测。希腊学者Zaphiro-poulos Yiorgos提出使用空间相关法对某山岭测风站实测风速建立“测量-相关计算-预测”模型，成功实现该点风速、风向预测。美国学者Li Shu hui提出使用神经网络法和卡尔曼滤波法建立预测改进模型，实现风速短时高精度预测。其建立改进模型的核心步骤是利用卡尔曼滤波法确定循环多层感知网络(RMLP网络)的初始连接权值。加拿大学者El-Fouly T.H.M.建立灰色GM(1，1)模型对风电场风速、风电功率实现短时预测。印度学者Basu S.提出使用遗传算法对北印度洋某测点海洋表面风速实现预测。美国学者Rife D.L.利用风土学逼近论成功对低水平风速建立中间比例预测模型，实现风速预测。

国内学者也对铁路沿线风速开展相关研究，如研建铁路沿线大风监测***。但所研建***仅涉及简单的风速采集监测、风速统计分析等功能，并没有开展实质性的铁路沿线风速预测研究工作。该方向到目前为止尚属空白。同时，国内对其他类型风速序列的预测研究也刚刚起步，其研究风速对象主要集中在风电场风速上。

就风电场风速预测领域而言，国内学者取得了不俗的成果，提出了利用时间序列法对我国某风电场风速开展建模预测，先对非平稳风速时间序列进行差分处理，获取平稳序列后再建立预测模型，最终依靠模型实现预测计算。还有的利用神经网络法对我国某风电场风速建立神经网络BP风速预测模型。另外一种基于时间序列法和神经网络法，提出可明显提高风速预测精度的混合建模算法。该方法先通过时间序列法对风速序列实现初步建模，获取直接反映风速特 征的时间序列模型方程，然后根据模型方程选择神经网络法BP模型的初始连接权值及神经元个数，使得混合算法获得了较高预测精度。除了风速预测之外，其他对象的预测研究也迅速发展，如对电力负荷量、桥梁振动信号、设备故障信号、铁路货运量等时间序列。

总结国内外参考文献可知，对铁路沿线风速开展预测研究属于前沿方向，其涉及交通运输工程、环境气象、***辨识、智能控制、信号处理等多个学科，研究难度较大。到目前为止，该方面涌现的研究成果非常稀少，几乎没有能够直接参考使用的现成算法。

发明内容

本发明提供了一种可用于不同预报场合、可获得高精度短期预测的优化算法的铁路沿线风速预测的方法，其基于成熟的小波分析法、时间序列分析法、卡尔曼滤波法三种智能算法，计算难度低，计算量小，能够进行超前一分钟以上的预测计算，并可获得显式预测方程式，进行超前多步高精度预测。

本发明所采用的技术方案如下：

一种铁路沿线风速预测的方法，用于铁路沿线风速的短期预测，其特征在于所述方法包括如下步骤：

(1)选择小波分析法的分解与重构算法对拟预测风速序列进行多层分解与重构计算，实现将原始非平稳风速序列分解为多层平稳的风速序列；

(2)利用改进时间序列分析法对分解重构后的各层风速序列分别建立相应的时序预测模型；

(3)利用各层所建的对应预测模型对分解层风速序列进行超前多步预测计算，然后对各层风速预测值加权计算，获取原始风速序列的相关预测值。

具体的讲，所述步骤(2)之后还包括如下步骤：

通过改进时间序列分析和卡尔曼滤波法混合建立预测精度模型，进一步提高改进时间序列法获得的超前一步预测精度。

所述原始风速序列是在铁路沿线上设置测风站，获得该测风站在某一短时间内的风速与时间关系序列。

所述步骤(1)包括将原始风速序列v(t)(t＝1，2，3，…)经小波分析法Mallat塔式算法进行风速信号分解与重构，分别获得高频风速序列vⁿ(t)(n＝1，2，3，…)和低频风速序列ω(t)。

所述步骤(2)包括：

对高频风速序列vⁿ(t)(n＝1，2，3，…)、低频风速序列ω(t)运用改进时间序列法建立时序预测模型，计算出各分解层的超前k步高频风速序列预测值 

低频风速序列预测值 

所述预测精度模型的建立包括：

对高频风速序列vⁿ(t)(n＝1，2，3，…)、低频风速序列ω(t)运用卡尔曼时间序列法建立预测精度模型，计算出各分解层的超前1步高频风速序列预测值 

低频风速序列预测值 

对各分解层的超前k步高频风速序列预测值 

低频风速序列预测值 

进行加权计算，输出最终风速预测值 

所述步骤(1)包括：选用具有时频紧支撑和高正则性的Daubechies 6小波对原始风速序列v(t)进行信号分解，分解深度n为3-6；选用Mallat塔式算法对分解后的风速序列在不同尺度上进行信号重构。

所述步骤(2)包括：

以原始风速序列v(t)进行信号分解后的某一层低频或高频分量序列为建模样本序列，取该层序列前多个数据完成滚动式时间序列法建模，并求出相对于采样点超前三步的预测风速值。

由于风速时间序列是一种变化异常复杂的非线性、非平稳性序列，通常难以使用单一的预测算法实现较高精度预测。同时，风速序列潜在的物理规律是不断变化的，预测模型只有能够实时跟踪这种变化规律，才有可能获取高精度的预测结果。铁路沿线风速预测的方法所采用的改进时间序列分析法具备强大的短期超前多步预测能力，对风速序列跳跃采样点的追踪能力表现异常出色。在对超前一步预测精度要求相对宽松的场合，运用改进时间序列分析法进行建模预测是一种理想的选择。卡尔曼时间序列分析法明显提高改进时间序列分析法获得超前一步预测精度，体现出在超前一步预测方面的优秀算法特性。该优化的算法体系(包括改进时间序列分析法和卡尔曼时间序列分析法)相对于其他算法而言，其混合建模过程并没有明显增加计算难度和计算量，通过借助相关高级语言容易编程实现其优化计算步骤，在工程中能够实时地获得超前风速预报值，可在相关风速预警***中推广使用。

为进一步说明本发明的有益效果，下面结合附图和具体实施方式对本发明作进一步的说明。

附图说明

图1是本发明具体实施方式中测定的原始风速序列图；

图2是本发明具体实施方式中v(t)序列小波分解与重构计算结果图；图2(a)是{X_1t}序列；图2(b)是{X_2t}序列；图2(c)是{X_3t}序列；图2(d)是{X_4t}序列；

图3是本发明具体实施方式中各层风速序列的后50个数据进行超前三步预测计算结果图；图3(a)是{X_1t}序列超前三步预测结果；图3(b)是{X_2t}序列超前三步预测结果；图3(c)是{X_3t}序列超前三步预测结果；图3(d)是{X_4t}序列超前三步预测结果；

图4是本发明具体实施方式中原始风速序列加权计算后的超前三步预测结果图的v(t)序列超前三步预测结果(1)；

图5是本发明具体实施方式中原始风速序列加权计算后的超前不同步的预测结果图；图5(a)是v(t)序列超前一步预测结果(1)；图5(b)是v(t)序列超前五步预测结果(1)；图5(c)是v(t)序列超前十步预测结果(1)；

图6是利用传统的时间序列分析法进行原始风速序列的超前不同步的预测结果图；图6(a)是v(t)序列超前一步预测结果(2)；图6(b)是v(t)序列超前三步预测结果(2)；图6(c)是v(t)序列超前五步预测结果(2)；图6(d)是v(t)序列超前十步预测结果(2)；

图7是本发明具体实施方式中引入卡尔曼滤波法的预测结果图。

具体实施方式

该铁路沿线风速预测的方法以青藏铁路格尔木-拉萨段18号测风站实测原始风速序列v(t)(每分钟一个采样点)进行建模与预测，原始序列v(t)如图1所示。取前150个数据建立预测模型，后50个数据检验模型。

首先，选择小波分析法的分解与重构算法对拟预测风速序列进行多层分解与重构计算，实现将原始非平稳风速序列分解为多层平稳的风速序列；

选用具有时频紧支撑和高正则性的Daubechies 6小波对原始序列v(t)进行信号分解，分解深度n取3。选用Mallat塔式算法对分解后的风速序列在不同尺度上进行信号重构。

vⁿ(t)(n＝1，2，3)分别表示在第1～3层上完成重构计算后获得的高频分量序列。ω(t)表示第3层重构处理后获得的低频分量序列。为了建模方便，将第1 层高频分量序列v¹(t)记为序列{X_4t}，第2层高频分量序列v²(t)记为序列{X_3t}，第3层高频分量序列v³(t)记为序列{X_2t}，第3层低频分量序列ω(t)记为序列{X_1t}。各序列如图2所示。

然后，利用改进时间序列分析法对分解重构后的各层风速序列分别建立相应的时序预测模型；

以第3层低频分量序列{X_1t}为建模样本序列，说明滚动式时间序列法的算法计算步骤，而其他分解层风速序列的建模步骤可参照完成。取序列{X_1t}前150个数据完成滚动式时间序列法建模，并求出相对于采样点X_1t(200)超前三步的预测风速值 

具体算法计算步骤如下所示：

(1)利用时间序列分析法对{X_1t(1)，X_1t(2)，...，X_1t(150)}序列进行BJ建模，选择最小信息准则AIC进行模型定阶，选择矩估计进行模型参数估计。确定序列{X_1t}的合适模型为AR IMA(6，1，0)，获得其对应的预测方程为：

X_1t(t)＝2.7423X_1t(t-1)-2.9532X_1t(t-2)

+1.8939X_1t(t-3)-1.3428X_1t(t-4)

+1.1826X_1t(t-5)-0.7913X_1t(t-6)

+0.2685X_1t(t-7)+a_t                           (1)

使用模型方程(1)实现风速序列{X_1t}超前一步预测计算，得到超前一步预测值 

(2)保持步骤(1)所辨识的合适模型类别，利用 

序列对模型参数进行重新估计，得到包含预测值 

信息特征的新模型方程，继续进行超前一步预测计算，得到相对于采样点X_1t(200)超前两步的预测值 

X_1t(t)＝2.4731X_1t(t-1)-2.1436X_1t(t-2)

+0.7840X_1t(t-3)-0.2292X_1t(t-4)

+0.2645X_1t(t-5)-0.3284X_1t(t-6)

+0.1796X_1t(t-7)+a_t                   (2)

(3)利用 

序列对模型重新参数估计，得到包含预测值 

信息特征的新模型方程，继续进行超前1步预测，得到相对于采样点X_1t(200)超前三步的预测值 

这就是滚动式时间序列法完成超前三步预测的一个完整计算周期。

X_1t(t)＝2.4271X_1t(t-1)-2.1279X_1t(t-2)

+0.9923X_1t(t-3)-0.4392X_1t(t-4)

+0.1067X_1t(t-5)-0.0212X_1t(t-6)

+0.0621X_1t(t-7)+a_t                        (3)

(4)当完成一个计算周期后，再根据最新的实测风速序列重新开始滚动式预测计算，获取采样点X_1t(201)超前三步的预测值 

再者，对各分解层进行超前多步预测计算，

参考滚动式时间序列法建模实例，分别对序列{X_1t}，{X_2t}，{X_3t}，{X_4t}分别建立预测模型，然后对各层风速序列的后50个数据进行超前三步预测计算，结果如图3所示。

当完成各分解层风速序列的预测计算后，根据下述公式加权计算序列v(t)的超前三步预测值 

${\hat{v}}_t\left(3\right)={\mathrm{\ρ}}_1{\hat{X}}_{1t}\left(3\right)+{\mathrm{\ρ}}_2{\hat{X}}_{2t}\left(3\right)+{\mathrm{\ρ}}_3{\hat{X}}_{3t}\left(3\right)+{\mathrm{\ρ}}_4{\hat{X}}_{4t}\left(3\right)$

取加权系数ρ₁＝ρ₂＝ρ₃＝1，预测结果如图4所示。此外，计算原始序列v(t)的超前一步预测值 

超前五步预测值 

超前十步预测值 

结果如图5所示。

比较例  为了进一步显示本发明的良好性能，再利用传统的时间序列分析法直接对原始序列v(t)的同一段风速样本进行建模，然后分别进行超前一步、三步、五步、十步预测计算，结果如图6所示。

引入通用的预测精度评价指标对图6所示预测结果进行评估，评价指标的相关计算公式如下所示：

(1)平均误差： $\mathrm{AE}=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}[X\left(i\right)-\hat{X}\left(i\right)]---\left(4\right)$

(2)平均绝对误差： $\mathrm{MAE}=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}|X\left(i\right)-\hat{X}\left(i\right)|---\left(5\right)$

(3)平均相对误差： $\mathrm{MRE}=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}\left|\frac{X\left(i\right)-\hat{X}\left(i\right)}{X\left(i\right)}\right|---\left(6\right)$

(4)均方根误差： $\mathrm{RMSE}=\sqrt{\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{[X\left(i\right)-\hat{X}\left(i\right)]}^2}---\left(7\right)$

相关指标的计算结果如表1所示。

表1预测结果分析表(1)

从表1可知：基于小波分析法和改进时间序列分析法的本发明具备出色算法特性，其超前一步、三步、五步、十步预测平均相对误差分别仅为0.30％，0.75％，1.15％，1.65％；比较图5-图6，在跳跃采样点预测上，本发明比传统时间序列分析法更优，传统时间序列分析法预测存在延时现象。

为进一步提高上述具体实施方式所得到的超前一步预测精度，引入卡尔曼滤波法与时间序列分析法混合建模来改善。

一般的线性离散***可表示为：

X(t+1)＝Φ(t+1，t)X(t)+Γ(t+1，t)w(t)     (8)

Z(t+1)＝H(t+1)X(t+1)+v(t+1)               (9)

式中：X(t)为n维状态向量；Z(t)为m维观测向量；w(t)为p维***噪声向量；v(t)为m维测量噪声向量；Φ(t+1，t)为从t时刻到t+1时刻的状态转移矩阵；Γ(t+1，k)为从t时刻到t+1时刻的激励转移矩阵；H(t+1)为t+1时刻的预测输出转移矩阵。式(8)称为状态方程，式(9)称为测量方程。

卡尔曼滤波法的核心优势是能够根据最新的测量值修正前一时刻的估计值，具有动态加权修正的特性。而要实现卡尔曼预测，必须推导出正确的状态方程和测量方程。

设X₁(t)＝X_1t(t)，X₂(t)＝X_1t(t+1)，…，X₇(t)＝X_1t(t+6)将式(5)改写成如下形式：

X_1t(t+1)＝2.7423X_1t(t)-2.9532X_1t(t-1)

+1.8939X_1t(t-2)-1.3428X_1t(t-3)           (10)

+1.1826X_1t(t-4)-0.7913X_1t(t-5)

+0.2685X_1t(t-6)+a_t

则有

X₁(t+1)＝2.7423X₁(t)-2.9532X₂(t)

+1.8939X₃(t)-1.3428X₄(t)

                            (11)

+1.1826X₅(t)-0.7913X₆(t)

+0.2685X₇(t)+a_t

由于X₂(t)＝X₁(t+1)，X₃(t)＝X₂(t+1)，…，X₇(t+1)＝X₆(t)，则可得：

$\left[\begin{array}{c}{X}_1(t+1)\\ X_2(t+1)\\ X_3(t+1)\\ X_4(t+1)\\ X_5(t+1)\\ X_6(t+1)\\ X_7(t+1)\end{array}\right]=\left[\begin{array}{ccccccc}2.7423& -2.9532& 1.8939& -1.3428& 1.1826& -0.7913& 0.2685\\ 1& 0& 0& 0& 0& 0& 0\\ 0& 1& 0& 0& 0& 0& 0\\ 0& 0& 1& 0& 0& 0& 0\\ 0& 0& 0& 1& 0& 0& 0\\ 0& 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 0& 1& 0\end{array}\right]\×\left[\begin{array}{c}{X}_1\left(t\right)\\ X_2\left(t\right)\\ X_3\left(t\right)\\ X_4\left(t\right)\\ X_5\left(t\right)\\ X_6\left(t\right)\\ X_7\left(t\right)\end{array}\right]+\left[\begin{array}{c}1\\ 0\\ 0\\ 0\\ 0\\ 0\\ 0\end{array}\right]a(t+1)---\left(12\right)$

易知，Z(k+1)＝X(k+1)+v(k+1)，其中，v(k+1)为测量附加噪声，为建模方便，可以假设为白噪声，则测量方程为：

Z(t+1)＝[1，0，0，0，0，0，0]×[X₁(t+1)，X₂(t+1)，X₃(t+1)，

                                                          (13)

X₄(t+1)，X₅(t+1)，X₆(t+1)，X₇(t+1)]^T+v(t+1)

将式(8)、(9)与式(12)、(13)比较，即可得到Φ(t+1，t)、Γ(t+1，t)、H(t+1)的具体取值。考虑到收敛的速度和参考工程习惯，取初始值如下：X(0|0)＝[0]、P(0|0)＝10I。取R(t)＝[1](t＝1，2，3，...)，

Q(t)＝[1](t＝1，2，3，...)。应用卡尔曼预测递推方程进行风速超前1步预测，如图7所示。

同理，运用方程(4)-(7)对图7预测结果进行评价，相关计算结果如表2所示。从表2可知，超前一步预测的效果得到了明显的提高，其优化效果非常显著。由此可知，在超前一步预测方面，本发明明显优于传统的时间序列分析法。

表2预测结果分析表(2)

Claims (8)

1.一种铁路沿线风速预测的方法，用于铁路沿线风速的短期预测，其特征在于所述方法包括如下步骤：

(1)选择小波分析法的分解与重构算法对拟预测风速序列进行多层分解与重构计算，实现将原始非平稳风速序列分解为多层平稳的风速序列；

(2)利用改进时间序列分析法对分解重构后的各层风速序列分别建立相应的时序预测模型；

(3)利用各层所建的对应预测模型对分解层风速序列进行超前多步预测计算，然后对各层风速预测值加权计算，获取原始风速序列的相关预测值。

2.根据权利要求1所述的铁路沿线风速预测的方法，其特征在于所述步骤(2)之后还包括如下步骤：

通过改进时间序列分析和卡尔曼滤波法混合建立预测精度模型，进一步提高改进时间序列法获得的超前一步预测精度。

3.根据权利要求1所述的铁路沿线风速预测的方法，其特征在于所述原始风速序列是在铁路沿线上设置测风站，获得该测风站在某一短时间内的风速与时间关系序列。

4.根据权利要求1所述的铁路沿线风速预测的方法，其特征在于所述步骤(1)包括将原始风速序列v(t)(t＝1，2，3，…)经小波分析法Mallat塔式算法进行风速信号分解与重构，分别获得高频风速序列vⁿ(t)(n＝1，2，3，…)和低频风速序列ω(t)。

5.根据权利要求4所述的铁路沿线风速预测的方法，其特征在于所述步骤(2)包括：

对高频风速序列vⁿ(t)(n＝1，2，3，…)、低频风速序列ω(t)运用改进时间序列法建立时序预测模型，计算出各分解层的超前k步高频风速序列预测值 

低频风速序列预测值 

6.根据权利要求2所述的铁路沿线风速预测的方法，其特征在于所述预测精度模型的建立包括：

将原始风速序列v(t)(t＝1，2，3，…)经小波分析法Mallat塔式算法进行风速信号分解与重构，分别获得高频风速序列vⁿ(t)(n＝1，2，3，…)和低频风速序列ω(t)；

对高频风速序列vⁿ(t)(n＝1，2，3，…)、低频风速序列ω(t)运用卡尔曼时间序列法建立预测精度模型，计算出各分解层的超前1步高频风速序列预测值 

低频风速序列预测值 

对各分解层的超前k步高频风速序列预测值 

低频风速序列预测值 进行加权计算，输出最终风速预测值 

7.根据权利要求1所述的铁路沿线风速预测的方法，其特征在于所述步骤(1)包括：选用具有时频紧支撑和高正则性的Daubechies 6小波对原始风速序列v(t)进行信号分解，分解深度n为3-6；选用Mallat塔式算法对分解后的风速序列在不同尺度上进行信号重构。

8.根据权利要求7所述的铁路沿线风速预测的方法，其特征在于所述步骤(2)包括：

以原始风速序列v(t)进行信号分解后的某一层低频或高频分量序列为建模样本序列，取该层序列前多个数据完成滚动式时间序列法建模，并求出相对于采样点超前三步的预测风速值。 

Images