基于GNSS TEC 的電離層自相關預報誤差分析及參數優化

熊雯，王博文，劉裔文,3，朱慶林

(1.中國電波傳播研究所,山東 青島 266107；2.南昌大學 信息工程學院,南昌 330031；3.上饒師范學院 物理與電子信息學院,江西 上饒 334001)

0 引言

電離層電子密度的變化是影響無線電系統工作性能的重要誤差源之一，對電離層總電子含量(TEC)進行預報是減緩電離層誤差的主要途徑.尤其對全球衛星導航系統(GNSS) 應用具有重要使用價值.例如，對GNSS 單頻用戶，TEC 預報是有效減緩電離層延遲誤差影響的重要途徑.另外，電離層預報也是提高基于GNSS 信標觀測的全球電離層地圖(GIM)精度的重要手段.

從預報時間提前量的角度看，電離層預報模型可分為長期預報模型和短期預報模型.長期預報模型一般是太陽活動指數或地磁活動指數驅動，能反應電離層周日、季節、11 年周期變化等大尺度變化的全球性模型，例如國際參考電離層模型、NeQuick 模型等.短期預報模型通常基于臨近的歷史觀測數據的建模來預測當地未來1～24 h 的結果.因短期預報具有較高的區域預測精度，且效率高，故常被用于一些對時效性和精度要求較高的導航、通信、雷達等系統的電離層延遲修正.

早在2000 年左右，就有學者提出了一些有效的短期預報方法，例如自相關方法、自回歸方法和人工神經網絡等[1-3].目前，工程上用的較多的是自相關或自回歸方法[4-8].制約電離層短期預報精度的主要因素主要是磁暴、亞暴等引起的電離層擾動和平靜期電離層存在的一些局地小尺度擾動.這些擾動的變化很復雜，至今依然是研究的熱點，目前還很難對其進行準確預報.總的來說，對于電離層短期預報，至今并沒有一種絕對優勢的方法，尋求對已有傳統方法進行優化仍然是提高預測精度的重要途徑之一[9-11].

自相關預報是一種模型簡單且精度較高的預報方法，本文針對自相關預報方法中的參數設置對誤差的影響進行深入分析，進而提出優化方案.首先介紹所使用的數據和處理方法，其次介紹自相關預報原理和方法，在再次比較自相關與自回歸方法的基礎上，開展誤差分析和參數優化分析；最后進行討論和總結.

1 數 據

本文研究所采用的數據主要來自麻省理工大學CEDAR Madrigal 數據庫，其數據來自于分布全球的地基衛星接收機網絡(Distributed Ground Based Satellite Receivers)中的GNSS 接收機網絡(World-wide GNSS Receiver Network) (http://cedar.openmadrigal.org/list).

Madrigal 收集的TEC 測量值均為垂直TEC，單位為TECU，1 TECU 表示每平方米有1016個電子.垂直TEC 是由斜TEC 通過電離層薄層假設(高度為350 km)轉換得到.最后通過格網化高精度TEC 地圖重構算法生成并提供1°×1°×5 min 時空分辨率的全球垂直TEC 地圖[12].

但是，Madrigal 數據庫并沒有直接提供單站垂直TEC 數據.為開展單站電離層預報方法研究，本文在全球TEC 地圖中從觀測站網密集區域任意選取了5 個格網點的TEC 數據，然后通過1 h 平均得到單站垂直TEC 日變化時間序列，以此作為本文的實測TEC 時間序列.圖1 給出了所選取的5 個網格點分布圖.圖1 中5 個網格點的坐標分別為：1(45°N，120°W)；2(30°N，105°W) ；3(40°N，90°W) ；4(30°S，25°E)；5(28°N，85°E).

圖1 選取的5 個網格點坐標示意圖

2 自相關預報原理和方法

在1 個月左右的時間尺度上，電離層參量的時間序列通常是平穩時間序列，在統計學上可將其視為具有非零自相關函數的隨機穩態過程的實現[1].假設某電離層參量小時值的時間序列為z(t)，其在t時刻的值可表示為前n個測量值的加權平均，如式(1)所示：

式(1)中的加權系數 λj通過以下方程組求解得到：

式中：μ 為拉格朗日系數；ρ(τ) 為z(t) 的自相關系數，其均值為0.對于1 h 分辨率的電離層時間序列，引入式(3)所示的擬合函數來進行加權系數的估計[1].

上述處理是從信號的角度將電離層參量的時間序列看成是一個周期信號疊加一個非周期性的擾動.對于周期性信號，我們熟知電離層受日照影響存在明顯的周日變化，根據周期信號的傅里葉級數理論，將自相關系數的周期成分用周期為24 h 的基波和2 次諧波的和來表示.而對于非周期的擾動信號，式(3)實際采用的是一個幅度受限的雙邊指數信號來模擬其變化，其中 τ0是時間常數.式(3)中的三個系數A24、A12、τ0通過最小二乘法來確定，然后，選取最近4 天的12 個最大的自相關系數代入式(2)求解，即可得到對應所選時間的權重系數λj[1].最后，將權重系數代入式(1)即可得到t時刻的預測值.

3 結果

3.1 自相關與自回歸預報方法比較

電離層TEC 觀測序列是一種典型的周期性時間序列，根據時間序列預測理論，對于數據量較少的情況，通常采用自相關或自回歸模型進行預報.本文首先對兩類模型進行了比較，自回歸模型選用了常用的ARIMA 模型.圖2 給出了一個結果示例，所用數據為網格點5 的數據，預報時間為2020 年2 月1 日至2 月10 日.從圖2 可以看出，自相關和ARIMA 預報結果均與觀測結果吻合較好，很好地再現了電離層TEC 的周日變化特征.

圖2 自相關與自回歸預報方法結果比較示例

進一步地，我們采用均方根偏差(RMSE)和歸一化均方根偏差(NRMSE)來衡量預報絕對誤差和相對誤差，計算公式為：

上兩式中：εRMSE和 εNRMSE分別為均方根偏差和歸一化均方根偏差；和zi分別為第i個點位的預測值和觀測值；N為總的點位數.

經估算，圖2 所示結果中，自相關預報的絕對誤差約為0.716 5 TECU，相對誤差約為10.87%；ARIMA預報絕對誤差和相對誤差分別為0.876 2 TECU 和13.29%.此外，我們還利用網格點1～4 的數據對兩種方法進行了比較分析，計算結果顯示，對于網格點1～4，自相關方法的絕對誤差分別為0.414 5 TECU、0.654 3 TECU、0.515 2 TECU、1.182 1 TECU，ARIMA方法絕對誤差分別為0.436 9 TECU、0.607 8 TECU、0.522 6 TECU、1.230 1 TECU；自相關方法的相對誤差分別為8.3%、9.52%、8.52%、13.28%，ARIMA 方法相對誤差分別為8.75%、8.84%、8.64%、13.82%.

從上述結果對比可以看出，總體上自相關方法的絕對誤差和相對誤差均略小于ARIMA 方法.此外，我們還對兩種預報方法需要花費的時間進行了測算.結果顯示，預報如圖2 所示的10 天結果，在相同的計算環境下，自相關預報時間約為10 s，ARIMA 預報時間約為50 s.還需指出的是，本文采用的ARIMA 模型使用之前已經定好了階數，后續的預測均使用固定階數，若采用動態自動定階，則預報時間將大大增加.

綜上分析，自相關預報方法預報誤差相對略優于ARIMA 方法，且自相關預報所需時間比ARIMA少得多.因此，本文后續將主要分析自相關方法.如前所述，自相關方法雖然已經取得了一些應用，但其依然有改進的空間，3.2 小節將在誤差分析的基礎上嘗試對自相關方法進行參數優化.

3.2 自相關誤差分析和參數優化

對于自相關預報方法，在預報參數選取方面，目前普遍的做法是，選取最近4 天相關系數最高的12 個時間點，然后通過式(2)解出相應的權重系數.由于這兩個參數是早期經驗值，其對預報誤差的影響目前還沒有文獻公開報道.為弄清楚這一影響，我們對具體預報過程和誤差產生原因進行了深入分析.

圖3 給出了參數選取4 天和12 點時的一個預報結果示例.圖3 中的TEC 數據來自網格點4 (30°S，25°E)，預測日期為2020 年2 月8 日，所用的歷史數據為2020 年1 月9 日至2 月7 日共30 天數據，分辨率為1 h.圖3 中黑色圓點表示實測值，2 月8 日的紅色圓點表示預測值.我們以第8 個預測點(如圖中黑色箭頭所指)為例，給出了用來計算該點預測值的歷史數據分布(如圖中藍色圓圈所示)，同時在每個歷史點位左側標出了權重系數.根據自相關預測理論，所選取的12 個采樣點是前4 天相關系數最大的采樣點，這些點主要是前4 天相同時刻和相鄰時刻的采樣點，且相同時刻的權重最大(為0.093)，而相鄰時刻采樣點的權重均為0.079.

圖3 2020 年2 月8 日“4+12”參數方案時的預報結果(紅色圓點)以及2 月4 至8 日實測結果(黑色圓點)

對比圖3 中2 月8 日的預報結果和實測結果可以看出，雖然兩條曲線的變化趨勢基本一致，但幾乎所有點位的預測值均要大于觀測值.對黑色箭頭所示第8 個點位，其觀測值為9.264 TECU，預測值為10.292 TECU，誤差達到1.029 TECU.從第8 個預測點所用歷史點位信息來看，距離預測點最遠的2 月4 日的兩個被選用的點位對誤差影響較大，這兩個采樣點數值明顯高于2 月5 日至2 月8 日其他點位的值，很可能是造成該點預測值明顯偏離實測值的主要原因.

圖3 結果反映出距離預測點較遠的歷史數據可能與當前時刻觀測值的偏差較大，從而影響預報精度.據此，我們推測，縮短參與預報的歷史數據天數，有可能能夠減少預報誤差.因此，我們將參與預測的天數減少為3 天，相應的點數減為9 個，對上述選例重新開展自相關預測，結果如圖4 所示.由圖4 可知，整體上，修改參數后的預測值與觀測值更接近.對黑色箭頭所示第8 個點位，其觀測值為9.264 TECU，預測值為9.863 TECU，誤差減小為0.599 TECU.

圖4 2020 年2 月8 日“3+9”參數方案時的預報結果(紅色圓點)以及2 月4 至8 日實測結果(黑色圓點)

此外，我們還分析了天數增加或減小時的結果，圖5 給出了不同參數設置下的所有24 個點位預測誤差曲線比較，圖5 中黑色圓圈代表參數設置為“4+12”時的誤差曲線(對應圖3 所示結果)，紅色圓圈代表參數設置為“3+9”時的誤差曲線(對應圖4 所示結果).對比上述兩條誤差曲線可以看出，除12 UT、13 UT、14 UT 三個點位外，其余點位的預測誤差都有所減小，尤其是10 UT 之前的點位減少更明顯.總體上，原方法“4+12”平均絕對預報誤差為0.8 TECU，參數調整為“3+9”后平均絕對預報誤差為0.63 TECU，平均絕對誤差減小了0.17 TECU.此外，圖5 還給出了參數設置分別為“5+12”(紫色圓圈)和“2+6”(綠色圓圈)時的預報誤差曲線，但是從對比結果看，兩種方案均存在誤差明顯增加的點位，而誤差減少的點位較少.經估計，圖5 中“5+12”案的平均絕對誤差為0.885 TECU，“2+6”方案的平均絕對誤差為0.892 TECU，相對于原方法“4+12”誤差還分別增加了0.085 TECU和0.092 TECU.

圖5 2020 年2 月8 日不同參數方案預報誤差比較

上述結果表明，在TEC 自相關短期預報方法中，參與預測的天數設置為3 天(相應的點位設置為9)可能是更優的方案.這從側面說明電離層當前的狀態可能主要取決于前3 天的狀態，或者說前3 天的電離層狀態對當前狀態的關聯最大.時間跨度越大，關聯性越小，引起誤差的可能性越大.

為進一步檢驗“3+9”參數方案的預報精度改善效果，我們對圖1 中所有5 個觀測點的數據進行了分析，且每個觀測點從2 月1 日至3 月1 日進行30 次1～24 h 短期預報.對每個觀測站，我們利用30×24 個時間點的觀測值和預測值來計算該站的均方根誤差(RMSE)和歸一化均方根偏差(NRMSE)，以進行兩種參數設置方案的誤差比較.兩種方案的綜合誤差比較結果如表1 所示.從表1 結果可以看出，對于本文隨機選取的5 個觀測站，“3+9”綜合預報誤差均要優于“4+12”方案.將5 個觀測站誤差值求平均值可得，所有站RMSE 平均值從0.825 TECU 降低至0.805 TECU，降低了0.02 TECU，NRMSE 平均值從11.756%降低至11.452%，降低了約0.3%.需要指出的是，由于本文所選用TEC 數據的取值水平較低，所以改進效果從絕對值上看較小，但本文經過一定規模的統計分析顯示所有站點均有一定改進，因此我們認為這一改進方案是依然有一定參考價值的.

表1 不同觀測站兩種方案預測誤差綜合比較

4 結果與討論

本文首先比較了自相關和自回歸滑動平均(ARIMA)方法，結果顯示兩種方法均能較好地預測TEC周日變化.對本文預報結果進行了RMSE 和NRMSE估算，結果表明總體上自相關預報誤差略低于ARIMA預報誤差，且在相同計算環境下，自相關預報所花費的時間明顯少于ARIMA 模型.

其次，本文對自相關預報方法的誤差進行了深入分析，重點分析了實際參與加權的觀測值覆蓋天數和觀測值數量這兩個參數的選取對預報誤差的影響.從Madrigal 數據庫中任意選取了5 個網格點的觀測數據開展了預報試驗，對不同參數設置方案的預報誤差進行了對比分析.2020 年2 月8 日網格點4 (30°S，25°E)的個例分析結果表明，“3+9”方案的預報誤差相對最小，天數增加或減少(例如“2+6”、“4+12”和“5+12”等方案) 都會使誤差變大.進一步地，我們對5 個觀測站分別進行了30 天的1～24 h 短期預報試驗，綜合分析了3 600 個預測點的RMSE 和NRMSE，結果顯示：相比于傳統的“4+12”方案，“3+9”方案在5 個觀測站中均取得了更優的預報性能.具體地，所有站RMSE 平均值從0.825 TECU 降低至0.805 TECU，降低了0.02 TECU，NRMSE 平均值從11.756%降低至11.452%，降低了約0.3%.

因此，我們推測，TEC 時間序列的當前狀態可能主要與前3 天的狀態有關.當參與加權的觀測天數增加(大于3 天)時，內在關聯性小的觀測值參與加權運算，可能會出現較大的狀態差異，從而造成誤差增大.當參與加權的觀測天數減少時，可能會丟失主要的內在關聯信息，從而也引起誤差增大.但需要指出的是，本文提出的參數優化方案只是對傳統參數方案的一個微小的改進，還無法明顯改善電離層預報性能，但可作為工程實現的參考方案.

電離層隨時間的變化非常復雜，除了熟知的周日變化、季節變化和太陽活動性11 年周期變化外，還存在很多小尺度的時間變化，目前我們對這些小尺度時間變化規律和控制因素的認識還不足.若要顯著提高電離層短期預報性能，必須深入研究電離層參量隨時間短期變化的規律和內在控制機理，這將是我們下一步研究工作的重點.