黃河流域降水時空演變規(guī)律研究

張金萍，張 航

（1.鄭州大學(xué)水利科學(xué)與工程學(xué)院，鄭州 450001；2.鄭州大學(xué)黃河生態(tài)保護與區(qū)域協(xié)調(diào)發(fā)展研究院，鄭州 450001；3.鄭州市水資源與水環(huán)境重點實驗室，鄭州 450001）

0 引 言

由于中上游水土流失、人類活動影響及氣候變化等因素導(dǎo)致黃河水生態(tài)愈發(fā)脆弱［1，2］。降水作為黃河流域水資源的主要制約因素，受地形地貌及季風(fēng)氣候的影響，在黃河流域內(nèi)部呈現(xiàn)明顯的空間異質(zhì)性。如何從已發(fā)生變異的降水變量中分析出水文特征規(guī)律［3-5］，對于準(zhǔn)確把握黃河流域水資源特征，促進流域內(nèi)工農(nóng)業(yè)生產(chǎn)，制定生態(tài)保護措施具有重要意義。

在黃河流域已有的降水研究中，楊特群［6］根據(jù)黃河流域1951-2007年降水資料認(rèn)為其汛期降水量和年降水量變化規(guī)律大致相同。王雁［7］則利用黃河流域1951-2008年月降水資料對流域上、中、下游的降水趨勢進行了分析。王國慶［8］、何金梅［9］分別從不同角度對蘭州上游區(qū)間的降水進行了研究分析。周帥［10］基于標(biāo)準(zhǔn)化降水指數(shù)（SPI）分析了黃河流域干旱時空演變規(guī)律。陳磊［11］基于黃河流域106 個代表氣象站1960-2010年逐日降水資料，對季節(jié)降水的變化趨勢和變異情況進行了相關(guān)研究。黃建平［12］則利用歐洲中期天氣預(yù)報中心ERA5再分析資料對黃河流域過去40年的降水變化情況進行了研究。紹曉梅［13］、李占杰［14］、魚京善［15］等人采用小波分析或最大熵譜分析的方法對黃河流域2003年之前的降水周期性規(guī)律進行了研究。

以上研究主要對黃河流域的降水趨勢進行了分析，并采用小波分析和最大熵譜分析對其年降水周期變化規(guī)律進行了揭示，但受限于小波窗口和最優(yōu)截止矩的確定，具有一定的局限性，而對其年降水和汛期降水相關(guān)關(guān)系的研究則更加少見。為了對黃河流域降水在不同時間尺度的演化規(guī)律進行更深層次挖掘，本研究在黃河流域年降水序列趨勢檢驗的基礎(chǔ)上，運用信息熵理論計算黃河流域年、汛期降水在不同階段的時空信息熵，從而對不同時間尺度上的降水復(fù)雜性進行識別；使用CEEMDAN 方法對流域年、汛期降水序列進行多時間尺度分解，進而提取原始序列蘊含的周期性規(guī)律；最后應(yīng)用Copula 理論建立年降水和汛期降水的聯(lián)合分布，從而更精確地把握兩者的相依關(guān)系。本研究不僅揭示了黃河流域降水的變化趨勢，而且從時間和空間維度對其演變規(guī)律進行了細致剖析，相關(guān)研究結(jié)果可為黃河流域生態(tài)保護和高質(zhì)量發(fā)展提供支撐作用。

1 數(shù)據(jù)與方法

1.1 數(shù)據(jù)來源與處理

黃河流域貫穿中國半干旱半濕潤地區(qū)，大部分區(qū)域?qū)贉貛Т箨懶约撅L(fēng)氣候。本研究將黃河流域分為8 個二級水資源區(qū)，分別為龍羊峽以上（簡稱龍上區(qū)間）、龍羊峽至蘭州（簡稱龍-蘭區(qū)間）、蘭州至河口鎮(zhèn)（簡稱蘭-河區(qū)間）、內(nèi)流區(qū)、河口至龍門（簡稱河-龍區(qū)間）、龍門至三門峽（簡稱龍-三區(qū)間）、三門峽至花園口（簡稱三-花區(qū)間）、花園口以下（簡稱花下區(qū)間）。

流域內(nèi)部分氣象站點由于遷移、新建原因?qū)е聰?shù)據(jù)資料序列較短且存在殘缺，若使用插補處理會對精度產(chǎn)生較大影響。因此以中國國家青藏高原數(shù)據(jù)中心（http：//data.tpdc.ac.cn/）空間分辨率為0.008 333 3°（約1 km）的中國逐月降水?dāng)?shù)據(jù)為基礎(chǔ)進行研究分析。該數(shù)據(jù)集是通過Delta 空間降尺度方案在中國地區(qū)降尺度生成，其結(jié)果已使用496 個獨立氣象觀測站點數(shù)據(jù)進行驗證，驗證結(jié)果可信。為了進行相關(guān)指標(biāo)的計算，在研究中以黃河流域內(nèi)321 個代表氣象站點基礎(chǔ)（如圖1），提取各個氣象站點所在柵格1960-2017年的月降水?dāng)?shù)據(jù)作為站點資料進行研究。

1.2 研究方法

1.2.1 信息熵

熵最初來源于熱力學(xué)，1997年，Singh［16］較為詳細地介紹了信息熵在水文水資源等方面的應(yīng)用，國內(nèi)也有學(xué)者［17，18］使用信息熵對部分區(qū)域的降雨規(guī)律進行了研究。熵的計算如公式（1）所示：

式中：pi為變量xi所對應(yīng)的概率。

Ebrahimi 等人［19］對信息熵的計算公式進行了改進，提出了在總體分布未知情況下的熵值計算公式，如公式（2）所示：

其中：

式中：yi是原始數(shù)據(jù)升序后的序列。當(dāng)i≤m時，yi-m=y1；i≥nm時，yi+m=yn。

1.2.2 CEEMDAN分解

CEEMDAN 算法通過在原始信號上加入有限方差約束的多組獨立同分布自適應(yīng)白噪聲實現(xiàn)對非平穩(wěn)原始信號不同頻率成分的提取［20］。設(shè)x(t)為待分解序列，w(t)表示均值為0，方差為1 的高斯白噪聲，r(t)為殘余差值，IMFk為第k階模態(tài)分量，CEEMDAN分解過程如下：

（1）向x(t)中加入高斯白噪聲，并進行I次實驗，第i次實驗信號可表示為：

（2）對第i次加入白噪聲后的信號xi(t)進行EMD 分解并進行I次實驗平均，得到第一階模態(tài)分量IMF1，然后計算x(t)與IMFi分量的殘余差值r1(t)。

（3）再次向殘余差值r1(t)中加入白噪聲，進行I次實驗平均，第i個殘余差值可表示為：

（4）對第i次加入白噪聲的信號ri1(t)進行EMD 分解，第i次為IMFi2分量，通過I次實驗平均，得到第二階IMF分量為：

（5）此時，分解得到的殘余差值為r2(t) =r1(t) -IMF2，返回步驟（4）、（5）計算下一階模態(tài)，則第k階模態(tài)為：

重復(fù)第3 步至第5 步直到差值rk(t)滿足以下條件之一時終止分解：①不能被EMD 進一步分解；②IMF分量滿足相應(yīng)條件；③局部極值點的個數(shù)小于3 個。最終，原始信號x(t)可分解為若干個IMF分量和一個趨勢項rK(t)：

1.2.3 Copula函數(shù)

Copula 函數(shù)具有將聯(lián)合分布和各自的邊緣分布聯(lián)結(jié)起來的作用，從而實現(xiàn)變量間隨機性和相依性的分離，廣泛應(yīng)用于降雨、洪水、干旱和枯水等多特征屬性頻率分析研究中［21-23］。應(yīng)用最廣泛的為Archimedean Copula，其類型主要包括Gumbel Copula、Clayton Copula 和Frank Copula 三種（見表1）。表中θ為參數(shù)，τ為Kendall秩相關(guān)系數(shù)。

表1 常用Archimedean Copula函數(shù)Tab.1 Three commonly used Archimedean Copula functions

2 結(jié)果與分析

2.1 降水變化趨勢

1960-2017年間，黃河流域年均降水量500.6 mm，龍上、龍-蘭、蘭-河、內(nèi)流區(qū)、河-龍、龍-三、三-花、花下區(qū)間各水資源區(qū)的年均降水量分別為543.4、458.5、250.3、349.1、443.3、563.6、619.6、654.5 mm。整體來看，黃河流域河源區(qū)及下游降水較為豐富，高于流域平均水平。

MK 檢驗結(jié)果（圖2）表明，黃河流域年降水序列自20 世紀(jì)70年代以來開始呈現(xiàn)一種明顯的下降趨勢，在20 世紀(jì)初，這種下降趨勢超過了0.05 的顯著性水平（u0.05= 1.96）。這和王飛等人［24］采用SPEI指數(shù)得到黃河流域干旱呈顯著增加趨勢的結(jié)論是一致的。河-龍、龍-三、三-花和花下4 個水資源區(qū)變化趨勢與流域整體比較接近；流域上游的龍上和龍-蘭水資源區(qū)年降水序列具有緩慢增加的趨勢，王遠見等人［25］通過研究同樣發(fā)現(xiàn)黃河源區(qū)年降雨量呈增加趨勢，只不過在其結(jié)論中認(rèn)為這種上升趨勢是顯著的；蘭-河、內(nèi)流區(qū)水資源區(qū)年降水序列主要呈下降趨勢，但相對于流域整體而言，其下降趨勢并不明顯。

此外，從圖2中UF統(tǒng)計量可以初步判斷黃河流域各水資源區(qū)年降水序列在2000年左右變化最為顯著，為進一步確定突變時間，結(jié)合滑動T檢驗綜合確定2002年為顯著突變時間點，并以此劃分突變前（1960-2001年）和突變后（2002-2017年）時段。

圖2 年降水序列MK檢驗結(jié)果Fig.2 MK test results of annual precipitation series

2.2 降水時空信息熵

2.2.1 降水時間熵

表2、3 分別為各水資源區(qū)年、汛期降水量突變前后兩階段的時間信息熵。由表2 可知，黃河流域年降水序列時間熵在突變后相比突變前降低5.8%，即2002年以后，黃河流域降水量在年際間的波動減弱。各水資源區(qū)中，龍-蘭水資源區(qū)時間熵減小幅度最大，降水量在年際間的波動相比突變前出現(xiàn)了顯著降低；三-花、花下水資源區(qū)年降水時間熵變化幅度較小，降水在年際間波動較為穩(wěn)定，其余水資源區(qū)降水在年際間的波動位于兩者之間。

表2 年降水序列突變前后時間信息熵Tab.2 Temporal information entropy value of annual rainfall in different stages

由表3 可知，蘭-河、黃河內(nèi)流區(qū)水資源區(qū)所對應(yīng)的汛期降水時間熵突變后相比突變前降幅均超過10%，說明2002年以后汛期降水量在年際間的波動顯著減弱；三-花、花下水資源區(qū)汛期時間熵變化較小，說明其汛期降水量在整個研究時段內(nèi)的波動比較平穩(wěn)。

表3 汛期降水突變前后時間信息熵Tab.3 Temporal information entropy value of flood season rainfall in different stages

此外，對比年和汛期時間信息熵可得到如下結(jié)論：

（1）無論年尺度還是汛期尺度，突變后的時間信息熵相比突變前均出現(xiàn)了降低，說明近些年的降水量雖然呈下降趨勢，但降水在黃河流域和各水資源區(qū)內(nèi)部沿時程的復(fù)雜程度也出現(xiàn)了下降，不確定性降低，有利于對降水規(guī)律進行研究。

（2）降水時間信息熵在不同水資源分區(qū)內(nèi)具有空間差異性。首先表現(xiàn)為突變前后各水資源區(qū)在汛期和年尺度上的時間信息熵值大小不同，且發(fā)生變異后下降幅度不一。其原因主要是不同水資源分區(qū)在地理位置、局部氣候、地勢地貌等方面不同所導(dǎo)致。其次，位于龍-三水資源分區(qū)內(nèi)的部分氣象站點汛期降水時間熵在突變后相比突變前出現(xiàn)了增加，但對龍-三水資源分區(qū)而言，突變后的信息熵相較于突變前是下降的，因此存在局部降水復(fù)雜性上升、不確定性增加的現(xiàn)象。

2.2.2 降水空間熵

由黃河流域年降水空間熵結(jié)果（見圖3）可以發(fā)現(xiàn)：

圖3 年降水空間信息熵Fig.3 Spatial information entropy value of annual rainfall

（1）2001年以后，黃河流域年降水空間熵出現(xiàn)小幅度降低，說明從流域整體角度而言，突變后年降水量在空間上的不確定性出現(xiàn)了下降，但這種趨勢并不明顯；

（2）河-龍及龍-三水資源分區(qū)空間熵出現(xiàn)了升高，其余水資源分區(qū)的空間熵下降，說明突變以后各水資源區(qū)降水復(fù)雜性出現(xiàn)了不同方向的改變，但這種改變以復(fù)雜性降低為主。

通過對空間熵進行MK 檢驗發(fā)現(xiàn)，在研究時段內(nèi)，黃河流域和龍上、龍-蘭、蘭-河、內(nèi)流區(qū)、河-龍、花下水資源區(qū)在1990年后空間熵呈現(xiàn)下降趨勢，龍-三、三-花水資源區(qū)的空間熵則呈現(xiàn)一種微弱的上升趨勢。說明黃河流域內(nèi)部降水在空間上的復(fù)雜性主要是降低的，即降水的空間分布更加均勻。

進一步對黃河流域年降水空間熵（Hs）和年降水變異系數(shù)（Cv）進行相關(guān)分析（圖4）可以看出，降水空間熵和變異系數(shù)之間存在明顯的正相關(guān)。變異系數(shù)代表著流域內(nèi)部各站點降雨的不均勻程度，變異系數(shù)越大，說明各流域內(nèi)部降水越不均勻。由此可見，空間差異性蘊藏在空間信息熵所蘊含的復(fù)雜性之中，即差異性是復(fù)雜性的一種體現(xiàn)。此外，黃河流域年降水空間熵＞汛期降水空間熵＞9月降水空間熵，說明時間跨度越長，所包含的不確定性越大，這與信息熵所代表的物理含義是一致的。就黃河流域而言，其汛期空間熵突變后是降低的，但9月份降水空間熵在突變后是增大的。由此可見，黃河流域9月份降水所包含的不確定性是增大的，汛期和年的變化相一致。因此在今后的降雨觀測過程中，應(yīng)密切注意9月份降水的變化情況。

圖4 Hs與Cv相關(guān)關(guān)系Fig.4 Correlation between Hs and Cv

2.3 多時間尺度分解

圖5（a）為黃河流域年降水序列，具有較強的非線性和非平穩(wěn)特性。圖5（b）、（c）為多時間尺度分解結(jié)果，突變前的原始序列最終被分解為5 層，包括4 個IMF分量和1 個趨勢項（RES余量）；突變后的原始序列被分解為4層，包括3個IMF分量和1個趨勢項。分解出的IMF分量較好地體現(xiàn)了原始信號序列在不同周期上的波動規(guī)律和變化特征。

圖5 黃河流域年降雨序列多時間尺度分解Fig.5 Multi-time scale decomposition results of annual rainfall series

突變前后各IMF分量周期波動和振幅變化特征如表4所示。

由表4可知，突變前：從IMF1分量到IMF4分量，周期增大，振幅減小；IMF1 分量的準(zhǔn)周期在2～3 a 范圍內(nèi)，說明其變化頻率較高，且其最大振幅和平均振幅都處于各分量的首位；IMF1～IMF4 和RES分量的方差貢獻率分別為74.41%、11.9%、2.64%、1.23%和9.81%，由此可以看出，具有短周期、高振幅的高頻分量包含了原始序列較多的信息量，在原始序列中起著主要作用；IMF1分量的平均振幅較大，IMF2分量的平均振幅位居其次，IMF3、IMF4 分量的平均振幅相差不大，這與最大振幅的變化規(guī)律是相一致的。從突變前的RES余量可以看出，黃河流域年降雨量在突變前呈現(xiàn)一種下降趨勢。

表4 年降水序列多時間尺度分解結(jié)果Tab.4 Multi-time scale decomposition results of annual rainfall series

突變后：原始序列分解出3 個IMF分量和一個趨勢項。IMF1～IMF3 分量周期逐漸增大，IMF1 分量包含2 a 和3 a 兩個準(zhǔn)周期，IMF2 分量為4 a 的準(zhǔn)周期，IMF3 分量為9 a 的準(zhǔn)周期。各分量的貢獻率分別為84.23%、7.08%、6.10%和2.59%；RES分量呈現(xiàn)先下降后上升的變化趨勢，前半程的下降趨勢對應(yīng)著突變前RES分量的下降趨勢，是相一致的，后半程的上升趨勢表明從全局而言，黃河流域降水在最近幾年逐漸開始恢復(fù)，這和原始序列的變化趨勢也是相一致的。

對各二級水資源區(qū)進行CEEMDAN分解，鑒于篇幅原因，圖6 僅列出龍上和花下區(qū)間年降水序列多周期分解結(jié)果，并得到以下結(jié)論：

圖6 各水資源區(qū)年降水多周期分解結(jié)果（以龍上和花下區(qū)間為例）Fig.6 Multi-time scale decomposition results of annual rainfall in secondary water resources area

（1）各水資源區(qū)IMF1 分量準(zhǔn)周期均為2～4 a，即占主要貢獻的高頻分量周期基本一致，表明黃河流域水資源區(qū)包含著相似的短周期變化規(guī)律。因此，無論是對黃河流域還是各水資源分區(qū)進行研究時都應(yīng)關(guān)注短周期尺度上的變化規(guī)律，即對2～4 a的周期進行重點研究。

（2）突變后序列分解出的IMF分量個數(shù)均減少，造成該結(jié)果的原因可能有兩點：①突變后黃河流域內(nèi)部降水的復(fù)雜程度和不確定性是降低的，包含的信息量減少，導(dǎo)致分解出的IMF分量減少；②突變后原始數(shù)據(jù)的時間序列較短。

（3）由RES分量可知，突變前各子水資源區(qū)降水量變化趨勢和黃河流域整體基本一致，呈下降趨勢；突變后，龍-蘭、蘭-河、內(nèi)流區(qū)、河-龍、龍-三5 個水資源區(qū)和流域整體在研究時段末呈現(xiàn)一種上升趨勢，龍上、三-花和花下3 個水資源區(qū)在研究時段后期則呈現(xiàn)一種下降趨勢。表明突變前各水資源分區(qū)和流域整體的變化趨勢是一致的，都呈現(xiàn)一種下降趨勢，但突變造成各水資源分區(qū)存在一定差異性。

2.4 多時間尺度熵

通過對各IMF分量進行歸一化處理，計算得到多時間尺度熵（表5），其代表了原始序列在不同時間尺度上的變化特征：

表5 黃河流域年降水序列突變前后多時間尺度熵Tab.5 Multi-time scale entropy of rainfall in different stages

（1）突變后各IMF分量熵值相比突變前均出現(xiàn)不同程度的下降，說明變異導(dǎo)致降水要素在不同尺度上的復(fù)雜程度下降，隨機性降低，這與原始序列時空信息熵的變化規(guī)律是一致的。

（2）無論突變前還是突變后，IMF1分量的熵值最大，表明高頻短周期分量所包含的信息量比較豐富，最能體現(xiàn)原始序列的變化規(guī)律。因此對于黃河流域而言，對降水的研究應(yīng)集中在2～4年的短周期上。

2.5 年-汛期降水聯(lián)合分布

黃河流域年降水和汛期降水關(guān)系密切，具有較強的相關(guān)關(guān)系。本研究使用參數(shù)法確定各變量邊緣分布，并利用Kolmogorov-Smirnov（KS）方法進行檢驗，檢驗結(jié)果顯示，突變前年降水、汛期降水分別以伽馬分布和韋伯分布進行擬合時效果最優(yōu)，見圖7。

圖7 突變前年降水序列經(jīng)驗頻率和理論頻率擬合圖Fig.7 Fitting results of empirical frequency and theoretical frequency of annual rainfall before mutation

采用相關(guān)性指標(biāo)法對Copula 函數(shù)進行參數(shù)估計，并以AIC信息準(zhǔn)則、OLS規(guī)劃原則作為擬合優(yōu)度評價指標(biāo)，結(jié)果（表6）表明Gumbel-Copula 函數(shù)在構(gòu)建突變前年降水量和汛期降水量之間的聯(lián)合分布時，擬合效果最優(yōu)。

表6 突變前擬合優(yōu)度檢驗結(jié)果Tab.6 Goodness of fit test results before mutation

突變前最終的Copula 函數(shù)如式（10）所示，其分布函數(shù)與等高線圖如圖8所示：

圖8 年雨量-汛期雨量聯(lián)合分布Fig.8 Joint distribution of annual and flood season precipitation

同理，求出突變后最優(yōu)Copula函數(shù)如式（11）所示：

以37.5%、62.5%分別作為豐枯分界，計算突變前后年、汛期降水豐枯遭遇概率，計算結(jié)果如表7、8所示，可以發(fā)現(xiàn)：

表7 突變前年-汛期降水聯(lián)合概率Tab.7 Joint probability of annual and flood season precipitation before mutation

（1）不同階段（突變前后），年和汛期降水聯(lián)合分布豐枯同步頻率大于豐枯異步概率，且同步頻率主要集中在同豐或同枯狀態(tài)；

（2）降水要素發(fā)生突變后，同步頻率增大，從突變前的0.655 增加到突變后的0.722。主要體現(xiàn)在同豐和同平狀態(tài)的增加。

上述兩個結(jié)論表明2002年以后黃河流域降水在年和汛期尺度上更加一致，即年內(nèi)降水更集中在汛期，這對于防洪是不利的，需要有關(guān)部門提前做好防汛措施。

3 結(jié) 論

以黃河流域1960-2017年月降水?dāng)?shù)據(jù)為基礎(chǔ)，基于信息熵、多時間尺度和Copula 函數(shù)分析了黃河流域降水在時間和空間上的變化規(guī)律，主要得出以下結(jié)論。

（1）過去50年間，黃河流域降水整體表現(xiàn)出一種下降趨勢，但源區(qū)降水在緩慢增加，且降水在黃河流域內(nèi)部呈現(xiàn)明顯的空間異質(zhì)性；

（2）突變后（2002年后），流域降水時空信息熵均出現(xiàn)了降低，其中龍羊峽-蘭州水資源區(qū)的時間信息熵減小幅度最大，超過10%，空間熵則在黃河流域內(nèi)部呈現(xiàn)不同方向的變化；

（3）通過CEEMDAN 方法將突變前后黃河流域年降水序列分別分解為4層（3個IMF分量和1個RES余量）和3層（2個IMF分量和1 個RES余量）。結(jié)果表明對黃河流域降水的觀測和研究應(yīng)主要集中在2～4 a 的短周期上，且RES余量反映了原始序列整體下降但在近幾年逐漸恢復(fù)的趨勢；

表8 突變后年-汛期聯(lián)合降水概率Tab.8 Joint probability of annual and flood season precipitation after mutation

（4）通過Gumbel-Copula 函數(shù)發(fā)現(xiàn)在不同階段黃河流域年、汛期降水聯(lián)合分布以豐枯同頻為主，同豐和同枯的重現(xiàn)期分別在3～4 a、4～5 a之間。此外，2002年以后，同步頻率增加，年內(nèi)降水更集中在汛期，對于防洪提出了更高要求。

本研究主要基于信息熵、CEEMDAN 和Copula 理論對黃河流域年、汛期降水的時空變化規(guī)律進行了研究，而黃河流域極值暴雨同樣是影響流域內(nèi)社會經(jīng)濟活動和生態(tài)穩(wěn)定的重要因素。受限于資料原因，在此并未對其進行研究分析，未來在獲取更短時間步長降水?dāng)?shù)據(jù)的基礎(chǔ)上，對極值暴雨事件開展研究，對于指導(dǎo)流域內(nèi)生產(chǎn)生活具有更加實際的意義。 □