基于ESMD的西江流域降水時(shí)空變化特征分析

李繼清，陳思雨

（華北電力大學(xué)水利與水電工程學(xué)院，北京102206）

0 引 言

氣候與人類社會(huì)生活密切相關(guān)，隨著全球變暖問題日益嚴(yán)重［1］，進(jìn)一步加劇了降水這一氣候基本因子的時(shí)空分布差異，進(jìn)而對(duì)水資源、生態(tài)環(huán)境系統(tǒng)以及人類的生產(chǎn)生活產(chǎn)生影響。西江流域擁有豐富的水資源，降水是其徑流的補(bǔ)給來源，決定徑流分布情況。研究西江流域降水中長期的時(shí)空變化特性，能夠掌握該流域水文資源的周期、趨勢和突變規(guī)律，為該地區(qū)的水資源利用和管理提供科學(xué)依據(jù)。

氣候由多個(gè)時(shí)空尺度構(gòu)成了多層次結(jié)構(gòu)，是一個(gè)具有非平穩(wěn)性的非線性復(fù)雜系統(tǒng)。氣候變化特性研究領(lǐng)域已有較為系統(tǒng)的分析方法，所采用的數(shù)學(xué)分析方法也在逐步改進(jìn)和完善。早期研究中多采用累積距平、等值線等方法［2-5］，后期逐漸引入數(shù)學(xué)分析法，國內(nèi)外在研究各地區(qū)降水變化時(shí)廣泛采用Mann-Kendall 檢驗(yàn)法［6，7］，并結(jié)合小波分析等多種方法進(jìn)行對(duì)比［8-12］。此外，有序聚類、全球氣候模型等方法也可應(yīng)用于氣候分析［13-17］。但傳統(tǒng)方法受其理論基礎(chǔ)限制，在實(shí)際應(yīng)用中存在一定的局限性，極點(diǎn)對(duì)稱模態(tài)分解（Extreme-point Symmetric Mode Decomposition，ESMD）法作為近年來發(fā)展起來的局部自適應(yīng)時(shí)間序列分析技術(shù)［18］，善于從非平穩(wěn)、非線性的觀測序列中尋找變化趨勢并對(duì)其做出異常診斷。王金良和李宗軍［19］以及房賢水［20］分別從ESMD 方法統(tǒng)計(jì)特征、工作原理等方面對(duì)其優(yōu)勢進(jìn)行了論證。近年來，該方法也逐漸應(yīng)用到氣候變化的實(shí)例研究［21-23］，進(jìn)一步驗(yàn)證了其優(yōu)勢和可靠性。本文采用ESMD 法，對(duì)西江流域的降水觀測資料進(jìn)行周期、突變及趨勢分析，并與線性趨勢分析、Mann-Kendall 檢驗(yàn)法和小波分析法的分析結(jié)果進(jìn)行對(duì)比驗(yàn)證，進(jìn)而詳細(xì)分析西江流域的降水時(shí)空變化特征，為流域水資源的開發(fā)利用提供決策依據(jù)。

1 研究方法

極點(diǎn)對(duì)稱模態(tài)分解法（ESMD）能夠通過模態(tài)分解和時(shí)-頻分析，全面分析時(shí)間序列的周期、突變和趨勢特征，在水文時(shí)間序列分析領(lǐng)域具有較好的優(yōu)勢。分別利用線性趨勢分析法判別數(shù)據(jù)序列的整體線性變化趨勢，Mann-Kendall 檢驗(yàn)法判別數(shù)據(jù)數(shù)列的變化趨勢和突變時(shí)間點(diǎn)，小波分析法提取數(shù)據(jù)序列的周期特征，與ESMD方法的計(jì)算結(jié)果進(jìn)行對(duì)比分析。

1.1 極點(diǎn)對(duì)稱模態(tài)分解法（ESMD）

極點(diǎn)對(duì)稱模態(tài)分解（ESMD）方法在希爾伯特-黃變換法之上發(fā)展而來，在保留了經(jīng)驗(yàn)?zāi)B(tài)分解法優(yōu)點(diǎn)的同時(shí)也解決了其中存在的“模態(tài)混疊問題”，該方法在氣候分析中的優(yōu)勢為：擅長尋找變化趨勢，能夠分離出數(shù)據(jù)序列的年際變化趨勢和整體趨勢；擅長異常診斷，能夠觀察數(shù)據(jù)序列各模態(tài)中的異常時(shí)頻段；擅長時(shí)-頻分析，采用直接插值法能夠更為直觀地分析各時(shí)間尺度上的頻率變化［24］。

第一部分模態(tài)分解具體計(jì)算步驟如下：

（1）找出序列X的全部極值點(diǎn)，包括極大值點(diǎn)和極小值點(diǎn)，并將它們依次記作Ei（i=1，2，3，…，n）；

（2）用線段連接所有相鄰的極點(diǎn)并將線段中點(diǎn)依次記為Fi（1≤i≤n-1）；并通過一定方式補(bǔ)充左、右兩端邊界中點(diǎn)F0和Fn；

（3）利用所獲取的n+1 個(gè)中點(diǎn)構(gòu)建p條內(nèi)插曲線L1，…，Lp（p≥1），并計(jì)算其均值曲線L*；

（4）對(duì)X-L*序列重復(fù)上述步驟，直到|L*|≤ε 或者篩選次數(shù)達(dá)到預(yù)設(shè)的最大值K，得到第一個(gè)經(jīng)驗(yàn)?zāi)1；

（5）對(duì)剩余序列X-M1重復(fù)上述4 個(gè)步驟，直到剩余序列R只剩余一定數(shù)量的極點(diǎn)，便可分別得到經(jīng)驗(yàn)?zāi)2，M3，…；

（6）讓最大篩選次數(shù)K在整數(shù)區(qū)間［Kmin，Kmax］內(nèi)變化并重復(fù)上述幾個(gè)步驟。計(jì)算方差比率，并畫出隨K的變化圖，找出方差比率最小值對(duì)應(yīng)的最大篩選次數(shù)K0，此時(shí)分解得到的趨勢余項(xiàng)R為數(shù)據(jù)的最佳擬合曲線。以K0作為限制條件再次重復(fù)上述步驟，得到最佳分解結(jié)果。最終利用ESMD 方法將序列X分解為一系列經(jīng)驗(yàn)?zāi)B(tài)和一個(gè)趨勢余項(xiàng)。

第二部分時(shí)-頻分析認(rèn)為應(yīng)當(dāng)尊重離散數(shù)據(jù)的離散性特征，提出了針對(duì)數(shù)據(jù)的“直接插值法”，不但可以直觀地體現(xiàn)各模態(tài)的振幅與頻率的時(shí)變性，還可明確地獲知總能量變化。直接插值法基本思路如下：

（1）尋找極值點(diǎn)，計(jì)算兩個(gè)相鄰極大值點(diǎn)和相鄰極小值點(diǎn)之間的時(shí)間差；

（2）將這些時(shí)間段視為局部周期賦給其中點(diǎn)，畫出時(shí)間-周期對(duì)應(yīng)點(diǎn)圖；

（3）將這些局部周期值取倒數(shù)得到局部頻率，再做3次樣條插值得到光滑的時(shí)間-頻率變化曲線［25］。

圖1 ESMD方法流程圖Fig.1 The step of ESMD

1.2 對(duì)比驗(yàn)證方法

1.2.1 線性趨勢法

線性趨勢法能夠?qū)?shù)據(jù)序列的線性變化趨勢進(jìn)行分析。用yi表示樣本量為n的氣候變量，ti為對(duì)應(yīng)時(shí)刻，建立一元線性回歸方程：

式中：a為回歸常數(shù)；b為回歸系數(shù)即氣候變量的傾向趨勢［26］。

1.2.2 Mann-Kendall檢驗(yàn)法

Mann-Kendall 檢驗(yàn)法是一種非參數(shù)統(tǒng)計(jì)檢驗(yàn)方法，廣泛應(yīng)用于水文序列的趨勢、突變特征分析。利用M-K 統(tǒng)計(jì)值Z進(jìn)行趨勢統(tǒng)計(jì)的顯著性檢驗(yàn)，先假設(shè)該序列無趨勢，通過兩尾檢驗(yàn)，在給定顯著性水平下，查得臨界值Z1-α/2，當(dāng)|Z|＜Z1-α/2時(shí)，接受原假設(shè)；若|Z|＞Z1-α/2，則拒絕原假設(shè)。進(jìn)行突變分析需要計(jì)算M-K統(tǒng)計(jì)值UFk和UBk：

通過分析UFk和UBk交點(diǎn)位置即可明確突變發(fā)生時(shí)間［27］。

1.2.3 小波分析法

小波分析的基本思想是用一簇小波函數(shù)系來表示或逼近某一信號(hào)或函數(shù)。在小波變換系數(shù)Wf（c，d）中，將小波變換系數(shù)的平方值在d域上積分可得到小波方差：

式中：c為伸縮尺度；d為平移參數(shù)。

小波方差隨尺度c的變化過程為小波方差圖，能夠反應(yīng)信號(hào)波動(dòng)的能量隨尺度c的分布，進(jìn)而獲得數(shù)據(jù)序列的周期變化規(guī)律［28］。

2 研究結(jié)果與分析

采用ESMD 法對(duì)于西江流域的降水資料進(jìn)行模態(tài)分解和時(shí)-頻分析，結(jié)合線性趨勢分析、Mann-Kendall 檢驗(yàn)法和小波分析法，對(duì)各降水系列的周期、突變和趨勢特征進(jìn)行分析討論。

2.1 研究區(qū)概況及數(shù)據(jù)資料

西江是珠江流域的主干河流，發(fā)源于云南省，流經(jīng)地區(qū)地勢西北高東南低，全長共2 075 km，主要河流由南盤江、紅水河等河段組成。西江流域?qū)儆趤啛釒夂騾^(qū)，春夏季降雨充沛，秋冬較為干旱，水資源量變化主要受到降水季節(jié)分配的影響，徑流年內(nèi)分配極不均勻，汛期徑流量比例占年徑流量的80%以上。降水量由東到西呈現(xiàn)遞減趨勢［29］。

分析數(shù)據(jù)來自于西江流域1954-2008年51 個(gè)氣象站的日降水資料，空間分布如圖2，基本能夠均勻覆蓋西江整個(gè)流域。對(duì)日降水資料進(jìn)行整理計(jì)算，得到西江流域年降水、年降水極大值、年降水極小值、汛期降水和非汛期降水?dāng)?shù)據(jù)系列。其中年降水極值為全年日降水極大或極小值。根據(jù)地理情況將西江流域劃分為500 m 以上的上游高原區(qū)以及500 m 以下的中下游平原盆地區(qū)。平原區(qū)年降水Cv值為0.121，高原區(qū)年降水Cv值為0.089，圖3為兩地區(qū)年降水量對(duì)比。高原區(qū)年降水量約為1 100 mm，平原區(qū)年降水量約為1 500 mm，兩地區(qū)降水趨勢基本相同。

圖2 西江流域氣象站分布圖Fig.2 Distribution of meteorological stations in the Xijiang River Basin

圖3 西江流域平原區(qū)和高原區(qū)年降水概況圖Fig.3 Annual precipitation in the plain and plateau of the Xijiang River Basin

2.2 周期性

ESMD 法具有自適應(yīng)性和基于信號(hào)的局部變化特性，能夠更好地提取非平穩(wěn)、非線性序列的變化規(guī)律。通過ESMD 法對(duì)西江流域各降水序列進(jìn)行模態(tài)分解，可得到3～5 個(gè)模態(tài)分量以及一個(gè)趨勢項(xiàng)余量R，玉溪和連縣年降水分解結(jié)果圖4，各模態(tài)分量能夠反應(yīng)原數(shù)據(jù)序列中固有的振蕩特征，即對(duì)應(yīng)為降水序列的周期規(guī)律。表1為各降水序列的周期結(jié)果。西江流域各站之間各降水序列的固定周期成分均出現(xiàn)接近或相同的情況，具有統(tǒng)一的空間周期特征。

圖4 年降水模態(tài)分解圖Fig.4 Exploded modal decomposition maps of annual precipitation

表1 西江流域內(nèi)各站各降水指標(biāo)周期統(tǒng)計(jì)表 aTab.1 Periodic statistics of precipitation indicators at stations in the Xijiang River Basin

對(duì)西江流域各站各降水指標(biāo)周期成分出現(xiàn)次數(shù)進(jìn)行統(tǒng)計(jì)整理，結(jié)果如圖5。其中橙色代表年降水，紫色代表年降水極大值，綠色代表年降水極小值，黃色代表汛期，藍(lán)色代表非汛期。年降水18 a 的周期規(guī)律在站中出現(xiàn)頻率最高，其次為28 a 的周期特征。年極大值和年極小值的周期規(guī)律出現(xiàn)頻次最高的均為28 a，其次為2 a、14 a 等。汛期降水以存在3 a 和28 a 的周期特征最為普遍。非汛期降水序列中84%的氣象站存在2 a 的周期變化特征，其次是28 a 和18 a。各降水指標(biāo)之間的周期特征具有較高一致性，主要存在2、14、18和28 a的周期變化規(guī)律。

圖5 各降水指標(biāo)周期成分出現(xiàn)次數(shù)統(tǒng)計(jì)圖Fig.5 Statistics of the number of occurrences of the periodic components of various precipitation index

采用小波分析法與ESMD 方法周期分析結(jié)果進(jìn)行對(duì)比，表2為貴陽、羅甸、通道和威寧年降水的周期對(duì)比結(jié)果。兩方法下存在部分一致或相近的周期結(jié)果，周期越短結(jié)果一致性越好，周期越長越容易出現(xiàn)差異。究其原因，小波分析是以傅里葉變換為理論依據(jù)，在小波基選取等方面存在一定局限性，限制了該方法的實(shí)際運(yùn)用效果。而ESMD方法具有較強(qiáng)的靈活性和自適應(yīng)性，可通過數(shù)據(jù)自身特點(diǎn)進(jìn)行分解，探索事物的內(nèi)在規(guī)律。

表2 ESMD法和小波分析法下年降水周期對(duì)比結(jié)果表 aTab.2 Comparison of sudden changes in precipitation in the next year by ESMD method and wavelet analysis

2.3 突變性

根據(jù)ESMD 法得到的數(shù)據(jù)序列模態(tài)分量，采用直接插值法得到其對(duì)應(yīng)的頻率-振幅時(shí)變圖，如圖6為玉溪和連縣站年降水時(shí)-頻分析結(jié)果圖，其中F代表頻率，A代表振幅，根據(jù)兩者對(duì)應(yīng)關(guān)系，可知圖像在出現(xiàn)低頻、大幅度振蕩和高頻、小幅度振蕩時(shí)均代表數(shù)據(jù)序列發(fā)生了降水突變。對(duì)比圖中F和A的變化情況，可知玉溪年降水在1983年和1993年發(fā)生，連縣年降水在1971年發(fā)生突變。

圖6 年降水頻率-振幅時(shí)變圖Fig.6 Frequency-amplitude time-varying data of annual precipitation

通過觀察頻率-振幅時(shí)變圖得到各降水序列突變年份，結(jié)果整理統(tǒng)計(jì)如表3。空間上，各降水指標(biāo)在站間的突變結(jié)果均不完全一致，在相鄰站間突變結(jié)果一致性稍高，圖7梧州、廣寧、高要3 個(gè)相鄰氣象站的年降水均在觀測期內(nèi)發(fā)生三次突變現(xiàn)象，且突變時(shí)間較為接近。但降水突變?cè)谖鹘饔蛘w空間范圍內(nèi)未顯示出明顯的空間分布規(guī)律和特征。

圖7 梧州、廣寧和高要站的年降水突變結(jié)果圖Fig.7 Abrupt changes in annual precipitation at Wuzhou，Guangning and Gaoyao stations

表3 西江流域內(nèi)各站各降水指標(biāo)突變年份統(tǒng)計(jì)表Tab.3 Statistics of abrupt change of precipitation index of each station in Xijiang River Basin

通過整理各降水指標(biāo)的突變時(shí)間點(diǎn)分布年代，分析突變的時(shí)間分布特征，結(jié)果如圖8。各降水指標(biāo)在各年代內(nèi)均有突變現(xiàn)象發(fā)生，其中年降水和汛期降水的突變結(jié)果相似，均在20世紀(jì)60年代出現(xiàn)突變現(xiàn)象最為頻繁，其次為90年代；非汛期降水在80年代出現(xiàn)突變現(xiàn)象的次數(shù)最多，60年代的突變次數(shù)與汛期降水持平。年極大值的突變現(xiàn)象頻發(fā)于60年代和90年代，年極小值則在80年代出現(xiàn)突變情況最多。各降水指標(biāo)突變現(xiàn)象多集中在60年代、80年代和90年代，在50年代和00年代較少發(fā)生突變。

圖8 各降水指標(biāo)突變時(shí)間出現(xiàn)年代統(tǒng)計(jì)圖Fig.8 The chronological statistics of the abrupt time of each precipitation index

將ESMD 法與Mann-Kendall 檢驗(yàn)法結(jié)果進(jìn)行對(duì)比分析，圖9為玉溪和連縣的年降水M-K 結(jié)果圖，其中玉溪在1956年、1977年、1993年、2001年和2008年發(fā)生降水突變，連縣站的突變年份為1960年、1963年、1970年和1986年。

圖9 年降水M-K圖Fig.9 Results of annual precipitation based on M-K method

表4為蒙自、獨(dú)山、欽州和百色站年降水M-K法和ESMD法的突變結(jié)果對(duì)比。雖然兩方法算得突變結(jié)果不完全一致，但均存在相近或相同的突變結(jié)果，可知ESMD 方法進(jìn)行西江流域各降水序列突變分析具有一定的可靠性。兩方法相比，M-K 法僅能對(duì)數(shù)據(jù)序列中的突變點(diǎn)進(jìn)行判別，而ESMD 方法能夠在通過各模態(tài)頻率與振幅的對(duì)應(yīng)關(guān)系判斷數(shù)據(jù)序列突變點(diǎn)的同時(shí)了解其局部變化情況，具有明顯優(yōu)勢。

表4 ESMD法和Mann-Kendall檢驗(yàn)法下部分站年降水突變年份對(duì)比結(jié)果表Tab.4 Results of abrupt change in precipitation of some stations under the ESMD method and the M-K test

2.4 趨勢性

線性趨勢分析能夠判別西江流域各降水指標(biāo)的多年線性變化特征，統(tǒng)計(jì)結(jié)果如表5。年降水、年極大值和汛期降水線性趨勢分布情況相似，正、負(fù)增長站數(shù)基本持平，3 個(gè)指標(biāo)的正增長率最大均為防城站，負(fù)增長率最大均為南寧站，其中年極大值的線性變化速率略小，多不超過1 mm/a，年降水和汛期降水線性變化速率多不超過4 mm/a。年極小值和非汛期降水的線性趨勢多為正增長，但線性變化趨勢均不顯著，年極小值正增長率最大為連縣站，負(fù)增長率最大為融安站，非汛期降水的正增長率最大為靖西站，負(fù)增長率最大為高要站。

表5 西江流域各站各降水指標(biāo)線性變化速率結(jié)果表Tab.5 Results of the linear change rate of precipitation indexes at various stations in the Xijiang River Basin

ESMD 方法模態(tài)分解得到的趨勢項(xiàng)R可描述降水序列的具體變化過程，進(jìn)而分析西江流域各降水指標(biāo)趨勢特征。圖10為年降水指標(biāo)的變化趨勢類型，主要變化趨勢為上升-下降-上升型和先降后升型，多分布于中下游地區(qū)，除此以外上游高原區(qū)域存在少量的先升后降型和下降-上升-下降型。年極大值指標(biāo)的變化趨勢以上升-下降-上升、先升后降以及先降后升型為主，其中前兩種變化類型均多分布于流域中下游，此外上游地區(qū)還存在少量持續(xù)上升和下降-上升-下降的變化特征。年極小值指標(biāo)中，持續(xù)上升型、先降后升型、先升后降型等多種變化類型數(shù)量較為均勻，在流域內(nèi)交錯(cuò)分布。汛期降水與年降水分布情況相似，出現(xiàn)了少量站點(diǎn)呈現(xiàn)持續(xù)下降的變化特征。非汛期降水指標(biāo)中近半數(shù)站呈現(xiàn)先升后降的變化特點(diǎn)，多分布于流域中下游。由此可見，西江流域中下游地區(qū)各降水指標(biāo)變化趨勢類型較為統(tǒng)一，而位于流域西部的上游地區(qū)由于山區(qū)地形因素影響，各降水指標(biāo)的趨勢變化存在較大差異。

圖10 ESMD方法下年降水的各類變化趨勢圖Fig.10 Various trends of annual precipitation based on the ESMD method

采用M-K 趨勢檢驗(yàn)法對(duì)ESMD 法的趨勢結(jié)果進(jìn)行對(duì)比分析，M-K 趨勢檢驗(yàn)法給定顯著水平α=0.05，臨界值為±1.96，兩方法計(jì)算得到的桂平站年降水趨勢結(jié)果對(duì)比如圖11。在M-K法分析中，桂平站年降水前期呈現(xiàn)較為平穩(wěn)的小幅波動(dòng)，在90年代后呈現(xiàn)增長趨勢，整體降水序列M-K 統(tǒng)計(jì)值Z為1.19，表現(xiàn)為不顯著上升趨勢，與ESMD 法的趨勢結(jié)果具有較高的一致性，證明采用ESMD方法分析數(shù)據(jù)序列的趨勢性具有可靠性。

圖11 桂平站年降水ESMD法與M-K法趨勢結(jié)果對(duì)比圖Fig.11 Comparison of trend results of ESMD and M-K methods of annual precipitation at Guiping station

3 結(jié) 論

氣候變化與水文資源有著密切關(guān)聯(lián)，通過研究降水變化特性能為流域水資源管理提供有效手段。在氣候分析領(lǐng)域，ESMD 方法的優(yōu)勢在于能夠更為全面有效地識(shí)別數(shù)據(jù)序列的趨勢、周期和突變特征。采用ESMD 方法結(jié)合小波分析等對(duì)比方法分析西江流域的各降水指標(biāo)變化特征，得到如下結(jié)論：①西江流域各站的降水指標(biāo)周期特征具有統(tǒng)一性，主要存在2、14、18、28 a等周期規(guī)律；②每個(gè)站的各降水序列均存在1～5個(gè)突變點(diǎn)不等，各降水指標(biāo)突變多集中在20世紀(jì)60年代、80年代和90年代；③年降水、年極大值和汛期降水的正、負(fù)增長線性趨勢站數(shù)量基本持平，年極小值和非汛期降水的線性趨勢多呈現(xiàn)出上升現(xiàn)象；④年降水和非汛期降水在中下游區(qū)域的ESMD 法趨勢變化具有較高的一致性，其余各降水指標(biāo)的變化趨勢多呈現(xiàn)散亂分布。

在全球變暖環(huán)境下，西江流域降水指標(biāo)年內(nèi)、年際周期變化規(guī)律相近，年極小值及非汛期降水量統(tǒng)一呈現(xiàn)微弱增長趨勢，中下游地區(qū)降水指標(biāo)趨勢特征相對(duì)一致，降水偏多和偏少時(shí)期交替出現(xiàn)，但由于地理位置等因素影響，各站點(diǎn)仍存在其獨(dú)有的突變和趨勢特征，需進(jìn)一步探尋各水文要素的變化特性及各要素之間的相關(guān)聯(lián)系和影響機(jī)制。