EEMD-ARIMA在干旱預測中的應用
——以新疆維吾爾自治區為例

許德合，丁 嚴，張 棋，黃會平

（華北水利水電大學，鄭州450000）

干旱災害是常見的自然災害之一，其發生的主要原因是一段時期內降水量的減少，若遇到高溫、大風和低溫等異常氣候，則會進一步加重干旱程度［1］。通常情況下，長時間的干旱將導致地區水資源供應不足，影響正常的農業生產和經濟發展，造成不可估量的經濟損失［2］。因此，如何對干旱的發生情況進行準確評估、監測和分析，成為了國內外學者研究的熱門議題［3，4］。在對干旱進行研究的過程中，常使用易于進行對比分析和計算的干旱指標來評估干旱程度、持續時間和影響范圍［5-7］。常用的干旱評價指標有標準化降水指數（Standard Precipitation Index，SPI）、標準化降水蒸散指數（Standard Precipitation Evaporation Index，SPEI）、帕默爾干旱指數（Plamer Drought Severity Index，PDSI）等［8-13］，其中SPI的運用最為廣泛。SPI是用來描述一段時間內降水量出現頻率多少的指標，可應用于不同時間尺度下的計算，干旱分級精度相對較高［14-16］。不同時間尺度的SPI能夠反映不同類型的干旱，1個月時間尺度的SPI適用于季節干旱的分析，6 個月時間尺度的SPI可用于分析農業干旱，長時間尺度的SPI適于分析水文干旱。

近幾年，國內干旱災害時常發生，其中以西北地區最為嚴重［17］。對新疆維吾爾自治區不同地區的干濕變化進行研究，探索不同時期的降水量對區域產生的影響，明確干旱的成因機制，可為相關部門在制定有效干旱應對措施時提供幫助，減少干旱災害所帶來的經濟損失［1］。目前，常用的干旱預測模型有差分自回歸移動平均模型（Autoregressive Integrated Moving Average Model，ARIMA）、支持向量機（Support Vector Machine，SVM）等。其中ARIMA 模型是最常用的數據驅動模型，常用于時間序列數據的預測。韓萍等利用ARIMA 模型對多個時間尺度的SPI進行預測，并對模型預測結果進行對比分析，得出ARIMA 模型較適合SPI3、SPI6、SPI9序列的短期預測，適合SPI12、SPI24序列長期預測的結論［18］；李佳佳等為識別長江流域月降水周期振蕩和長期趨勢的顯著影響因子，應用集合經驗模態分解（Ensemble Empirical Mode Decomposition，EEMD）方法，分別對各站點的月降水序列進行分解［19］；王佳等基于EEMD-ANN（Artificial Neural Network，ANN）預報模型預測入庫徑流量，結果表明使用EEMD-ANN 可以較為精準地進行徑流量預測［20］；李輝等為提取機組故障特征，將EEMD 與SOM 神經網絡結合進行故障自動識別，結果表明該方法可以準確提取機組故障特征且計算速度快［21］；李勃旭等為確定地鐵門傳動系統的退化狀態，采用EMD-ARIMA 模型對夾緊力峰值的均值和標準差進行預測，結果表明EMD-ARIMA 模型能夠較好的預測夾緊力的峰值［22］。利用氣象站點監測到的降雨量值計算得到的SPI序列屬于典型的非平穩時間序列，而ARIMA模型能夠更多的提取原始序列的信息，常用于處理非平穩的時間序列，因此可以利用ARIMA模型預測SPI序列［23］。EEMD 能夠提取出原始信號在不同尺度的局部特征信號，精準地反映出原時間序列信號的物理特性，從而為模型預測提供穩定的前提。EEMD-ARIMA 模型集合了EEMD 和ARIMA 模型的優點，通過EEMD 分解得到時間序列平穩的局部特征，有利于ARIMA 模型的預測，可以提高模型的預測精度。因此，本文利用EEMD-ARIMA 組合模型對多尺度SPI序列進行預測，分析模型在干旱預測中的應用及優勢，以期為干旱防控工作提供幫助。

本文通過對新疆維吾爾自治區32 個氣象站點于1960-2019年間收集到的降水量數據進行計算，得到不同時間尺度的SPI，通過ARIMA 模型和EEMD-ARIMA 組合模型，分別對1、3、6、9、12、24，6個時間尺度的SPI進行預測，并采用平均絕對誤差（Mean Absolute Error，MAE）、均方根誤差（Root Mean Square Error，RMSE）、決定系數（R-Square，R2）對兩種模型的預測精度進行評價。由于計算處理的數據量大，受論文篇幅限制無法展示所有數據的處理結果，加之本文側重于講述組合模型在干旱預測中的優勢，故在文中僅以分布在整個新疆地區3 個不同方位上的福海站、巴音布魯克站、莎車站為例，分析組合模型的預測結果，并進行兩種模型預測結果的對比和分析。兩種模型對32 個站點多尺度SPI的預測結果，結合ArcGIS 的經驗貝葉斯克里金插值法在文末進行可視化展示和分析。

1 研究區概況

新疆維吾爾自治區地理坐標介于34°25'N～48°10'N、73°40'E～96°18'E 之間。省域內平均海拔約1 000 m，區域北部有阿爾泰山，南部有昆侖山系，省域中部的天山以北為準噶爾盆地，南部為塔里木盆地。新疆維吾爾自治區遠離海洋，區域地形復雜，四周高山環繞，海洋氣流不易到達，致使區域內氣候干燥，是典型的干旱、半干旱地區。區域年平均降水量在150 mm 左右，區域北部的降水量高于南部，各地降水量相差很大。受氣候和地理環境的影響，新疆維吾爾自治區生態環境脆弱，各種自然災害頻繁發生，在全球氣候變暖的形勢下，干旱災害產生的影響呈現出擴大化的趨勢［24，25］。

圖1 研究區域及氣象站點分布Fig.1 Distribution of meteorological stations in study area

2 數據來源與方法

2.1 數據來源

本文所用逐日降水量數據來源于國家氣象科學數據中心（http：//data.cma.cn/）中新疆維吾爾自治區氣象站觀測數據。所用新疆維吾爾自治區地理高程數據來源于地理空間數據云（http：//www.gscloud.cn/search）。表1 為所選3 個示例站點的信息。

表1 示例氣象站點信息Tab.1 Information of sample meteorological stations

2.2 研究方法

2.2.1 標準化降水指數（SPI）

近年來，國內外學者常使用SPI指數進行氣象干旱研究。SPI指數是通過計算出某時段內降水量的Γ 分布概率，經正態標準化處理后，根據處理結果劃分干旱等級。它可以用來反映一段時間內降水量出現頻率的多少，易于計算且能直觀反映氣象干旱的程度［26］。SPI的具體計算過程參見氣象干旱等級（GB/T20481-2017）。依據國家標準氣象干旱等級（GB/T20481-2017）規定的干旱分級標準（表2），根據計算得到的SPI，確定區域的干旱類型。

表2 標準化降水指數干旱分級Tab.2 Drought classification based on SPI

SPI1、SPI3、SPI6、SPI9、SPI12、SPI24分別對應1、3、6、9、12、24個月時間尺度的SPI。短時間尺度的SPI常用來反映短期內的降雨量變化。SPI1適用于氣象干旱，可用于反映短期降水變化情況。3個月時間尺度下的SPI用來描述降雨量的季節變化，計算出的SPI3中5、8、11、2月數據可分別用于描述春、夏、秋、冬的干旱情況。6 個月時間尺度下的SPI適用于農業干旱。長時間尺度的SPI常用來分析長期的干旱趨勢。SPI9可用于表征較長時間內的地下水位變化。SPI12的時間周期較長，能較清楚的描述年際降雨量變化情況，24 個月時間尺度的SPI序列則用于分析長期的降水變化所引起的干旱［27-29］。

2.2.2 ARIMA模型

自回歸移動平均模型（Auto-Regressive and Moving Average Model，ARMA）的建模思想是將預測值假定為一組隨機序列，確定能夠近似描述這組序列的模型，之后根據該序列的過去值和現在值，對未來值進行預測［30］。ARMA 模型分為自回歸模型（Auto-Regressive，AR）、滑動平均模型（Moving-Average，MA）以及ARMA。只要是平穩且非白噪聲的時間序列皆可通過建立ARMA 模型進行預測，但大多數時間序列都是非平穩時間序列。因此，在時間序列具有某種趨勢時，需要對時間序列進行d次差分，使其成為平穩序列。進行d次差分后的ARMA 模型即為ARIMA模型。ARIMA（p，d，q）模型的一般式為：

式中：Ht為時間序列值；φi(i= 1，2，…，p)和θj(j= 1，2，…，q)分別為自回歸系數和滑動平均系數；ut為白噪聲序列，且ut～N(0，σ2)。

ARIMA模型建模流程為：

（1）平穩性檢驗。本文使用單位根檢驗（Augmented Dickey-Fuller Test，ADF）法進行判斷。在ADF 檢驗中，原假設為非平穩時間序列且存在單位根，給定顯著水平α=0.05，若檢驗統計量對應的概率值（P）小于0.05，則拒絕原假設［31，32］。對于非平穩的時間序列需要進行差分，得到平穩序列。

（2）p、q取值范圍確定。根據數據的自相關函數和偏自相關函數來確定模型階數p、q的取值范圍。

（3）模型定階。采用赤池信息準則（Akaike Information Criterion，AIC）、貝葉斯信息準則（Bayesian Information Criterion，BIC）進行模型定階，AIC、BIC公式如下：

式中：N為參數個數。

在不同p和q組合的模型中，選擇AIC、BIC的最小值所對應的參數，從而得到最優ARIMA模型。

2.2.3 EEMD分解

1998年N E Huang 等人提出了經驗模態分解（Empirical Mode Decomposition，EMD），EMD 能夠完美地適應于全部的非線性、非平穩信號的處理，并且經過該方法處理后的結果具有相當高的信噪比［33］。原始序列輸入EMD 后得到有限個固有模態函數（Intrinsic Mode Function，IMF）和余量，各IMF分量包含了原始序列不同時間尺度的局部特征，盡可能地保留原始數據的 特 性［22］。EEMD 是EMD 的 一種改進方法，與EMD相比，EEMD在信號中加入了高斯白噪聲，高斯白噪聲以其均勻分布的特性補償了IMF分量的損失，其步驟算法如下［34］：

（1）在原始數據Y(t)中添加正態分布的白噪聲序列ωi(t)，試驗次數為i；

（2）進行EMD 分解，將含有白噪聲的原始數據分解為IMF的組合；

（3）每次都加入服從同一分布的不同白噪聲序列，重復步驟1和2，得到一組不同的IMF成分和殘差；

（4）以所有IMF的均值作為最終的IMF組。

最終的分解結果為：

2.2.4 EEMD-ARIMA模型

通過Python編程語言，將EEMD 與ARIMA 模型結合為新的預測方法EEMD-ARIMA 模型。組合模型先利用EEMD 將有非平穩特征的時間序列SPI分解為N項含有原序列局部特征的序列，再利用ARIMA 模型對這些序列分別進行預測，最后把各項預測結果求和得到該序列的最終預測結果，這樣的預測結果比直接用ARIMA 模型預測的非平穩時間序列SPI的結果具有更高的精度。使用EEMD-ARIMA模型進行預測的步驟如下：

（1）EEMD 處理。將降水量數據導入EEMD 進行分解，原始序列分解為從高頻到低頻的IMF1、IMF2、…、IMFn以及殘差量。殘差量即趨勢，將其記為IMFn+1。

（2）ARIMA 模型處理。將IMF1、IMF2、…、IMFn+1，分別導入ARIMA模型進行預測，對各個分量進行平穩性檢驗和模型定階，并將預測結果輸出。預測結果記為P1、P2、…、Pn+1。

（3）對P1、P2、…、Pn+1進行相加求和。

EEMD-ARIMA模型建立流程見圖2。

圖2 EEMD-ARIMA模型建立流程Fig.2 EEMD-ARIMA combined model forecast flow chart

2.3 評價指標

常見的回歸預測評估指標有MAE、均方誤差（Mean Square Error，MSE）、RMSE，其中MAE和MSE是基礎的評估指標，RMSE是MSE指標的擴展，相較MSE指標更加準確。為了比較各模型預測精度的高低，本文選定了RMSE、MAE和R23 種評價指標對ARIMA 模型和EEMD-ARIMA 組合模型的預測結果進行評價。RMSE是用來衡量觀測值與真實值之間的偏差，MAE是絕對誤差的平均值，能夠更好地反映預測值誤差的實際情況，RMSE和MAE的值越小，模型效果越好。R2是將預測值與均值進行對比，R2越大，表示擬合效果越好，最大值為1。

式中：xi是觀測值；yi是真實值；是yi的平均值為預測值；N為樣本數。

3 結果與分析

3.1 ARIMA模型建模和預測

依據國家氣象信息中心提供的中國地面氣候資料月值數據集，選擇1960-2019年新疆維吾爾自治區境內32 個氣象站點持續測定的逐月降水量數據，進行多尺度SPI的計算。由于SPI在不同時間尺度適用于不同種類的干旱，因此選取了1、3、6、9、12、24 共6 個時間尺度［18］。將計算得到的SPI序列，分別記為SPI1、SPI3、SPI6、SPI9、SPI12、SPI24。

通過Python3.6 對ARIMA 建模。經過ADF 檢驗，SPI1、SPI3、SPI6、SPI9、SPI12、SPI24的P值均小于0.05，故為平穩時間序列。通過ACF、PACF進行模型定階，選擇當AIC、BIC值最小時對應的p、q值［35］，各序列的模型定階結果見表3，即為各時間尺度的最優模型。選取1960-2007年數據作為訓練集，2008-2019年數據作為測試集。應用32個氣象站點各時間尺度的ARIMA 最優模型對SPI序列進行預測。示例站點的預測結果分別見圖4、圖5和圖6。

表3 六尺度SPI序列的ARIMA模型定階Tab.3 Order the ARIMA model based on six scales SPI values

3.2 EEMD-ARIMA模型建立和預測

EEMD能夠將原始信號逐級分解并提取出其在不同尺度的局部特征信號，從而準確反映出原時間序列信號的物理特性。在分解過程中，需要對Nstd和NE進行設置。Nstd用于設置添加高斯白噪聲的標準差，取值范圍一般為0.01～0.4，具體設置數值根據原信號中的噪聲干擾大小，視具體情況確定。NE是添加噪聲的次數，通常取值為50或100。利用EEMD對多尺度SPI進行分解，經過多次參數修改及其分解結果的對比，最終將Nstd設置為0.2，NE設置為100。將加入白噪聲序列的原始序列進行分解，得到8 個IMF分量和1 個殘差數據，如圖3 所示為對福海站SPI3序列進行EEMD 分解得到的結果。從圖3 中可以看出，分解出的時間序列均圍繞0 值上下波動，趨于平穩，說明通過EEMD進行分解能夠降低原始序列的非平穩性。

圖3 EEMD分解結果Fig.3 EEMD decomposition results

選取1960-2007年數據作為訓練集，2008-2019年數據作為測試集。將EEMD 分解得到的1960-2007年數據輸入ARIMA 模型進行預測，并以ARIMA 模型預測結果之和作為預測的最終結果。將6 個時間尺度的SPI實際計算值、ARIMA 模型預測值、EEMD-ARIMA 模型預測值進行對比展示。3 個示例站點的預測對比圖分別為圖4、圖5和圖6。

圖4 基于ARIMA模型與EEMD-ARIMA組合模型對福海站多時間尺度SPI值的預測結果與觀測值計算結果對比（2008-2019）Fig.4 Comparison of predicted and observed value of multi-time scale SPI of ARIMA model and EEMD-ARIMA combined model in Fuhai Station（2008-2019）

圖5 基于ARIMA模型與EEMD-ARIMA組合模型對巴音布魯克站多時間尺度SPI值的預測結果與觀測值計算結果對比（2008-2019）Fig.5 Comparison of predicted and observed value of multi-time scale SPI of ARIMA model and EEMD-ARIMA combined model in Bayinbuluke Station（2008-2019）

圖6 基于ARIMA模型與EEMD-ARIMA組合模型對莎車站多時間尺度SPI值的預測結果與觀測值計算結果對比（2008-2019）Fig.6 Comparison of predicted and observed value of multi-time scale SPI of ARIMA model and EEMD-ARIMA combined model in Shache Station（2008-2019）

由圖4（a）可知，2010年、2012年、2018年，福海站旱情最為嚴重，EEMD-ARIMA 組合模型的計算結果與實際情況最為接近，較為精確的預測到了嚴重干旱的發生，在這些時間段中組合模型的預測值比ARIMA 模型的預測值更接近SPI實際計算值。在圖5（a）中2018年，組合模型的預測結果更符合實際情況，這表明組合模型對極端天氣預測結果的準確度要高于單一模型。根據組合模型的預測結果可知，3 個站點在2008年和2014年有嚴重干旱發生，這與新疆年鑒中對大范圍干旱的記錄一致。圖4（b）2010-2011年，圖5（a）2016-2017年SPI值異常，在快速增長后快速下降，對于這種異常變化，通過組合模型預測得到的干旱情況更接近實際的變化情況。根據圖4、圖5 和圖6 的（c）、（d）、（e）、（f）圖可以看到，兩種模型的預測值均與實際計算得到的SPI值接近，但隨著時間尺度的增大，單一模型預測結果延遲越發明顯，在24 個月時間尺度，與實際計算結果相比有近一個月的延遲。而組合模型隨著時間尺度的增大，預測結果與實際計算結果更為接近。從3個站點的多尺度SPI對比中可以明顯看到，兩種模型在1 個月時間尺度下的預測值與實際值有明顯差異。在各時間尺度下，組合模型的預測值與SPI實際值之間的差值均小于單一模型預測值與SPI實際值的差值，表明組合模型在6 個時間尺度下均能更好的預測SPI序列。

為驗證EEMD-ARIMA 組合模型和ARIMA 模型在6個時間尺度下的預測精度，通過MAE、RMSE、R2對兩種模型在6個時間尺度下的預測結果進行了評價，結果見表4。RMSE和MAE的取值范圍為［0，+∞］，取值越小，模型效果越好，當值為0 時，說明模型完美。R2取值越大，表示擬合效果越好，最大值為1，當值為1 時表明模型預測效果最佳。對比3 個站點在不同時間尺度下的組合模型和單一模型的MAE、RMSE、R2值發現，組合模型在任一時間尺度下的MAE、RMSE值均小于單一模型的值，且在任一時間尺度下的R2值均大于單一模型的值，說明EEMDARIMA 組合模型在各個時間尺度下的預測精度均高于單一模型，表明組合模型更適用于多時間尺度下的時間序列預測。組合模型在1個月時間尺度時，預測精度遠高于單一模型，隨著時間尺度的增加，兩種模型預測精度之間的差距逐漸縮小，在24個月時間尺度時達到最小，表明組合模型在1、3 個月時間尺度下的預測結果遠優于單一模型，在6、9個月時間尺度時，優于單一模型，在12、24 個月時間尺度下略優于單一模型。隨著時間尺度的增加，SPI序列特征趨于平穩，ARIMA模型對序列的擬合效果逐漸提升。

表4 ARIMA模型和EEMD-ARIMA組合模型的MAE、RMSE、R2值Tab.4 MAE、RMSE and R2 values for ARIMA and EEMD-ARIMA models

結合ArcGIS的空間分析功能，利用經驗貝葉斯克里金插值法對2019年32 個站點的SPI實際觀測計算值和預測值進行可視化（如圖7）。由于不同時間尺度的SPI適用情況不同，且新疆的春旱威脅最大，其次為夏旱和秋旱［36］，故此處選擇能夠進行降雨量季節變化分析的SPI3，對新疆維吾爾自治區春、夏、秋、冬的干旱情況進行插值可視化。從圖7 中可以看到EEMD-ARIMA 組合模型在四季的預測結果與實際情況接近，較ARIMA 模型的預測結果更符合實際。

圖7 基于SPI實際值、ARIMA模型預測值和EEMD-ARIMA模型預測值的季尺度干旱空間分布Fig.7 Spatial distributions of seasonal drought levels based on the actual values and the predicted results of models

4 討 論

ARIMA 模型和EEMD-ARIMA 組合模型對SPI的預測精度隨時間尺度的增加而提高，并在24 個月時間尺度時達到最高。從SPI的計算公式來看，隨著時間尺度的增加，SPI序列中集合了原始數據中更多的信息，因此，預測值對實際計算值的擬合更為充分。從結果來看，EEMD-ARIMA 組合模型的預測結果比單一模型的預測結果更接近實際情況，這是因為通過ARIMA模型對多尺度SPI進行預測存在特征不平穩的情況，而EEMD可以通過提取出原始時序在不同尺度的局部特征，從而將SPI數據序列平穩化，為ARIMA 預測提供穩定的前提。因此，EEMD-ARIMA 組合模型優于單一模型，能夠用于新疆地區的干旱預測。

新疆四季降水不均，四季的干旱情況不同，以春旱最為嚴重，秋旱和夏旱次之。因此，文中主要考慮了降水量對區域的影響，選取了3 個月尺度的SPI進行可視化展示。在對預測結果進行分析時，也側重于對3 個月尺度的計算結果的分析。但河流、植物覆蓋等都會對干旱情況產生影響，1978-2008年新疆耕地面積總體呈增加的趨勢，農作物的增加對地區的水利設施造成負擔，減弱了地區的抗旱能力，造成干旱程度的進一步加劇。在之后的干旱研究中，將通過對長時間尺度的SPI值的分析，進一步研究地下水位的變化，植被的覆蓋等對干旱產生的影響。從不同尺度的SPI值進行研究得到更加準確詳細的結果。

對極端干旱的預測是目前研究的熱點和干旱預測問題的研究重點。本文利用SPI對新疆地區的干旱時空變化情況進行分析，基于計算出的數據，對EEMD-ARIMA 模型在新疆地區的干旱預測中的適用性進行研究發現，組合模型可以較好的捕獲2008-2019年的干旱事件，適用于新疆地區的干旱預測。北疆和南疆的干旱情況有所不同，南疆大部分地區以農業為主，水利設施投入不足；而北疆的抗旱澆灌面積持續增長、農業和水利設施投入程度增加，因此南疆的受旱面積增長趨勢大于北疆［37］。在后續的研究中，針對不同地區的實際情況，進一步細化干旱評價指標的選取，找出不同地區最為適合的干旱指數，提高模型在新疆的適用性。

5 結 論

本文基于新疆維吾爾自治區32個氣象站點的降水量數據，計算出不同時間尺度下的SPI值，通過對SPI序列的分析、建模及預測，得到如下結論。

（1）單一的ARIMA 模型，在1 個月尺度下的預測精度最低，在24個月尺度的預測精度最高。主要是由兩個原因造成的，其一，ARIMA 模型是一個整體線性自回歸模型，隨著預測數據的增多，精度會下降。在1 個月尺度下，通過ARIMA 模型得到的預測數據多于其他時間尺度，故精度相較其他時間尺度而言較低。其二，長時間尺度的SPI序列集合了更多原始數據的信息，有利于模型的預測。在不同時間尺度下，組合模型的預測精度均高于單一模型。

（2）EEMD-ARIMA 組合模型預測出的干旱發生年份與實際情況近似一致，能夠用于對新疆地區干旱的預測。通過EEMD對數據進行分解，精簡了模型對原始數據的讀取過程，提高了預測精度。利用EEMD 對各尺度的SPI進行分解，得到一組比原始序列特征更加平穩的分量，將各分量導入ARIMA進行預測，并將各分量的預測結果累加得到每一時間尺度下的預測數據，得到的預測數據的精度均高于單一模型。從預測結果來看，EEMD-ARIMA 組合模型優于ARIMA 模型，即EEMDARIMA組合模型能夠更好的擬合不同時間尺度下的SPI。

氣象干旱能夠為區域的干旱情況提供預警，氣象站點監測到的氣溫、降水數據只能近似的反映區域的情況，粗略的估計區域的干旱變化情況，需要研究更加優化的干旱評價指標使其更精確的描述區域的變化情況。