王雪妮，郝舒哲，常建波，鄧顯羽

（1.山西省地質(zhì)調(diào)查院,太原 030006;2.太原理工大學(xué)水利科學(xué)與工程學(xué)院,太原 030024;3.中水東北勘測(cè)設(shè)計(jì)研究有限責(zé)任公司,長(zhǎng)春 130021）

干旱對(duì)農(nóng)業(yè)發(fā)展具有廣泛而持續(xù)的影響，其被認(rèn)為是最嚴(yán)重的自然災(zāi)害之一，并引起眾多學(xué)者的關(guān)注。干旱成因較為復(fù)雜，導(dǎo)致其在持續(xù)時(shí)間、強(qiáng)度和空間延伸等方面難以被量化[1]。鑒于此，國(guó)內(nèi)外通常采用帕默爾干旱指數(shù)、標(biāo)準(zhǔn)化降水指數(shù)和綜合氣象干旱指數(shù)等作為度量標(biāo)準(zhǔn)監(jiān)測(cè)干旱事件[2]。其中，帕默爾干旱指數(shù)及自適應(yīng)帕默爾干旱指數(shù)(self-calibrated palmer drought severity index,scPDSI)為應(yīng)用最為廣泛的干旱指標(biāo)之一。YANG 等[3]通過(guò)對(duì)7 個(gè)氣象干旱指數(shù)的分析，指出scPDSI 在中國(guó)區(qū)域內(nèi)更為適用。此后，scPDSI 作為較成熟的干旱監(jiān)測(cè)指標(biāo)被應(yīng)用于我國(guó)各地的氣象干旱研究[4-5]。而對(duì)于氣象干旱的預(yù)測(cè)方法，各種機(jī)器學(xué)習(xí)算法則成為目前研究的熱點(diǎn)。吳晶等[6]將隨機(jī)森林(random forest,RF)模型應(yīng)用于干旱等級(jí)預(yù)測(cè)研究中，并取得較為準(zhǔn)確的預(yù)測(cè)結(jié)果。措姆等[7]在基于機(jī)器學(xué)習(xí)模型開(kāi)展干旱趨勢(shì)預(yù)測(cè)的研究中發(fā)現(xiàn)，支持向量機(jī)(support vector machine,SVM)模型具有更高的預(yù)測(cè)精度。隨著“先分解，后重構(gòu)”思想的提出，耦合模態(tài)分解的組合機(jī)器學(xué)習(xí)模型開(kāi)始被應(yīng)用于氣象干旱預(yù)測(cè)研究當(dāng)中。丁嚴(yán)等[8]建立了基于互補(bǔ)集合經(jīng)驗(yàn)?zāi)B(tài)分解(complementary ensemble empirical mode decomposition,CEEMD)的長(zhǎng)短期記憶(long short term memory,LSTM)網(wǎng)絡(luò)模型和差分自回歸移動(dòng)平均(autoregressive integrated moving average,ARIMA)模型，對(duì)新疆多時(shí)間尺度標(biāo)準(zhǔn)化降水指數(shù)(standardized precipitation index,SPI)序列進(jìn)行干旱預(yù)測(cè)。李子陽(yáng)等[9]將變分模態(tài)分解(variational mode decomposition,VMD)和門(mén)控循環(huán)單元網(wǎng)絡(luò)(gated recurrent unit,GRU)相結(jié)合，建立了VMD-GRU 干旱預(yù)測(cè)模型，據(jù)此對(duì)趙口大型灌區(qū)逐月SPI值進(jìn)行了預(yù)測(cè)。此后，為進(jìn)一步提高氣象干旱預(yù)測(cè)精度，智能優(yōu)化算法逐漸被融入至耦合模態(tài)分解的組合機(jī)器學(xué)習(xí)模型中。劉雪梅等[10]利用改進(jìn)的混沌量子粒子群優(yōu)化算法(chaotic quantum particle swarm optimization,CQPSO)優(yōu)化了GRU 神經(jīng)網(wǎng)絡(luò)超參數(shù)，構(gòu)建了VMD-CQPSO-GRU 智能預(yù)測(cè)模型，通過(guò)實(shí)例分析發(fā)現(xiàn)，相比于單一循環(huán)神經(jīng)網(wǎng)絡(luò)模型，其所建模型的預(yù)測(cè)精度顯著提升。其他相關(guān)領(lǐng)域研究表明，自適性噪聲完備經(jīng)驗(yàn)?zāi)B(tài)分解法(complete ensemble empirical mode decomposition with adaptive noise, CEEMDAN)的效果通常優(yōu)于CEEMD[11]，而采用灰狼優(yōu)化算法(grey wolf optimization,GWO)優(yōu)選參數(shù)比傳統(tǒng)的粒子群優(yōu)化算法(particle swarm optimization,PSO)得到的SVM模型預(yù)測(cè)結(jié)果更為準(zhǔn)確[12]。

根據(jù)上述分析，本研究提出將CEEMDAN、GWO 和SVM 相結(jié)合，構(gòu)建CEEMDAN-GWO-SVM 模型，以北方典型的半干旱區(qū)——山西省為研究對(duì)象，基于其1956-2020年逐月scPDSI序列，采用CEEMDAN 對(duì)序列進(jìn)行分解，應(yīng)用GWO-SVM模型實(shí)現(xiàn)序列的預(yù)測(cè)重構(gòu)。研究結(jié)果可為山西省氣象干旱預(yù)測(cè)提供參考。

1 研究區(qū)概況

山西省位于黃土高原地區(qū)，介于34°34~40°44'N、110°14~114°33'E 之間，地跨黃河、海河兩大水系。山西省屬溫帶大陸性季風(fēng)氣候，四季分明、雨熱同季，年平均氣溫為9.1°C，年平均降水量為486.9 mm。山西省降水量時(shí)空分布不均，山多地瘠，水土流失嚴(yán)重，歷史上旱災(zāi)頻發(fā)[13]。近年來(lái)，山西省極端氣象災(zāi)害更是多發(fā)重發(fā)，多項(xiàng)指標(biāo)突破歷史極值[14]。

2 數(shù)據(jù)與研究方法

2.1 數(shù)據(jù)來(lái)源與處理

本研究所采用scPDSI 數(shù)據(jù)來(lái)源于英國(guó)環(huán)境數(shù)據(jù)分析中心(http://www.badc.rl.ac.uk/)，時(shí)間范圍為1956 年1月-2020年12月，空間分辨率為0.5°×0.5°。

在數(shù)據(jù)處理方面，首先對(duì)原始序列集進(jìn)行重采樣處理，形成空間分辨率為0.25°×0.25°的數(shù)據(jù)集；其次，提取山西省范圍內(nèi)空間逐月scPDSI數(shù)據(jù)；最后，計(jì)算出山西省1956年1月-2020年12月逐月平均scPDSI值，并以1956年1月-2007年12月scPDSI數(shù)據(jù)作為預(yù)報(bào)訓(xùn)練集，2008年1月-2020年12月數(shù)據(jù)作為測(cè)試集。

2.2 研究方法

2.2.1 CEEMDAN算法 CEEMDAN算法是由TORRES等[11]基于EMD[15]、EEMD[16]、CEEMD[17]提出的一種信號(hào)分解算法，該方法通過(guò)加入高斯白噪聲和迭代分解余量信號(hào)，可將原始時(shí)間序列分解為多個(gè)固有模態(tài)函數(shù)(intrinsic mode function,IMF)，據(jù)此較好地解決此前信號(hào)分解方法存在的模態(tài)混疊和白噪聲殘留問(wèn)題。CEEMDAN分解的主要步驟如下[11]。

(1)向原始信號(hào)x(t)中添加i(i=1,2,…,m)次服從標(biāo)準(zhǔn)正態(tài)分布的高斯白噪聲wi(t)，得到新的信號(hào)xi(t)＝x(t)+wi(t)。對(duì)xi(t)實(shí)施EMD分解，且只保留第一個(gè)分解分量，取均值得到原始信號(hào)的第一個(gè)IMF分量IIMF1：

對(duì)應(yīng)的殘余分量r1(t)為：

(2)向殘余分量r1(t)中添加I(i=1,2,…,m)次服從標(biāo)準(zhǔn)正態(tài)分布的高斯白噪聲wi(t)得到新分量ri(t)＝r(t)+wi(t)。對(duì)(t)進(jìn)行EMD分解，且僅保留第1個(gè)模態(tài)分量，則原始信號(hào)的第2個(gè)IMF分量IIMF2為：

(3)重復(fù)分解j(j=1,2,…,n)次直至殘差余量無(wú)法被分解時(shí)，可依次得到IIMF1，IIMF2，…，IIMFn及對(duì)應(yīng)的殘余分量。最終，原始信號(hào)可表示為：

基于CEEMDAN 方法，本研究發(fā)現(xiàn)將山西省逐月平均scPDSI 序列分解為10 個(gè)模態(tài)分量時(shí)，剩余殘余分量無(wú)法被再次分解，因此取子模態(tài)數(shù)量為10個(gè)。

2.2.2 GWO-SVM模型 SVM模型[18]是一種利用拉格朗日乘子法和對(duì)偶性質(zhì)求解優(yōu)化目標(biāo)的機(jī)器學(xué)習(xí)方法，可通過(guò)高斯核函數(shù)將非線性問(wèn)題轉(zhuǎn)化為線性問(wèn)題，被廣泛應(yīng)用于時(shí)間序列預(yù)測(cè)。相較于常見(jiàn)的BP神經(jīng)網(wǎng)絡(luò)模型，SVM模型的預(yù)測(cè)效果通常更好[19]。然而，采用SVM模型進(jìn)行預(yù)測(cè)時(shí)其參數(shù)取值對(duì)預(yù)測(cè)結(jié)果影響較大，選取合適參數(shù)對(duì)模擬效果至關(guān)重要。GWO是MIRJALILI等[20]基于自然界4種灰狼（α狼、β狼、δ狼、ω狼）捕食獵物活動(dòng)提出的一種群智能優(yōu)化算法，具有結(jié)構(gòu)簡(jiǎn)單、參數(shù)較少、收斂性能好等特點(diǎn)。采用GWO對(duì)SVM模型的兩個(gè)重要參數(shù)（懲罰系數(shù)C和高斯徑向基核函數(shù)參數(shù)σ）進(jìn)行優(yōu)選，可以提高SVM模型的預(yù)測(cè)精度和速度。

采用GWO選擇SVM 最優(yōu)參數(shù)C和σ時(shí)，首先，輸入樣本數(shù)據(jù)，設(shè)置C與σ取值范圍，初始化狼群數(shù)量和迭代次數(shù)，并計(jì)算每只灰狼的個(gè)體適應(yīng)度。然后，設(shè)置目標(biāo)函數(shù)，進(jìn)行SVM模型訓(xùn)練，在GWO模型中以C和σ為獵物進(jìn)行優(yōu)化，當(dāng)目標(biāo)函數(shù)值小于灰狼的個(gè)體適應(yīng)度時(shí)，將個(gè)體適應(yīng)度更新為當(dāng)前最優(yōu)目標(biāo)函數(shù)值。最后，隨著迭代次數(shù)不斷增加，精度逐漸提高，達(dá)到終止條件時(shí)輸出GWO 全局最優(yōu)值。采用最優(yōu)參數(shù)C與σ建立SVM 模型，并對(duì)測(cè)試樣本進(jìn)行預(yù)測(cè)、分析。GWO-SVM模型參數(shù)尋優(yōu)流程如圖1。

圖1 GWO-SVM模型參數(shù)尋優(yōu)流程圖Figure 1 GWO-SVM model parameters optimization flow diagram

2.2.3 基于CEEMDAN-GWO-SVM的scPDSI預(yù)測(cè)模型 基于CEEMDAN-GWO-SVM組合模型的山西省氣象干旱預(yù)測(cè)流程如圖2。本研究在每個(gè)IMF模態(tài)分量的預(yù)測(cè)中，使用前一個(gè)月的值對(duì)下一個(gè)月進(jìn)行預(yù)測(cè)，將預(yù)測(cè)所得值補(bǔ)充至原序列末尾構(gòu)成新序列，運(yùn)用新序列完成后續(xù)遞歸預(yù)測(cè)。

圖2 CEEMDAN-GWO-SVM模型流程圖Figure 2 CEEMDAN-GWO-SVM model flow diagram

2.2.4 RF模型 RF是由BREIMAN[21-22]在2001年提出的一種統(tǒng)計(jì)學(xué)習(xí)理論，RF模型基于自助重采樣技術(shù)，從總體訓(xùn)練樣本集S中有放回等概率地重復(fù)抽樣生成K個(gè)新的訓(xùn)練樣本集1，2，…，K，其中每個(gè)訓(xùn)練樣本集對(duì)應(yīng)一棵決策樹(shù)。在每棵樹(shù)的結(jié)點(diǎn)，隨機(jī)選取若干個(gè)特征因素進(jìn)行節(jié)點(diǎn)分裂，并按照“節(jié)點(diǎn)不純度最小”原則選擇一個(gè)特征因素對(duì)該節(jié)點(diǎn)進(jìn)行分裂。對(duì)于新輸入數(shù)據(jù)，對(duì)所有決策樹(shù)的預(yù)測(cè)結(jié)果取平均值作為最終預(yù)測(cè)結(jié)果。RF模型算法流程如圖3。

圖3 RF模型算法流程圖Figure 3 Algorithm flow diagram of RF model

RF 模型參數(shù)包括決策樹(shù)參數(shù)和Bagging 框架參數(shù)，其中各重要參數(shù)及其含義如表1。

表1 RF模型重要參數(shù)說(shuō)明表Table 1 Table of important parameters for RF model

2.2.5 模型預(yù)測(cè)效果評(píng)價(jià)指標(biāo) 在預(yù)測(cè)模型當(dāng)中，均方誤差(mean square error,MSE)和平均絕對(duì)誤差(mean absolute error, MAE)是常用的基礎(chǔ)預(yù)測(cè)評(píng)估指標(biāo)，均方根誤差(root mean square error, RMSE)是MSE的擴(kuò)展，其和決定系數(shù)R2也是預(yù)測(cè)模型中最常用的性能度量指標(biāo)[23]。因此，本研究采用RMSE、MAE和R2作為模型預(yù)測(cè)效果的評(píng)價(jià)指標(biāo)。其中，R2的取值范圍為[0,1]，R2越接近1，RMSE和MAE值越小，說(shuō)明模型的預(yù)測(cè)效果越好。三種指標(biāo)的計(jì)算公式分別如下[23]。

式中：yi和分別表示實(shí)測(cè)值和預(yù)測(cè)值；是yi的平均值；n為樣本數(shù)量。

3 結(jié)果與分析

3.1 scPDSI時(shí)間序列的CEEMDAN分解

本研究采用CEEMDAN 算法進(jìn)行scPDSI 序列分解時(shí)取噪聲系數(shù)為0.2，噪聲添加次數(shù)為500，EMD 內(nèi)部最大包絡(luò)次數(shù)為5 000。scPDSI 序列分解結(jié)果如圖4。由圖4 可知，原始scPDSI 序列的波動(dòng)幅度較大，而10 個(gè)IMF 分量的波動(dòng)幅度隨著分解的進(jìn)行逐漸減小，說(shuō)明CEEMDAN 算法有效降低了scPDSI序列的非平穩(wěn)性，有利于提升序列預(yù)測(cè)的準(zhǔn)確性。

圖4 基于CEEMDAN算法的scPDSI序列分解波形圖Figure 4 Decomposition waveform of the scPDSI series based on CEEMDAN

3.2 GWO-SVM模型參數(shù)優(yōu)化

GWO-SVM 模型參數(shù)率定時(shí)，取GWO 狼群數(shù)量為30，最大迭代次數(shù)為300，各子模態(tài)的SVM 模型的懲罰參數(shù)C和高斯核函數(shù)參數(shù)σ取值如表2。

表2 各子模態(tài)的SVM模型參數(shù)取值Table 2 Values of SVM model parameters for each IMF

3.3 RF模型的參數(shù)確定

RF 模型參數(shù)優(yōu)選時(shí)首先對(duì)外層的Bagging 框架進(jìn)行參數(shù)擇優(yōu)，然后再對(duì)內(nèi)層的決策樹(shù)模型進(jìn)行參數(shù)擇優(yōu)。在優(yōu)化RF 模型某一參數(shù)時(shí)，需要把其他參數(shù)設(shè)置為常數(shù)。經(jīng)多次調(diào)參，選取出山西省逐月scPDSI 序列的RF模型各參數(shù)最優(yōu)值(表3)。

表3 山西省逐月scPDSI的RF模型最優(yōu)參數(shù)Table 3 Optimal parameters of random forest model for scPDSI series in Shanxi Province

3.4 各模型預(yù)測(cè)效果比較

利用訓(xùn)練好的CEEMDAN-GWO-SVM 模型對(duì)測(cè)試集scPDSI 各IMF 分量進(jìn)行預(yù)測(cè)，并將預(yù)測(cè)結(jié)果進(jìn)行重構(gòu)，得到該模型最終預(yù)測(cè)結(jié)果。同時(shí)采用GWO-SVM 和RF 模型對(duì)預(yù)測(cè)結(jié)果進(jìn)行對(duì)比分析。以RMSE、MAE 和R2共3種評(píng)價(jià)指標(biāo)對(duì)各模型預(yù)測(cè)結(jié)果進(jìn)行評(píng)價(jià)(表4)。

表4 各模型評(píng)價(jià)指標(biāo)計(jì)算結(jié)果Table 4 Calculation results of prediction evaluation indicators of each model

由表4 可知，在測(cè)試集中，CEEMDAN-GWO-SVM 模型預(yù)測(cè)結(jié)果R2為0.762，RMSE 和MAE 分別為0.890 和0.557。與GWO-SVM 模型和RF 模型相比，CEEMDAN-GWO-SVM 模型預(yù)測(cè)結(jié)果R2分別提高5.6%和20.3%，RMSE 分別下降7.6%和19.5%，MAE 則分別降低19.9%和27.4%。其原因在于，CEEMDAN-GWO-SVM 模型進(jìn)行scPDSI序列預(yù)測(cè)時(shí)增加了基于CEEMDAN 算法的scPDSI序列分解環(huán)節(jié)，較好地緩解了模態(tài)混疊和噪聲的干擾，進(jìn)而提高了干旱特征因子R2預(yù)報(bào)水平，降低了預(yù)報(bào)誤差，預(yù)測(cè)效果優(yōu)于GWO-SVM模型和RF模型。

3.5 基于CEEMDAN-GWO-SVM預(yù)測(cè)模型的不同干旱等級(jí)預(yù)測(cè)效果評(píng)價(jià)

為進(jìn)一步明確CEEMDAN-GWO-SVM 模型對(duì)山西省氣象干旱預(yù)測(cè)的適用性，基于文獻(xiàn)[24]提出的氣象干旱等級(jí)劃分標(biāo)準(zhǔn)（表5），對(duì)scPDSI 測(cè)試集的預(yù)測(cè)值和原始值進(jìn)行干旱等級(jí)劃分(圖5)。各干旱等級(jí)月份數(shù)量統(tǒng)計(jì)值則見(jiàn)表6。

表5 scPDSI等級(jí)劃分Table 5 scPDSI ranks

表6 各干旱等級(jí)月份數(shù)量統(tǒng)計(jì)值Table 6 Statistical values of months for each drought grade

圖5 測(cè)試集預(yù)測(cè)值與原始值干旱等級(jí)劃分Figure 5 The drought classification for predicted values and original values in test set

通過(guò)分析表6可知，測(cè)試集中存在110個(gè)月預(yù)測(cè)氣象干旱等級(jí)與實(shí)際氣象干旱等級(jí)完全一致，預(yù)測(cè)準(zhǔn)確率為70.5%，未能準(zhǔn)確預(yù)測(cè)的級(jí)別中，有67.4%的月份預(yù)測(cè)等級(jí)差距為1 個(gè)級(jí)別，93.5%的月份預(yù)測(cè)等級(jí)差距為2個(gè)級(jí)別以?xún)?nèi)。表明基于CEEMDAN-GWO-SVM 模型的氣象干旱等級(jí)預(yù)測(cè)結(jié)果和實(shí)際情況接近，符合山西省氣象干旱趨勢(shì)。

4 討論與結(jié)論

目前，常用于氣象干旱預(yù)測(cè)的研究方法包括LSTM、ARIMA、SVM 和RF 模型等單一機(jī)器學(xué)習(xí)模型，對(duì)干旱指數(shù)進(jìn)行評(píng)估和預(yù)測(cè)的相關(guān)研究也取得了一定成果。吳晶等[6]基于RF 模型對(duì)流域內(nèi)空間和時(shí)間尺度的干旱進(jìn)行了預(yù)測(cè)；方秀琴等[25]側(cè)重于分析不同時(shí)間尺度和干旱等級(jí)劃分標(biāo)準(zhǔn)對(duì)RF模型性能的影響；PARK等[26]則研究了16 個(gè)干旱因素對(duì)氣象和農(nóng)業(yè)干旱的相對(duì)重要性。已有研究在預(yù)測(cè)時(shí)多選用SPI 作為干旱指數(shù)，而關(guān)于scPDSI的氣象干旱預(yù)測(cè)研究，特別是基于組合機(jī)器學(xué)習(xí)模型的scPDSI氣象干旱預(yù)測(cè)研究則較為罕見(jiàn)。本研究通過(guò)CEEMDAN和GWO對(duì)單一機(jī)器學(xué)習(xí)模型SVM進(jìn)行了優(yōu)化，相較于GWO-SVM和RF模型，本研究所提出的模型預(yù)測(cè)精度有較明顯的改善。然而，由于采用CEEMDAN 算法分解得到的前幾個(gè)模態(tài)分量仍然包含有噪聲，一定程度上干擾了預(yù)測(cè)效果，因此，降低了預(yù)測(cè)的精度。未來(lái)，本研究將進(jìn)一步考慮對(duì)分解后的各模態(tài)分量進(jìn)行篩選和聚類(lèi)，從而盡可能降低模態(tài)分量噪聲干擾對(duì)預(yù)測(cè)結(jié)果的影響，提高模型預(yù)報(bào)精度，為各地區(qū)氣象災(zāi)害預(yù)警預(yù)報(bào)提供技術(shù)支撐。

干旱預(yù)測(cè)是提高防旱抗旱能力的重要非工程措施。本研究提出將組合機(jī)器學(xué)習(xí)模型CEEMDAN-GWOSVM 應(yīng)用于山西省氣象干旱預(yù)測(cè)，通過(guò)與GWO-SVM 和RF 模型對(duì)比分析，得到主要結(jié)論：（1）基于“先分解，后重構(gòu)”思想的CEEMDAN-GWO-SVM 模型，將原始序列分解為10 個(gè)IMF 分量，對(duì)各IMF 分量進(jìn)行預(yù)測(cè)，并將預(yù)測(cè)結(jié)果進(jìn)行重構(gòu)，與GWO-SVM 相比較，預(yù)測(cè)結(jié)果R2提高了5.6%，RMSE 下降7.6%，MAE 則降低19.9%。和RF模型相比較，相應(yīng)值分別為20.3%、19.5%和27.4%，說(shuō)明CEEMDAN-GWO-SVM 模型在預(yù)測(cè)山西省逐月scPDSI序列上具有一定優(yōu)勢(shì)。（2）從干旱等級(jí)預(yù)測(cè)結(jié)果可發(fā)現(xiàn)，CEEMDAN-GWO-SVM 模型預(yù)測(cè)結(jié)果有70.5%的月份與實(shí)際完全一致，而未能準(zhǔn)確預(yù)測(cè)的月份中有93.5%差距為2 個(gè)級(jí)別以?xún)?nèi)，說(shuō)明本研究所提出的CEEMDANGWO-SVM模型適用于山西省氣象干旱預(yù)測(cè)研究，可為該省氣象干旱監(jiān)測(cè)提供一定的技術(shù)支持。