高佳佳， 杜 軍

（1.中國氣象局成都高原氣象研究所拉薩分部，西藏 拉薩 850000； 2.西藏高原大氣環(huán)境科學(xué)研究所，西藏 拉薩 850000；3.西藏高原大氣環(huán)境研究重點實驗室，西藏 拉薩 850000）

0 引言

極端降水是全球最受關(guān)注、影響最大的自然災(zāi)害之一，是短期氣候預(yù)測研究的重點［1］。IPCC［2］曾指出，隨著全球氣候變暖，內(nèi)陸地區(qū)的極端降水事件頻率呈現(xiàn)出增加趨勢。近年來，區(qū)域性洪水、干旱、高溫、雨雪冰凍等極端事件頻發(fā)，尤其是20世紀(jì)80 年代以后，頻繁的極端事件給生態(tài)環(huán)境、經(jīng)濟(jì)發(fā)展和人民生活造成了嚴(yán)重影響［3-4］。因此，研究區(qū)域性極端降水事件對科學(xué)認(rèn)識氣候變化背景下水循環(huán)的時空演變，把握氣候異常對極端降水的影響規(guī)律，評估水資源管理及區(qū)域水資源安全具有十分重要的理論和實踐意義。

極端氣候事件歸根到底是氣候極值問題，氣候極值是極端事件產(chǎn)生的必要條件，極端事件發(fā)生發(fā)展的預(yù)測首先要考慮極值的分布規(guī)律。國內(nèi)外部分學(xué)者研究表明，中國西部的極端降水天數(shù)呈增加趨勢，尤其是西北地區(qū)［5-8］。Fischer 等［9］分析了珠江的極端降水分布特征，并估算了極端降水指數(shù)。Hong 等［10］認(rèn)為海河流域的極端降水主要發(fā)生在38°N，大部分站點的降水呈現(xiàn)出減少趨勢。程炳巖等［11］、江志紅等［12］研究認(rèn)為廣義帕累托分布函數(shù)（generalized Pareto distribution，GPD）在重慶、中國東部的日降水模擬中具有更高的擬合度；李占玲等［13］基于GPD 函數(shù)分析了黑河流域的極端降水頻率特征，得出該流域20 世紀(jì)60 年代發(fā)生的極端降水次數(shù)最多，90 年代以后次數(shù)較少。Eylon 等［14］運用極值理論分析了巴拿馬運河的極端降水分布特征，并估算了極端降水的重現(xiàn)期及相應(yīng)的置信區(qū)間。劉彩紅等［15］運用CMIP5 模式指出，青海高原的降水有極端化的趨勢，極端降水頻次增加，強(qiáng)度增大。韓國軍［16］、游慶龍等［17］運用統(tǒng)計方法計算出青藏高原極端降水大部分呈增加趨勢，且逐年平均降水強(qiáng)度和逐年連續(xù)降水天數(shù)均有所增加，90 年代以來增加明顯。

雅魯藏布江發(fā)源并流經(jīng)西藏高原，地理位置特殊，是世界上海拔最高的大河之一，平均海拔4 000 m 以上，是全球氣候主要變化區(qū)與敏感區(qū)。流域沿岸為西藏主要農(nóng)、牧業(yè)生產(chǎn)區(qū)，其洪澇和干旱的頻繁發(fā)生導(dǎo)致了水資源分布不均，進(jìn)而影響了流域的用水矛盾和生態(tài)環(huán)境的惡化，而極端事件的發(fā)生是對區(qū)域氣候、環(huán)境變化的重要響應(yīng)。目前對氣候極值進(jìn)行定量評估的方法以氣候動力模式為主，從概率論角度對極端氣候事件及可預(yù)測性研究并不多，尤其是預(yù)測方法。IPCC報告中特別強(qiáng)調(diào)統(tǒng)計方法對極端氣候事件的重要性［2］，因此，本文基于廣義帕累托方法（GPD）的分布參數(shù)模型，針對超出閾值的數(shù)據(jù)作為樣本數(shù)據(jù)來建模，從氣候極值的分布規(guī)律出發(fā)，揭示極端降水的發(fā)生發(fā)展規(guī)律，探索極端降水的可預(yù)測性，從而更好地預(yù)估極端事件，為提高防災(zāi)減災(zāi)能力提供科學(xué)依據(jù)。

1 研究區(qū)與數(shù)據(jù)

1.1 研究區(qū)域

雅魯藏布江（簡稱雅江）全部在中國境內(nèi)，橫貫西藏高原南部，干流全長約2.1×103km，流域面積2.4×105km2，雅魯藏布江干流河谷沿東西向的斷裂帶發(fā)育，流域呈東西向的狹長帶，支流多而短小，較大支流有拉薩河、雅魯藏布江帕隆藏布、易貢藏布、拉喀藏布、尼澤曲、年楚河等。干流在仲巴縣里孜以上為上游，是高寒河谷地帶。里孜到米林縣派區(qū)為中游，支流眾多，流量增大，河谷展寬，氣候溫和，水利條件較好，是西藏農(nóng)業(yè)最發(fā)達(dá)的地區(qū)。下游位于林芝一帶。截止國境線，年徑流總量為1.1×1010m3，洪水由強(qiáng)降水形成，持續(xù)時間較長。流域東部地區(qū)降水量充足，年平均降水量超過900 mm，達(dá)到半濕潤地區(qū)水平，西部地區(qū)降水量少，年均降水量不足100 mm，為干旱地區(qū)，整個流域的降水量從東至西呈現(xiàn)遞減的趨勢。

圖1 研究區(qū)域位置示意Fig.1 Location of study area

1.2 數(shù)據(jù)來源

選取流域內(nèi)9 個氣象站近50 年（1959—2017年）5—9 月逐日降水資料作為研究對象。部分站點開始于1978年和1979年。

2 研究方法

2.1 GPD分布函數(shù)和密度函數(shù)

GPD 可以直接利用歷年的原始數(shù)據(jù)，人為設(shè)置閾值，在設(shè)置好閾值后，以此為標(biāo)準(zhǔn)來抽取每一年超過此閾值的極大或極小值，即“超門限峰值POT”（peaks over threshold），可以提高估算精度［18］。具體為：

GPD的分布函數(shù)為：

其相應(yīng)的密度函數(shù)為：

式中：ξ為門限值；a為尺度參數(shù)；k為形狀參數(shù)，自變量x的取值，取決于k的值，當(dāng)k< 0時，ξ 0 時，GPD 為指數(shù)分布；當(dāng)k=1 時，GPD 在區(qū)間[ξ，ξ+a]上均勻分布。

表1 雅魯藏布江流域內(nèi)各站點資料長度Table 1 Basic information of each stations in Yarlung Zangbo River basin

2.1.1 閾值的選取

GPD 模型的核心在于確定閾值，它是正確估計參數(shù)的前提。如果閾值選取的過高，會使得超額數(shù)據(jù)量太少，導(dǎo)致估計出來的參數(shù)方差很大；如果閾值選取的過低，則不能保證函數(shù)的收斂性，所估計參數(shù)有較大偏差。本文主要使用Hill 圖估計、百分位法和年交叉率法來確定各站點日降水量閾值。

Hill 圖法是基于Hill 估計量的一種閾值圖形法［19］，由點（k，1/H（k，n））構(gòu)成的曲線，通過觀察圖中尾部指數(shù)穩(wěn)定的區(qū)域來選擇閾值。其定義為：

百分位法是將該站點的日降水量從小到大排序，并計算相應(yīng)的累計百分位，某一百分位所對應(yīng)數(shù)據(jù)的值就稱為這一百分位的百分位數(shù)，文中分別計算了90、95、97、98、99百分位。

年交叉率法是假定超閾值降水極值出現(xiàn)次數(shù)服從泊松分布，以一年為時段所量度，極值超過閾值的次數(shù)，λ=n/T，其中，λ為年交叉率，T為資料總年數(shù)。

2.1.2 參數(shù)估計

對GPD 分布進(jìn)行參數(shù)估計的方法有很多，本文主要通過極大似然估計法對參數(shù)進(jìn)行估計。該方法具有很強(qiáng)的靈活性，可以適應(yīng)不同模型的需求，且統(tǒng)計特性良好，能夠綜合各種有關(guān)信息到統(tǒng)計推斷中去。模型估計時，樣本變異可能會導(dǎo)致模型參數(shù)的不確定性，由于極大似然方法具有漸進(jìn)正態(tài)性，容易給出估計值及其標(biāo)準(zhǔn)誤差（標(biāo)準(zhǔn)誤差是參數(shù)不確定性或變異性的度量之一）。該方法唯一的缺點是計算時迭代繁瑣。具體方法見文獻(xiàn)［20］。

2.1.3 重現(xiàn)期計算

某一指定重現(xiàn)期T時間的降水量分位計算公式：

2.2 獨立性和穩(wěn)定性檢驗

GPD 擬合需要超閾值數(shù)據(jù)序列滿足平穩(wěn)性的條件，因此，擬合之前需對超閾值序列進(jìn)行穩(wěn)定性檢驗。文章中使用Mann-Kendall（M-K 檢驗）對序列的變化趨勢和突變點進(jìn)行檢驗。M-K 檢驗［21］是氣象學(xué)、氣候?qū)W中經(jīng)常用來進(jìn)行突變檢驗的一種非參數(shù)檢驗方法，它不要求樣本符合一定的分布。即給定顯著性水平α=0.05，則統(tǒng)計量的臨界值為±1.96。統(tǒng)計量大于0，表示序列呈上升趨勢；反之，表明呈下降趨勢，大于或小于±1.96，表示上升或下降趨勢明顯。該方法能有效區(qū)分某一自然過程是處于自然波動還是存在確定的變化趨勢，常用于氣候變化影響下的降水、干旱頻次趨勢檢測。

2.3 擬合度檢驗

通過對超閾值序列進(jìn)行平穩(wěn)性檢驗后，即可對GPD 進(jìn)行擬合，擬合結(jié)果需通過Kolmogorov-Simimov（K-S 檢驗）。一般在K-S 檢驗中，先計算需要做比較的兩組觀察數(shù)據(jù)的累積分布函數(shù)，然后求這兩個累積分布函數(shù)的差的絕對值中的最大值D。最后通過查表以確定D值是否落在所要求對應(yīng)的置信區(qū)間內(nèi)。若D值落在了對應(yīng)的置信區(qū)間內(nèi)，說明被檢測的數(shù)據(jù)滿足要求。

3 結(jié)果與分析

3.1 閾值的選取

通過Hill 圖尾部特征的穩(wěn)定性來選取閾值（圖2），由于Hill圖具有較強(qiáng)的主觀性，從而會導(dǎo)致選取的閾值不同。以中游地區(qū)的拉孜、拉薩和墨竹工卡站為例，拉孜站點的Hill 圖尾部指數(shù)趨于穩(wěn)定大致位于65～70 附近，對應(yīng)的降水量為7.9～15.4 mm；拉薩和墨竹工卡站的尾部指數(shù)均位于160左右時趨于穩(wěn)定，對應(yīng)的降水量分別為8.7～18.5 mm 和9.9～20.3 mm。同樣的，下游地區(qū)的林芝站在80～120位置時趨于穩(wěn)定，對應(yīng)的降水量為14.4～17.6 mm。為更好的確定閾值，在此基礎(chǔ)上，我們結(jié)合百分位法和年交叉率法進(jìn)行閾值確定（表2）。

圖2 各站點Hill圖指數(shù)分布Fig.2 The distribution of index of Hill plot for each stations

表2 各站點的百分位閾值（mm）和年交叉率Table 2 Precipitation threshold selection in GPD and the average annual occurrence number for each stations

與Hill 圖相比，百分位法可以更精確的確定閾值。以拉孜站為例，Hill 圖顯示的閾值為7.9～15.4 mm，而對應(yīng)的是該站點的95 百分位；同樣拉薩站、墨竹工卡站的Hill 圖確定的閾值均位于93 百分位和95 百分位。總體而言，Hill 圖確定的閾值要小于百分位確定的閾值。根據(jù)前人研究結(jié)果得出，年交叉率為1～2 時對應(yīng)的閾值可作為GPD 擬合分析時的參考閾值。結(jié)合表1 得出，當(dāng)各站點日降水量達(dá)到99 百分位時，年交叉率均穩(wěn)定在1.5 附近，因此我們確定99 百分位時的閾值為最佳閾值。

3.2 超閾值序列的平穩(wěn)性檢驗

根據(jù)閾值的平穩(wěn)性和穩(wěn)定性要求，我們用M-K法在顯著性水平α=0.05 條件下檢測超閾值降水序列的變化趨勢和突變。結(jié)果顯示，下游地區(qū)的林芝和米林站的統(tǒng)計量略大于1.96，分別為2.3 和2.8，其余大部分站點均位于臨界區(qū)域內(nèi)，通過了顯著性水平檢驗。通過時間序列曲線（UFk曲線）可以看出（圖3），米林、墨竹工卡、南木林站呈顯著增長趨勢，說明三個站點的超閾值序列的日降水量呈逐漸增加趨勢；其余站點呈下降趨勢（此部分只給出拉孜、墨竹工卡、日喀則、米林站的趨勢變化圖，其余圖表省略）。各站點超閾值序列沒有明顯的突變性。

圖3 流域內(nèi)四個站點的時間序列統(tǒng)計量變化圖Fig.3 Time series statistics of four sites in the basin

3.3 模型擬合效果檢驗

通過對流域內(nèi)各站點進(jìn)行GPD 擬合（圖4），并使用K-S檢驗看其是否符合已知理論分布函數(shù)。大部分站點的統(tǒng)計值小于0.01顯著性水平，接受原假設(shè)。說明雅江流域各站點之間雖存在差異，但由GPD 擬合曲線可知，理論頻數(shù)和實測頻數(shù)基本相符。另外，雖然西藏地區(qū)大部分臺站的觀測記錄起始年代不一致，如拉薩站、拉孜站，資料開始時間分別為1955 年、1977 年，但從圖4 中可看出，數(shù)據(jù)滿足方程需求，擬合結(jié)果表明資料長度并不影響降水極值的統(tǒng)計推斷，且資料年限越長擬合結(jié)果越好。

圖4 各站點累積頻率和實際頻率分布曲線對比Fig.4 The distribution of cumulative frequency and empirical frequency over the observation stations

3.4 雅江流域極端降水的分布特征

從各站點的極端降水閾值分布來看，林芝地區(qū)的閾值最大，閾值最小的站點為拉孜站。這與我們的觀測事實一致，林芝地區(qū)日降水量大，連續(xù)降雨日數(shù)長，得到的閾值就大。

尺度參數(shù)主要是描述極值分布的變率，尺度參數(shù)越大，極值波動范圍越大，表明打破極端降水的記錄值也越大。整體而言，雅江流域的尺度參數(shù)由下游向上游是逐漸減小的，平均值為5.95。由表3可知，下游地區(qū)的尺度參數(shù)最大，約為7.00，表明這一區(qū)域的極端降水變化幅度很大。從氣候背景來看，該地區(qū)位于高原季風(fēng)區(qū)，受印緬槽和西風(fēng)帶影響，季節(jié)性降水較大［22-23］，5—9 月的降水總量可達(dá)600 mm，是西藏地區(qū)夏季降水量最大的區(qū)域，因此可能出現(xiàn)的破極端降水記錄值要高于其他地區(qū)。尺度參數(shù)最小的區(qū)域位于流域中上游地區(qū)，受地理位置和大氣環(huán)流影響，雨期短，降水量少，且連續(xù)降水日數(shù)也少，降水極值的范圍比較小，區(qū)域打破極端降水的記錄值要比下游地區(qū)低。這與前人研究結(jié)果一致［24］，我國的干旱地區(qū)大部分位于非季風(fēng)區(qū)，降水極值范圍較小，破紀(jì)錄的可能性較季風(fēng)區(qū)小。

表3 流域內(nèi)各站點的GPD模型參數(shù)估計及檢驗Table 3 Estimation and validation of parameters in GPD model

形狀參數(shù)作為模型的第二個重要參數(shù)，不同的形狀有不同的尾部分布特征，它表示該區(qū)域極端降水的破紀(jì)錄率。由表3 看出，形狀參數(shù)正值區(qū)主要位于拉孜地區(qū)，這些地區(qū)發(fā)生破紀(jì)錄降水事件的可能性比其他地區(qū)大。正是由于該區(qū)域降水日數(shù)少，所以一旦有降水過程，就可能會打破降水極值。而下游地區(qū)的形狀參數(shù)為負(fù)值，說明這些區(qū)域的降水發(fā)生破紀(jì)錄的概率偏小。因為夏季，降雨越頻繁的區(qū)域，極值變率大，較均值離散程度大，則破紀(jì)錄的概率較小。有降水的地區(qū)，不是有較大的形狀參數(shù)就是有較大的尺度參數(shù)，不可能兩個參數(shù)都大。形狀參數(shù)較大的地區(qū)，稱為“形狀參數(shù)主導(dǎo)區(qū)”，該區(qū)域多持續(xù)降水，極端氣候事件較少；尺度參數(shù)較大的區(qū)域稱為“尺度參數(shù)主導(dǎo)區(qū)”，該區(qū)域降水較少，且多變，極端氣候事件較多［23］。

3.5 雅江流域極值重現(xiàn)期預(yù)估

極值模型建立最重要的目的之一就是預(yù)測極端事件的重現(xiàn)期或重現(xiàn)水平。如表4 所示，根據(jù)重現(xiàn)期公式可以得出，從5 年一遇和10 年一遇的極端降水值來看，雅江流域除拉孜站外，其他地區(qū)降雨極值均超過30 mm，日喀則地區(qū)的降水極值達(dá)50 mm，其中拉薩、澤當(dāng)、墨竹工卡和南木林站5 年一遇和10年一遇的極端降水量分別在40 mm 左右。相關(guān)研究指出［26-27］，當(dāng)溫度上升為2 ℃時，青藏高原的強(qiáng)降水距平百分率平均增多44.5%～59.5%，大值區(qū)出現(xiàn)在山南附近，這與我們的研究結(jié)論相符合。各站點在15年一遇的極端降水值之后，極值水平的增長變得非常緩慢，其中林芝地區(qū)的增長最緩慢，以0.32 mm·a-1的速率增長；日喀則地區(qū)的降雨極值增長率最快，約0.72 mm·a-1。西藏地區(qū)極端降水頻率一般為每年4.3 次，強(qiáng)度在20 mm·d-1以上，林芝地區(qū)為極端降雨的高值區(qū)，且沿雅江一線極端降水的頻次呈增加趨勢［28-29］。

表4 雅江流域各站點日降水量極大值重現(xiàn)水平（單位：mm）Table 4 The maximum daily precipitation of flood season reappeared in Yarlung Zangbo River（unit：mm）

為檢驗各站點不同重現(xiàn)期水平的合理性，將其帶入各站點逐日降水序列中進(jìn)行驗證。以拉薩和拉孜站為例，拉薩站5 年重現(xiàn)期水平值為39 mm，在1967—2017 年所有逐日降水中，共有10 次大于39 mm的降水過程，平均每5年一次。拉孜站5年重現(xiàn)期水平值為35.1 mm，在1977—2017 年期間，共有9 次大于35.1 mm 的降水過程，平均每4.4 年一次。由GPD 擬合計算出的極端降水重現(xiàn)期水平基本符合實際。值得注意的是這里的“重現(xiàn)期”并不意味著經(jīng)過T年之后一定會出現(xiàn)的“周期”，它是概率意義上的“統(tǒng)計周期”。

4 結(jié)論

雅魯藏布江作為高原河流，由于強(qiáng)降水的時空分布不均而引起洪澇和干旱，并對流域內(nèi)的農(nóng)牧業(yè)產(chǎn)生重要威脅。因此在氣候日益增暖趨勢下，評估極端降水規(guī)律及其發(fā)生概率十分必要。本文通過引進(jìn)GPD 概率分布模型，對西藏地區(qū)汛期強(qiáng)降水規(guī)律進(jìn)行模擬。結(jié)果表明：

（1）通過Hill 圖法選取的流域內(nèi)各站點的閾值序列小于百分位法選取的閾值序列，綜合考慮Hill圖法、百分位法及年交叉率法最終確定99百分位時的閾值為最佳閾值。

（2）各站點閾值序列在M-K 顯著性水平檢驗下，無明顯突變。擬合效果通過K-S檢驗，各站點擬合的理論頻數(shù)和實測頻數(shù)基本相符，且資料長度并不影響降水極值的統(tǒng)計推斷。

（3）通過分析流域內(nèi)各站點擬合的極端降水特征可知，尺度參數(shù)的大值區(qū)位于流域下游，即林芝、米林地區(qū)，表明該地區(qū)的極值波動大；相反地，小值區(qū)位于流域中上游的拉孜站附近，表明極值波動小。形狀參數(shù)正值區(qū)位于流域中上游地區(qū)，說明發(fā)生破紀(jì)錄的降水事件概率較大，擬合結(jié)果與實際觀測一致。

（4）從5 年一遇和10 年一遇的極值水平看，雅江流域除拉孜站外，其他地區(qū)降雨極值均超過30 mm，日喀則地區(qū)的降水極值達(dá)50 mm；各站點在15 年一遇的極端降水值之后，極值水平的增長變得非常緩慢，由GPD 擬合計算出的降水極值具有一定的合理性。