陳雨田,劉立龍,黎峻宇,田祥雨
(桂林理工大學 a.測繪地理信息學院;b.廣西空間信息與測繪重點實驗室,廣西 桂林 541006)
電離層是地球大氣受太陽輻射和宇宙射線影響而形成的一個電離區域。 在衛星通訊、 導航和定位等應用領域, 由于電磁波穿過電離層時發生折射、 散射、 吸收等現象[1], 產生延遲造成誤差, 而誤差大小與傳播路徑上的總電子含量(total electron content,TEC)成正比[2]。 一些學者通過研究地震多發區上空的電離層電子濃度, 發現震前震區上空的電離層會產生異常擾動, 其TEC較平時發生明顯增加或減少[3-4]。 因此,對電離層總電子含量的高精度預報既是導航定位的迫切需要, 又是對地震等自然災難分析研究的重要手段。
目前, 常用的電離層總電子預報模型主要有經典電離層模型[5-6](Bent模型和Klobuchar模型等)和線性、 非線性預測模型[7-8](時間序列和神經網絡等)。 經典電離層模型在本領域經過較長時間的發展和改進具有易于操作、 計算簡單等優點, 但其固有的缺陷并沒有得到很好的彌補, 模型預報精度較低, 一般情況下預報值精度僅約50%~60%, 且受地域影響較大[9-10]。 自回歸積分滑動平均模型(autoregressive integrated moving average model,ARIMA)是時間序列中預測精度較高、 實際應用廣泛的一種模型, 能夠對非平穩、 非線性的時間序列進行有效預報[11-12], 但電離層TEC的時空分布差異較大, 單一ARIMA模型的預報精度較低。 國內一些學者通過利用小波分解或經驗模態分解(empirical mode decomposition, EMD)對電離層的TEC變化進行探索,并改進了ARIMA預報模型[13-14], 獲得了一些有益的結論, 然而改進模型的精度仍有提升的空間。 本文在過往研究的基礎上進一步探索, 將小波分解與EMD方法相結合, 共同納入ARIMA模型,構建WEARIMA預報模型, 并采用IGS提供的GIM數據驗證改進模型的預報精度。
小波分析是一種對信號的局部化分析,通過小波分析可以將信號序列分解為一個低頻序列和多個高頻序列。低頻序列表現為高頻率分辨率和低時間分辨率,高頻序列表現為低頻率分辨率和高時間分辨率。在進行小波分析時,定義平方可積函數ψ(t)為基本小波函數,滿足條件:
(1)
將基本小波函數ψ(t)伸縮和平移后即可得到一個小波序列:
(2)
式中:a和b分別為伸縮因子和平移因子,a、b∈R,a≠0, 則任意函數x(t)的連續小波變換及其逆變換可表示為
(3)
(4)
經驗模態分解(EMD)是一種能夠有效處理非線性、非平穩信號的新型信號處理方法。根據信號的特征,通過EMD分解可以將復雜信號序列x(t)從高頻到低頻進行多尺度分解, 從而得到頻率由高到低排列的若干個IMF分量和殘余函數R(n), 其實質是將信號序列由高頻至低頻進行平穩化處理[15],即
(5)
EMD分解要求分解得到的每個IMF分量滿足兩個條件[16]:第一,在整個信號序列中極值點與零點數目相同或最多相差一個;第二,由極大值和極小值構成的上下包絡關于時間軸局部對稱。由此可以看出,對于幅度與頻率相對較為對稱的信號,EMD分解會取得較好的效果,而對時域波形發生畸變的信號分解效果則較差。
自回歸積分滑動平均模型(ARIMA)是對自回歸滑動平均模型(ARMA)的擴展。ARMA(p,q)模型的一般形式為
xt=φ1xt-1+φ2xt-2+…+φpxt-p+εt+θ1εt-1+θ2εt-2+…+θqεt-q,
(6)
式中:εt為白噪聲序列, 并且與t時刻前的原始序列x(k)(k ARIMA(p,d,q)×(P,D,Q)考慮了季節因素影響模型[11], 具有建模簡單、 操作方便等特點,在對具有季節性的非平穩時間序列的預報中具有廣泛的應用。 模型中的p、d、q分別為自回歸階數、 差分階數和移動平均數; 而參數P為季節性自回歸階數,D為常規差分階數,Q為季節性移動平均數。 WEARIMA模型的建模流程如圖1所示,建模步驟如下: ① 提取IGS中心發布的30天電離層TEC格網數據,將前25天作為樣本序列,利用db4正交小波對樣本序列進行一級分解(用以提高序列的平穩性及相對于時間軸的對稱性),對分解出的參數選取合適的閾值進行分析,通過逆小波變換對分解結果重構,得到高頻序列d1和低頻序列a1; ② 將小波分解得到的低頻序列a1進行EMD分解,得到若干個IMF分量和單調殘余函數R(n)(EMD分解后的趨勢余量); ③ 將高頻序列d1和EMD分解后的各個分量分別代入ARIMA模型中,根據不同分量表現出的自相關與偏相關特性,確定各模型的參數; ④ 對各分量的預報值進行求和,得到最終的預報結果。 圖1 小波-EMD改進ARIMA模型(WEARIMA)建模流程Fig.1 Modeling process of improved ARIMA model based on wavelet-EMD 由EMD分解后IMF分量的特性可知,當原始信號序列的幅度與頻率較為對稱時,分解效果較好,利用ARIMA模型預報的精度也較高。太陽輻射對大氣分子的作用是形成電離層的主要因素[17],在太陽活動的活躍期往往會伴隨磁暴等現象,造成電離層總電子含量(TEC)異常擾動,不利于EMD分解后的預報,且使得預報結果具有偶然性,故選取太陽活動較為平靜的2009年年積日122~151共30天的GIM電離層格網數據進行實驗,利用前25天的TEC觀測值作為樣本序列,預報后5天的TEC值。 圖2給出了預報期內表征地磁活動強度的Kp指數與表征全球環電流強度的Dst指數。Kp指數可以反映3 h的地磁活動強度,一般認為,Kp值小于5時不會發生不同等級的地磁風暴;每日Kp指數之和小于30表示地磁活動平靜。而Dst指數的時間分辨率為1 h,當Dst指數不低于-30 nT時不會發生任何等級的磁暴。預報期的地磁指數均滿足上述條件,故可排除磁暴與地磁活動對預報結果的影響。 圖3是對年積日122~146的數據在格網點(45°N, 110°E)處的小波-EMD分解, 采用db4正交小波對原始序列進行一級分解, 得到低頻序列a1與高頻序列d1, 高頻序列d1的極值介于-2~2 TECu, 使用ARIMA模型直接進行預報。 低頻序列a1極值介于5~15 TECu, 利用EMD方法對其進一步分解得到5個IMF分量和1個趨勢余量R(6), 將得到的各分量分別使用ARIMA模型對后5天的TEC進行預報,結果如圖4所示。除IMF1與IMF2以外,其余各分量的預報值與真實值之間相差很小,其中IMF5與趨勢余量R(6)和真實值基本一致。 實驗中為了方便表示,將僅用小波分解改進的ARIMA模型記為WARIMA模型,小波-EMD改進的ARIMA模型記為WEARIMA模型,通過ARIMA、WARIMA和WEARIMA分別對上述格網點(45°N, 110°E)進行5天預報,得到各模型預報結果與IGS觀測值的對比情況,如圖5所示。可以看出,僅采用ARIMA模型得到的預報結果在每日TEC的極值附近與觀測值有較大偏差,而WARIMA與WEARIMA相對于原模型能夠有效削弱極值點附近預報值的誤差,在預報性能上有較大的提升,且WEARIMA模型的改進效果更好。對預報值殘差的絕對值進行統計分析,得到表1。 從后5天的總體預報來看,ARIMA模型在改進前僅有60%左右的殘差絕對值小于1 TECu,通過兩種方法改進后誤差絕對值小于1 TECu的可增大約10%;改進后誤差絕對值大于3 TECu的較原模型明顯減少,且WEARIMA模型的預報值在3天中沒有出現大于3 TECu的誤差。 圖2 預報期Dst與Kp指數(年積日122~151)Fig.2 Dst and Kp index during forecast period(DOY 122-151) 圖3 樣本序列小波-EMD分解圖(年積日122~146)Fig.3 Wavelet-EMD decomposition of sample sequence(DOY 122-146) 僅通過單個格網點數據與預報值殘差的絕對值比例并不能很好地體現出改進模型的優勢,為進一步說明小波-EMD改進模型在電離層TEC預報上較原模型與WARIMA模型有更高的預報精度,本文將IGS中心發布的觀測值視為真實值,采用平均相對精度(P)和均方根誤差(RMSE)進一步評定3種模型預報結果的優劣。P、RMSE指標為 (7) (8) 式中:IPRE和IIGS分別為模型的預報值與IGS觀測值;m和n分別為當天觀測的起始歷元和結束歷元。 由于電離層TEC值受太陽直射影響較大,其周日變化幅度隨地理經度的變化較小,而隨緯度的升高呈現明顯的降低趨勢,表2給出了北半球低緯(15°N, 110°E)、 中緯(45°N, 110°E)和高緯(75°N, 110°E)3個緯度地區ARIMA、WARIMA和WEARIMA模型預報值的日平均相對精度。從緯度變化上看,3種模型預報精度隨緯度的升高而增加。 在同一格網點上, 除低緯度WEARIMA模型的相對精度低于WARIMA以外,其余緯度的平均相對精度均表現為WEARIMA模型高于WARIMA模型,高于ARIMA模型。各模型5 d預報值均方根誤差的計算結果見表3,均方根誤差隨緯度的升高而減小,在同一格網點上WEARIMA模型預報值的均方根誤差要明顯低于WARIMA和ARIMA模型,具有最為可靠的預報性能。 圖4 分解后各分量的5 d預報值(點線)與觀測值(實線)對比Fig.4 Comparison between the 5 d prediction value and the observed value of each component obtained after decomposition 圖5 各模型預報結果對比Fig.5 Comparison of prediction results of each model 圖6為中國及周邊地區上空TEC預報值的均方根誤差(RMSE)分布情況。 從空間上看,隨緯度的降低,RMSE呈現明顯的增大趨勢, 且在中國西南部地區達到最大值, 即表現為“駝峰”; 從時間上分析, 除第一天的預報結果在中國南部地區有較大誤差, 其余預報值隨時間推移誤差逐漸增大, 在最后一天達到極大值。 通過對比每天的RMSE值以及5天的RMSE平均值, 可以發現WEARIMA模型在該區域的預報精度相較于ARIMA模型有明顯提升, 在預報誤差較大的中國西南部地區, WEARIMA模型能夠有效削弱RMSE峰值,提高預報精度。 表1 5 d預報值殘差絕對值(Δ)統計 Table 1 Absolute value of the residuals(Δ)of the predicted values over the next 5 days % 預報時間/dARIMA模型/WARIMA模型/WEARIMA模型Δ<1TECu1TECu≤Δ<2TECu2TECu≤Δ<3TECuΔ≥3TECu163.89/66.67/66.6719.44/13.89/19.448.33/11.11/5.568.33/8.33/8.33275.00/86.11/69.4419.44/8.33/19.440.00/0.00/11.115.56/5.56/0.00358.33/72.22/75.0027.78/13.89/13.895.56/8.33/8.338.33/5.56/2.78461.11/72.22/69.4422.22/13.89/16.6755.56/5.56/13.8911.11/8.33/0.00552.78/66.67/77.7819.44/16.67/16.6711.11/8.33/5.5616.67/8.33/0.00 表2 不同模型5 d預報值相對精度對比 Table 2 Comparison of relative accuracy of predicted values of different models over the next 5 days % 經緯度模 型預報天數/d12345平均相對精度ARIMA82.9791.0084.9485.8078.46846315°N,110°EWARIMA87.9693.4990.3789.6185.2189.33WEARIMA80.3083.5989.9379.4789.9984.66ARIMA91.3992.1590.7789.3689.2690.5945°N,110°EWARIMA88.8492.3391.8793.0094.0292.01WEARIMA95.0394.9193.3195.1195.4094.75ARIMA95.0394.1391.3393.5295.2193.8475°N,110°EWARIMA95.5594.8692.1393.0995.1994.16WEARIMA95.5793.8493.3094.6996.3194.74 表3 不同模型5 d預報值均方根誤差(RMSE) 圖6 RMSE空間分布情況對比Fig.6 Comparison of RMSE spatial distribution 采用小波與EMD分解相結合的方法改進ARIMA預報模型,根據IGS提供的GIM電離層格網數據,分析改進模型在中國及周邊地區的預報精度,并與單一小波分解改進的模型及原模型進行對比分析,得到如下結論: (1)對小波分解后的低頻序列進行EMD分解可得到若干個IMF分量與一個趨勢余量,對其分別預報得到的預報結果表現為隨頻率降低與真值的擬合度升高,并且較單一小波分解具有一定的提升,實驗表明該方法在提高電離層TEC預報精度上有一定作用。 (2)分析各模型預報值的相對精度與均方根誤差發現,總體而言,3種模型的預報精度隨緯度的升高而增加,改進模型的預報精度較原模型有較大提高,3種模型預報精度表現為WEARIMA模型>WARIMA模型>ARIMA模型。 (3)對中國及周邊地區上空進行預報并分析區域均方根誤差分布情況,得到WEARIMA模型能有效削弱中國西南部地區預報值均方根誤差的峰值,與原模型相比精度有較大的提升。 本文基于小波-EMD分解改進ARIMA模型,為電離層總電子預報提供了新的思路。但是,電離層各項參數受太陽活動及其他因素的影響較大,如何顧及這些因素進行預報以及如何選取小波基、確定ARIMA模型階數等問題仍有待進一步研究。
2 實驗分析
2.1 實驗數據準備
2.1 實驗與精度分析








3 結 論