朱曉樂 王 華 符菊梅 陳景鵬 徐臘萍

（1裝備指揮技術(shù)學(xué)院 2中國西昌衛(wèi)星發(fā)射中心）

1 引言

航天科研試驗中，飛行數(shù)據(jù)可以作為故障診斷的一種數(shù)值依據(jù)，需要對運載火箭飛行數(shù)據(jù)進行監(jiān)測，以便及時了解運載火箭的工作狀態(tài)和工作環(huán)境，為指揮控制中心提供決策支持[1]。動力系統(tǒng)作為運載火箭的重要組成部分，其系統(tǒng)參數(shù)變化是否正常，能夠直接反映運載火箭飛行情況，目前，國內(nèi)衛(wèi)星發(fā)射中心的飛行狀態(tài)快速評估系統(tǒng)能夠?qū)\載火箭的飛行結(jié)果進行快速、科學(xué)的評估，如果能夠?qū)︼w行數(shù)據(jù)，尤其動力系統(tǒng)緩變數(shù)據(jù)變化趨勢進行實時預(yù)測，將可提前預(yù)測出潛在的故障趨勢，為指揮控制中心進行決策提供更有力的支持。

在眾多預(yù)測方法中，回歸分析預(yù)測法主要基于因果關(guān)系分析，并且需要大量現(xiàn)有數(shù)據(jù)進行分析，比較適合事后處理，實時性不夠[2]?；疑A(yù)測法要求數(shù)據(jù)資料具有確定性趨勢，靈活性不夠[3]。神經(jīng)網(wǎng)絡(luò)模型預(yù)測法需要大量現(xiàn)有數(shù)據(jù)來訓(xùn)練模型，模型訓(xùn)練辨識完成后才能進行預(yù)測，一旦數(shù)據(jù)序列特征發(fā)生較大變化，就要重新學(xué)習(xí)與建模，靈活性不夠，實時性也不夠[4]。指數(shù)平滑法建模簡單，適用范圍廣，對歷史數(shù)據(jù)的依賴程度比較低，實時性也比較好，但是對數(shù)據(jù)突變趨勢的跟蹤能力比較弱，容易導(dǎo)致預(yù)測結(jié)果滯后[5]。本文采用的ARMA模型是一種時間序列預(yù)測方法，它將預(yù)測對象隨時間變化形成的序列，看作是一個隨機時間序列。其基本思想是：一串隨時間變化而又相互關(guān)聯(lián)的數(shù)字序列，可以用相應(yīng)的模型加以近似描述，通過對相應(yīng)數(shù)學(xué)模型的分析研究，能更本質(zhì)的認識這些動態(tài)數(shù)據(jù)內(nèi)結(jié)構(gòu)和復(fù)雜性，從而達到在最小方差意義下的最佳預(yù)測。

2 ARMA模型的定義

假設(shè)｛xi｝為隨機時間序列，Box-Jenkins模型理論認為xi的取值不僅與其前p步的各值xi-1，xi-2，…，xi-p有關(guān)，而且同前q步的隨機干擾ai-1，ai-2，…，ai-q也有關(guān)，且均為線性關(guān)系，從而得到自回歸移動平均模型 ARMA（p，q）模型如下[6]：

令

則上式可簡記為

其中B為線性推移算子且有如下性質(zhì)：

這一模型就稱作p階自回歸-q階滑動平均混合模型，記為 ARMA（p，q）模型，特殊地，若p=0，稱作純滑動平均模型，記為MA（q）;若q=0，稱作純自回歸模型，記為 AR（p）；若p=q=0，模型退化為Xt=at，即{Xt}為白噪聲序列。

3 ARMA模型辨識和預(yù)測

3.1 辨識模型類別與確定階次。

通過計算時間序列的自相關(guān)函數(shù)和偏相關(guān)函數(shù)，根據(jù)截尾性和拖尾性來確定是采用AR（p）、MA（q）模型或者 ARMA（p，q）模型。屬于 AR（p）過程的時間序列，它的自相關(guān)函數(shù)隨著滯后期k的增加呈現(xiàn)幾何衰減形式（即具有拖尾性），而偏相關(guān)函數(shù)則應(yīng)在k>p時截止為0；而屬于MA（q）過程的時間序列的自相關(guān)函數(shù)在k>p后為0（即具有截尾性），而偏相關(guān)函數(shù)或者呈指數(shù)衰減，或者呈正弦函數(shù)衰減；如果時間序列的樣本自相關(guān)函數(shù)和偏相關(guān)函數(shù)均不截止，但較快收斂到0，則很可能屬于ARMA過程；此外，當自相關(guān)函數(shù)的尾部隨k增加不是趨近于0，而是呈周期起伏，則表明時間序列中含有周期分量。當樣本數(shù)量較小時，樣本自相關(guān)函數(shù)將會偏離實際的自相關(guān)函數(shù)，用它們來識別模型，工作量大且效果不好，因此，采用赤池弘次提出的信息量準則（稱為AIC準則）來判斷模型階次p和q。

定義模型的AIC統(tǒng)計量：

在模型階次確定的情形下，進行模型參數(shù)的估計。模型參數(shù)系數(shù)確定一般分兩步完成：先用矩估計或逆函數(shù)法粗估計，再用粗估計的值作為疊代初值進行最小二乘法精估計。兩步估計完成后，則得到具體ARMA模型。由于ARMA（p，q）模型參數(shù)的估計非常復(fù)雜，論文在進行數(shù)據(jù)仿真的時候直接調(diào)用MATLAB自帶函數(shù)，在此不再介紹參數(shù)估計的數(shù)學(xué)推導(dǎo)[7]。

3.2 模型檢驗和數(shù)據(jù)預(yù)測

數(shù)據(jù)序列通過平穩(wěn)性檢驗，并建立了相應(yīng)的ARMA模型之后，為考核所建模型的優(yōu)劣，一般還需檢驗ARMA模型殘量e1，e2，…，en是不是白噪聲。也就是說，如果殘量經(jīng)檢驗確定是白噪聲序列，則認為模型是合理的，否則應(yīng)當進一步改進模型。

獲得較為滿意的時間序列模型后，采用最小方差線性估計的原則對飛行數(shù)據(jù)進行預(yù)報。用記號表示為時間序列…，Z-1，Z0，Z1，…，Zk，…，Zk+L，其中k≥1，L≥1。若已觀測到Z1，…，Zk的數(shù)值，要估計Z k+L的數(shù)值，稱為在k時刻作L步預(yù)報，Zk+L的估計值記為。

3.3 預(yù)測結(jié)果分析

預(yù)測結(jié)果分析包括評定預(yù)測的精度，以及評價預(yù)測模型的合理性。

目前，對預(yù)測的精度評定主要基于誤差理論，即用平均絕對誤差（Mean Absolute Error，MAE）和均方誤差（Mean Squared Error，MSE）來衡量[8]。

首先，計算預(yù)測誤差et

然后分別計算平均絕對誤差和均方誤差，其計算公式如下。

對預(yù)測模型的合理性評價主要從兩個方面來衡量：預(yù)測精度和預(yù)測跟蹤速度。

4 趨勢預(yù)測實例

4.1 預(yù)測對象

運載火箭飛行數(shù)據(jù)多達幾百甚至上千路，這些數(shù)據(jù)按其頻率初步可以分為緩變參數(shù)（頻率在10Hz以內(nèi)）和速變參數(shù)（頻率遠大于10Hz）兩種。速變參數(shù)變化速度快，不同時刻的數(shù)據(jù)相關(guān)性小，對其進行趨勢預(yù)測意義不大。當然也不是所有緩變參數(shù)都適合進行趨勢預(yù)測，只有那些連續(xù)緩變型參數(shù)，例如供電電壓參數(shù)，壓力參數(shù)，過載參數(shù)等，由于數(shù)據(jù)變化的連續(xù)性，前后時刻的數(shù)據(jù)間存在一定的相關(guān)性，對它們進行預(yù)測有工程實踐價值[1]。實時飛行數(shù)據(jù)中，動力系統(tǒng)壓力參數(shù)的數(shù)據(jù)曲線特征很有代表性，下面以這種類型的參數(shù)為例，進行趨勢預(yù)測仿真。

圖1是一次飛行試驗任務(wù)中采集壓力數(shù)據(jù)繪制成的曲線，圖2是引入一個陡然下降的故障趨勢后的數(shù)據(jù)曲線。

圖1 原始數(shù)據(jù)曲線

圖2 含故障趨勢的數(shù)據(jù)曲線

4.2 預(yù)測步驟

預(yù)測實現(xiàn)步驟：

（1）讀入數(shù)據(jù)；

（2）對原始數(shù)據(jù)進行零均值化，平穩(wěn)化處理，然后進行分析；

（3）分別計算置信度為95%的自相關(guān)函數(shù)和偏相關(guān)函數(shù)，并畫出其自相關(guān)函數(shù)和偏相關(guān)函數(shù)曲線；

（4）由自相關(guān)函數(shù)拖尾性和偏相關(guān)函數(shù)的截尾性，初步判斷ARMA模型的階次；

（5）在步驟（4）的基礎(chǔ)上建立一系列模型，然后根據(jù)AIC準則判別最優(yōu)模型；

（6）對數(shù)據(jù)進行一步預(yù)測，并計算預(yù)測誤差，進行擬和誤差的自相關(guān)檢驗。

4.3 預(yù)測結(jié)果

理論上，當p，q的階次越高時，平均絕對誤差和均方誤差越小，也就意味著預(yù)測精度越好。實際情況是當p，q達到一定階次時預(yù)測精度不再顯著的提高，但計算量卻以指數(shù)上升，這樣必定會影響預(yù)測的跟蹤速度。為了克服這一問題，模型階次的確立采用了AIC準則。

按照4.2中預(yù)測實現(xiàn)步驟對引入故障趨勢后的壓力參數(shù)進行趨勢預(yù)測，由圖4自相關(guān)函數(shù)、圖5偏相關(guān)函數(shù)確定了基本模型參數(shù)ARMA（41，8），然后建立了九個模型參數(shù) ARMA（40，7），ARMA（40，8），ARMA（40，9），ARMA（41，7），ARMA（41，8），ARMA（41，9），ARMA（42，7），ARMA（42，8），ARMA（42，9）。分別計算九個模型參數(shù)所對應(yīng)的AIC值，如表1所示。

表1 不同階次對應(yīng)的AIC值

根據(jù)AIC準則選取了最優(yōu)模型參數(shù)ARMA（41，7）建立了模型：

最終采用上述模型對參數(shù)進行了趨勢預(yù)測，得到預(yù)測曲線圖6，誤差函數(shù)曲線圖7，誤差自相關(guān)系數(shù)曲線圖8。預(yù)測模型的平均絕對誤差MAE=7.7137×e-4，均方誤差 MSE=1.3992×e-6。

4.4 預(yù)測結(jié)果分析

對一步預(yù)測結(jié)果圖6、圖7、圖8進行分析可以看出：

圖3 去除均值后的數(shù)據(jù)曲線

圖4 自相關(guān)函數(shù)

圖5 偏相關(guān)函數(shù)

圖6 預(yù)測曲線

圖7 誤差函數(shù)

圖8 誤差自相關(guān)系數(shù)

（1）模型的跟蹤速度較好，在引入劇烈變化趨勢時，能夠快速響應(yīng)時間序列數(shù)據(jù)的急劇變化，快速跟蹤與自修正能力比較強，使得預(yù)測值曲線與實測值曲線擬合得很好，達到95.77%；

（2）模型的預(yù)測精度比較高，平均絕對誤差MAE=7.7137×e-4，均方誤差 MSE=1.3992×e-6；

（3）模型的預(yù)測誤差函數(shù)和預(yù)測自相關(guān)函數(shù)顯示誤差序列符合白噪聲序列的要求，表明預(yù)測模型是適合的；

（4）當預(yù)測步數(shù)增大時，預(yù)測時間無大影響，但是預(yù)測的平均誤差和均方誤差都變大，不同預(yù)測步數(shù)的誤差結(jié)果如表2所示。

表2 不同預(yù)測步數(shù)的誤差分析表

5 結(jié)束語

論文結(jié)合航天科研試驗，利用ARMA模型對運載火箭飛行數(shù)據(jù)中的動力系統(tǒng)緩變數(shù)據(jù)進行了實時趨勢預(yù)測，數(shù)據(jù)仿真分析表明，ARMA模型預(yù)測方法在預(yù)測精度與預(yù)測跟蹤速度方面都有不錯的效果，對數(shù)據(jù)表現(xiàn)出的故障趨勢，在保證較高預(yù)測精度的同時能夠快速跟蹤預(yù)測，為計算機輔助決策和飛行數(shù)據(jù)的自動判讀提供了一種解決方法。

［1］劉蘊才.導(dǎo)彈衛(wèi)星測控系統(tǒng)工程（下）［M］.北京:國防工業(yè)出版社，1996:138-157.

［2］周紀薌.實用回歸分析方法［M］.上海:上?？茖W(xué)技術(shù)出版社，1990.

［3］袁嘉祖.灰色系統(tǒng)理論及其應(yīng)用［M］.北京:科學(xué)出版社，1991.

［4］汪成亮，宋軍，胡炳權(quán)，張勤.智能神經(jīng)網(wǎng)絡(luò)在時序信號預(yù)測上的應(yīng)用［J］.重慶大學(xué)學(xué)報，2003，26（1）:34-36.

［5］余國浩，蔡遠文.自適應(yīng)指數(shù)平滑法用于遙測數(shù)據(jù)實時趨勢預(yù)測研究［J］.裝備指揮技術(shù)學(xué)院學(xué)報，2007，08:48-51.

［6］韓路躍，杜行檢.基于MATLAB的時間序列建模與預(yù)測［J］.計算機仿真，2005.04:105-107

［7］王沫然.MATLAB與科學(xué)計算［M］.北京:電子工業(yè)出版社，2004.

［8］張有為.預(yù)測的數(shù)學(xué)方法[M].北京:國防工業(yè)出版社，1991:12-14.