電離層TEC預(yù)報(bào)模型的卡爾曼濾波改進(jìn)方法

王建敏，黃佳鵬，祝會(huì)忠，馬天明

(遼寧工程技術(shù)大學(xué) 測繪與地理科學(xué)學(xué)院，遼寧 阜新 123000)

0 引言

由于電離層自身的不穩(wěn)定性，對于穿過該層的全球衛(wèi)星導(dǎo)航系統(tǒng)(global navigation satellite system，GNSS)的電波會(huì)產(chǎn)生折射、散射等物理現(xiàn)象。正常情況下電離層的反射有助于GNSS無線電信號(hào)傳播，但是E層和F層劇烈變化對GNSS信號(hào)會(huì)產(chǎn)生嚴(yán)重的影響。電離層總電子數(shù)(total electron content，TEC)是電離層的重要特性之一，測量電離層中的總電子含量、研究其變化規(guī)律并對其進(jìn)行預(yù)報(bào)是一項(xiàng)重大的研究課題。因此研究預(yù)報(bào)模型，對電離層TEC的變化情況進(jìn)行高精度預(yù)報(bào)，進(jìn)而揭示電離層的物理機(jī)制，具有重要的科研及實(shí)用意義[1-8]；文獻(xiàn)[2]使用國際全球衛(wèi)星導(dǎo)航系統(tǒng)服務(wù)組織(International GNSS Service，IGS)公布的電離層TEC進(jìn)行自回歸滑動(dòng)平均模型(auto regressive and moving average model，ARMA)的建模與預(yù)測，以電離層平靜期和活躍期為研究對象，預(yù)報(bào)6 d相對精度達(dá)80 %以上；文獻(xiàn)[3]使用ARMA模型對電離層TEC進(jìn)行預(yù)報(bào)，以24 h為預(yù)報(bào)時(shí)間，精度可達(dá)90 %；文獻(xiàn)[4]利用夾角余弦和聚類分析對電離層TEC進(jìn)行混沌預(yù)測，得到的預(yù)測值在標(biāo)準(zhǔn)差和均方根(root mean square，RMS)2個(gè)指標(biāo)方面得到較好結(jié)果；文獻(xiàn)[5]利用神經(jīng)網(wǎng)絡(luò)模型完成1 d的電離層TEC預(yù)報(bào)，表明訓(xùn)練好的神經(jīng)網(wǎng)絡(luò)模型能夠充分反映出不同季節(jié)電離層TEC的變化；文獻(xiàn)[6]通過分析TEC的各個(gè)趨勢項(xiàng)，確定模型相應(yīng)參數(shù)，完成時(shí)間序列模型的預(yù)測；文獻(xiàn)[7]將奇異譜分析方法(singular spectrum analysis，SSA)與ARMA模型結(jié)合，使用SSA方法分解TEC數(shù)據(jù)，再利用ARMA模型完成電離層TEC預(yù)報(bào)，結(jié)果表明組合模型的預(yù)報(bào)精度相對僅利用ARMA模型預(yù)報(bào)提高了4 %；文獻(xiàn)[8]利用自回歸模型(auto regressive model，AR)模型與誤差反向傳播(back propagation，BP)神經(jīng)網(wǎng)絡(luò)模型等組合模型完成電離層TEC預(yù)報(bào)，得出組合模型相對單一模型完成預(yù)報(bào)可得到更好的預(yù)報(bào)結(jié)果的結(jié)論。目前所有研究均是直接基于IGS公布的數(shù)據(jù)完成預(yù)報(bào)，對于預(yù)處理方式大多以差分方式使數(shù)據(jù)平穩(wěn)化；但差分方式多是以簡單數(shù)學(xué)模型完成，并不完全遵循數(shù)據(jù)的變化規(guī)律，只能完成數(shù)據(jù)平穩(wěn)化處理。

卡爾曼濾波是以最小無偏差性為準(zhǔn)則的數(shù)據(jù)處理模型。本文將卡爾曼濾波引入電離層TEC預(yù)報(bào)的預(yù)處理，以濾波算法對IGS網(wǎng)站公布的電離層原始數(shù)據(jù)進(jìn)行降噪處理，減小噪聲對于組合模型建立和預(yù)報(bào)的影響，再使用AR模型和BP神經(jīng)網(wǎng)絡(luò)模型的組合模型進(jìn)行預(yù)測，分析比較濾波前后數(shù)據(jù)預(yù)測精度及模型適用條件。

1 模型基本原理

1.1 卡爾曼濾波

卡爾曼濾波是一種利用線性系統(tǒng)狀態(tài)方程對系統(tǒng)狀態(tài)進(jìn)行最優(yōu)估計(jì)的算法，它能減少隨機(jī)噪聲對電離層TEC觀測量的影響。它具有最小無偏差性，是當(dāng)前應(yīng)用最廣的一種數(shù)據(jù)預(yù)處理方法。它處理數(shù)據(jù)不僅會(huì)去除突變的數(shù)據(jù)，而且會(huì)保持?jǐn)?shù)據(jù)原有的變化趨勢，這樣有利于預(yù)報(bào)模型的建立。離散線性系統(tǒng)的卡爾曼濾波方程包括狀態(tài)方程和觀測方程[9]，即：

(1)

L(k)=B(k)X(k)+Δ(k)。

(2)

通過最小二乘原理可推得卡爾曼濾波方程為：

(3)

(4)

Dx(k)=(E-J(k)×B(k))×Dx(k-1)]。

(5)

1.2 AR模型

AR模型全稱為自回歸模型，是使用過去時(shí)間段的數(shù)據(jù)建立模型獲得權(quán)值，利用權(quán)值及過去時(shí)間段的數(shù)據(jù)進(jìn)行加權(quán)生成預(yù)測值。因此通過建立AR模型，可以對電離層TEC數(shù)據(jù)進(jìn)行預(yù)測。AR模型是使用一段隨時(shí)間變化而變化的數(shù)據(jù)序列，用相應(yīng)的模型加以描述，基于原始觀測量可以充分反映相關(guān)因素對預(yù)測數(shù)據(jù)的影響和作用，不受模型各個(gè)變量相互獨(dú)立的假設(shè)條件約束，所生成的描述模型可以消除一般預(yù)測方法中由于參數(shù)不明確、自變量選擇等造成的困難。數(shù)學(xué)模型的建立更多基于原始數(shù)據(jù)，更能發(fā)現(xiàn)這些數(shù)據(jù)本質(zhì)的、內(nèi)在的規(guī)律。

自回歸AR模型為

(6)

(7)

式中：

1.3 BP神經(jīng)網(wǎng)絡(luò)模型

BP神經(jīng)網(wǎng)絡(luò)是神經(jīng)網(wǎng)絡(luò)中被廣泛使用的一種模型。BP神經(jīng)網(wǎng)絡(luò)的特點(diǎn)是在運(yùn)算過程中可以進(jìn)行學(xué)習(xí)修正功能，包括自我學(xué)習(xí)、訓(xùn)練以及泛化等功能，可以實(shí)現(xiàn)對于數(shù)據(jù)反復(fù)運(yùn)算，不斷修正相應(yīng)參數(shù)和權(quán)值，直至達(dá)到預(yù)期數(shù)值。BP神經(jīng)網(wǎng)絡(luò)是一種多層反復(fù)運(yùn)算神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)。如圖1所示，BP神經(jīng)網(wǎng)絡(luò)最基本的拓?fù)浣Y(jié)構(gòu)是3層，包括輸入層、輸出層和一個(gè)隱層，上一層各個(gè)神經(jīng)元與下一層所有的神經(jīng)元連接，但是同一層各個(gè)神經(jīng)元之間沒有聯(lián)系，導(dǎo)致數(shù)據(jù)在運(yùn)算過程中，每層的神經(jīng)元只會(huì)影響到下一層的權(quán)值，對其他層沒有影響。若運(yùn)算過程并沒有得到理想的數(shù)值輸出，會(huì)計(jì)算誤差的值，隨后進(jìn)行學(xué)習(xí)過程，這一過程是反向傳播的，通過不斷學(xué)習(xí)訓(xùn)練誤差值，調(diào)節(jié)參數(shù)和權(quán)值，直到輸出值小于期望值，則輸出運(yùn)算結(jié)果值[13-18]。以電離層TEC濾波值作為BP神經(jīng)網(wǎng)絡(luò)的輸入值完成訓(xùn)練、建模和預(yù)報(bào)。

1.4 組合模型

由于AR模型和BP神經(jīng)網(wǎng)絡(luò)模型在完成電離層TEC預(yù)測時(shí)都留有殘差，為減少誤差較大的模型對組合模型的影響，組合過程中選擇改進(jìn)的加權(quán)組合方式，其表達(dá)式為

C(n)=αA(n)+βB(n)。

(8)

式中；A(n)為AR模型對電離層TEC的預(yù)測值；B(n)為BP神經(jīng)網(wǎng)絡(luò)模型對TEC的預(yù)測值；C(n)為組合模型的預(yù)測值，其中n表示組合模型的第n項(xiàng)；α與β分別為組合模型中2個(gè)單一模型的權(quán)值。為充分利用2個(gè)模型預(yù)測值，減少模型在組合過程中的運(yùn)算量，需要限制α與β的取值，即當(dāng)AR模型殘差優(yōu)于BP神經(jīng)網(wǎng)絡(luò)模型，取AR模型權(quán)值α為1，BP神經(jīng)網(wǎng)絡(luò)模型權(quán)值β為0，以此方式減少誤差較大模型對于組合模型的影響。下文簡稱使用IGS公布數(shù)據(jù)的BP神經(jīng)網(wǎng)絡(luò)與AR模型的組合模型(combined model of BP neural network and AR model using IGS data)為BA模型，使用經(jīng)過卡爾曼濾波數(shù)據(jù)的BP神經(jīng)網(wǎng)絡(luò)與AR模型的組合模型(combined model of BP neural network and AR model based on Kalman filter data)為KBA模型[19]。圖2為這2種模型的流程圖。

選擇的精度驗(yàn)證方法是殘差V分布情況和均方根(root mean square，RMS)，即

V=pi-qi。

(9)

式中：pi為不同預(yù)測模型的預(yù)測數(shù)據(jù)；qi為IGS網(wǎng)站公布的實(shí)際TEC數(shù)據(jù)。

(10)

式中：TECi(pre)為預(yù)測數(shù)據(jù)；TECi(means)為實(shí)際數(shù)據(jù)；n為預(yù)測TEC數(shù)據(jù)的個(gè)數(shù)。

2 實(shí)驗(yàn)與結(jié)果分析

使用IGS網(wǎng)站下載的2015年電離層網(wǎng)格點(diǎn)TEC數(shù)據(jù)對該方法進(jìn)行分析?？紤]自2015年全球網(wǎng)格點(diǎn)數(shù)據(jù)以1 h為時(shí)間間隔(1 d共24個(gè)觀測值)，限于篇幅，選用2015年4個(gè)月份前10 d數(shù)據(jù)，這些數(shù)據(jù)的可預(yù)測性處于平均水平。為驗(yàn)證該方法是否與數(shù)據(jù)的時(shí)間和緯度有關(guān)，取年積日為第1天至第10天、第90天至第100天、第181天至第191天、第273天至第283天，經(jīng)緯度為(87.5°N，125°E)、(45°N，125°E)、(0°N，125°E)。利用BA模型和KBA模型對上述基礎(chǔ)數(shù)據(jù)進(jìn)行提前1 d的預(yù)報(bào)。圖3為年積日為第181天至第190天電離層TEC原始數(shù)據(jù)與卡爾曼濾波后的數(shù)據(jù)對比，可以看出通過卡爾曼濾波后的電離層TEC數(shù)據(jù)仍然保持原始數(shù)據(jù)的變化趨勢，但是相對原始數(shù)據(jù)更加平緩(電離層TEC的計(jì)算單位為TECU，1個(gè)TECU表示每平方米有1016個(gè)電子)。

使用第1天至第9天數(shù)據(jù)進(jìn)行濾波前后模型的對比，圖4至圖6分別表示 87.5°N、 45°N、赤道數(shù)據(jù)在濾波前后建模預(yù)報(bào)的情況。

統(tǒng)計(jì)預(yù)測值的殘差分布和RMS值結(jié)果如表1所示。

表1 不同緯度的殘差(V)分布與均方根(RMS)比較情況(1月)

由表1可知：比較1月份各個(gè)緯度的預(yù)測情況，使用BA模型完成預(yù)測，殘差在±2個(gè)TECU以內(nèi)的預(yù)測值比例占36 %，殘差值小于±1個(gè)TECU的預(yù)測值比例占25 %，殘差平均值為-1.725 9個(gè)TECU；使用KBA模型完成預(yù)測，殘差在±2個(gè)TECU以內(nèi)的預(yù)測值比例占41.67 %，殘差值小于±1個(gè)TECU的預(yù)測值比例占25 %，殘差平均值為-0.654 9個(gè)TECU。結(jié)合圖4～圖 6，組合模型對于數(shù)據(jù)具有較好的預(yù)測作用，使用經(jīng)過卡爾曼濾波的數(shù)據(jù)與直接使用原始數(shù)據(jù)進(jìn)行預(yù)測相比，預(yù)測值曲線變化的趨勢非常接近，使用經(jīng)過卡爾曼濾波的數(shù)據(jù)完成建模和預(yù)報(bào)，殘差分布優(yōu)于直接使用原始數(shù)據(jù)進(jìn)行建模和預(yù)報(bào)，且使用經(jīng)過卡爾曼濾波后的數(shù)值進(jìn)行預(yù)測，預(yù)測結(jié)果的RMS值更小，擬合程度更好。

使用第90天至第99天的數(shù)據(jù)進(jìn)行濾波前后組合模型的對比，圖7～圖9分別表示4月87.5°N、 45°N、赤道數(shù)據(jù)在濾波前后建模預(yù)報(bào)的情況。

統(tǒng)計(jì)預(yù)測值的殘差分布和RMS值結(jié)果如表2所示。

表2 不同緯度的殘差(V)分布與均方根(RMS)比較情況(4月)

由表2可知：使用BA模型預(yù)測，殘差在±2個(gè)TECU以內(nèi)的預(yù)測值占總預(yù)測值的18.05 %，殘差在±1個(gè)TECU以內(nèi)的預(yù)測值占總預(yù)測值的9.72 %；使用KBA模型完成預(yù)報(bào)，預(yù)報(bào)殘差值小于±2個(gè)TECU的預(yù)測值比例占23.61 %，殘差在±1個(gè)TECU以內(nèi)的預(yù)測值占總預(yù)測值的18.05 %以上。結(jié)合圖7～圖9，BA模型和KBA模型對于數(shù)據(jù)都具有較好的預(yù)測作用，而KBA模型預(yù)測值相對BA模型預(yù)測值在殘差平均值、殘差分布和RMS值等方面占優(yōu)，表示KBA模型預(yù)測值擬合程度更好。

使用第181天至190天數(shù)據(jù)進(jìn)行濾波前后組合模型的對比，圖10～圖12分別表示7月87.5°N、45°N、赤道數(shù)據(jù)在濾波前后建模預(yù)報(bào)的情況。

統(tǒng)計(jì)預(yù)測值的殘差分布和RMS值結(jié)果如表3所示。

表3 不同緯度的殘差(V)分布與均方根(RMS)比較情況(7月)

由表3可知：使用BA模型完成預(yù)測，殘差在±1個(gè)TECU以內(nèi)的預(yù)測值占總預(yù)測值的 27.78 %，RMS值為3.968 6個(gè)TECU，殘差平均值為-1.439 4個(gè)TECU；使用KBA模型完成預(yù)報(bào)，預(yù)報(bào)殘差值小于±1個(gè)TECU的預(yù)測值比例占30.55 %，RMS值為3.193 3個(gè)TECU，殘差平均值為-0.938 7個(gè)TECU。結(jié)合圖10～圖12，組合模型對于數(shù)據(jù)具有較好的預(yù)測作用，BA模型和KBA模型預(yù)報(bào)數(shù)據(jù)曲線變化趨勢與實(shí)際數(shù)據(jù)非常接近，且KBA模型殘差統(tǒng)計(jì)方面和RMS值均優(yōu)于BA模型，因此KBA模型相對BA模型具有更好的預(yù)測能力。

使用第273天至第282天數(shù)據(jù)進(jìn)行濾波前后組合模型預(yù)測，圖13～圖15分別表示10月87.5°N、45°N、赤道建模預(yù)報(bào)的情況。

統(tǒng)計(jì)預(yù)測值的殘差分布和RMS值結(jié)果如表4所示。

表4 不同緯度的殘差(V)分布與均方根(RMS)比較情況(10 月)

由表4可知：BA模型殘差在±2個(gè)TECU以內(nèi)的預(yù)測值比例占72.22 %，殘差值小于±1個(gè)TECU的預(yù)測值比例占38.89 %，殘差平均值為-0.944 2個(gè)TECU；KBA模型殘差在±2個(gè)TECU以內(nèi)的預(yù)測值比例占77.78 %，殘差值小于±1個(gè)TECU的預(yù)測值比例占45.83 %，殘差平均值為-0.298 1個(gè)TECU。KBA模型RMS的平均值為2.609 0個(gè)TECU，BA模型RMS的平均值為3.002 1個(gè)TECU。RMS值越小精度越高，表示使用經(jīng)過濾波后的數(shù)據(jù)生成模型具有更高的預(yù)報(bào)精度。

3 結(jié)束語

本文采用卡爾曼濾波對IGS公布的電離層TEC值進(jìn)行預(yù)處理，再使用AR模型和BP神經(jīng)網(wǎng)絡(luò)模型得到的組合模型進(jìn)行建模和預(yù)測。原始數(shù)據(jù)通過濾波能夠降低原始狀態(tài)偏差的影響，得到更高精度的狀態(tài)估值，再利用高精度狀態(tài)值進(jìn)行建模和預(yù)報(bào)。通過實(shí)驗(yàn)證明，使用卡爾曼濾波處理原始數(shù)據(jù)不僅保持了數(shù)據(jù)趨勢性，而且使用處理后的數(shù)據(jù)進(jìn)行建模，模型預(yù)測精準(zhǔn)度更高。KBA模型在殘差值大小、殘差分布和RMS值方面均優(yōu)于BA模型，表示KBA模型相對BA模型具有更高的預(yù)報(bào)精度；KBA模型在中低緯度預(yù)測更加精準(zhǔn)，且10月是實(shí)驗(yàn)部分中預(yù)測精度最高的月份。將卡爾曼濾波應(yīng)用于預(yù)測模型預(yù)處理為模型的建立提供了新的思路。

[1] RATNAM V D，DINESH V B，TEJASWI B，et al.TEC prediction model using neural networks over a low latitude GPS station[J].International Journal of Soft Computing and Engineering，2012，2(2)：2231-2307.

[2] 張小紅，任曉東，吳風(fēng)波，等.自回歸移動(dòng)平均模型的電離層總電子數(shù)含量短期預(yù)報(bào)[J].測繪學(xué)報(bào)，2014，43(2)：118-124.

[3] 李磊，張淑芳，王軼卓，等.利用ARMA模型進(jìn)行電離層TEC預(yù)報(bào)[J].應(yīng)用基礎(chǔ)與工程科學(xué)學(xué)報(bào)，2013，21(5)：814-821.

[4] 孫佳龍，郭金云，郭淑艷.利用夾角余弦和聚類分析的電離層TEC混沌預(yù)測[J].武漢大學(xué)學(xué)報(bào)(信息科學(xué)版)，2014，39(4)：441-456.

[5] 李淑慧，彭軍還，徐偉超，等.利用神經(jīng)網(wǎng)絡(luò)預(yù)報(bào)短期電離層TEC變化[J].測繪科學(xué)，2013，38(1)：8-12.

[6] 陳鵬，姚宜斌，吳寒.利用時(shí)間序列分析預(yù)報(bào)電離層TEC[J].武漢大學(xué)學(xué)報(bào)(信息科學(xué)版)，2011，36(3)：267-270.

[7] 盧辰龍，匡翠林，張晉升，等.組合SSA與ARMA模型預(yù)報(bào)電離層TEC[J].大地測量與地球動(dòng)力學(xué)，2014，34(6)：44-49.

[8] 王建敏，黃佳鵬，祝會(huì)忠，等.電離層總電子數(shù)預(yù)報(bào)方法研究[J].測繪科學(xué)，2016，41(12)：47-53.

[9] 李奕.自適應(yīng)卡爾曼濾波在變形監(jiān)測數(shù)據(jù)處理中的應(yīng)用研究：以太原萬達(dá)廣場A2區(qū)基坑監(jiān)測為例[D].成都：成都理工大學(xué)，2012.

[10] 宋雪航，史良勝，楊金忠.基于集合卡爾曼濾波的潛水動(dòng)態(tài)預(yù)測方法[J].武漢大學(xué)學(xué)報(bào)(工學(xué)版)，2015，47(3)：324-331.

[11] 魏二虎，張帥，WEI Jianan.自適應(yīng)濾波在環(huán)火探測器脈沖星自主導(dǎo)航中的應(yīng)用[J].武漢大學(xué)學(xué)報(bào)(信息科學(xué)版)，2015，40(5)：701-705.

[12] 李秀海，郭達(dá)志.應(yīng)用半?yún)?shù)AR模型的電離層TEC建模與預(yù)測[J].測繪科學(xué)，2011，36(2)：149-151.

[13] 李捷斌，孔令杰.基于Kalman濾波的BP神經(jīng)網(wǎng)絡(luò)方法在大壩變形預(yù)測中的應(yīng)用[J].大地測量與地球動(dòng)力學(xué)，2009，29(4)：124-126.

[14] HABARULEMA J B，MCKINNELL L A，CILLIERS P J.Prediction of global positioning system total electron content using neural networks over South Africa[J].Journal of Atmospheric and Solar-Terrestrial Physics，2007，69：1842-1850.

[15] 田間.一種訓(xùn)練BP神經(jīng)網(wǎng)絡(luò)的融合算法[D].長春：吉林大學(xué)，2011.

[16] 王德明，王莉，張廣明.基于BP神經(jīng)網(wǎng)絡(luò)的短期風(fēng)速預(yù)報(bào)模型[J].浙江大學(xué)學(xué)報(bào)(工學(xué)版)，2012，46(5)：837-841.

[17] 李聰.基于BP神經(jīng)網(wǎng)絡(luò)的股票指數(shù)期貨價(jià)格預(yù)測[D].青島：青島大學(xué)，2015.

[18] 閃麗潔，張利平，劉戀，等.基于多方法優(yōu)選因子和人工神經(jīng)網(wǎng)絡(luò)耦合模型的枯水期徑流預(yù)報(bào)[J].武漢大學(xué)學(xué)報(bào)(工學(xué)版)，2015，48(6)：758-763.

[19] 趙傳華，黨亞民，秘金鐘.組合方法在電離層TEC短期預(yù)報(bào)中的應(yīng)用[J].大地測量與地球動(dòng)力學(xué)，2012，32(3)：76-79.