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

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

1 引言

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

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

2 ARMA模型的定義

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

令

則上式可簡記為

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

4 趨勢預(yù)測實(shí)例

4.1 預(yù)測對象

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

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

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

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

4.2 預(yù)測步驟

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

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

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

（3）分別計(jì)算置信度為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準(zhǔn)則判別最優(yōu)模型；

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

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

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

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

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

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

最終采用上述模型對參數(shù)進(jìn)行了趨勢預(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進(jìn)行分析可以看出：

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

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

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

圖6 預(yù)測曲線

圖7 誤差函數(shù)

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

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

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

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

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

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

5 結(jié)束語

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

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

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

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

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

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

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

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

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