基于二阶自适应Volterra级数的船舶运动极短期预报研究

翁震平，顾 民，刘长德

（中国船舶科学研究中心，江苏 无锡 214082）

基于二阶自适应Volterra级数的船舶运动极短期预报研究

翁震平，顾 民，刘长德

（中国船舶科学研究中心，江苏 无锡 214082）

针对随机海浪作用下船舶运动的非平稳、非线性特性，文章提出了基于卡尔曼（Kalman）滤波原理的非线性二阶Volterra级数自适应预报模型。通过把Volterra级数核向量作为状态向量，利用随机游动模型建立系统的状态方程，一步Volterra级数预报模型作为系统的观测方程，从而进一步提高了Volterra级数模型的核估计的收敛速度。同时验证了利用AIC准则对Volterra级数预报模型定阶的可行性，通过迭代法实现了自适应多步预报。仿真结果表明文中提出的基于Kalman滤波算法的自适应预报模型应用于船舶运动极短期预报是可行的，该方法在理论和工程应用方面具有重要的意义。

卡尔曼滤波；二阶Volterra级数；AIC准则；船舶运动预报

1 引 言

船舶在航行过程中，由于受海浪、海风等环境因素影响，不可避免地会产生摇摆，尤其是在恶劣的海况下对船舶的海上作业会造成很大的隐患。如果能提前几秒或十几秒预报出船舶的运动姿态或其运动趋势，可以为重要设备的作业提供指导，从而可有效提高作业安全性，并在特定条件下可改善船舶的耐波性和适航性。

目前，国内外对船舶运动姿态建模预报都非常重视并展开了许多研究。文献[1]对国内外研究现状进行了归纳总结，其中利用时间序列分析法[2-3]对船舶运动姿态进行极短期预报越来越受到重视，这种方法的最大优点是无需知道海浪的任何先验信息和船舶航行姿态的状态方程，仅仅利用历史数据寻求规律进行预报。这种方法假设船舶在海浪中的运动姿态为一平稳的窄带随机过程，从而利用线性回归（AR）模型拟合这一过程，但是实际情况下的船舶运动是非线性的，文献[4-7]利用神经网络进行非线性建模取得了一定成效，但对输入层以及隐层节点的选取往往通过试凑法，不宜处理时变系统，难以做到对海况实时变化的自适应。因此本文提出采用非线性二阶Volterra级数预报模型，结合AIC准则实现了船舶运动的自适应多步预报。

在非线性自适应预报中，基于Volterra级数的模型已经得到了广泛的应用。Volterra级数理论之所以具有如此大的吸引力，原因在于它与幂级数有着天然的密切联系和相似之处且具有鲜明的物理意义，为解决非线性实际问题提供了强有力的方法和工具。张家树、肖先赐等人[8-9]提出了多种非线性自适应预报算法来提高预报性能。其中归一化LMS算法[8]和递推最小二乘（RLS）算法[9-10]应用最为广泛，文献[8]中归一化LMS（NLMS）算法的特点是运算简单，无须计算相关函数，但收敛速度较慢，且在对时间序列递推过程中，收敛参数μ的选取不容易掌握，对数据的相关阵特征值分布也比较敏感，不利于实现自适应在线预报。文献[9-10]递推最小二乘（RLS）算法能实现快速收敛，该算法旨在使期望输出信号与模型输出之差的平方和最小，有效地提高了算法的鲁棒性，其收敛速度明显高于NLMS算法，但RLS算法得到的数据自适应滤波器是对每一组输入数据而言的，未能体现信息的完整性。本文提出利用自适应Kalman滤波算法对Volterra级数核进行估计，并结合AIC准则实现了自适应船舶运动极短期预报，自适应Kalman滤波算法的特点类似于RLS算法，但Kalman滤波算法是从总体平均导出的在均方意义上最优的滤波器，且其收敛速度更快，对数据相关阵的特征值分布具有鲁棒性，能较好地实现多步预报。

2 自适应Volterra级数预报模型

2.1 Volterra级数模型

设非线性离散动力系统的输入为XT(n)=｛x( n- 1 ),x( n- 2 ),…,x( n- p )｝(n=p+1,p+2,…,N )，输出为(n+ 1 )，则非线性系统函数的Volterra级数展开式为

式中：p为滤波器长度，hp(m1,m2,…,mp)为p阶Volterra核。

在实际应用中，这种无穷级数展开式难以实现，必须采用有限次求和的方式。最常用的是二阶截断求和形式，则用于时间序列预报的Volterra自适应滤波器的形式为

2.2 多步迭代预报

定义Volterra级数滤波器的状态扩展的输入矢量为系数向量为则（2）式所示非线性Volterra自适应滤波器可描述成如图 1 所示结构。 （2）式由（3）、（4）式可以表示为

图1 Volterra自适应滤波器结构Fig.1 Volterra adaptive filter

2.3 预报原理

假设所用建模数据个数为N，预报数据的个数为M。则实时递推预报原理如图2所示。

Step1.利用1～N个数据建立模型，预报未来N+1～N+M个数据。

Step2.获得第N+1个数据后，利用2～N+1个数据重新建立模型，预报N+2～N+M+1个数据。

此后一直重复这个过程，即获得新数据后，建模数据个数始终保持为N，预报未来M个数据。

图2 递推预报原理示意图Fig.2 Recursive prediction principle

2.4 滤波器长度p的确定

实际应用中，我们采用AIC准则[2]（Akaike information criterion）来确定模型的阶数。该准则是一种最佳准则函数法，即确定出一个准则函数，该函数既要考虑用某一个模型拟和对原始数据的接近程度，同时又考虑模型中所含待定参数的个数，具体算法如下：

设一个最大阶数M，一般我们取N1=N2，M=N/3（N为建模时所用数据个数），Sp()N 为预报模型（6）-（8）式的残差平方和,则

3 基于自适应Kalman滤波算法的Volterra级数核估计

3.1 Kalman滤波原理及状态估计[11]

设以z(t)表示t时刻离散线性系统状态m维向量，以y(t)表示t时刻系统n维观测数据向量，则系统的状态空间模型如下：

状态方程：

观测方程：

其中，φ( t+ 1 ,t)为系统的m×m维状态转移阵，v1(t)为m×1过程噪声向量，它是零均值的白噪声过程，v1(t)的相关矩阵为Q1(t)。c(t)为已知的n×m维观测矩阵，v2(t)为n×1维观测噪声向量，v2(t)的相关矩阵为Q2(t)。其中v1(t)与v2(t)是统计独立的。

新息过程：

Kalman增益：

一步预报：

状态估计：

状态滤波误差相关阵：

状态预报误差相关阵：

3.2 Volterra级数核的状态空间模型的建立

利用自适应Kalman滤波算法对Volterra级数核进行递推估计时，因为模型参数随时间推移而发生变化，因此我们把参数估计过程看作非平稳过程，引入过程噪声v()n，假定Volterra级数核的状态方程是随机游动模型[12]：

其中，v(n)是零均值平稳随机过程，相关阵为Q(n)=E[v( n ) vT(n)]=qI（q 为标量，是 v(n)的方差）。显然这里的转移矩阵为单位阵。

观测方程为

其中，U(n)=[1,x(n),x( n- 1 ),…,x( n- p +1),x2(n),x(n) x( n- 1 ),…,x2(n- p +1)]T为量测矩阵，ε(n)均值为零，方差为，其中 ｛v( n )｝和 ｛ε (n )｝相互独立。

由（21）式可以看出，观测方程实际为一步Volterra级数预报模型。

3.3 基于自适应Kalman滤波算法的Volterra级数核估计

设（20）、（21）式确立的状态空间模型中状态滤波误差定义为e(n)，状态预报误差相关阵为K（n+1,n），则由Kalman滤波算法递推公式可得：

方程（22）-（26）式构成了自适应Kalman滤波算法对Volterra级数核的状态估计过程。

则对应输入向量XT(p+ 1 ),XT(p+ 2 ),…,XT(N)的希望输出为x( p+ 1 ),x( p+ 2 ),…,x(N ) ；对应的Volterra级数的状态扩展为U(p),U( p+ 1 ),…,U(N)。

基于Kalman滤波算法的Volterra级数核估计的具体推导过程如下：

令p=1,2,…,M/3；计算t=p,p+2,…,N-1时的

由（27）、（31）式我们可以看到在该算法过程中需预先假定:

（2） v(t)的方差q。

4 实例预报结果与分析

4.1 数据预处理

数据预处理的目的是使时间序列变为零均值，使其变为近似平稳序列。

4.2 预报结果误差分析

本文分别采用浪向为150°和180°工况下的模型试验数据转化为实船后的数据进行验证分析，采样频率2Hz，建模数据个数为200，采用滑动窗法，连续递推预报次数为300（即连续预报时间长度为150s），并对这300次多步预报在未来tk=2s，4s，6s，8s时刻的预报误差进行统计分析。

图3、图4分别给出了顶浪状态下航速为18kns、24kns的纵摇数据在未来tk=6s时刻的预报曲线与实际曲线比较。

由图3、图4可以观察到，在顶浪状态高航速和低航速的纵摇预报值曲线与实测值曲线都比较接近。由表1可得，对应于tk=2s，4s，6s，8s时刻，低航速下的纵摇运动预报精度EPA要高于高航速状态下的运动预报精度EPA，且随着预报时间的增长，低航速下的纵摇预报精度明显高于高航速状态下的运动预报精度。

表1 纵摇预报误差统计表（顶浪）Tab.1 Prediction errors statistics of pitch motion(head sea)

表示预报精度，它反映了预报曲线与真实曲线之间的形状相似程度。EPA越接近1，表明预报曲线与实际数据曲线相似程度越高。其中，n为预报的数据个数，x(t)、(t)分别为t时刻的真实值与预报值，(t)、xm(t)分别为n个预报值和真实值的平均期望；σ、σx分别代表预报序列和期望序列的标准差。

为进一步验证本文自适应预报模型的可行性，下面利用浪向150°、航速为18kns工况下的横摇、纵摇数据对其进一步验证分析。

对图5、图6预报值曲线与实测值曲线比较以及表2中对横摇、纵摇误差统计可以看出，本文方法对该船在浪向150°、航速18kns工况下的横摇、纵摇数据同样取得了良好的预报效果。纵摇预报精度EPA要明显好于横摇预报精度EPA。

4.3 自适应性与收敛性

利用本文方法，以顶浪、航速18kns的船舶运动纵摇数据考察模型自适应性与收敛性。

（1） 自适应性

在预报过程中，本文采用滑动窗法建模，通过引入AIC准则，结合Kalman滤波算法，模型根据输入数据自动定阶同时完成Volterra核估计，完全实现了自适应建模与预报。随着新数据的到来，模型阶数与参数自适应变化，使得一步预报模型始终保持最佳，从而极大地保证了预报精度。

表2 横摇、纵摇误差统计表（浪向150°，航速18kns）Tab.2 Prediction errors statistics on pitch motion and roll motion(wave direction 150°,forward speed 18kns)

（2） 收敛性

由图7可以得出，对于连续300次滑动建模预报，利用AIC准则对每次建模进行定阶，模型阶数变动较大，因此，需要对每次建模进行重新定阶。图8-10分别为利用第11-210个采样点、第151-250个采样点、第221-420个采样点建模时的Volterra核收敛曲线，通过图7可以得出，对应于这三次建模的阶数分别为3、6、7，其核向量收敛曲线分别如图8-10所示。

图8中h1(0)，h1(1)，h1(2)所示曲线为Volterra级数模型中对应一次项的核收敛曲线，其余为二次项对应核收敛曲线，从图8中可以看出，在建模过程中一次项核起主要作用，二次项基本趋于零，当迭代次数为200时，Volterra核趋于收敛。图9、图10只画出了一次项核收敛曲线，同理也可以看出当迭代次数为200时，Volterra核趋于收敛，因此验证了本文选取建模数据为200是可行的。

5 结 论

本文基于二阶Volterra级数模型，提出了利用Kalman滤波算法对Volterra级数核估计的一种新的辨识方法，并将其应用到船舶运动极短期预报当中，得出了以下结论:

（1）在顶浪状态下高航速和低航速的纵摇预报效果都比较理想，但低航速下的纵摇运动预报精度要高于高航速状态下的纵摇运动预报精度；斜浪状态下的纵摇运动预报效果明显好于横摇运动预报效果。

（2）利用AIC准则对非线性Volterra级数定阶是可行的。通过采用滑动窗法与迭代法，实现了船舶运动自适应多步预报。

（3）通过Volterra核收敛曲线验证了Kalman滤波算法的快速收敛性，同时通过收敛性能验证了所选建模数据个数的可行性，并作为建模数据个数选取的规范，从而在Volterra核估计中，避免了不必要的计算，同时又保证了预报精度，为今后实时在线预报提供了可行性依据。

[1]赵希人,彭秀艳.舰船运动极短期预报建模的研究现状[J].船舶工程,2002(3):4-8.

[2]彭秀艳,赵希人,魏纳新等.大型舰船姿态运动极短期预报的一种AR算法[J].船舶工程,2001(5):6-8.

[3]彭秀艳,赵希人,高奇峰.船舶姿态运动实时预报算法研究[J].系统仿真学报,2007(1):267-269.

[4]徐 培,金鸿章,王科俊等.一种新型的船舶横摇运动实时预报方法[J].中国造船,2002,43(1):71-73.

[5]谢美萍,赵希人.基于小波网的船舶运动极短期建模预报[J].船舶力学,2002,6(2):16-20.

[6]Richards R A.Application of multipl artificial intelligence techniques for an aircraft carrier landing decision support tool[C].IEEE International Conference on Plasma Science,2002,1:7-11.

[7]Shen Yan,Xie Meiping.The extreme short prediction of ship pitch based on diagonal recurrent neural network[J].Journal of Marine Science and Application,2005,4(2):56-60.

[8]张家树,肖先赐.用于混沌时间序列的Volterra自适应预测[J].物理学报,2000,49(3):403-408.

[9]张家树,肖先赐.用于混沌时间序列自适应预测的一种少参数二阶Volterra滤波器[J].物理学报,2001,50(07):1248-1254.

[10]孔祥玉,韩崇昭,马红光等.一种总体最小二乘算法及在Volterra滤波器中的应用[J].西安交通大学学报,2004,38(4):340-342.

[11]Altunkaynak A.Adaptive estimation of wave parameters by Geno-Kalman filtering[J].Ocean Engineering(2008),doi:10.1016/j.oceaneng.2008.03.010:1-7.

[12]Louka P,Galanis G,Siebert N,et al.Improvements in wind speed forecasts for wind power prediction purposes using Kalman filtering[J].Wind Eng Ind.Aero-dyn.(2008).doi.:10.1016/j.jwe-ia.2008.03.013:5-15.

Extreme short-term prediction of ship motion based on second-order adaptive Volterra series

WENG Zhen-ping,GU Min,LIU Chang-de
(China Ship Scientific Research Center,Wuxi 214082,China)

Aiming at non-stationary and non-linear nature of the ship motion,the application of Kalman filtering in non-linear second-order Volterra series model is proposed.The system’s state space model is established by taking the Volterra series kernals vector as the state vector,the observation equation can be obtained by Volterra series model.Then the convergence rate is improved further and carried out adaptive ship motion prediction combining AIC criterion.The multi-step simulation results obtained show that non-linear second-order Volterra series model based Kalman filtering is feasible in the predicting skill.The approach has vital significance in the theory and the engineering application.

Kalman filtering;second-order Volterra series;AIC criterion;ship motion prediction

U661.32

A

1007-7294（2010）07-0732-09

2010-05-05

翁震平（1958-），男，博士生，中国船舶科学研究中心研究员。