艦船運動姿態極短期預報算法研究

吳爽，焦淑紅，任慧龍

1. 哈爾濱工程大學 信息與通信工程學院，黑龍江 哈爾濱 150001 2. 哈爾濱工程大學 船舶工程學院，黑龍江 哈爾濱 150001

艦船在航行過程中，因受到外加作用會產生晃動，在海浪較大情況下，表現更加明顯。艦船的晃動會對需要在艦船進行一系列特殊作業(如：艦載機起降作業中指導艦載機安全起降、導彈發射、風浪中航行的控制等[1])產生巨大影響，甚至危害生命財產安全。若提前對艦船運動姿態進行預測，可以為決策者提供決策依據，減少安全事故的發生，對航海事業具有重要意義。

目前,國內外對艦船運動姿態極短期預報非常重視并進行了許多研究。近年來，隨著混沌系統理論的快速發展，混沌預報方法被引入到船舶運動姿態的預報中來。文獻[2]對于艦船運動姿態序列的混沌特性進行了分析，說明艦船運動姿態序列屬于低維混沌序列。

已有研究表明，很多低維混沌時間序列可用二階 Volterra 自適應濾波器進行較為精確的自適應預報[3,4]。文獻[5]利用最大 Lyapunov 指數法對船舶運動姿態時間序列樣本進行了混沌性判別，針對船舶運動姿態的非線性和不確定性與混沌特性的緊密關系，給出了船舶運動姿態混沌時間序列的二階 Volterra 自適應預報模型。文獻[6]提出基于Volterra級數的RLS[7-10]自適應算法，驗證其有較好的收斂性并在預測低維混沌系統取得很好的效果。但是，該方法只是對數據進行一次性處理，沒有討論其對于數據流的預測效果如何。本文在Volterra 濾波器與 RLS 自適應算法的基礎上，提出一種基于小波變換的RLS的Volterra級數核估計自適應算法，應用于艦船運動姿態數據流的極短期實時預報研究中。

1 理論基礎

1.1 RLS算法原理

RLS算法全稱是遞推最小二乘算法，它是一種最小二乘自適應橫向濾波器的遞推算法，它是由n-1時刻濾波器抽頭權向量的最小二乘軌跡來遞推n時刻權向量的最新估計。

(1)

式中：

φ(t)=[φ1(t),φ2(t),…,φn(t)]T

X(i)=[x(i),x(i-1),…,x(i-n+1)]T

估計誤差為：

設i<0數據為零， 則有：

最小二乘的性能指標是使ξ(t)最小，式中β(t,i)代表遺忘因子，0<β(t,i)<1,i=1,2,…,t,遺忘因子使建模過程對非平穩情況下數據統計特性的變動更具有適應性。通常取指數式因子作為遺忘因子，即

β(t,i)=λt-i(i=1,2,…,t)

式中λ為趨近于1的正常數，在取指數式因子為遺忘因子時有：

(2)

式中：n×n相關陣N(t)和n×1互相關M(t)分別為：

(3)

將式(3)中i=t的項單獨分開，則有：

λN(t-1)+X(t)XT(t)

同理有：

M(t)=λM(t-1)+X(t)XT(t)

引理若A、B為兩個n×n維正定陣，且滿足：

A=B-1+CD-1CT

式中：D為另一個m×m正定陣;C為n×m陣，則有：

A-1=B-BC(D+CTBC)-1CTB

令N(t)=A,λN(t-1)=B-1,X(t)=C,1=D，經推導可得：

其中a(t)稱作新息，定義為

RLS算法如下：

設定初始值為

計算t=1,2,…時

Z(t)=λ-1P(t-1)X(t)

K(t)=[1+XT(t)Z(t)]-1Z(t)

P(t)=λ-1P(t-1)-K(t)ZT(t)

(4)

定義Fn(t)=λFn(t-1)+ηn(t)fn(t)，后驗估計誤差定義為ε(t)=y(t)-φT(t)X(t)，則最小二乘估計的加權誤差平方和為：

(5)

由式(4)以及式(5)可得:

λξmin(t-1)+a(t)ε(t)

1.2 RLS算法應用于Volterra核估計原理

根據RLS算法的收斂性可知，當輸入向量相關矩陣的特征值擴展很大， RLS算法使得收斂速度也得以保證，并且RLS算法通過選取合適的自適應濾波器的系數，來保證輸出信號y(t)與期望信號盡可能地匹配。

RLS算法中確定目標函數為

(6)

式中：e(i)表示i時刻的輸出誤差，d(t)表示輸入信號為U(i)向量時的期望輸出；并且H(t)、U(t)分別表示輸入向量及自適應濾波器的系數向量；λ為指數加權因子，應該在0<λ<1范圍內進行選擇。

定義系統的輸入矩陣為

P(t)=[U(m),U(m+1),…,U(t-1)]T

式中，U(t)與式(6)中的U(t)相同，m代表二階Volterra濾波器的輸入項數。則輸出為

Y(t)=[y(m),y(m+1),…,y(t-1)]T

式中y(i)(i=m,m+1,…,t-1)代表系統的實際輸出，且滿足y(i)=HT(i)U(i)。

利用傳統最小二乘法求解Volterra級數核可得

直接利用Volterra 濾波器對艦船運動姿態進行實時在線預報時，隨著新數據的不斷獲取，P(t)、Y(t)維數將不斷增大，勢必會耗費大量存儲空間，為此，使用遞推最小二乘法對Volterra級數核進行估計。這樣在預報過程中P(t)維數是確定的，可以減少了數據對存儲空間的占用，其具體步驟總結如下：

設t+1時刻系統的輸入矩陣為P(t+1)，對應輸出向量為Y(t+1)，輸入輸出滿足：

P(t+1)=[U(m),U(m+1),…,U(t+1)]T

Y(t+1)=[y(m),y(m+1),…,y(t+1)]T

記：

φt=[PT(t)P(t)]-1

若增加輸入數據x(t+1)后可得：

φt+1=[PT(t+1)P(t+1)]-1=

[P-1(t+1)+U(t+1)UT(t+1)]-1

令

可得：

φt+1=[P-1(t+1)+U(t+1)UT(t+1)]-1

綜上，可對Volterra RLS濾波器的核估計算法總結如下：

依次令t=m+1,m+2,…,m+N,N為樣本數據個數。則有：

K(t+1)=φtU(t+1)/[1+UT(t+1)φtU(t+1)]

φt+1=φt-K(t+1)UT(t+1)φt

φm=[PT(m)P(m)]-1

1.3 小波去噪細節選擇

文獻[11]中對小波的定義及小波降噪過程已進行較詳細的闡述，本文不再贅述。小波降噪中母小波的選擇較為重要。而母小波的選擇需要考慮到母小波的數學特性及實際應用情景，以保證精度需求。

其中，母小波的數學特性主要包括正交性、正則性、對稱性、高階消失矩、緊支性[12]。考慮姿態預測需要母小波具備上述數學特性，故最終選擇可以很好滿足這些數學特性的Daubechies小波函數，即dbN系列小波，經試驗驗證，該系列中db6小波在分解層數為6時降噪效果最優。

1.4 基于小波去噪的Volterra-RLS實時預報方法

對于Volterra +RLS自適應預報模型，其預報函數是在擬合混沌序列的基礎上獲得的，其中，運動姿態序列本身就不可避免的噪聲，故經過一系列變換之后，使其得到的預報值包含更多的噪聲，使得預報精度不僅要受到模型的影響，還受到噪聲影響。同時，Volterra +RLS算法只能對”靜”的姿態序列進行單次預報，還不具備對運動姿態數據流進行實時預報的能力。

運動姿態時間序列屬于非平穩時間序列，對于分平穩時間序列的降噪方法有很多，但是較為常用的經驗模態分解(empirical mode decomposition，EMD)、EMD的優化方法集合經驗模態分解(ensemble empirical mode decomposition，EEMD)及小波變換。由于EEMD是對EMD的優化形式，故本文這里主要討論EEMD和小波降噪方法。EEMD是對序列添加白噪聲并且進行多次EMD分解，所以其降噪需要時間較多，本文需要進行實時預測，對于算法的效率要求較高，故對于降噪本文最終選取小波閾值降噪方法。

同時考慮Volterra +RLS自適應預報模型存在的問題，本文最終采用的運動姿態時間序列數據流預測的整體框架如圖1所示。傳感器采集到的運動姿態數據數據流進入工作站，首先利用滑動窗口技術進行處理，生成概要數據結構，然后概要數據結構經小波分解分解獲得各尺度系數，經閾值處理去除噪聲，再經Volterra+RLS算法進行預報，得到最終預報值。

圖1 運動姿態時間序列數據流預測框架

對于概要數據結構通過滑動窗口獲得，對于窗口長度的選擇，綜合考慮算法時間復雜度、預報精度及計算機資源占用等問題后，選用窗口長度為400，預報時長為15 s。

2 仿真

對于艦船運動來說，橫搖，又稱側滾角，指浸于水中的物體繞最長延伸方向或波浪入射方向的水平軸的旋轉振蕩運動，即以船舶重心所在的前后軸線(縱軸線)為中心的回轉搖晃。船舶在海上最容易發生橫搖且搖擺幅度最大，故本文主要研究橫搖運動姿態。

本文將某船模總長12.48 m、船寬1.568 m、設計吃水0.404 m、排水量4.672 t、航速6 kn，在浪向45°、浪向90°工況中20 min橫搖運動姿態作為試驗數據。

每次選用400個數據作為建模樣本，對Volterra+RLS、EEMD+Volterra+RLS及小波+Volterra+RLS等3種預報模型進行預報仿真分析。時序序列實時預報中常用均方根誤差RMSE對預報結果進行分析，故本文也以均方根誤差RMSE對精度進行評估。

1)預報精度

先以浪向45°橫搖為例，利用3種預報模型進行預報起始點為2 182 s，建模數據長度為400，預報時長為15 s的姿態運動預報，預報結果如圖2。從圖2中可以看出在使用Volterra+RLS對運動姿態序列預報前使用EEMD、小波閾值降噪使得預報精度較Volterra+RLS模型預報精度均有所提高，并且小波閾值降噪效果更好。

為了驗證在不同工況下，小波閾值+Volterra+RLS的適用性，再選用90°工況下部分橫搖試驗數據進行預報仿真，對預報誤差記錄如表1所示。

圖2 3種模型預報15 s預報結果

預報模型預報起始點/sRMSE耗時/sVolterra+RLSEEMD+Volterra+RLS小波+Volterra+RLS6970.0690.0238970.2310.0241 2970.4220.0241 4240.1270.02321820.2940.0246970.00911.4468970.22111.4781 2970.19911.3811 4240.03611.5072 1820.16611.2106970.0030.1518970.0630.1601 2970.1590.1551 4240.0170.1522 1820.1540.149

從表1中可以看出，EEMD+Volterra+RLS和小波+Volterra+RLS對于Volterra+RLS來說預報精度都有所提高，但小波+Volterra+RLS的預報精度更高，這樣再次對圖2進行驗證。

2)時間復雜度

根據表1所示，可以分析得出：

a)對于相同的橫搖運動姿態序列，預報時長相同時Volterra+RLS算法多次進行預報耗時基本穩定，而后兩種預報方法耗時會產生一定范圍的波動。

b)由于EEMD需要進行多次EMD分解，并尋找符合條件的IMF分量，所以利用EEMD進行降噪的Volterra+RLS預報算法耗時最長。

c)通過小波降噪優化的Volterra+RLS預報算法耗時可以發現，小波分解需要的時間遠小于EEMD分解所需要的時間，盡管小波降噪耗時比Volterra+RLS要長，但預報15 s需要時間穩定在0.2 s，還是完全可以滿足預報需求的。

3 結論

本文提出一種基于小波去噪的Volterra+RLS實時預報方法應用到數據量大、流速快的運動姿態數據流 的極短期實時預報中。首先對模型中涉及到的基本理論進行介紹，考慮到噪聲對精度影響及基于RLS的Volterra核估計算法只能對橫搖運動姿態進行單次預報，不能很好地處理動態的橫搖運動姿態數據流等問題，提出使用對于時序數據流常用的滑動窗口方法，對實時數據流進行概要數據結構的獲取，然后在概要數據結構上進行小波閾值處理及結合Volterra+RLS方法進行預測。對于算法驗證工作，則采用不同工況下橫搖運動姿態試驗數據對本文提出的模型與EEMD+Volterra+RLS及Volterra+RLS模型進行仿真，并從預報精度和時間復雜度兩個方面進行分析討論，經討論得出，本文提出的基于小波去噪的Volterra+RLS實時預報方法在預測精度和耗時方面較其他兩種方法都具有明顯優勢。