基于深度學(xué)習(xí)的多步預(yù)測方法在太陽風(fēng)速度預(yù)測上的應(yīng)用

毛凱舟，謝宗霞，孫彥茹，劉梓璇

（天津大學(xué)智能與計算學(xué)部，天津 300350）

0 引言

惡劣空間天氣對依賴于電磁環(huán)境的高科技系統(tǒng)具有重要影響。這些空間天氣包括太陽耀斑、日冕物質(zhì)拋射、太陽質(zhì)子事件以及電子暴、地磁暴等由于太陽的劇烈活動引起的短時間尺度的變化。在太陽系內(nèi)，空間天氣受到太陽風(fēng)的巨大影響。

太陽風(fēng)是一種連續(xù)存在、來自太陽并以200～800 km/s的速度運(yùn)動的高速帶電粒子流，對人類的生產(chǎn)生活有著重要影響：當(dāng)太陽風(fēng)掠過地球時，會使電磁場發(fā)生變化，引起地磁暴、電離層暴，并影響通信，特別是短波通信；對地面的電力網(wǎng)、管道和其他大型結(jié)構(gòu)發(fā)送強(qiáng)大元電荷，影響輸電、輸油、輸氣管線系統(tǒng)的安全；對衛(wèi)星的安全運(yùn)行也會產(chǎn)生影響；還可能引發(fā)火山爆發(fā)、地震災(zāi)害，等等。對近地環(huán)境太陽風(fēng)狀況進(jìn)行可靠的預(yù)測，旨在減輕空間天氣對人類生產(chǎn)生活的影響，保證近地空間區(qū)域的航天活動安全，導(dǎo)航與無線電系統(tǒng)的正常工作，以及地面重要電力網(wǎng)絡(luò)系統(tǒng)的安全。

目前，對太陽風(fēng)的各項參數(shù)進(jìn)行預(yù)測的技術(shù)和方法主要可以分為3類：1）基于物理的磁流體動力學(xué)（Magnetohydro-dynamics,MHD）模型的方法，通過檢測太陽表面日冕面積的大小，判斷當(dāng)前的太陽狀態(tài)[1-4]；2）基于人工觀測的方法，由大量的觀測員對太陽的活動進(jìn)行實時監(jiān)測，通過他們的經(jīng)驗知識對太陽每個時刻的狀態(tài)進(jìn)行判斷[5-7]；3）機(jī)器學(xué)習(xí)方法，采用支持向量機(jī)（support vector machine,SVM），由人工神經(jīng)網(wǎng)絡(luò)（artificial neural network,ANN）對數(shù)據(jù)進(jìn)行訓(xùn)練來判斷太陽狀態(tài)。

由于太陽活動的復(fù)雜性和災(zāi)害性空間天氣事件預(yù)報時效性的迫切需求，需要尋找不同于傳統(tǒng)MHD模型與經(jīng)驗?zāi)Ｐ偷奶幚硗緩?，ANN方法正是能夠適應(yīng)新需求的一種技術(shù)手段。Liu D D等[8]使用SVM算法預(yù)測了太陽風(fēng)速度，但他們的預(yù)測僅能提前3 h。Wintoft 等[9]開發(fā)了一種結(jié)合半經(jīng)驗PFSS模型和ANN 的混合智能系統(tǒng)，可以預(yù)測太陽風(fēng)的日平均速度。Jones等[10]研究發(fā)現(xiàn)，使用基于圖像的信息可以改善全球日冕磁場模型。2018年，楊易和沈芳等[11]開發(fā)了一種基于ANN技術(shù)的、利用更多觀測數(shù)據(jù)驅(qū)動的近地太陽風(fēng)速度預(yù)報模型，稱為Hybrid Intelligent Source Surface（HISS）模型；他們對該模型2007年—2016年（包含1個完整太陽活動周的極小、上升、極大、下降全部4 個活動相）太陽風(fēng)速度的計算結(jié)果與近地觀測數(shù)據(jù)進(jìn)行了比較以及誤差分析和高速流事件分析，結(jié)果顯示模型計算結(jié)果與觀測數(shù)據(jù)較為一致。Upendran 等[12]使用預(yù)訓(xùn)練的ConvNet 網(wǎng)絡(luò)結(jié)合LSTM網(wǎng)絡(luò)開發(fā)出優(yōu)于基準(zhǔn)持久模型和自回歸模型的ANN模型。

本文希望基于ANN 方法，利用多步預(yù)測思想進(jìn)行太陽風(fēng)速度預(yù)測，擬使用多元時間序列數(shù)據(jù)，并基于長短時記憶（long short-term memory networks,LSTM）網(wǎng)絡(luò)開發(fā)單一LSTM模型、編解碼器LSTM模型和CNN-LSTM 編解碼器模型3種深度學(xué)習(xí)預(yù)測模型，用以對未來24、48、72、96 h 的太陽風(fēng)速度進(jìn)行預(yù)測。對于每種模型，除了直接對目標(biāo)時間點進(jìn)行預(yù)測，還引入多步預(yù)測的方法對包含目標(biāo)時間點的連續(xù)多個時間步進(jìn)行預(yù)測。

1 多元時間序列數(shù)據(jù)及評價指標(biāo)

1.1 數(shù)據(jù)選擇

1.1.1 時間序列數(shù)據(jù)

NASA 發(fā)布的OMNI數(shù)據(jù)集（https://omniweb.gsfc.nasa.gov/index.html）包括：

1）低分辨率OMNI（low resolution OMNI, LRO）數(shù)據(jù)，包括每小時的近地太陽風(fēng)磁場、等離子體數(shù)據(jù)，高能質(zhì)子（能量超過1～60 MeV）通量，以及地磁和太陽活動指數(shù)等。除每小時分辨率外，LRO數(shù)據(jù)還有每日和每27天2種分辨率。

2）高分辨率OMNI（high resolution OMNI）數(shù)據(jù)，包括1 min 和5 min 的地球弓形激波鼻處的太陽風(fēng)磁場和等離子體數(shù)據(jù)，以及地磁活動指數(shù)等。

本文選用的數(shù)據(jù)集為2011年—2017年、分辨率為1 h 的LRO數(shù)據(jù)集，包含小時分辨率的太陽風(fēng)磁場和等離子體數(shù)據(jù)。這些數(shù)據(jù)來自日地拉格朗日點1（L1）附近軌道上的許多航天器觀測數(shù)據(jù)。

預(yù)測模型的目標(biāo)輸出是在L1處的每小時平均太陽風(fēng)速度。由于每小時的風(fēng)速變化不大，所以使用每小時平均值，并且平均值的變化表明了風(fēng)速值的不確定性。在OMNIWEB檔案庫中（可在https://omniweb上在線獲得），將變化量（或標(biāo)準(zhǔn)偏差σ）計算為1天中每小時測量的變化量。

1.1.2 圖像數(shù)據(jù)

太陽動力學(xué)天文臺（Solar Dynamics Observatory,SDO）是NASA 的“與星共棲”（LWS）計劃發(fā)起的第一個任務(wù)，旨在了解太陽變化的原因及其對地球的影響。SDO通過在小范圍的時空和同時在多個波長下研究太陽大氣來了解太陽對地球和近地空間的影響，其目標(biāo)是通過確定太陽磁場的產(chǎn)生和結(jié)構(gòu)方式以及如何將存儲的磁能轉(zhuǎn)換并釋放到日光層和地球空間中，來理解太陽活動（包括太陽風(fēng)、高能粒子和太陽輻照度等變化）及其對空間天氣的影響，并推動對空間天氣的預(yù)測。

SDO任務(wù)搭載了大氣成像組件（AIA）、EUV 可變性實驗（EVE）、地震和磁成像儀（HMI）3個科學(xué)載荷。這3臺儀器將同時觀測太陽，覆蓋了解太陽變化所需的整個測量范圍，以期測量表征太陽如何變化以及為什么變化。

AIA 是一個四望遠(yuǎn)鏡陣列，可觀測可見光、UV 和EUV 波長下的整個太陽，對多種波長的太陽大氣成像，以將太陽表面變化與其內(nèi)部變化聯(lián)系起來。由SDO整理的機(jī)器學(xué)習(xí)數(shù)據(jù)集（以下稱為SDOML）已公開可用，其中的AIA 圖像已被重新采樣到512×512像素的網(wǎng)格上，并以Python NumPy格式存儲為二進(jìn)制數(shù)組，數(shù)據(jù)每10 s包含10 個波長的太陽圖像。本文的訓(xùn)練和測試數(shù)據(jù)包括每小時UTC拍攝的AIA 太陽圖像（見圖1），所選圖像構(gòu)成了整個觀察時間的代表；如果不存在該小時處的圖像（很多天都是這種情況），則將當(dāng)天最接近此時間點的圖像用作該小時圖像的代替。

圖1 AIA 太陽圖像數(shù)據(jù)Fig.1 AIA solar image data

1.2 數(shù)據(jù)預(yù)處理

1.2.1 輸入特征

1）引入冕洞特征——將0131波段的AIA 太陽表面圖像轉(zhuǎn)換為灰度圖，計算圖中像素值小于15的像素值之和，代表冕洞的影響，記為area。

2）引入行星際日冕物質(zhì)拋射（interplanetary coronal massejection,ICME）因素——將ICME列表二值化為一個新特征，有日冕物質(zhì)拋射事件對應(yīng)時刻的特征值為1，其他時刻特征值為0。

3）除了area、ICME 事件之外，選擇5個影響因素一起作為輸入——bulk speed（歷史太陽風(fēng)速度）、proton density（質(zhì)子密度）、proton temperature（質(zhì)子溫度）、flow pressure（流壓力）和σB（場矢量的均方根標(biāo)準(zhǔn)偏差）。其中σB是由分量平均值的標(biāo)準(zhǔn)偏差形成的向量的長度，即

1.2.2 異常值處理

1）對于數(shù)據(jù)中存在的異常值，采用三階B樣條曲線插值的方法進(jìn)行填充。B樣條曲線的總方程為

式中：Pi為控制曲線的特征點；Fi,k(t)為k階B樣條基函數(shù)。

三階B樣條曲線方程中的基函數(shù)為

2）圖像的異常值處理——圖像的異常值多為缺失或者拍攝不清晰，可采用高頻率數(shù)據(jù)集中臨近時間的圖像作為代替補(bǔ)充。

1.3 數(shù)據(jù)標(biāo)準(zhǔn)化

不同評價指標(biāo)往往具有不同的量綱，不利于數(shù)據(jù)分析。為此，需要進(jìn)行數(shù)據(jù)標(biāo)準(zhǔn)化處理，使不同數(shù)據(jù)指標(biāo)間具有可比性；且原始數(shù)據(jù)經(jīng)過標(biāo)準(zhǔn)化處理后，可使各指標(biāo)處于同一數(shù)量級，便于進(jìn)行綜合對比評價。

本文的數(shù)據(jù)標(biāo)準(zhǔn)化采用了標(biāo)準(zhǔn)差標(biāo)準(zhǔn)化，根據(jù)原始數(shù)據(jù)的均值和標(biāo)準(zhǔn)差進(jìn)行數(shù)據(jù)標(biāo)準(zhǔn)化，經(jīng)過處理的數(shù)據(jù)符合標(biāo)準(zhǔn)正態(tài)分布，即均值為0，方差為1。標(biāo)準(zhǔn)化處理公式為

式中：x*和x分別為標(biāo)準(zhǔn)化的數(shù)據(jù)和原始數(shù)據(jù)；μ為所有樣本數(shù)據(jù)的均值；σ為所有樣本數(shù)據(jù)的標(biāo)準(zhǔn)差，計算時對每個屬性/每列分別進(jìn)行。

1.4 數(shù)據(jù)集劃分

本文實驗采用2011年—2017年的數(shù)據(jù)，其中，2011年—2015年的數(shù)據(jù)作為訓(xùn)練集，共有43 824個；2016年的數(shù)據(jù)作為驗證集，共有8784個；2017年的數(shù)據(jù)作為測試集，共有8760個。數(shù)據(jù)構(gòu)建完成后的結(jié)構(gòu)如表1所示。

表1 太陽風(fēng)時間序列數(shù)據(jù)集劃分Table 1 Partitioning of the solar wind time series data sets

1.5 評價指標(biāo)

為明確量化擬合的優(yōu)劣，定義了3個指標(biāo)來評估擬合優(yōu)度——均方根誤差（RMSE）、平均絕對誤差（MAE）和相關(guān)系數(shù)（CC）。它們的定義如下：

這些指標(biāo)從不同的角度評估模型，RMSE 和MAE表征測量誤差，CC表征模型預(yù)測速度和觀測速度之間的擬合程度。因此，RMSE 和MAE 越小，CC越大，則表明預(yù)測效果越好。

2 預(yù)測方法

2.1 多步預(yù)測

多步預(yù)測又稱為長期預(yù)測，其含義是運(yùn)用歷史的時間序列[y1,···,yn]預(yù)測之后的t步序列[ym+1,···,ym+t]，其中：n表示已觀測數(shù)據(jù)的個數(shù)；m是≥n的一個時間點；t＞1，表示預(yù)測范圍。多步預(yù)測方法的有效性已在機(jī)器學(xué)習(xí)領(lǐng)域得到驗證以及快速發(fā)展和完善[13-18]。

本文將多步預(yù)測的思想應(yīng)用于太陽風(fēng)速度預(yù)測，并采用一種較好的策略——多目標(biāo)輸出的多步預(yù)測策略，進(jìn)行后續(xù)的建模和實驗，遍歷時間步并將數(shù)據(jù)劃分為重疊的窗口，每次迭代都前進(jìn)1個時間步，并預(yù)測包含目標(biāo)時間點的連續(xù)多個時間步的太陽風(fēng)速度。該方法的單步預(yù)測和多步預(yù)測模型的輸入和輸出可以表示為

相比于單步預(yù)測，采用多步預(yù)測的方法在模型訓(xùn)練過程中添加了輸出端太陽風(fēng)速度的連續(xù)變化信息，在預(yù)測過程中攜帶了目標(biāo)時間步之前的預(yù)測信息。因此，我們預(yù)期多步預(yù)測能夠?qū)崿F(xiàn)更好的預(yù)測效果。

2.2 基于LSTM神經(jīng)網(wǎng)絡(luò)的太陽風(fēng)速度預(yù)測模型

循環(huán)神經(jīng)網(wǎng)絡(luò)（recurrent neural networks,RNN）用前一個時間步的輸出作為其后續(xù)時間步的輸入，允許模型根據(jù)當(dāng)前時間步的輸入和前一個時間步的輸出的直接知識來共同決定模型下一步的操作。

LSTM神經(jīng)網(wǎng)絡(luò)可能是一種最成功、應(yīng)用最廣泛的RNN。LSTM神經(jīng)網(wǎng)絡(luò)克服了RNN 所面臨的長依賴關(guān)系，得到了穩(wěn)定可靠的模型；除了利用前一個時間步的輸出的循環(huán)連接之外，LSTM 神經(jīng)網(wǎng)絡(luò)還有一個內(nèi)部存儲，其操作類似于一個局部變量，允許它們在輸入序列上累積狀態(tài)。因此，當(dāng)涉及多步時間序列預(yù)測時，LSTM神經(jīng)網(wǎng)絡(luò)具備很多優(yōu)勢：1）LSTM 神經(jīng)網(wǎng)絡(luò)是一種循環(huán)網(wǎng)絡(luò)，被設(shè)計成以序列數(shù)據(jù)作為輸入，不同于其他必須以滯后觀測作為輸入特征的模型。2）LSTM神經(jīng)網(wǎng)絡(luò)直接支持多元輸入的多個并行輸入序列，而不像其他模型中多元輸入以扁平結(jié)構(gòu)表示。3）與其他神經(jīng)網(wǎng)絡(luò)一樣，LSTM 神經(jīng)網(wǎng)絡(luò)能夠?qū)⑤斎霐?shù)據(jù)直接映射到一個可能表示多個輸出時間步長的輸出向量。

下面介紹本文基于LSTM神經(jīng)網(wǎng)絡(luò)構(gòu)建的3種預(yù)測模型。

2.2.1 單一LSTM模型

單一LSTM 模型使用1個普通的LSTM層，以20個LSTM 神經(jīng)元作為參數(shù)；LSTM 層之后經(jīng)過1 個dropout 層減小模型的過擬合；然后是1個具有10個節(jié)點的全連接層（dense_1），用于解釋LSTM層習(xí)得的特性；最后，以另1個全連接層（dense_2）作為輸出層，用于直接輸出1個預(yù)測向量，該向量的維度與訓(xùn)練輸出的維度一致，在輸出序列中每個時間步對應(yīng)1個元素。圖2展示了24 h 多步預(yù)測實驗中連續(xù)5步預(yù)測所用的單一LSTM 模型結(jié)構(gòu)。

圖2 單一LSTM 模型結(jié)構(gòu)Fig.2 Structure of the single LSTM model

2.2.2 編解碼器LSTM 模型

編解碼器LSTM模型由2個子模型組成——用于讀取和編碼輸入序列的編碼器；讀取編碼的輸入序列并對輸出序列中的每個元素進(jìn)行1步預(yù)測的解碼器。該模型首先定義1個包含20個單元的LSTM 層（lstm_1）作為編碼器模型，用于讀取輸入序列并輸出1個包含20個元素的向量（每個單元1 個輸出），以便從輸入序列中獲取特征信息；然后，對輸入序列的內(nèi)部表示進(jìn)行多次重復(fù)，對于輸出序列中的每個時間步重復(fù)1次，這個向量序列將被呈現(xiàn)給LSTM解碼器；之后，使用repeat vector 將編碼器的輸出序列的最后1個作為內(nèi)部表示復(fù)制N份，作為LSTM 解碼器的N次輸入（N代表N步預(yù)測）；定義另1個包含20個單元的LSTM層（lstm_2）作為解碼器，輸出整個序列而不僅僅是輸出序列的末尾；經(jīng)過dropout 層后，使用全連接層（dense_1和dense_2）來解釋解碼器輸出序列中的每個時間步長，并輸出預(yù)測。這意味著須對輸出序列中的每個步驟使用相同的層。為了實現(xiàn)這一點，將解釋層和輸出層封裝在1個time_distributed 包裝器中，該包裝器允許解碼器每次執(zhí)行操作時都使用所封裝的層。圖3展示了24 h 多步預(yù)測實驗中連續(xù)5步預(yù)測所用的編解碼器LSTM 模型結(jié)構(gòu)。

圖3 編解碼器LSTM 模型結(jié)構(gòu)Fig.3 Structure of the codec LSTM model

2.2.3 CNN-LSTM 編解碼器模型

本文將使用CNN 和LSTM的混合模型稱為CNN-LSTM編解碼器模型。該模型以一維CNN作為編解碼器結(jié)構(gòu)中的編碼器，能夠直接支持序列輸入并自動學(xué)習(xí)顯著特征；使用LSTM解碼器解釋這些特征。CNN-LSTM編解碼器模型將編碼器定義為一個簡單而有效的CNN 架構(gòu)，由2個卷積層（convld_1和convld_2）以及1個最大池化層（max_pooling ld）組成，然后將結(jié)果扁平化（flatten）。第1個卷積層讀取輸入序列并將結(jié)果投射到特征圖上；第2個卷積層對第1層創(chuàng)建的特征圖執(zhí)行相同的操作，以期放大所有顯著的特征。每個卷積層使用16個特征圖，讀取卷積核大小為3個時間步長的輸入序列。最大池化層通過保留最大信號值的1/4來簡化特征映射。在池化層之后提取的特征映射被扁平化為一個長向量，可以用作解碼器的輸入。解碼器與編解碼器LSTM模型定義的解碼器相同。圖4展示了24 h 多步預(yù)測實驗中連續(xù)5步預(yù)測所用的CNN-LSTM 編解碼器模型結(jié)構(gòu)。

圖4 CNN-LSTM 編解碼器模型結(jié)構(gòu)Fig.4 Structure of the CNN-LSTM codec model

3 預(yù)測實驗結(jié)果分析

本文使用2.2節(jié)所述的3種模型對未來24、48、72、96 h 的太陽風(fēng)速度進(jìn)行預(yù)測；并且對于每個預(yù)測時間點和每種模型，分別進(jìn)行單步預(yù)測和多步預(yù)測實驗。在單步預(yù)測中，一系列訓(xùn)練輸入對應(yīng)未來T時間步后的1個訓(xùn)練輸出；在N步預(yù)測中，一系列訓(xùn)練輸入對應(yīng)未來連續(xù)N個訓(xùn)練輸出VT-N+1,VT-N+2,…,VT。每個預(yù)測時間點和每種模型的實驗重復(fù)3次，記錄下最優(yōu)結(jié)果。

24 h 預(yù)測的實驗結(jié)果如表2所示，多步預(yù)測為5步預(yù)測。為了說明模型的有效性，以24 h 持久模型作為基準(zhǔn)模型——預(yù)測結(jié)果為24 h 之前的太陽風(fēng)速度。從實驗結(jié)果可以看出，3種模型都顯示了較好的預(yù)測性能，且均表現(xiàn)為多步預(yù)測效果優(yōu)于單步預(yù)測效果。這可以歸因于多步預(yù)測攜帶了目標(biāo)時間步之前的預(yù)測信息。3種模型中，最簡單的單一LSTM 模型的預(yù)測效果最佳，其多步預(yù)測結(jié)果與觀測數(shù)據(jù)的相關(guān)系數(shù)為0.789，RMSE 為62.469 km/s，MAE為46.476 km/s。

表2 24 h預(yù)測實驗結(jié)果Table 2 Predicted experimental results at 24 h

48 h 預(yù)測的實驗結(jié)果如表3所示，多步預(yù)測為5步預(yù)測，3種模型均表現(xiàn)為多步預(yù)測效果優(yōu)于單步預(yù)測效果，其中編解碼器LSTM 模型的預(yù)測效果最佳，其多步預(yù)測結(jié)果與觀測數(shù)據(jù)的相關(guān)系數(shù)為0.522，RMSE為87.310 km/s，MAE 為67.976 km/s。

表3 48 h預(yù)測實驗結(jié)果Table 3 Predicted experimental resultsat 48 h

72 h 和96 h 預(yù)測的實驗中，為了避免預(yù)測步長過大，對數(shù)據(jù)的尺度進(jìn)行放大，由原先每小時1個數(shù)據(jù)改為每天（24 h）1個數(shù)據(jù)，取每天00：00時的數(shù)據(jù)進(jìn)行訓(xùn)練和預(yù)測。這種情況下，多步預(yù)測為3步預(yù)測。72 h 和96 h 預(yù)測的實驗結(jié)果分別如表4和表5所示。從實驗結(jié)果可以看出：3種模型均表現(xiàn)為多步預(yù)測效果優(yōu)于單步預(yù)測效果；對于遠(yuǎn)期的預(yù)測，CNN-LSTM 編解碼器模型依然有效，且預(yù)測效果最佳，而單一LSTM模型和編解碼器LSTM模型的預(yù)測效果急劇下降。其原因是，在遠(yuǎn)期的預(yù)測中CNN 對于數(shù)據(jù)特征的提取依然可靠有效。CNNLSTM編解碼器模型72 h 的多步預(yù)測結(jié)果與觀測數(shù)據(jù)的相關(guān)系數(shù)為0.519，RMSE 為89.493 km/s，MAE 為72.556 km/s；96 h 的多步預(yù)測結(jié)果與觀測數(shù)據(jù)的相關(guān)系數(shù)為0.428，RMSE 為93.924 km/s，MAE為75.605 km/s。

表4 72 h預(yù)測實驗結(jié)果Table 4 Predicted experimental resultsat 72 h

表5 96 h預(yù)測實驗結(jié)果Table 5 Predicted experimental results at 96 h

4 結(jié)束語

本文采用深度學(xué)習(xí)的方法，根據(jù)多元時間序列數(shù)據(jù)來預(yù)測在太陽與地球之間的拉格朗日點1（L1）處的太陽風(fēng)速度；基于LSTM 開發(fā)了3種深度學(xué)習(xí)預(yù)測模型——單一LSTM 模型、編解碼器LSTM 模型和CNN-LSTM編解碼器模型；并使用這3種模型對未來24、48、72、96 h 的太陽風(fēng)速度進(jìn)行預(yù)測。對于每種模型，除了直接對目標(biāo)時間點進(jìn)行預(yù)測，還引入多步預(yù)測的方法，對包含目標(biāo)時間點的連續(xù)多個時間步進(jìn)行預(yù)測。實驗結(jié)果表明，每種模型的多步預(yù)測效果均優(yōu)于單步預(yù)測效果。

將來，深度學(xué)習(xí)和多步預(yù)測的方法可能會幫助我們發(fā)現(xiàn)太陽物理學(xué)中的新關(guān)系，或者進(jìn)一步提高太陽風(fēng)速度的預(yù)測精度，支持建立有效的空間天氣預(yù)報體系。