基于時(shí)變矩模型的湄公河上游非一致性洪水頻率分析

龍?jiān)剖妫?葉長(zhǎng)青， 李運(yùn)剛， 熊合勇， 趙承遠(yuǎn)

(1.云南大學(xué) 云南省國(guó)際河流與跨境生態(tài)安全重點(diǎn)實(shí)驗(yàn)室， 云南 昆明 650504； 2.海南大學(xué) 生態(tài)與環(huán)境學(xué)院， 海南 海口 570228； 3.華能瀾滄江水電股份有限公司， 云南 昆明 650214)

1 研究背景

在水利水電工程設(shè)計(jì)中，水文頻率分析是設(shè)計(jì)洪水計(jì)算及確定重現(xiàn)期的關(guān)鍵[1-2]。傳統(tǒng)水文頻率計(jì)算建立在水文時(shí)間序列一致的基礎(chǔ)上，近些年來(lái)，在氣候變化和人類(lèi)活動(dòng)的影響下許多河流的水文序列一致性遭到破壞，傳統(tǒng)水文頻率分析方法的適用性受到質(zhì)疑[3-5]。變化環(huán)境下非一致性水文頻率分析已成為水工設(shè)計(jì)規(guī)劃、洪水分析和風(fēng)險(xiǎn)管理中的關(guān)鍵技術(shù)環(huán)節(jié)，也是水文科學(xué)研究中的熱點(diǎn)問(wèn)題[6-7]。

非一致性洪水頻率研究主要集中在水文序列非一致性診斷和非一致性洪水頻率分析兩個(gè)方面。國(guó)內(nèi)外代表性的非一致性檢驗(yàn)方法有Mann-Kendall、Pettitt、Spearman、Brown-Forsythe檢驗(yàn)和貝葉斯方法等[7-8]。在非一致性水文頻率分析方面，常用水文極值序列還原/還現(xiàn)的方法，如降雨徑流關(guān)系法、時(shí)間序列的分解與合成法及水文模型法等[9-10]。這類(lèi)方法對(duì)水文序列進(jìn)行還原/還現(xiàn)處理后，采用一致性假設(shè)的方法對(duì)水文序列進(jìn)行頻率計(jì)算，然而還原/還現(xiàn)方法的一致性修正成果往往存在不確定性[6-7]。因此，基于非平穩(wěn)極值序列的直接水文頻率分析逐漸受到重視，如混合分布法[11]、條件概率分布法[12]、時(shí)變矩模型[13]和GAMLSS(generalized additive models for location, scale and shape)模型[14]等。其中，時(shí)變矩模型考慮均值和方差的趨勢(shì)性，將趨勢(shì)性成分嵌入到了概率分布的一、二階矩中(時(shí)變矩)，可得到設(shè)計(jì)值隨時(shí)間的變化關(guān)系。時(shí)變矩模型是描述單變量水文序列非一致性的有效工具，已經(jīng)被應(yīng)用到許多地區(qū)[15-17]。

瀾滄江-湄公河流經(jīng)中國(guó)、緬甸、老撾、泰國(guó)、柬埔寨和越南6個(gè)國(guó)家，是東南亞最重要的國(guó)際河流。在氣候變化和人類(lèi)活動(dòng)影響尤其是大規(guī)模梯級(jí)水電開(kāi)發(fā)背景下，瀾滄江-湄公河干流徑流特征已經(jīng)改變，水文序列一致性遭到破壞[8,18]。本研究以干流泰國(guó)清盛水文站為例，利用Mann-Kendall趨勢(shì)突變檢測(cè)方法對(duì)年最大日流量序列進(jìn)行非一致性分析，基于時(shí)變矩模型選擇6種概率分布和7種趨勢(shì)模型進(jìn)行組合，產(chǎn)生42種競(jìng)爭(zhēng)模型進(jìn)行比較擇優(yōu)，分析變化環(huán)境下洪水設(shè)計(jì)值的響應(yīng)規(guī)律，以期為流域水工規(guī)劃、防洪安全和跨境水資源利用管理等提供科學(xué)依據(jù)。

2 數(shù)據(jù)來(lái)源與研究方法

2.1 研究區(qū)概況及數(shù)據(jù)來(lái)源

瀾滄江-湄公河河長(zhǎng)4 880 km, 流域面積為81×104km2，多年平均徑流量為4 750×108m3(圖1)。流域主要受東南季風(fēng)和西南季風(fēng)的影響，降水量豐富但季節(jié)分配不均，有著較為明顯的干濕季。2008年全流域水庫(kù)有效庫(kù)容(86×108m3)僅占多年平均徑流量的2%，而預(yù)計(jì)到2025年水庫(kù)有效庫(kù)容將增加到868×108m3，達(dá)到多年平均徑流量的19%[19]。截止目前中國(guó)境內(nèi)瀾滄江干流已建成11個(gè)梯級(jí)大壩，總庫(kù)容約444.88×108m3(表1)。位于泰國(guó)境內(nèi)的清盛水文站是湄公河干流上離中國(guó)最近的水文站，控制面積為18.9×104km2，多年平均流量為2 653 m3/s。

表1 瀾滄江干流已建梯級(jí)電站情況

圖1 瀾滄江-湄公河上游流域水系及已建水電站分布

本研究中采用的清盛站1960-2019年逐日流量數(shù)據(jù)來(lái)源于湄公河委員會(huì)(https://www.mrcmekong.org/)。

2.2 研究方法

首先對(duì)水文序列采用Mann-Kendall法[20]進(jìn)行趨勢(shì)分析和突變檢驗(yàn)。若序列的一致性遭到破壞，則選擇時(shí)變矩方法進(jìn)行分析。時(shí)變矩方法主要分析水文頻率曲線(xiàn)特征參數(shù)隨時(shí)間變化的影響，認(rèn)為水文序列的第一階矩均值(m)和第二階矩標(biāo)準(zhǔn)差(σ)隨時(shí)間具有線(xiàn)性或拋物線(xiàn)性趨勢(shì)變化特征，并引入時(shí)間t對(duì)水文頻率曲線(xiàn)進(jìn)行函數(shù)描述[12]。時(shí)變矩模型可由不同的概率分布線(xiàn)型和水文序列的第一、二階矩的趨勢(shì)模型相互組合得到。本研究選擇P-Ⅲ型分布(P-Ⅲ)、廣義極值分布(GEV)、廣義邏輯斯諦分布(GLO)、正態(tài)分布(NORM)、二參數(shù)對(duì)數(shù)正態(tài)分布(LN2)和耿貝爾分布(GMB)共6種概率分布函數(shù)。為分析前兩階矩隨時(shí)間的變化關(guān)系，選取合適的趨勢(shì)模型嵌入到分布模型中。考慮到曲線(xiàn)的外延性，采用線(xiàn)性(L)和拋物線(xiàn)性(P)趨勢(shì)。將趨勢(shì)嵌入頻率曲線(xiàn)的第一、二階矩中可得到7類(lèi)趨勢(shì)模型，具體描述參考文獻(xiàn)[13]、[16]。采用極大似然法進(jìn)行參數(shù)估計(jì)，并使用AIC(Akaike information criterion)準(zhǔn)則選擇最優(yōu)模型作為時(shí)變矩模型進(jìn)行相應(yīng)的頻率計(jì)算。時(shí)變矩方法具體的計(jì)算流程框架如圖2所示。

圖2 時(shí)變矩方法計(jì)算流程框架

3 結(jié)果與分析

3.1 洪水系列非一致性診斷

圖3為1960-2019年清盛水文站年最大日流量序列變化過(guò)程及Mann-Kendall突變檢驗(yàn)。由圖3(a)可以看出，1960-2019年年最大日流量呈“顯著下降-平穩(wěn)-下降”的變化過(guò)程。其中，1966和2006年分別發(fā)生了洪量為23 500和29 300 m3/s的特大洪水。自2008年后，年最大日流量雖有波動(dòng)，但均小于10 000 m3/s。圖3(b)中趨勢(shì)和突變分析結(jié)果表明，清盛水文站的年最大日流量呈顯著的下降趨勢(shì)(p<0.01)，突變點(diǎn)出現(xiàn)在2008年。

圖3 1960-2019年清盛水文站年最大日流量變化及Mann-Kendall突變檢驗(yàn)

3.2 時(shí)變矩模型優(yōu)選

基于時(shí)變矩方法，選擇6種概率分布線(xiàn)型和7種趨勢(shì)模型相互組合，共得到42種備選模型。采用極大似然法給出最優(yōu)線(xiàn)型的參數(shù)估計(jì)值，并計(jì)算各種模型組合的AIC準(zhǔn)則值，結(jié)果見(jiàn)表2。根據(jù)AIC準(zhǔn)則，選擇AIC值最小的模型組合作為最優(yōu)模型組合，則清盛水文站年最大日流量序列最優(yōu)擬合分布為GLO分布，最優(yōu)擬合趨勢(shì)模型為AL趨勢(shì)模型，時(shí)變矩最優(yōu)模型為GLO-AL模型，即均值具有線(xiàn)性趨勢(shì)而標(biāo)準(zhǔn)差不考慮趨勢(shì)變化。

表2 清盛水文站時(shí)變矩模型AIC擬合檢驗(yàn)值

3.3 水庫(kù)修建不同階段洪水線(xiàn)型響應(yīng)規(guī)律及重現(xiàn)期變化

3.3.1 不同時(shí)間基點(diǎn)下洪水線(xiàn)型響應(yīng)規(guī)律 結(jié)合Mann-Kendall法突變點(diǎn)檢驗(yàn)并考慮瀾滄江梯級(jí)水電工程的建設(shè)時(shí)間(表1)，選擇時(shí)間基準(zhǔn)點(diǎn)分析水電站工程建設(shè)前后洪水線(xiàn)型的變化情況。瀾滄江干流第1個(gè)梯級(jí)水電站漫灣水電站于1986年開(kāi)工，1992年一期工程結(jié)束；第2個(gè)梯級(jí)水電站大朝山水電站于2003年完工投產(chǎn)，其他電站于2008年以后陸續(xù)投產(chǎn)運(yùn)行。因此選擇水庫(kù)修建前的1970、1985年以及1994、2005和2019年作為時(shí)變矩模型的基準(zhǔn)時(shí)間點(diǎn)。

圖4為清盛水文站不同基準(zhǔn)時(shí)間點(diǎn)的洪水?dāng)M合線(xiàn)型變化情況。在尚未進(jìn)行梯級(jí)水電開(kāi)發(fā)時(shí)期即1960-1970年， GLO-AL線(xiàn)型高水尾端位于所有線(xiàn)型的最上方，同一量級(jí)洪水出現(xiàn)的概率偏大。1971-1985年GLO-AL曲線(xiàn)下移，但下移的幅度不大。在梯級(jí)水電站建成影響期，1986-1994年(漫灣水電站運(yùn)行)GLO-AL曲線(xiàn)有所下降；1995-2005年(漫灣和大朝山水電站聯(lián)合運(yùn)行)GLO-AL曲線(xiàn)相較于1994年又有所下降；2006-2019年(11個(gè)水電站陸續(xù)運(yùn)行)GLO-AL曲線(xiàn)位于所有時(shí)間基點(diǎn)線(xiàn)型的下方，同量級(jí)洪水發(fā)生概率最低。瀾滄江中下游以小灣、糯扎渡水電站為核心的梯級(jí)水電開(kāi)發(fā)對(duì)洪水線(xiàn)型的影響較小，但水庫(kù)聯(lián)合運(yùn)行以后清盛站同量級(jí)洪水發(fā)生的概率降低，設(shè)計(jì)洪水量級(jí)也明顯減小。這說(shuō)明瀾滄江水電站的建設(shè)對(duì)清盛水文站洪水具有明顯的調(diào)控作用。

圖4 清盛水文站不同基準(zhǔn)時(shí)間點(diǎn)的洪水?dāng)M合線(xiàn)型變化

3.3.2 設(shè)計(jì)洪峰流量與重現(xiàn)期的變化過(guò)程 用傳統(tǒng)頻率分析方法計(jì)算得到清盛水文站100年一遇洪水的設(shè)計(jì)值為22 881 m3/s，用時(shí)變矩模型分析該設(shè)計(jì)值下重現(xiàn)期的變化，如圖5所示。由圖5可見(jiàn)，與傳統(tǒng)頻率分析下某一指定流量標(biāo)準(zhǔn)值重現(xiàn)期不變相比，采用時(shí)變矩模型得到的重現(xiàn)期隨時(shí)間呈上升趨勢(shì)。如1960-1965年為小于100年一遇，1966-2011年為大于100年一遇且小于200年一遇，2012年后為大于200年一遇。選取指定重現(xiàn)期T=100 a，采用時(shí)變矩模型分析指定重現(xiàn)期下設(shè)計(jì)流量的變化，如圖5所示。

圖5 清盛水文站設(shè)計(jì)流量下重現(xiàn)期及指定重現(xiàn)期下設(shè)計(jì)流量變化

由圖5可見(jiàn)，與傳統(tǒng)頻率分析下指定重現(xiàn)期標(biāo)準(zhǔn)下的設(shè)計(jì)流量值不變相比，采用時(shí)變矩方法計(jì)算得到的設(shè)計(jì)流量隨時(shí)間呈減小趨勢(shì)。清盛水文站100年一遇洪水量級(jí)從1960年的23 274 m3/s減小至2019年的19 085 m3/s。如果不考慮序列非一致性處理，得到的100年一遇設(shè)計(jì)值為22 881 m3/s，該值比時(shí)變矩模型求得的2014-2019年糯扎渡水電站運(yùn)行后的設(shè)計(jì)值高出18%～20%，即傳統(tǒng)頻率分析方法會(huì)高估洪水量級(jí)。從時(shí)變矩方法計(jì)算得到的設(shè)計(jì)流量來(lái)看，瀾滄江-湄公河上游梯級(jí)水庫(kù)尤其是糯扎渡水庫(kù)運(yùn)行后明顯減小了清盛水文站100年一遇的洪峰流量，對(duì)下游防洪具有積極作用。

4 討 論

氣候和人類(lèi)活動(dòng)是河川徑流變化的兩大主要驅(qū)動(dòng)因素?；谡军c(diǎn)觀測(cè)數(shù)據(jù)分析表明，近50年來(lái)瀾滄江-湄公河上游地區(qū)降水量減少趨勢(shì)不明顯，氣溫升高和潛在蒸散發(fā)增加是徑流量顯著減少的重要因素[21-22]。相對(duì)于氣候變化，人類(lèi)活動(dòng)對(duì)徑流的影響在特定時(shí)期內(nèi)則更為明顯。研究表明，瀾滄江允景洪水文站1987-2014年多年年均徑流量較1980-1986年減少了6%，人類(lèi)活動(dòng)對(duì)徑流量變化的貢獻(xiàn)率從1987-2007年的43%增加到2008-2014年的95%，人類(lèi)活動(dòng)對(duì)徑流量變化的影響逐漸加劇[23]。20世紀(jì)80年代以來(lái)，瀾滄江梯級(jí)水庫(kù)建設(shè)導(dǎo)致徑流年內(nèi)分配呈現(xiàn)均勻化，下游清盛站年最大日流量呈顯著的減少趨勢(shì)[24-25]。總的來(lái)說(shuō)，變化環(huán)境尤其是梯級(jí)水庫(kù)的建設(shè)對(duì)允景洪-萬(wàn)象區(qū)間徑流年內(nèi)分配的影響顯著[26]。

本研究結(jié)果表明，1960-2019年來(lái)清盛水文站年最大日流量序列呈現(xiàn)出顯著減小趨勢(shì)，水文極值序列突變年份(2008年)與景洪水電站蓄水、小灣水電站投產(chǎn)時(shí)間吻合，說(shuō)明變化環(huán)境下水文極值序列發(fā)生了變異。未來(lái)湄公河上游地區(qū)極端降水量將趨于增多且降水集中程度將增大，流域洪水風(fēng)險(xiǎn)勢(shì)必增加[27]。流域內(nèi)水庫(kù)聯(lián)合調(diào)度可以在一定程度上緩解洪水的風(fēng)險(xiǎn)[28-29]，如全流域水庫(kù)聯(lián)合調(diào)度可以將湄公河干流洪水從200年一遇減至20～50年一遇[30]。在未來(lái)氣候變化和人類(lèi)活動(dòng)加劇的背景下，使用傳統(tǒng)頻率分析方法得到的重現(xiàn)期和洪水設(shè)計(jì)值將會(huì)“失真”。本研究采用時(shí)變矩方法分析指定設(shè)計(jì)流量下重現(xiàn)期和指定重現(xiàn)期下設(shè)計(jì)流量，發(fā)現(xiàn)重現(xiàn)期和設(shè)計(jì)流量均隨時(shí)間變化。在非一致背景下，如果采用傳統(tǒng)水文頻率分析方法而不對(duì)序列的非一致性進(jìn)行考慮，將會(huì)錯(cuò)估洪水重現(xiàn)期及設(shè)計(jì)洪水量級(jí)，因此傳統(tǒng)洪水頻率計(jì)算方法的使用應(yīng)該慎重。

5 結(jié) 論

(1)1960-2019年清盛水文站年最大日流量序列呈顯著減小的趨勢(shì)并在2008年發(fā)生突變，水文極值序列呈現(xiàn)非一致性。

(2)采用時(shí)變矩方法選擇6種概率分布和7種趨勢(shì)模型進(jìn)行組合，得到42種競(jìng)爭(zhēng)模型進(jìn)行比較，結(jié)果表明清盛水文站廣義邏輯斯諦分布搭配均值具有線(xiàn)性趨勢(shì)的模型擬合最優(yōu)。

(3)基于傳統(tǒng)頻率分析方法求得的100年一遇洪水設(shè)計(jì)值，在時(shí)變矩模型下重現(xiàn)期由水電開(kāi)發(fā)前的小于100年一遇變化到2012年后的大于200年一遇；100年一遇洪水設(shè)計(jì)值也隨時(shí)間呈減小趨勢(shì)。

(4)湄公河上游梯級(jí)水庫(kù)尤其是糯扎渡水庫(kù)運(yùn)行后明顯減小了清盛站100年一遇的洪峰流量，對(duì)下游防洪具有積極作用。