基于深度學習的全球電離層TEC預測

張富彬 周晨 王成 趙家奇 劉祎 夏國臻 趙正予

（1. 武漢大學電離層實驗室，武漢 430072；2. 北京航空航天大學前沿科學技術創(chuàng)新研究院，北京 100191；3. 哈爾濱工業(yè)大學（深圳）空間科學與應用技術研究院，深圳 518055）

引 言

電離層總電子含量(total electron content，TEC)是電離層研究過程中的一個非常重要的參數[1]. 從理論研究意義上來說，分析電離層的活動以及電離層的不同時空尺度的變化規(guī)律，包括電離層的逐日變化、電離層擾動等現象，這些研究都可以基于TEC參數的獲取和研究[2]；從應用角度來看，TEC的建模和預報對于衛(wèi)星定位、導航、測控、時頻傳遞等星地鏈路電波傳播有著至關重要的意義[3-4].

電離層TEC的建模研究最早采用的是經驗模型方法. Rawer等人首次提出了國際電離層參考(international reference ionosphere，IRI)模型，它是統(tǒng)計預報模型[5]，被空間研究委員會(Committee on Space Research,COSPAR)和國際無線電科學聯(lián)盟(International Union of Radio Science, URSI)共同主導，而且此項目一直在不斷更新. 由于眾多學者的共同努力和新的建模技術的使用，IRI模型已經得到了很大的完善和改進[6-7]. 2015年，Aa等人提出使用卡爾曼濾波的數據同化算法對中國及其鄰近地區(qū)的TEC進行建模，極大地改進了經驗模型的預報效果[8].

人工神經網絡的出現，為TEC建模提供了新的研究方向. 人工神經網絡能夠學習不同變量之間的復雜性，通過大量的神經元把輸入的向量與輸出結果建立函數映射關系. 正是因為這個特性，神經網絡逐漸地被應用到電離層TEC的建模中[9-11]. 1997年，Cander等人開始利用神經網絡來分析預測非線性時間序列，使用多層感知機的方法對不同電離層參數進行建模[12]. 2011年，Habarulema等人針對南非地區(qū)的電離層TEC建模，提出了一個基于神經網絡的預測模型，該模型具有一個輸入層、兩個隱含層和一個輸出層[13]. 2018年，Song等人提出一種遺傳算法改進的神經網絡進行中國區(qū)域TEC的預測，它克服了傳統(tǒng)神經網絡容易陷入局部最優(yōu)的問題，預測誤差在5～8 TECU左右[14]. 2020年，Cesaroni等人使用神經網絡的方法對全球范圍內不同地理位置的臺站TEC進行預測，誤差在3～5 TECU左右[15].

深度學習是一種特定類型的機器學習，它是在傳統(tǒng)人工神經網絡基礎上的升級，可以使用訓練網絡自動提取數據的特征. 它的網絡層數相比傳統(tǒng)的神經網絡更多，從理論上可以實現任意函數的擬合.2017年，Sun等人采用長短時記憶(long short-term memory, LSTM)網絡的方法對北京地區(qū)的TEC進行預測，誤差在3～4 TECU左右[16]. 2019年，Chen等人將深度卷積生成對抗網絡(deep convolutional generative adversarial networks, DCGAN)應用到空間物理領域的TEC圖片補全，通過添加一個額外的判決器以及使用IGS (International GNSS Service)提供的TEC地圖作為先驗知識，提出了一種基于正則化的DCGAN[17]. 文獻[18]對電離層TEC預測的滑動平均自相關法、LSTM以及encode-to-decode(即編碼器-解碼器)[19-20]三種方法進行比較，結果顯示在地磁活動頻繁時期模型的預測效果都會變差，但是LSTM在磁暴時期會取得更好的效果.

本文采用深度學習的基本方法來進行電離層TEC的預測，時間精度相比以往的大部分論文來說更高，為1 h. 預測的TEC空間尺度是全球，空間經度分辨率為5°，緯度分辨率為2.5°，即每個TEC map的形狀尺寸為71行、73列. 神經網絡采用的是編碼器-解碼器結構[19]，其中，使用卷積層并且配合池化層來進行特征提取和降維便于加快網絡訓練和減小數據規(guī)模，中間層使用基于卷積操作的LSTM網絡——ConvLSTM網絡來實現時間序列預測的功能[21].

1 數據和實驗方法

1.1 數據來源和預處理

電離層TEC數據集為IGS提供的標準IONEX(IONosphere map EXchange format)文件格式生成的全球電離層TEC數據. 時間分辨率是1 h；空間經度范圍從西經180°到東經180°，分辨率是5°；緯度范圍從北緯87.5°到南緯87.5°，分辨率是2.5°. 數據集的時間范圍從1998-01-01—2018-12-31，每一天的數據構成一個文件. 根據以上信息即可得到全球TEC map網格點的規(guī)模是71×73.

另外，根據電離層TEC數據存在周期性的日變化規(guī)律，我們可以把一天的TEC map數據，即24個按時間順序的TEC map合并作為一個單個樣本的輸入，該樣本的輸出為下一天的24個TEC map. 按照這種結構對所有的原始IONEX格式文件進行提取、拼接和重新設置矩陣維度信息等操作，我們可以得到一個五維張量，它的各個維度的大小是(number, 24,71, 73, 1). 其中：number代表了數據集的總樣本個數；24代表每個樣本包含24個小時，即24張如圖1的TEC map；71和73則是每個整點時刻的TEC map的維度；1表示該map的顏色通道數目，因為每個經緯度的TEC是一個一維的數值.

圖1 2018-07-07T18:00的全球TEC mapFig. 1 The global TEC map of 2018-07-07T18:00

我們可以很容易地使用Python代碼把一個時間序列數據轉變成一組有監(jiān)督學習的數據樣本，數據序列的變化如圖2所示，其中的每個數字代表了一整天24個整數時刻的TEC數據，是一個(24,71, 73, 1)形狀的四維張量，且這些數據是按日期順序遞增排列的. 例如，圖中的數字1,2,3可以表示為2017-01-01—03的TEC數據.

圖2 時間序列數據轉換為有監(jiān)督數據的過程Fig. 2 Transform a time series data into a set of supervised learning data samples

本文選取1998-01-01T0:00—2017-12-31T23:00的數據作為訓練集(占據總數的95%)和驗證集(占據總數的5%)，2018-01-01T0:00—12-31T23:00的數據作為測試集，數據集劃分和每個樣本的輸入及輸出維度具體如表1所示.

表1 數據集劃分及樣本維度Tab. 1 The scale and division of data set

1.2 神經網絡模型的選取

本文采用深度學習的方法來對全球電離層的TEC map進行預測，使用經典的編碼器-解碼器框架來搭建深度神經網絡[19]. 這個框架也被稱為自編解碼器，它主要的目的就是將輸入網絡的向量通過編碼器轉換為維度更低的中間向量，再利用解碼器把這個中間變量轉換為輸出變量. 在深度學習中，編碼器實現了特征降維的功能，相當于一個特征提取器；解碼器會把編碼的結果還原為需要的數據. 在圖像處理領域，如手寫數字的識別，自編解碼器有著非常廣泛的應用.

編碼器-解碼器的典型結構如圖3所示. 其中：xi表示輸入此框架的原始數據展開成一維后的第i個分量，pi表示經過編碼器結構后得到的向量的第i個分量，yi表示通過此框架后輸出的向量的第i個分量.

圖3 編碼器-解碼器模型結構Fig. 3 The structure of encoder-decoder model

目前深度學習常用的網絡模型是卷積神經網絡(convolutional neural networks, CNN)和循環(huán)神經網絡(recurrent neural network, RNN)，從它們的實現原理來看，前者更擅長提取空間結構信息，后者比較適合處理時間序列數據. RNN模型的眾多實現中包括獨立循環(huán)神經網絡(independently recurrent neural networks，IndRNN)、門控循環(huán)單元(gate drecurrent unit, GRU)、LSTM網絡等結構. LSTM網絡憑借它的性能和效率成為了使用最廣泛的一個網絡模型. 但是單純的LSTM網絡對于同時存在時間和空間特征的電離層TEC數據就存在一定的局限性，所以本文采用了施行健等人[21]在2015年提出的ConvLSTM網絡模型，其原理如圖4所示. 它通過將矩陣乘法更改為元素卷積的方式實現了對空間局部特征的提取.

圖4 ConvLSTM網絡原理框圖Fig. 4 The functional block diagram of ConvLSTM網絡

LSTM網絡使用輸入門、輸出門、遺忘門以及與隱藏狀態(tài)相同的記憶細胞來記錄額外的信息. 所有門的輸入均由當前時刻t的輸入Xt與上一時刻t?1的隱藏狀態(tài)Ht?1組成，用W與b表示各部分的權重和偏置系數，Ct表示t時刻的狀態(tài)，Ft為經過遺忘門后的輸出結果，It為經過輸入門后的輸出結果，Ot為經過輸出門后的結果，圖4中的×表示元素對應位置點乘，公式中使用符號⊙表示，?表示卷積，根據原理框圖有以下式子成立：

1.3 網絡模型的構建和網絡參數的選取

本文采用基于深度學習的神經網絡結構，主要可以分為三個部分：第一部分為編碼器，也就是卷積層和池化層的成對組合，它們可以將輸入的電離層TEC map進行下采樣，有效提取出數據的空間關鍵特征，同時池化操作也可以降低數據的維度，此操作可以減少神經網絡的訓練參數來加快訓練過程；第二部分是一層ConvLSTM網絡和一層批次標準化Batch Normalization，前者是實現時間序列預測的關鍵部分，可以提取出TEC map的時序特征，后者是對數據進行規(guī)范化，還可以防止過擬合；最后一個部分是根據編碼器-解碼器框架創(chuàng)建的與第一部分對應的卷積層和上采樣層的成對組合，即解碼器. 圖5是構建好的網絡模型結構.

圖5 本文搭建的模型結構Fig. 5 The structure of deep learning model in this paper

每一層的激活函數均采用elu函數，然后進行編譯. 編譯主要涉及損失函數、優(yōu)化器和指示列表三個參數. 本文調整參數多次試驗最終確定損失函數采用均方誤差(mean-square error, MSE)，優(yōu)化器使用改進后的自適應調整學習率的Adamax算法[22]. 訓練輪數為100，批次大小為32.

2 結果和分析

按照圖5中的網絡結構框圖，使用Python語言基于高度抽象的Keras庫編寫深度學習代碼(后端采用TensorFlow)，搭建神經網絡. 通過不斷地調整每一層的內核大小以及濾波器數目最終得到訓練好的模型，圖6為訓練過程中損失函數隨輪數的變化曲線，可以看出驗證集誤差和訓練集誤差都隨訓練輪數的增加而逐漸下降并最終趨于水平穩(wěn)定.

圖6 訓練誤差隨輪數的變化Fig. 6 The relationship between loss and epoch

本文采用每個TEC map的均方根誤差(root mean square error, RMSE)來對模型的性能進行評價，RMSE公式為

式中：N表示每一個TEC map的數據個數，即71×73；TECp表示預測得到的TEC map；TECr表示根據實測數據得到的TEC map；i表示此TEC map的第i個數據點. 另外，由于模型的輸入和輸出為某一天的TEC數據，即24個TEC map，所以我們定義模型的RMSE公式為

式中，ERMS,map,n表示第n個TEC map的RMSE.

下面討論在不同情況下本文模型對全球TEC分布圖的預測性能.

文中使用地磁活動指數Dst來判斷是否處于磁暴時期[23-24]，如表2所示.

表2 磁暴類型的劃分Tab. 2 The thresholds for storms

1)模型在訓練集上，且處于磁平靜時期(無磁暴)，提前1天的預測結果(圖7). 圖7(a)表示訓練集上，當輸入樣本為2017-12-13的TEC數據時，即輸入24×71×73×1的四維張量，模型預測得到的下一天的數據. 它與輸入數據的維度是相同的，其中包含24個時刻的TEC map，即從當天的0:00到23:00，圖中我們選擇10:00的TEC map進行展示. 2017-12-14T10:00的實際臺站觀測并處理得到的TEC數據如圖7(b)所示，使用公式(2)計算得到預測的RMSE為1.366 8，這表明預測的結果與實際TEC的數據誤差很小.

圖7 磁平靜時期在訓練集的提前1天的預測結果(RMSE=1.366 8)Fig. 7 The one-day-prediction of model during the geomagnetic quiet time on the training set(RMSE=1.366 8)

2)模型在訓練集上，且處于磁平靜時期，提前2天的預測結果(圖8). 模型在訓練過程中的輸入是前一天的24張TEC map數據，輸出是神經網絡預測得到的當天的24張TEC map數據. 在這里，我們把預測得到的這一天內的TEC map當作模型的輸入，再次使用模型來進行預測，此時我們就可以得到提前2天的預測結果. 當模型的輸入數據仍然是2017-12-13的TEC數據，根據這一方法，我們可以預測得到2017-12-15的數據. 圖8(a)是這個數據在18:00的TEC map，圖8(b)是實際觀測的TEC數據.

圖8 磁平靜時期在訓練集的提前2天的預測結果(RMSE=1.564 4)Fig. 8 The two-day-prediction of model during the geomagnetic quiet time on the training set(RMSE=1.564 4)

3)模型在訓練集上，且處于磁平靜時期，提前7天的預測結果(圖9). 其中輸入樣本仍為2017-12-13的TEC數據，模型預測得到2017-12-20的數據，圖9(a)是這個數據在22:00的TEC map，圖9(b)是實際觀測的TEC數據.

圖9 磁平靜時期在訓練集的提前7天的預測結果(RMSE=2.129 1)Fig. 9 The seven-day-prediction of model during the geomagnetic quiet time on the training set(RMSE=2.129 1)

使用已經訓練好的模型對整個測試集的數據進行提前1天的電離層TEC預測誤差分析，利用式(3)得到每次預測的RMSE，誤差的分布直方圖如圖10所示.

圖10 模型在測試集上的預測誤差分布直方圖Fig. 10 the histogram of prediction error on the test set

從圖10可以看出，模型在測試集上的誤差的均值為1.679 9 TECU，而且誤差的標準差為0.362 0 TECU，這表明模型的預測效果很穩(wěn)定. 圖中誤差在3.5 TECU的兩次預測是在中等磁暴和強磁暴的情況下得到的預測結果. 模型在弱磁暴和磁平靜時期的預測效果非常穩(wěn)定，而且誤差很小，下面給出了在不同地磁活動水平下模型的預測效果.

1)模型在測試集上，且處于磁平靜時期，提前2天的預測結果(圖11). 其中輸入為2018-12-15的TEC，預測得到2018-12-17的數據，圖11(a)是這個數據在5:00的TEC map，圖11(b)是實際觀測的TEC數據.

圖11 磁平靜時期在測試集的提前2天的預測結果(RMSE=1.426 0)Fig. 11 The two-day-prediction of model during the geomagnetic quiet time on the test set(RMSE=1.426 0)

2)模型在測試集上，且處于弱磁暴時期，提前1天的預測結果(圖12)，其中輸入樣本為2018-09-10的TEC數據，模型預測得到2018-09-11的數據，圖12(a)是這個數據在5:00的TEC map，圖12(b)是實際觀測的TEC數據.

圖12 弱磁暴時期在測試集的提前1天的預測結果(RMSE=2.204 7)Fig. 12 The one-day-prediction of model during the geomagnetic small storm on the test set(RMSE=2.204 7)

3)模型在測試集上，且處于中等磁暴時期，提前1天的預測結果(圖13). 其中輸入樣本為2018-04-19的TEC數據，模型預測得到2018-04-20的數據，圖13(a)是這個數據在5:00的TEC map，圖13(b)是實際觀測的TEC數據.

圖13 中等磁暴時期在測試集的提前1天的預測結果(RMSE=2.441 7)Fig. 13 The one-day-prediction of model during the geomagnetic moderate storm on the test set, RMSE=2.441 7

4)模型在測試集上，且處于強磁暴時期，提前1天的預測結果(圖14). 模型預測得到2018-08-26T5:00的TEC map為圖14(a)，圖14(b)是實際觀測的TEC數據.

圖14 強磁暴時期在測試集的提前1天的預測結果(RMSE=4.549 1)Fig. 14 The one-day-prediction of model during the geomagnetic intense storm on the test set(RMSE=4.549 1)

圖7～14表明模型在訓練集和測試集取得相近的效果. 另外，模型對于提前1天的TEC預測效果要好于提前2天和7天的預測效果，即隨著天數的增加，模型的預測效果會逐漸變差；模型在磁暴時期的預測相對于磁平靜時期誤差更大. 圖11～14是在不同的地磁活動條件下，模型對于提前1天預測得到的24個TEC map的其中一張與實際進行對比并計算誤差. 圖15是對這24張TEC map的每個時刻分別計算誤差，可以看出隨著磁暴的劇烈程度增加，模型的預測效果也會變差.

圖15 不同磁暴條件下提前1天預測的不同時刻的RMSEFig. 15 One-day prediction of RMSE at each hour in different storms

下面研究同一天不同經緯度位置下模型提前預測1～7天的結果，如圖16所示. 可以看出，模型在高緯度的預測效果要優(yōu)于低緯度的預測效果.

圖16 模型在不同緯度連續(xù)提前7天預測的RMSEFig. 16 Seven-day prediction of RMSE at different latitudes

接下來研究模型進行提前1天TEC預測時的誤差與當天不同時刻的關系. 圖17中的每條曲線是模型提前1天預測的RMSE，曲線上的點表示在特定時刻的誤差，mean_rmse表示圖中5次預測誤差的平均值. 這表明模型對于TEC數據的預測效果基本不會受到一天不同時間段的影響.

圖17 提前1天預測的不同時刻的RMSEFig. 17 One-day prediction of RMSE at every hour

3 結 論

針對全球電離層TEC的時空預測問題，本文提出了一種基于深度學習方法的全球電離層TEC預測模型，采用編-解碼器結構配合卷積優(yōu)化的長短時記憶網絡(ConvLSTM)實現了全球TEC數據的高時間精度和空間分辨率的預測. Song等人使用遺傳算法改進的神經網絡對中國地區(qū)TEC的預測誤差在5～8 TECU，Cesaroni等人使用神經網絡的預測誤差在3～5 TECU. 本文將全球范圍的TEC預測誤差下降到1.5 TECU左右，這一方面得益于CNN能夠代替人工對數據進行特征工程，適合處理具有局部關聯(lián)性和統(tǒng)計平穩(wěn)性的數據，可以提取得到樣本的空間特征；另一方面是運用ConvLSTM將卷積結合到LSTM網絡實現了對數據的時序預測功能. 兩者結合并選用合適的參數，最終得到了更好的預測效果.

通過研究模型在訓練集和測試集上以及在不同磁暴條件和不同緯度下的效果，我們可以得到以下結論：

1)模型的預測誤差隨著預測時間的增加而逐漸變大. 對于提前1天的TEC預測，RMSE基本穩(wěn)定在1.5 TECU以下；對于提前2天甚至7天的TEC預測也具有相當令人滿意的結果.

2)隨著磁暴的劇烈程度增加，模型的預測效果也會變差.

3)模型在高緯度下的預測效果要優(yōu)于低緯度的預測效果，而且對于TEC數據的預測效果基本不會受到一天不同時間段的影響.

但是，模型對于磁暴時期的TEC數據預測還存在較大誤差，需要進一步的研究工作：一方面，電離層與太陽風、磁層以及電離層不同緯度區(qū)之間有著復雜的耦合過程[25]，可以考慮將太陽風等參數加入模型來進行訓練；另一方面，目前除了ConvLSTM網絡以外還有Attention機制、GRU以及Transform等內容可在深度學習算法中使用，可以借助相關技術實現更高精度的TEC預測.