張軍杰

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

0.引言

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

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

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

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

1.LS+AR模型

1.1 LS模型

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

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

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

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

1.2 AR模型

AR模型是依據(jù)自身過去的時間序列的變化規(guī)律來預(yù)測下一時刻的值，要求建模的時間序列是零均值的平穩(wěn)隨機(jī)序列。所以建模之前需要對LS模型擬合后的殘差序列進(jìn)行預(yù)處理，通常是先進(jìn)行一階差分后再進(jìn)行建模。最終的極移預(yù)報值為LS模型的外推值與AR模型預(yù)報值之和。AR模型的數(shù)學(xué)表達(dá)式可表示為隨機(jī)序列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種，分別是最終預(yù)測誤差準(zhǔn)則、信息論準(zhǔn)則、傳遞函數(shù)準(zhǔn)則。理論上這三種方法是等效的，本文采用最終預(yù)測誤差準(zhǔn)則來確定AR模型的階。

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

1.3 精度評定

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

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

2.改進(jìn)的LS+AR模型

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

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

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

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

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

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

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

本次預(yù)報所用到的數(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）起開始預(yù)報，每7天預(yù)報一次，一共預(yù)報了350期，預(yù)報長度為50天。為進(jìn)一步削弱一次擬合后的殘差項里隱含周期項對預(yù)報精度的影響，首先用LS模型對極移參數(shù)進(jìn)行擬合，采用傅里葉分析法對一次擬合后的殘差項時間序列進(jìn)行頻譜分析，目的是將時間域中的隱藏波形轉(zhuǎn)變?yōu)轭l率域的頻譜，顯示其隱藏周期。每140天對極移X和Y分量的殘差項的時間序列去趨勢項和功率譜分析，對極移參數(shù)第一次擬合后的殘差結(jié)果的頻譜分析值（如表1所示），在18次的頻譜分析值中，每次頻譜分析的周期值變化差異較大，說明了十分有必要對其進(jìn)行二次擬合，削弱明顯周期項的影響。

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

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

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

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

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

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

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

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

圖8 精度提高百分比

5.結(jié)束語

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