李秀海,郭達(dá)志,張隨甲

(1.中國礦業(yè)大學(xué) 資源與安全學(xué)院,北京 100083;2.黑龍江工程學(xué)院 測繪工程系,黑龍江 哈爾濱 150050;3.中國人民解放軍武裝警察黃金部隊第三支隊,黑龍江 哈爾濱 150086)

電離層 TEC的預(yù)測模型

李秀海1,2,郭達(dá)志1,張隨甲3

(1.中國礦業(yè)大學(xué) 資源與安全學(xué)院,北京 100083;2.黑龍江工程學(xué)院 測繪工程系,黑龍江 哈爾濱 150050;3.中國人民解放軍武裝警察黃金部隊第三支隊,黑龍江 哈爾濱 150086)

電離層總電子含量(TEC)的精確預(yù)報對提高 GNSS導(dǎo)航精度,保障無線電空間遠(yuǎn)程通訊具有重要作用。分析了IGS發(fā)布的電離層格網(wǎng)點總電子含量(TEC)的時間序列特點,基于時間序列分析理論,以AR模型對格網(wǎng)點TEC隨機(jī)時間序列平穩(wěn)化后建模和預(yù)報。實例分析表明,研究的預(yù)報技術(shù)和方法是可行的。

電離層;時間序列分析;總電子含量(TEC);AR模型;預(yù)報

自2000-05-01SA政策取消后,電離層延遲便成為單頻接收機(jī)定位和導(dǎo)航的最大誤差源,在太陽活動強烈的年份,電離層延遲誤差可達(dá)到150m以上,對導(dǎo)航和定位的影響不能忽視[1]。單頻 GPS接收機(jī)通常是通過建立電離層延遲誤差改正模型來提高定位精度,電離層延遲誤差的大小與電離層總電子含量(Total Electron Content,TEC)有關(guān) ,因此 ,研究精密測量與預(yù)報 TEC的變化規(guī)律對于建立精確的電離層改正模型,提高 GNSS定位和導(dǎo)航的精度具有十分重要的意義。此外,在航空航天、通訊、國防等領(lǐng)域利用電磁波進(jìn)行遠(yuǎn)程空間通訊時,研究電離層對電磁波傳播的影響也需要精密測定與預(yù)報電離層總電子含量 TEC[2]。

國內(nèi)外的一些學(xué)者對電離層 TEC的預(yù)報進(jìn)行了一些有益的探索[2-6]。國際 GPS服務(wù)(International GPS Service,IGS)通過對全球數(shù)百個 GPS跟蹤站的觀測數(shù)據(jù)處理,每兩小時公布地面上空固定高度的5°(經(jīng)度間隔)×2.5°(緯度間隔)的電離層格網(wǎng)點電子含量,推算的電子含量 TEC具有較高的精度,但由于不是實時發(fā)布,所以一些實時應(yīng)用尚需根據(jù)公布的數(shù)據(jù)進(jìn)行推估。電離層格網(wǎng)點的 TEC可看作是一種動態(tài)數(shù)據(jù)序列,這可利用時間序列分析理論對電離層格網(wǎng)點的 TEC進(jìn)行分析。本文利用IGS發(fā)布的電離層格網(wǎng)點 TEC時間序列數(shù)據(jù),基于時間序列分析理論,建立電離層格網(wǎng)點 TEC的預(yù)報模型,實現(xiàn)電離層 TEC的預(yù)報,為實時應(yīng)用提供 TEC數(shù)據(jù)。

1 時間序列建模及預(yù)報

1.1 時間序列模型簡介

電離層格網(wǎng)點 TEC的時間序列 Xt通常具有周期性的變化,引起這種周期性的主要原因是季節(jié)變化或其他物理因素所致,故有時也把這種序列稱為季節(jié)性時間序列。對這類時間序列數(shù)據(jù)進(jìn)行平穩(wěn)化處理時應(yīng)根據(jù)時序觀測圖、觀測序列的自相關(guān)函數(shù)等來判別其周期S,而后實行季節(jié)差分▽S=(1-BS)Xt(B為后移算子,▽為差分算子),差分后的時間序列如仍有趨勢,則再進(jìn)行正常差分,直到成為平穩(wěn)時間序列 z(t)[7]。

平穩(wěn)時間序列 z(t)的建模方法可利用ARMA(p,q)模型

式中:p、q為模型的階,φi、θi為模型的參數(shù),a(t)是均值為0,方差為σ的白噪聲。

當(dāng)θi=0時,式(1)可改化為自回歸AR(p)模型;當(dāng)φi=0時,式(1)則轉(zhuǎn)化為滑動平均MA(q)模型。

通過考察樣本自相關(guān)函數(shù)和偏相關(guān)函數(shù)以及A IC、B IC準(zhǔn)則等對模型進(jìn)行辨識和定階。吳賢銘和Pandit從分析系統(tǒng)特性出發(fā),提出了動態(tài)數(shù)據(jù)系統(tǒng)(Dynam ic Data System,DDS)建模方法。基于如下的ARMA(n,n-1)模型對時間序列統(tǒng)一建模

建模時依次由低階到高階,并用 F檢驗法判別模型的適用性,避免了定階的困難。

對于平穩(wěn)可逆的ARM A模型,如下關(guān)系成立:

當(dāng) n充分大時,AR(n)可以任意精度逼近ARMA(p,q)模型。因此,可基于有限階AR(n)模型對時間序列統(tǒng)一建模,優(yōu)點是對AR(n)模型可用遞推線性最小二乘法估計其參數(shù),其算法簡單,計算量和存儲量小,從而極大改進(jìn)了傳統(tǒng)的建模方法。對電離層格網(wǎng)點 TEC的預(yù)報而言,建立AR模型也很方便。

1.2 AR(n)模型參數(shù)估計的遞推最小二乘(RLS)

為便于分析，假設(shè)儲罐內(nèi)LNG經(jīng)BOG換熱后完全變?yōu)轱柡蛻B(tài)LNG，分析所需要的LNG過冷度。其中，LNG儲罐內(nèi)壓力為0.3 MPa，對應(yīng)飽和溫度為-146.6℃。

算法

式(3)改寫成向量形式為

定義向量AR(n)模型的遞推最小二乘估值為

P(t)為未知數(shù)的協(xié)因數(shù)陣,計算時可取 P(0)=a In,a為很大的實數(shù),例如,取 a=105時,Φ^(0)=Φ0,σ^2a(1)=^a2(1)[8]。

當(dāng)模型參數(shù)隨時間變化時,隨著觀測數(shù)據(jù)的增加,用以上遞推最小二乘法估計參數(shù)的誤差有時較大,這是因為新數(shù)據(jù)被舊數(shù)據(jù)所淹沒的原因。為了體現(xiàn)過程的時變性,應(yīng)強調(diào)新數(shù)據(jù)的作用,逐漸淡化陳舊數(shù)據(jù)的影響。為此,可引入遺忘因子,并采用加入遺忘因子的遞推最小二乘解法。對于AR模型,對殘差平方和的每項采用一個加權(quán)系數(shù)λt-i(λ為遺忘因子,0＜λ≤1),則AR(n)模型的遞推最小二乘估值為

當(dāng)λ＜1時,體現(xiàn)了新近數(shù)據(jù)的作用,當(dāng)λ=1時化為普通最小二乘法。一般情況下,對慢時變參數(shù),λ應(yīng)取較大的值,而對快時變參數(shù),λ應(yīng)取較小的值。關(guān)于λ的取值可參考文獻(xiàn)[4]。

1.3 AR(n)模型階的 F檢驗法

模型定階是基于已知的 N個觀測數(shù)據(jù){z(i),i=1,2,…,N},由低階到高階遞增地對系統(tǒng)擬合AR模型,并依次對相鄰的兩個AR模型用以下檢驗判定模型的合適性[9]。

若記高、低階模型階次分別為 n1、n2,構(gòu)造統(tǒng)計量

式中:S1、S2分別為低、高階模型的殘差平方和,N為觀測序列長度。對于給定的顯著水平α,一般選α=0.05～0.10(對應(yīng)的置信度為90%～95%)。當(dāng)按式(12)算得的 F＞Fα(n2-n1,N-n2)時,表明S1、S2有顯著差異,低階模型不適用,可升階建模;反之,則認(rèn)為低階模型適用。

經(jīng)F檢驗后的模型具有AR(n)的形式,但并不一定為最終的合適模型。AR(n)模型中某些參數(shù)φ1,φ2,…,φn的估值可能接近于 0,即實際上可能某些參數(shù)為0。上述情況歸結(jié)為統(tǒng)計假設(shè)檢驗所得到的AR(n)中某些系數(shù)為0,這可通過模型參數(shù)估計的置信區(qū)間和F檢驗進(jìn)行模型修正。先刪去一個95%置信區(qū)間含零點的參數(shù),重新建立模型,并用F檢驗法與原來擬合的AR(n)進(jìn)行比較。如不顯著,正式刪去這個參數(shù);如顯著,則保留這個參數(shù)。對其余95%置信區(qū)間含零點的參數(shù)作類似的處理,則可得到修正后的模型。

2 實例分析與應(yīng)用

在IGS網(wǎng)站下載2008-01-01—2008-01-09的東經(jīng)125°,北緯 45°地區(qū)的電離層格網(wǎng)點的 TEC,共108個數(shù)據(jù)(每2 h 1個數(shù)據(jù),TEC單位以 TECU表示,1TECU=1016個電子m-2),對該點的 TEC進(jìn)行分析,以前7 d觀測數(shù)據(jù)建立數(shù)學(xué)模型,并對后2 d的電離層格網(wǎng)點TEC進(jìn)行預(yù)測,同時與觀測值進(jìn)行比較分析,驗證模型預(yù)測的正確性和可靠性。

首先,對該時間序列數(shù)據(jù)進(jìn)行平穩(wěn)化預(yù)處理。建立時間序列模型之前,首先對觀測數(shù)據(jù)進(jìn)行零均值、平穩(wěn)化處理,使其變成平穩(wěn)時間序列。圖1為原始時間序列數(shù)據(jù)隨時間的變化情況,可以看出原始數(shù)據(jù)大致有時間間隔為12的周期,不是平穩(wěn)時間序列。進(jìn)一步求出原始數(shù)據(jù)的自相關(guān)系數(shù)如圖2所示,原序列周期為 12(在ρ(12)=0.757 9,ρ(24)=0.624 1,ρ(36)=0.508 4,ρ(48)=0.411 3 等處出現(xiàn)峰值)。對原始數(shù)據(jù)進(jìn)行季節(jié)差分和一次正常差分,則原序列變?yōu)閦(t)=(1-B)(1-BS)Xt。經(jīng)統(tǒng)計檢驗 z(t)序列的均值 ˉz ＜ 2σˉz,故可以認(rèn)為 z(t)是零均值平穩(wěn)時間序列。用AR(n)模型經(jīng)程序辨識、F檢驗,獲得模型合適的階為3。

圖1 觀測數(shù)據(jù)分布曲線

圖2 時間序列自相關(guān)系數(shù)圖

模 型 參 數(shù) Φ =(φ^1,φ^2,φ^3)T=(-0.470 3,-0.232 3,-0.163 7)T。電離層格網(wǎng)點 TEC的預(yù)報模型為

應(yīng)用建立的AR(3)模型,按最小方差預(yù)報方法進(jìn)行步長為24的預(yù)報。對時序還原后,得到預(yù)測數(shù)據(jù)的時序曲線如圖3所示。通過試算遺忘因子λ,表明本例λ=1較合適。另外,也利用前6 d的數(shù)據(jù)建模,得到AR(1)模型,利用該模型預(yù)報后3 d的結(jié)果見圖4。

由AR(3)和 AR(1)的預(yù)報標(biāo)準(zhǔn)差分別為0.73TECU和0.94TECU。利用AR(1)模型預(yù)報的2008-01-07—2008-01-09的精度較差,尤其是預(yù)測最后一天的TEC誤差較大,說明短時間的預(yù)報精度比長時間的預(yù)報精度高。另外,從圖3、圖4可看出AR(n)模型預(yù)測的結(jié)果很好地反映了電離層格網(wǎng)點 TEC的變化趨勢。因此,用AR(n)模型對電離層格網(wǎng)點 TEC進(jìn)行預(yù)報時,時間越短效果越好。值得注意的是,個別時刻預(yù)測的誤差較大,需要做進(jìn)一步研究。本例中兩種預(yù)報模型的預(yù)報精度都比較高,用時間序列分析理論對電離層格網(wǎng)點TEC的預(yù)報取得了初步成果。

在電離層格網(wǎng)點TEC預(yù)報過程中,用已取得的觀測時間序列建立AR模型,對后面的電離層格網(wǎng)點 TEC進(jìn)行預(yù)測;當(dāng)再一次取得新值時,把它加入觀測時間序列,重新建立AR模型,模型得到修正,再去進(jìn)行后面的預(yù)測;如此反復(fù)建模、預(yù)測、修正…,可達(dá)到動態(tài)預(yù)測電離層格網(wǎng)點TEC的目的。

3 結(jié) 論

本文應(yīng)用時間序列分析理論和方法對電離層格網(wǎng)點 TEC進(jìn)行預(yù)報,推出了簡單的AR模型,利用F檢驗法對模型定階,采用遞推最小二乘解法進(jìn)行模型參數(shù)估計,便于計算機(jī)上實現(xiàn)建模、預(yù)報。通過上述研究和分析,得到以下幾點結(jié)論:

1)利用時間序列分析法建立的AR模型應(yīng)用于電離層格網(wǎng)點 TEC預(yù)報,具有建模容易、計算簡單和預(yù)報快速的特點。

2)在建立AR模型過程中,由于需要反復(fù)利用歷史數(shù)據(jù),本文推出了遞推最小二乘解法,使計算得到簡化,計算量和存儲量較小,便于上機(jī)實現(xiàn)。

3)時間序列分析方法對電離層格網(wǎng)點 TEC具有較高的模型擬合及預(yù)報精度,尤其是短期預(yù)報,效果更佳。

[1]蔡昌盛,高井祥,李征航.利用 GPS監(jiān)測電離層總電子含量的季節(jié)變化[J].武漢大學(xué)學(xué)報:信息科學(xué)版,2006,31(5):451-453.

[2]李志剛,程宗頤,馮初剛,等.電離層預(yù)報模型研究[J].地球物理學(xué)報,2007,50(2):327-337.

[3]劉瑞源,劉順林,徐中華,等.自相關(guān)分析法在中國電離層短期預(yù)報中的應(yīng)用[J].科學(xué)通報,2005,50(24):2781-2784.

[4]KN IGHT,M.M.CERVERA,A.FINN,“A Comparison of Predicted and Measured GPS Perfo rmance in an I-onospheric scintillation Environment,”Proceedings of ION GPS-99,1999:1437-1450.

[5]P.K.BHUVAN,RASHM IREKHA BORAH,“TEC derived from GPSnetwo rk in India and comparison w ith the IRI,”Advances in Space Research,2007:830-840.

[6]ZH IZHAO L IU,Ionosphere Tomographic Modeling and Applications U sing Global Positioning System (GPS)Measurements[D].Calgary:Diversity of Calgary,2004.

[7]楊叔子,吳 雅.時間序列分析的工程應(yīng)用(上、下冊)[M].武漢:華中理工大學(xué)出版,2007:170-175.

[8]鄧自立.建模于與估計[M].北京:科學(xué)出版社,2007:62-67.

[9]潘國榮.基于時間序列分析的動態(tài)變形預(yù)測模型研究[J].武漢大學(xué)學(xué)報:信息科學(xué)版,2005,30(6):483-487.

The ionospheric TEC forecastmodel

L IXiu-hai1,2,GUO Da-zhi1,ZHANG Sui-jia3
(1.School of Safety and Resource Engineering,China University of M ining&Technology,Beijing 100083,China;2.Department of Surveying and Mapping,Heilongjiang Institute of Technology,Harbin 150050,China;3.The Third Branch M ilitary of the Chinese Peop le’s Liberation Gold Police A rmy,Harbin 150086,China)

Accurately fo recasting the ionospheric total electron content p lays an impo rtant role in imp roving GNSS navigation and p rotecting remote radio comm unication.The characteristic of the TEC time series values fo r the ionosphere grid point observed by International GPSServicer is analyzed in this paper.Based on time seriesanalysis,ARmodelsare built fo r the stationary time series valuesat the ionospheric grid and used for forecasting ionospheric TEC.The p reliminary result show s the p roposed methodology and technique is feasible by illustration.

ionosphere;time series analysis;total electron content(TEC);AR model;forecast

P352

A

1006-7949(2010)01-0005-04

2009-03-13

李秀海(1962-),男,副教授,博士研究生.

[責(zé)任編輯劉文霞]