地球自轉(zhuǎn)參數(shù)短期預(yù)報方法及其實現(xiàn)*

嚴(yán) 鳳 姚宜斌

(武漢大學(xué)測繪學(xué)院，武漢 430079)

地球自轉(zhuǎn)參數(shù)短期預(yù)報方法及其實現(xiàn)*

嚴(yán) 鳳 姚宜斌

(武漢大學(xué)測繪學(xué)院，武漢 430079)

根據(jù)極移及日長變化的特性，采用LS與AR模型相結(jié)合的方法，建立了適合地球自轉(zhuǎn)參數(shù)(ERP)的預(yù)報模型。為減少相鄰殘差存在的強(qiáng)相關(guān)性，先對殘差序列進(jìn)行差分預(yù)處理，然后利用AR模型進(jìn)行ERP的預(yù)報。對ERP進(jìn)行了1—10天的短期預(yù)報，其預(yù)報精度接近國際最好水平。

地球自轉(zhuǎn)參數(shù);三角多項式模型;自回歸模型;預(yù)報精度;平穩(wěn)性檢驗

1 引言

地球自轉(zhuǎn)參數(shù)是指由極移、日長變化以及歲差和章動構(gòu)成的描述地球自轉(zhuǎn)變化的參數(shù)，簡稱為ERP。ERP表達(dá)了地固參考系相對于瞬時天球參考系的運動，是實現(xiàn)天球坐標(biāo)系與地球坐標(biāo)系相互轉(zhuǎn)換的必要參數(shù)。極移分量表示天球歷書軸在地球參考系中的運動。世界時與協(xié)調(diào)時之差或者日長變化則反映了地球自轉(zhuǎn)的不均勻性。地球定向參數(shù)包含有豐富的地球動力學(xué)信息，而且還是導(dǎo)航與衛(wèi)星精密定軌中必需的參數(shù)［1］，對其進(jìn)行實時高精度的預(yù)報具有重要的科學(xué)意義和應(yīng)用價值。

近幾十年來，隨著甚長基線干涉測量、人衛(wèi)激光測距和全球定位系統(tǒng)等現(xiàn)代測量技術(shù)的飛速發(fā)展，極大地推動了地球自轉(zhuǎn)變化的研究。然而，由于測量模型和資料處理的復(fù)雜性，目前VLBI和SLR等技術(shù)獲取的ERP參數(shù)往往需要延遲2～5天，難以滿足導(dǎo)航與衛(wèi)星精密定軌對ERP參數(shù)的時效性需求。因此，尋求高精度ERP預(yù)報方法成為一項值得深入研究的課題。目前國際上常用預(yù)報方法有最小二乘外推法、協(xié)方差法、自回歸移動平均法、神經(jīng)網(wǎng)絡(luò)法、小波分析、卡爾曼濾波等方法，其預(yù)報精度已經(jīng)達(dá)到一定水平。國內(nèi)對ERP參數(shù)進(jìn)行預(yù)報起步較晚，但隨著探月工程以及北斗二代衛(wèi)星導(dǎo)航系統(tǒng)等空間計劃的逐步拓展，提高ERP參數(shù)預(yù)報的精度和時效性已變得十分迫切，目前多家科研機(jī)構(gòu)對此展開了研究。1977年湯家豪教授［2］在AR模型基礎(chǔ)上首次提出門限序列自回歸模型;2007年王琪潔將大氣角動量的影響引入地球自轉(zhuǎn)參數(shù)預(yù)報中，建立了基于神經(jīng)網(wǎng)絡(luò)技術(shù)的地球自轉(zhuǎn)變化預(yù)報;2008年楊元喜等人［3］提出為了帶多適應(yīng)因子的自適應(yīng)抗差濾波的動態(tài)數(shù)據(jù)預(yù)報方法，有效改善了動態(tài)預(yù)報精度;2010年許雪晴［4］等人比較LS+AR模型和人工神經(jīng)網(wǎng)絡(luò)模型的預(yù)報精度，發(fā)現(xiàn)AR模型在短期的預(yù)報中精度較高，而ANN模型則在中長期預(yù)報中有優(yōu)勢。為了研究不同技術(shù)和方法對ERP預(yù)報精度的影響，2005年10月—2008年2月，國際上舉行了一次大型的地球定向參數(shù)預(yù)報競賽。這次比賽總共由13組來自不同國家及不同學(xué)院的專家參與，經(jīng)過兩年多預(yù)報與比較，發(fā)現(xiàn)采用LS+AR模型對ERP參數(shù)進(jìn)行短期預(yù)報的精度很高。因此本文采用LS+AR模型對地球自轉(zhuǎn)參數(shù)進(jìn)行短期預(yù)報，并與EOP PCC的結(jié)果進(jìn)行比較，以驗證預(yù)報模型的有效性。

2 LS+AR模型

2.1 LS模型

地球自轉(zhuǎn)參數(shù)刻畫了極移、日長變化、歲差、章動等自然現(xiàn)象的變化，這些變化包括周期性、趨勢性等復(fù)雜過程，時間序列中主要的部分是周期項和趨勢項，本文采用三角多項式模型，并用最小二乘方法來擬合這兩項。極移參數(shù)計算公式為:

式中a、b為趨勢項的擬合系數(shù)，c、d為周期項的擬合系數(shù)，P1、P2分別為 Chandler項和周年項的周期［5］，ω1、ω2分別為Chandler項和周年項的初始相位。當(dāng)方程個數(shù)大于未知數(shù)個數(shù)時，可根據(jù)最小二乘原理求解參數(shù)。

UT1-UTC、Δlod計算公式為:

式中:a、b、c為趨勢項的擬合系數(shù)，d、e、f為周期項的擬合系數(shù)，P1、P2、P3分別為Chandler項、周年項和半周年項的周期［6］，φ1、φ2、φ3為這3個周期項的初始相位。

2.2 AR模型

1)預(yù)報模型

剔除周期項和趨勢項后的殘差序列，沒有明顯的規(guī)律特征。經(jīng)分析，發(fā)現(xiàn)該時間序列的后一項與前幾項具有相關(guān)性，根據(jù)這一特點采用自回歸(AR)模型來擬合殘差序列。自回歸模型就是根據(jù)變量自身過去的規(guī)律來建立預(yù)測模型，在動態(tài)數(shù)據(jù)處理中存在廣泛的應(yīng)用。本文采用的自回歸模型為:

式中:xt為模型變量;b1，b2，…，bp為模型的回歸系數(shù);εt為模型的隨機(jī)誤差;p為模型階數(shù)［7］。

自回歸模型是時間序列分析的主要模型之一。在建立自回歸模型時，難點在于模型階數(shù)的確定，而這個問題可以利用Matlab高階統(tǒng)計量工具箱中的函數(shù)來解決。確定階數(shù)后，即可用最小二乘方法估計出模型參數(shù)。Matlab高階統(tǒng)計工具箱也提供了相關(guān)函數(shù)求解模型參數(shù)，計算十分簡便。

2)外推預(yù)報

當(dāng)回歸系數(shù)b1，b2，…，bp確定時，可根據(jù)方程進(jìn)行預(yù)報。

第一步預(yù)報值為:

第二步預(yù)報值為:

一般地，第l步預(yù)報值為:

每一個基礎(chǔ)時間序列確定一個模型階數(shù)，然后根據(jù)模型階數(shù)求出回歸系數(shù)，再用公式(4)即可進(jìn)行預(yù)報。

2.3 地球自轉(zhuǎn)參數(shù)的預(yù)報

極移參數(shù)預(yù)報流程如圖1，UT1-UTC的預(yù)報流程如圖2。

從www.iers.org下載eopc04.62-now文件，抽取極移參數(shù)序列、UT1-UTC序列和Δlod序列(圖3)，采用從1995年至今的數(shù)據(jù)進(jìn)行建模分析。

3 差分處理

圖1 極移參數(shù)的預(yù)報流程Fig.1 Prediction process of polar motion parameters

AR模型要求建模的時間序列是平穩(wěn)隨機(jī)序列，即滿足平穩(wěn)、正態(tài)、零均值的條件［8］。剔除了周期項和趨勢項影響的殘差序列，其平穩(wěn)性還不能滿足平穩(wěn)隨機(jī)序列的條件。本文提出先對殘差序列進(jìn)行一次差分，再對差分后的殘差序列建立AR模型的處理方法。一次差分后的殘差序列平穩(wěn)性大大增強(qiáng)，以極移數(shù)據(jù)為例，首先利用時間序列的散點圖來判斷殘差序列以及差分后殘差序列的平穩(wěn)性(圖4)，可以看出，差分前時間序列都繞其均值上下波動，但差分后的時間序列波動更趨于平穩(wěn)，更接近零均值。

計算以250天作為子時間序列的極移參數(shù)的殘差和差分殘差的均值統(tǒng)計量以及方差統(tǒng)計量，結(jié)果如圖5、6。從圖5、6可以看出差分后時間序列的均值統(tǒng)計量及方差統(tǒng)計量更接近零，并且每段子時間序列的均值統(tǒng)計量及方差統(tǒng)計量基本不變，說明差分后的時間序列平穩(wěn)性更強(qiáng)。

極移參數(shù)以6年觀測數(shù)據(jù)為基礎(chǔ)時間序列，采用對殘差序列進(jìn)行差分處理和不做差分處理兩種方案分別進(jìn)行解算，圖7為兩種方案解算精度的比較，從圖7可以看出進(jìn)行一次差分處理后，大大提高了預(yù)報精度。

Δlod可以與UT1-UTC相互轉(zhuǎn)換，UT1-UTC去掉秒跳后，再作差分，結(jié)果的相反數(shù)即為Δlod。因此，對Δlod的預(yù)報可以轉(zhuǎn)換為對UT1-UTC的預(yù)報，再將UT1-UTC的預(yù)報值換算為Δlod的值即可。

圖2 UT1-UTC的預(yù)報流程Fig.2 Flow chart of UT1-UTC prediction

圖3 極移參數(shù)、UT1-UTC、Δlod基礎(chǔ)序列Fig.3 Basic sequences of polar motion parameters、UT1-UTC and Δlod

4 精度檢驗及結(jié)果分析

4.1 精度檢驗

圖4 極移參數(shù)的殘差序列和差分后殘差序列Fig.4 Residual sequences of the polar motion parameters before and after the difference

圖5 殘差序列和差分后殘差序列的均值統(tǒng)計量Fig.5 Mean statistics of residual sequences before and after the difference

圖6 殘差序列和差分后殘差序列的方差統(tǒng)計量Fig.6 Variance statistics of residual sequences before and after the difference

為與國際結(jié)果進(jìn)行比較，用預(yù)報誤差絕對值的平均值(MAE)來評定預(yù)報精度，計算式為:

圖7 做差分與不做差分預(yù)報精度比較Fig.7 Comparison of prediction accuracy between doing and not doing the difference

式中，εi，n是預(yù)報值與它在IERS C04中相應(yīng)值的差值，N為預(yù)報次數(shù)，I為預(yù)報長度。

本文采用與EOP PCC相同的預(yù)報間隔和期數(shù)，并預(yù)報相同時間點的ERP參數(shù)。EOP PCC從2005年10月1日起開始預(yù)報，每隔7天預(yù)報10天的值，總共預(yù)報了126期，取預(yù)報誤差絕對值的平均值作為預(yù)報精度(圖8(a)［9］)。本文分別采用2005年10月1日前3年、前6年、前9年的數(shù)據(jù)作為基礎(chǔ)序列，建立預(yù)報模型，其預(yù)報精度如圖8(b)。圖8反映了本文方法與EOP PCC所用方法的精度比較，其具體分析見表1。

表1 ERP絕對誤差(MAE)比較Tab.1 Comparison between ERP mean absolute prediction error

由表1可以看出Xp、Yp方向1～5天的預(yù)報精度高于EOP PCC絕大部分預(yù)報方法，6～10天預(yù)報精度達(dá)到EOP PCC領(lǐng)先水平。UT1-UTC 1～10天預(yù)報精度均達(dá)到EOP PCC領(lǐng)先水平。對于Δlod的預(yù)報，EOP PCC中大部分方法的預(yù)報精度都不高，精度較高的預(yù)報方法1～10天的精度為0.13～0.38 ms/day;本文直接由UT1-UTC的預(yù)報結(jié)果轉(zhuǎn)換而來的Δlod，其預(yù)報精度也相當(dāng)高，1～10天的預(yù)報精度為0.17～0.25 ms/day。綜上可見，本文采用三角多項式模型與自回歸模型相結(jié)合的方法，其預(yù)報精度已經(jīng)達(dá)到國際領(lǐng)先水平。

圖8 ERP預(yù)報精度Fig.8 ERP prediction accuracy

4.2 基礎(chǔ)序列長度對預(yù)報精度的影響

從上一節(jié)可以看出，不同長度的基礎(chǔ)序列其預(yù)報結(jié)果是不一樣的。目前，國內(nèi)外對ERP的預(yù)報絕大多數(shù)沒有考慮基礎(chǔ)序列長度對ERP預(yù)報精度的影響。本文分別以3年、6年、9年、12年作為基礎(chǔ)序列長度對地球自轉(zhuǎn)參數(shù)進(jìn)行預(yù)報，尋找基礎(chǔ)序列長度與ERP預(yù)報精度的聯(lián)系(圖9)。圖9中反映，基礎(chǔ)序列長度越長，其預(yù)報精度反而越差。這主要是因為相鄰數(shù)據(jù)之間具有明顯的相關(guān)性，而隨著時間間隔的增加，相關(guān)性逐漸下降。同時ERP基礎(chǔ)序列過長的話，其數(shù)據(jù)的精度不一致，從而影響預(yù)報精度。因此，在建立預(yù)報模型時，基礎(chǔ)序列長度不宜太長，這樣一方面可以保證基礎(chǔ)序列精度的一致，提高預(yù)報精度;另一方面，有效縮短了預(yù)報時間。另外，6年接近Chandler周期與周年周期的最小公倍數(shù)［9］，選取6年作為基礎(chǔ)序列長度，更能體現(xiàn)地球自轉(zhuǎn)參數(shù)的周期特征。

圖9 3、6、9、12年的預(yù)報精度Fig.9 Prediction accuracy of 3，6，9 and 12 years

5 結(jié)語

本文采用LS+AR模型對地球自轉(zhuǎn)參數(shù)進(jìn)行預(yù)報，并與EOP PCC的結(jié)果進(jìn)行比較，以驗證預(yù)報模型的有效性。提出先對殘差序列進(jìn)行一次差分處理的方法，求一次差分后序列的平穩(wěn)性大大增強(qiáng)。本文通過LS+AR模型預(yù)報精度達(dá)到國際領(lǐng)先水平。國內(nèi)外對EOP的預(yù)報絕大多數(shù)沒有考慮基礎(chǔ)序列長度對EOP預(yù)報精度的影響，本文通過試驗發(fā)現(xiàn)基礎(chǔ)序列長度越長，其預(yù)報精度反而越差，以6年作為基礎(chǔ)序列長度比較合適。

1 黨亞民，秘金鐘，成英燕.全球?qū)Ш叫l(wèi)星系統(tǒng)原理與應(yīng)用［M］.北京:測繪出版社，2007.(Dang Yamin，Mi Jinzhong and Cheng Yingyan.Principles and applications of global navigation satellite system［M］.Beijing:Chinese Surveying and Mapping Press，2007)

2 王琪潔.基于神經(jīng)網(wǎng)絡(luò)技術(shù)的地球自轉(zhuǎn)變化預(yù)報［D］.中國科學(xué)院上海天文臺，2007.(Wang Qijie.The earth rotation changes prediction based on neural network technology［D］.Shanghai Astronomical Observatory，2007)

3 Yang Y X and Cui X Q.Adaptively robust filter with multi adaptive factors［J］.Survey Review，2008，40(309):260 -270.

4 許雪晴，周永宏.地球定向參數(shù)高精度預(yù)報方法研究［J］.飛行控制學(xué)報，2010，29(2):70-76(Xu Xueqing and Zhou Yonghong.The research of the earth orientation parameters high accuracy prediction method［J］.Journal of Flight Control，2010，29(2):70-76)

5 錢昌夏.極移速度的周期分析［J］.中國科學(xué)院上海天文臺年刊，1995，16:35-40.(Qian Changxia.Periodic analysis of velocity of polar motion［J］.Annals of Shanghai Observatory Academia Sinica，1995，16:35-40)

6 韓永剛，李志安.地球自轉(zhuǎn)速率變化主要中短周期的分析［J］.地球物理學(xué)進(jìn)展，2002，17(2):349-352.(Han Yonggang and Li Zhian.Analysis of main middle and short periods in the variation of the Earth’s rotation［J］.Progress in geophysics，2002，17(2):349-352)

7 張善文，雷英杰，馮有前.MATLAB在時間序列分析中的應(yīng)用［M］.西安:西安電子科技大學(xué)出版社，2007.(Zhang Shanwen，Lei Yingjie and Feng Youqian.Application of MATLAB in time series analysis［M］.Xi’an:Xi’an University of Electronic Science and Technology Press，2007)

8 王新洲，等.高等測量平差［M］.北京:測繪出版社，2006.(Wang Xinzhou，et al.Advanced surveying adjustment［M］.Beijing:Chinese Surveying and Mapping Press，2006)

9 李征航，等.空間大地測量學(xué)［M］.武漢:武漢大學(xué)出版社，2010.(Li Zhenghang，et al.Space geodesy［M］.Wuhan:Wuhan University Press，2010)

SHORT-TERM PREDICTION METHODS AND REALIZATION OF EARTH ROTATION PARAMETERS

Yan Feng and Yao Yibin
(School of Geodesy and Geomatics，Wuhan University，Wuhan 430079)

According to the respective characteristics of polar motion and length of day，combining LS and AR model，a suitable prediction model for the Earth Rotation Parameter(ERP)was builded.In order to reduce the strong correlation，the difference before establishing AR model was down and ERP 1-10 days were predicted.The results show that the prediction accuracy has reached an international leading level.

the Earth rotation parameters(ERP);trigonometric polynomial model;auto-regressive model;prediction accuracy;stationary test

1671-5942(2012)04-0071-05

2012-01-24

國家自然科學(xué)基金創(chuàng)新研究群體項目(41021061);國家自然科學(xué)基金(41174012)

嚴(yán)鳳，女，碩士生，主要研究方向為地球自轉(zhuǎn)參數(shù)短期預(yù)報.E-mail:642304401@qq.com

P227;P207

A