張軍杰

（上海市巖土地質(zhì)研究院有限公司，上海 200072）

0.引言

極移是指地球的極點在地球表面發(fā)生有規(guī)律的緩慢變化的現(xiàn)象，極移參數(shù)作為地固坐標系和天球坐標系的轉(zhuǎn)換參數(shù)之一，在深空探測、衛(wèi)星定軌和激光測月等領(lǐng)域中均有著廣泛地應用[1]。但是由于復雜的數(shù)據(jù)處理過程，得到的最終極移參數(shù)存在著時間上的滯后性，無法滿足某些要求實時或準實時極移參數(shù)的領(lǐng)域[2]。尤其是隨著北斗二代衛(wèi)星導航系統(tǒng)的逐步拓展，對極移參數(shù)的時效性和精度的需求日益迫切[3]，因此，如何提高極移參數(shù)的短期預報精度成為亟需解決的問題。

近年來，國內(nèi)外學者對于極移參數(shù)預報做了大量有益的研究工作，提出了很多非常實用的模型。神經(jīng)網(wǎng)絡(luò)模型（Artificial Neural Network,ANN）是常見的模型之一。人工神經(jīng)網(wǎng)絡(luò)是由大量處理單元互聯(lián)組成的非線性、自適應信息處理系統(tǒng)，王琪潔等人研究和探索了應用非線性的神經(jīng)網(wǎng)絡(luò)技術(shù)，將大氣角動量時間序列引入到地球自轉(zhuǎn)變化預報中，改善地球自轉(zhuǎn)參數(shù)（ERP）的預報精度，以及應用神經(jīng)網(wǎng)絡(luò)技術(shù)預報El Nino南方濤動（ENSO）事件，結(jié)果證實了神經(jīng)網(wǎng)絡(luò)具有良好的非線性預報能力。張志、廖瑛等人利用傅里葉分析法，獲得了插值基礎(chǔ)序列的周期，驗證了極移參數(shù)基礎(chǔ)序列重采樣的可行性，然后提取插值后的基礎(chǔ)序列數(shù)據(jù)的趨勢項，利用多輸入—單輸出反向傳播（Back Propagation，BP）神經(jīng)網(wǎng)絡(luò)進行建模，預報不同跨度的殘差序列項,最后合并趨勢項和殘差序列得到最終的極移參數(shù)預報。

最小二乘和自回歸預報模型（Least squares and autoregressive,LS+AR）[4,5]也較常用。雷雨等人分別利用遞推、迭代和間隔3種預報方法對地球自轉(zhuǎn)參數(shù)進行預報。預報結(jié)果顯示，這3種方式對日長變化（length of day，LOD）參數(shù)的所有跨度預報的精度相當,遞推方式在極移所有跨度的預報結(jié)果上精度最高，迭代方式預報結(jié)果最差。姚宜斌等人通過對LS+AR模型短期預報殘差的時間序列統(tǒng)計分析，發(fā)現(xiàn)臨期的模型預報殘差具有極高的相關(guān)性，由此提出利用上一期模型預報殘差和經(jīng)驗調(diào)節(jié)矩陣對下一期預報結(jié)果進行修正，從而提高了極移參數(shù)的超短期和短期預報精度。張昊、王琪潔等人利用最小二乘外推與ARIMA（P，1，0）組合模型對極移參數(shù)進行預報，ARIMA的全稱為求和自回歸移動模型，是d階差分運算（對原始時間序列進行差分的階數(shù)）和ARIMA（p,q）模型的結(jié)合，p為AR模型的階數(shù)，q為WA模型的階數(shù)。ARIMA（p,d,q）進行預報的思路是首先對原始極移參數(shù)時間序列進行d階差分，然后運用ARMA（p,q）對差分項進行預報，最后逆運算得到預報值結(jié)果。最小二乘外推與ARIMA(P,1,0)組合模型對于提高超短期極移參數(shù)預報精度有一定幫助。王小輝等人利用LS-SVM模型，將赤道向大氣角動量（AAM）和赤道向海洋角動量（OAM）時間序列引入到極移參數(shù)序列的預報中，改善了極移參數(shù)的預報精度。張昊根據(jù)最小二乘模型中錢德勒周期項的時變特性，提出了基于錢德勒參數(shù)時變修正的CLS模型,并利用CLS+AR模型對極移參數(shù)序列進行預測，CLS+AR模型在極移參數(shù)預報精度上較LS+AR模型有所改善。趙丹寧等人利用經(jīng)驗模態(tài)分解對LS+AR模型進行改進，首先利用經(jīng)驗模態(tài)分解方法對極移序列進行分解,獲得極移的高頻分量和低頻分量；然后采用最小二乘外推模型對極移低頻分量進行擬合，獲得最小二乘擬合殘差；之后利用自回歸模型對極移高頻分量和最小二乘擬合殘差之和進行建模預報；最后將最小二乘模型和自回歸模型外推值相加獲得極移參數(shù)的預報值。結(jié)果表明可以有效改善極移參數(shù)的預報精度。為檢驗各種模型的極移預報精度，維也納理工大學舉行了全球性的EOP預報對比運動（Earth Orientation Parameters Prediction Comparison Campaign，EOP PCC），確定了目前預報精度最高的模型是諧波最小二乘外推+自回歸預報模型，即LS+AR模型。

本文將基于國際上公認為最優(yōu)的LS+AR模型，對其進行改進，以期進一步提高極移參數(shù)的預報精度，達到國際先進水平。

1.LS+AR模型

1.1 LS模型

自1900年國際緯度局利用目視天頂儀觀測極移以來，天文工作者對極移已進行了近百年的研究。根據(jù)長期的觀測和研究發(fā)現(xiàn)，極移參數(shù)存在多個基本的周期項，如，錢德勒項、周年項和半周年項等周期項。在利用最小二乘方法對極移參數(shù)序列進行建模時，首先要確定其固定周期項和趨勢項，在本文中包含趨勢項、錢德勒項、周年項和半周年項等四項。LS模型的具體數(shù)學表達式如式（1）、式（2）所示：

式（1）中，a0為極移長趨勢項；t為極移參數(shù)對應的時間，（單位為年）；a1為趨勢項系數(shù)；a2、a3為錢德勒項系數(shù)，對應的T1為錢德勒周期項的周期時間1.183a；a4、a5為周年項系數(shù)；對應的T2為周年項周期1a；a6、a7為半周年項的系數(shù)，對應的T3為半周年項周期0.5a。極移Y分量對應參數(shù)含義與式（1）中的一致。

LS模型的參數(shù)具體解算方法如式（3）～式（6）：

式（3）～式（6）中，X為最小二乘的模型參數(shù)矩陣；B為模型參數(shù)的系數(shù)矩陣；L為X分量的觀測值序列組成的矩陣。極移Y分量的計算方法與X分量的算法一致。

1.2 AR模型

AR模型是依據(jù)自身過去的時間序列的變化規(guī)律來預測下一時刻的值，要求建模的時間序列是零均值的平穩(wěn)隨機序列。所以建模之前需要對LS模型擬合后的殘差序列進行預處理，通常是先進行一階差分后再進行建模。最終的極移預報值為LS模型的外推值與AR模型預報值之和。AR模型的數(shù)學表達式可表示為隨機序列zt（t＝1，2，…，n）與t時刻以前的規(guī)律性變化和t時刻的白噪聲的關(guān)系：

式（7）中，η1，η2，…，ηp為模型參數(shù)；at為白噪聲序列；p為模型階數(shù)；式（7）中zt稱為p階自回歸模型，簡記為AR（P）。

建立AR模型就必須首先確定模型的階，然后再確定AR模型參數(shù)。定階方法主要有3種，分別是最終預測誤差準則、信息論準則、傳遞函數(shù)準則。理論上這三種方法是等效的，本文采用最終預測誤差準則來確定AR模型的階。

式（9）中，PM為用AR模型擬合zt（t=1，2，3，…，n）序列殘差的均方差，當FPE（M）為最小值時取M作為AR模型的階。定階也可以利用Matlab高階統(tǒng)計量工具箱中的函數(shù)來解決，計算十分方便。

1.3 精度評定

為了便于評估預測精度，本文采用國際上通用的平均絕對誤差（Mean Absolute Error,MAE）作為精度評定標準，數(shù)學公式如式（10）所示：

式（10）中，n為預報次數(shù)；k為預報長度；εi，j為預報值與真值的差值。

2.改進的LS+AR模型

極移是用來描述地極在地球表面位置隨時間而變化的現(xiàn)象，其自身激發(fā)機制十分復雜，包括太陽輻射、月球引力、海洋潮汐和地下水分布，還有季節(jié)的變化引起的降雪變化和大氣流動變化等因素。極移現(xiàn)象激發(fā)因素的不確定性，導致極移分量的周期項和趨勢項具有時變性的特征。為了提高極移參數(shù)的預報精度，本文在獲得基本的極移參數(shù)后，利用基本的趨勢項和德勒項、周年項和半周年項等周期項進行初步擬合，在獲得一次擬合后的殘差項后，利用頻譜分析探測殘差項中隱含的殘差周期項，此時殘差項依然包含比較明顯的未知趨勢項和周期項，利用LS模型建模對殘差項再次進行擬合，獲得第一次殘差項的趨勢項和周期項，然后對二次殘差項進行一階差分，對差分結(jié)果進行AR建模，然后利用AR模型對差分項進行預報，最后根據(jù)LS模型得到的兩次趨勢項和周期項以及AR模型預報的殘差差分項獲得最終的預報極移參數(shù)。具體流程（如圖1所示）：

圖1 極移參數(shù)預報流程

3.基礎(chǔ)序列的選取

極移運動成因復雜，不同長度的基礎(chǔ)序列對短期預報結(jié)果的精度也是不同的，本文首先以2a、4a、6a、8a和10a作為基礎(chǔ)序列長度對極移參數(shù)進行預報，尋找更加合適的基礎(chǔ)序列長度?；A(chǔ)序列長度并非越長越好，主要是因為隨著時間間隔的增加，地球季節(jié)性變化、海洋潮汐、月球引力、太陽輻射等各種偶然物理影響因素也會隨之增加，極移序列之間的相關(guān)性會逐步下降。因此在利用LS+AR模型對極移參數(shù)進行建模時，基礎(chǔ)序列長度不宜過長，這樣不僅可以提高預報效率，還可以保障預報精度的穩(wěn)定性。對比基礎(chǔ)序列長度與預報精度的關(guān)系（如圖2、圖3所示），結(jié)果顯示當基礎(chǔ)序列為4a時預報精度最佳，所以本次預報采用4a長度作為預報極移參數(shù)的基礎(chǔ)序列。

圖2 Xpole參數(shù)2、4、6、8年的預報精度

圖3 Ypole參數(shù)2、4、6、8年的預報精度

4.預報結(jié)果與精度分析

本次預報所用到的數(shù)據(jù)來自國際地球自轉(zhuǎn)和參考系服務(wù)組織（IERS）所發(fā)布的EOP 08C04序列（https://www.iers.org/IERS/EN/DataProducts/EarthOrientationData/eop.html）。該序列包含從1962年1月1日至今的極移參數(shù)的Xpole、Ypole分量，時間間隔為1天。本文自2005年1月1日（MJD：55317）起開始預報，每7天預報一次，一共預報了350期，預報長度為50天。為進一步削弱一次擬合后的殘差項里隱含周期項對預報精度的影響，首先用LS模型對極移參數(shù)進行擬合，采用傅里葉分析法對一次擬合后的殘差項時間序列進行頻譜分析，目的是將時間域中的隱藏波形轉(zhuǎn)變?yōu)轭l率域的頻譜，顯示其隱藏周期。每140天對極移X和Y分量的殘差項的時間序列去趨勢項和功率譜分析，對極移參數(shù)第一次擬合后的殘差結(jié)果的頻譜分析值（如表1所示），在18次的頻譜分析值中，每次頻譜分析的周期值變化差異較大，說明了十分有必要對其進行二次擬合，削弱明顯周期項的影響。

表1 擬合殘差的頻譜分析周期表

利用頻譜分析得到的殘差周期再次對極移殘差序列進行擬合，圖4和圖5中藍色為一次擬合后的殘差序列，紅色為兩次擬合后的殘差序列?？梢悦黠@看出極移分量兩次擬合后的殘差時間序列更加趨向于平穩(wěn)、零均值，不包含明顯的周期項。然后對得到的兩次擬合后的殘差項進行一階差分，利用AR模型對差分后的極移殘差序列進行建模，得到一階差分項的預報值，最后將兩次擬合值和一階差分值預報值逆運算結(jié)果相加，獲得最終極移參數(shù)預報值。

圖4 Xpole參數(shù)第二次擬合前后殘差對比

圖5 Ypole參數(shù)第二次擬合前后殘差對比

為便于比較最小二乘和自回歸預報（LS+AR）模型改進后的預報精度，本文將采用三種預報方式進行短期預報。方案一：最小二乘和自回歸預報模型（LS+AR）；方案二：首先利用最小二乘進行擬合，然后對殘差項進行一階差分，最后對差分項進行自回歸預報（LS+一階差分+AR）；方案三：首先利用最小二乘進行擬合，然后對殘差項進行頻譜分析，獲得一次擬合后的殘差周期項，再次進行擬合，最后對二次擬合后的殘差項進行一階差分，利用自回歸模型對差分項進行預報（LS+LS（頻譜分析）+一階差分+AR）。最終試驗的預報結(jié)果（如圖6、圖7所示）：

圖6 Xpole參數(shù)預報精度比較

圖7 Ypole參數(shù)預報精度比較

由圖6和圖7可以看出：前10天三種預報方案預報精度差別不大，10～50天之間方案2和方案3的極移參數(shù)預報精度比方案1的極移參數(shù)預報精度有著明顯提高，再次驗證了AR模型要求建模的時間序列必須是平穩(wěn)的隨機序列。方案2和方案3的極移參數(shù)預報精度整體比較接近，圖8是方案3較之方案2的精度提高百分比，可以明顯看出極移X分量的預報精度提高較為明顯，極移Y分量除了在3～6日有少許下降外，其他時間段預報精度均有明顯提升。由此得出改進后的LS+AR模型是三種方案里預報精度最高的，驗證了本文所提出改進方法的可行性。

圖8 精度提高百分比

5.結(jié)束語

本文考慮到極運動周期時變性的特征，對LS+AR預報模型進行改進，對首次擬合后的極移參數(shù)殘差序列進行頻譜分析，利用頻譜分析得到的周期再次進行擬合，然后再對殘差項進行一階差分，最后對一階差分項進行自回歸預報，有效地提高了極移參數(shù)的預報精度。該方法僅從數(shù)學角度進行極移參數(shù)短期預報，并未分析頻譜分析得到周期的實質(zhì)的物理意義。研究極移序列潛在周期的具體激發(fā)源，并進一步改進將是下一步的工作重點。此外，本文推薦以4年作為極移基礎(chǔ)序列進行短期預報較為合適。