基于LS+AR模型的極移預(yù)報及其周期項影響分析

孫張振，徐天河，范朋飛

(1．長安大學地質(zhì)工程與測繪工程學院，陜西西安710054;2．地理信息工程國家重點實驗室，陜西西安710054;3．西安測繪研究所，陜西西安710054)

一、引 言

地球的自轉(zhuǎn)運動表征著地球整體的運動狀況，由地球各圈層地核、地幔地殼和大氣之間的相互作用引起的，可用地球定向參數(shù)(earth orientation para meters，EOP)描述[1]。EOP 包括 3 部分：① 章動和歲差;② 極移;③ 日長變化。其中，極移是表征地球運動的重要參數(shù)，它用來描述地球瞬時自轉(zhuǎn)軸在地球本體的運動而導(dǎo)致地球極點在地球表面上的位置發(fā)生緩慢變化的現(xiàn)象。高精度的極移觀測和預(yù)測數(shù)據(jù)對于天球參考框架和地球參考框架具有重要意義，是提高空間導(dǎo)航定位技術(shù)精度的必要部分。針對極移的預(yù)報，國內(nèi)外許多學者建立了各種預(yù)報模型對其進行研究、預(yù)報和監(jiān)測[2-8]?？偟膩碚f，這些模型大部分是將極移分為有規(guī)律的確定性部分和無規(guī)律的隨機性部分，即首先利用最小二乘方法(least squares，LS)來擬合極移時間序列中的趨勢項，這些有規(guī)律的確定性部分包括極移的長期趨勢(LT)、錢德勒擺動(CW)、周年擺動(AW)和半周年擺動(SAW)，并對其進行外推預(yù)測;然后利用AR 模型[6]、

ARMA模型或者人工神經(jīng)網(wǎng)絡(luò)[9]等確定性或者隨機性模型對其殘差部分進行預(yù)測，最后得到預(yù)報的極移結(jié)果。在上述預(yù)測方法中，利用最小二乘方法擬合其趨勢項時，國內(nèi)外學者在周期項的具體數(shù)值上雖大體一致但普遍不同，國外學者在錢德勒周期中選擇在433～435 d之間較多。國內(nèi)學者如許雪晴在極移預(yù)報時選取錢德勒周期約為435 d[10]，張昊等錢德勒周期選取為 1．183 a(432．08 d)[11]，在周年項和半周年項上國內(nèi)外選取比較一致，都為1 a(365．24 d)和0．5 a(182．62 d)。由于在極移觀測數(shù)據(jù)中其周期項和趨勢項都有時變性的特點[12]，為了更加精確地進行極移預(yù)報，本次在LS+AR預(yù)報模型的基礎(chǔ)上，在利用最小二乘擬合數(shù)據(jù)趨勢項之前，首先對62 a來的極移數(shù)據(jù)做頻譜分析[13]，由計算結(jié)果顯示其周期項中錢德勒周期為435．17 d，周年項為365．54 d，半周年項峰值對于X、Y分量有兩個峰值，卻有微小不同，鑒于其波峰值較小，故仍統(tǒng)一采用182．62 d。因此，此次在進行LS擬合時，應(yīng)用頻譜分析得出的周期項參數(shù)，聯(lián)合AR模型預(yù)報極移并與其他學者所用周期項預(yù)報結(jié)果相比較。計算結(jié)果表明，應(yīng)用本文所選取的周期項進行預(yù)報時，其預(yù)報精度最高。

二、LS+AR模型

LS+AR模型為一種組合的模型，在對極移進行預(yù)報時，首先應(yīng)用最小二乘擬合趨勢項，扣除后得到殘差部分，再利用AR模型對殘差部分進行預(yù)測，結(jié)合最小二乘外推部分即可得到預(yù)測值。

1．LS模型

在本文中，利用最小二乘外推模型中包含有線性項、半周年項、周年項和錢德勒項，模型公式為

式中，a0為常數(shù)項;a1為趨勢項;B1、B2、C1、C2、D1、D2…為各周期項擬合系數(shù);R1、R2、R3…分別為周期項中的周期;t為UTC時間。

2．AR模型

式(2)中，φ1，φ2，…，φp為模型參數(shù);at為白噪聲序列;p為模型階數(shù);式(2)稱為p階自回歸模型，簡記為AR(p)。

AR模型要求建模的時間序列是平穩(wěn)隨機序列，即滿足平穩(wěn)、正態(tài)、零均值的條件。運用AR模型的關(guān)鍵是確定AR模型的階數(shù)p，一般采用最終誤差準則、信息論準則、和傳遞函數(shù)準則，實際中這3種方法是等效的。本文采取最終誤差準則來確定AR模型的階數(shù)p：

最終誤差準則

其中

式(4)中使FPE(p)達到極小值的p就取作AR模型的階數(shù)。φ1，φ2，…，φp等模型參數(shù)可以通過求解Yule-Walker方程來確定。

3．精度評定

為了便于評估預(yù)測精度，采用平均絕對誤差(mean absolute error，MAE)作為精度評定標準，其計算公式為

式中，Pi為第i期的預(yù)報值;Oi為相應(yīng)的觀測值;Ei為跨度第i期預(yù)報的真誤差(假定觀測值為真值);n為總預(yù)報期數(shù);MAEj為預(yù)報跨度為j的平均絕對誤差。

三、數(shù)據(jù)計算和結(jié)果分析

1．數(shù)據(jù)說明

本次計算所用到的數(shù)據(jù)來自國際地球自轉(zhuǎn)和參考系服務(wù)IERS發(fā)布的EOP 05C04序列(http：∥hpiers．obspm．fr/eoppc/eop/)，該序列包含1962 年1 月1日至今的極移分量Xp和Yp、日長變化LOD、地球自轉(zhuǎn)速率UT1-UTC等，時間間隔為1 d。首先對極移分量Xp和Yp做頻譜分析，圖1為極移頻譜周期分布圖，從圖上可明顯看出極移中錢德勒周期約在435 d和周年周期項約為365 d。結(jié)合實際計算數(shù)據(jù)，可知其中錢德勒周期為435．17 d，周年項為365．54 d。

圖1 極移頻譜分析圖

2．數(shù)據(jù)處理

為了確定極移中具體周期對極移預(yù)報精度的影響，本文此次在進行最小二乘擬合時根據(jù)以往所采用的周期經(jīng)驗值和本文中所確定的周期值分別對極移進行預(yù)報，設(shè)計方案如下。

方案一：周期項選取錢德勒周期項432．08 d，周年項365．24 d，半周年項182．62 d。

方案二：周期項選取錢德勒周期項433 d，周年項365．24 d，半周年項182．62 d。

方案三：周期項選取錢德勒周期項435 d，周年項365．24 d，半周年項182．62 d。

方案四：周期項選取錢德勒周期項435．17 d，周年項365．24 d，半周年項182．62 d。

方案五：周期項選取錢德勒周期項435．17 d，周年項365．54 d，半周年項182．62 d。

在對極移(PM)進行預(yù)報時，其觀測基礎(chǔ)序列長度對預(yù)報精度有很大的影響，為了提高預(yù)報精度，本次預(yù)報所應(yīng)用的都是最優(yōu)基礎(chǔ)序列長度，極移分量X和Y基礎(chǔ)序列長度都為10 a觀測序列。預(yù)報精度結(jié)果統(tǒng)計是從2004年1月1日至2012年4月8日，分別對其進行跨度為1～360 d的預(yù)報。圖2～4給出了以上5種方案的預(yù)報精度對比圖。

從圖2可以看出，錢德勒周期的具體數(shù)據(jù)選取對極移的預(yù)報精度有很大關(guān)系，雖然其變化只在3 d內(nèi)，但越是接近錢德勒周期真值，其極移預(yù)報就會越精確。從圖3中可知，方案三和方案四預(yù)報精度非常接近，方案四更精確一些，也說明本文所用的錢德勒周期數(shù)據(jù)更加接近于真值。而從圖4可以看出，方案五的預(yù)報精度要高于方案四，這也可以說明在極移數(shù)據(jù)中，周年項的具體真值并不是只局限于365．24 d，本文中所應(yīng)用的周年項數(shù)值365．54 d更接近與真值。

圖2 方案一、二、三預(yù)報精度對比

圖3 方案三、四預(yù)報精度對比

圖4 方案四、五預(yù)報精度對比

四、結(jié)束語

本文基于LS+AR預(yù)測模型對極移進行高精度的預(yù)報，由于在極移數(shù)據(jù)中其周年項、錢德勒項等周期項具有時變性的特點，并不是只局限于某些具體的數(shù)值，故在預(yù)報之前首先對極移數(shù)據(jù)做頻譜分析，確定出所應(yīng)用觀測數(shù)據(jù)中的最優(yōu)周年項和錢德勒項數(shù)值，然后應(yīng)用所確定的具體數(shù)值來進行最小二乘的擬合進而預(yù)報未來極移數(shù)據(jù)，通過此種方法預(yù)報的極移精度要比以往國內(nèi)外學者應(yīng)用同種方法的具體周期數(shù)所得的預(yù)報精度都要高，因此，建議在進行地球定向參數(shù)預(yù)報之前，應(yīng)首先確定所用數(shù)據(jù)的具體周期頻譜值，在極移預(yù)報中加以考慮，可進一步提高預(yù)報的精度和可靠性。

[1]鄭大偉，虞南華．地球自轉(zhuǎn)及其地球物理現(xiàn)象的聯(lián)系：I日長變化[J]．地球物理學進展，1996，11(2)：81-101．

[2]AKYILMAZ O，KUTTERER H．Prediction of Earth Rotation Parameters by Fuzzy Inference Systems[J]．Journal of Geodesy，2004(78)：82-93．

[3]SCHUH H，ULRICH M，EGGER D，et al．Prediction of Earth Orientation Parameters by Artificial Neural Networks[J]．Journal of Geodesy，2002(76)：247-258．

[4]KOSEKW，KALARUSM，NIEDZIELSKIT．Forecasting of the Earth Orientation Parameters—Comparison of Different Algorithms[C]∥ Proceedings of the Journèes 2007，SystèMes De Référence Spatio-temporels“The Celestial Reference Frame for the Future”．Paris：[s．n．]，2007：155-158．

[5]KOSEKW．Future Improvements in EOP Prediction[C]∥Proceedings of the IAG 2009，“Geodesy for Planet Earth”．Argentina：[s．n．]，2009．

[6]KOSEK W．Causes of Prediction Errors of Pole Coordinates Data[C]∥Proceedings of the 6th Orlov’s Conference，“The Study of the Earth As a Planet by Methods of Geophysics，Geodesy and Astronomy”Ukraine：[s．n．]，2009：96-103．

[7]KALARUSM，SCHCH H，KOZEK W，et al．Achievements of Earth Orientation Parameters Prediction Comparison Campaign[J]．Journal of Geodesy，2010(84)：587-596．

[8]SCHUH H，NAGEL S，SEITZ T．Linear Drift and Periodic Variations Observed in Long Time Series of Polar Motion[J]．Journal of Geodesy，2001(74)：701-710．

[9]王琪潔．基于神經(jīng)網(wǎng)絡(luò)技術(shù)的地球自轉(zhuǎn)變化預(yù)報[D]．上海：中國科學院上海天文臺，2007．

[10]XU X Q，ZHOU Y H，LIAO X H．Short-term Earth O-rientation Parameters Predictions by Combination of the Least-squares，AR Model and Kalman Filter[J]．Journal of Geodynamics，2011(1107)：4．

[11]張昊，王琪潔，朱建軍，等．加權(quán)最小二乘法與AR組合模型在極移預(yù)測中的應(yīng)用研究[J]．天文學進展，2011，29(3)：343-352．

[12]KOSEK W，MCCARTHY D D，LUZUM B J．El Nino Impact on Polar Motion Prediction Errors[J]．Stadia Geophysica et Geodaetica，2001(45)：347-361．

[13]DENNISDM，BRIAN JL．Prediction of Earth Orientation[J]．Journal of Geodesy，1991(65)：18-21．