雅魯藏布江流域極端降水模擬及預估
2021-05-24 02:31:04高佳佳
冰川凍土 2021年2期
高佳佳， 杜 軍
（1.中國氣象局成都高原氣象研究所拉薩分部，西藏 拉薩 850000； 2.西藏高原大氣環境科學研究所，西藏 拉薩 850000；3.西藏高原大氣環境研究重點實驗室，西藏 拉薩 850000）
0 引言
極端降水是全球最受關注、影響最大的自然災害之一，是短期氣候預測研究的重點［1］。IPCC［2］曾指出，隨著全球氣候變暖，內陸地區的極端降水事件頻率呈現出增加趨勢。近年來，區域性洪水、干旱、高溫、雨雪冰凍等極端事件頻發，尤其是20世紀80 年代以后，頻繁的極端事件給生態環境、經濟發展和人民生活造成了嚴重影響［3-4］。因此，研究區域性極端降水事件對科學認識氣候變化背景下水循環的時空演變，把握氣候異常對極端降水的影響規律，評估水資源管理及區域水資源安全具有十分重要的理論和實踐意義。
極端氣候事件歸根到底是氣候極值問題，氣候極值是極端事件產生的必要條件，極端事件發生發展的預測首先要考慮極值的分布規律。國內外部分學者研究表明，中國西部的極端降水天數呈增加趨勢，尤其是西北地區［5-8］。Fischer 等［9］分析了珠江的極端降水分布特征，并估算了極端降水指數。Hong 等［10］認為海河流域的極端降水主要發生在38°N，大部分站點的降水呈現出減少趨勢。程炳巖等［11］、江志紅等［12］研究認為廣義帕累托分布函數（generalized Pareto distribution，GPD）在重慶、中國東部的日降水模擬中具有更高的擬合度；李占玲等［13］基于GPD 函數分析了黑河流域的極端降水頻率特征，得出該流域20 世紀60 年代發生的極端降水次數最多，90 年代以后次數較少。Eylon 等［14］運用極值理論分析了巴拿馬運河的極端降水分布特征，并估算了極端降水的重現期及相應的置信區間。劉彩紅等［15］運用CMIP5 模式指出，青海高原的降水有極端化的趨勢，極端降水頻次增加，強度增大。韓國軍［16］、游慶龍等［17］運用統計方法計算出青藏高原極端降水大部分呈增加趨勢，且逐年平均降水強度和逐年連續降水天數均有所增加，90 年代以來增加明顯。
雅魯藏布江發源并流經西藏高原，地理位置特殊，是世界上海拔最高的大河之一，平均海拔4 000 m 以上，是全球氣候主要變化區與敏感區。流域沿岸為西藏主要農、牧業生產區，其洪澇和干旱的頻繁發生導致了水資源分布不均，進而影響了流域的用水矛盾和生態環境的惡化，而極端事件的發生是對區域氣候、環境變化的重要響應。目前對氣候極值進行定量評估的方法以氣候動力模式為主，從概率論角度對極端氣候事件及可預測性研究并不多，尤其是預測方法。IPCC報告中特別強調統計方法對極端氣候事件的重要性［2］，因此，本文基于廣義帕累托方法（GPD）的分布參數模型，針對超出閾值的數據作為樣本數據來建模，從氣候極值的分布規律出發，揭示極端降水的發生發展規律，探索極端降水的可預測性，從而更好地預估極端事件，為提高防災減災能力提供科學依據。
1 研究區與數據
1.1 研究區域
雅魯藏布江（簡稱雅江）全部在中國境內，橫貫西藏高原南部，干流全長約2.1×103km，流域面積2.4×105km2，雅魯藏布江干流河谷沿東西向的斷裂帶發育，流域呈東西向的狹長帶，支流多而短小，較大支流有拉薩河、雅魯藏布江帕隆藏布、易貢藏布、拉喀藏布、尼澤曲、年楚河等。干流在仲巴縣里孜以上為上游，是高寒河谷地帶。里孜到米林縣派區為中游，支流眾多，流量增大，河谷展寬，氣候溫和，水利條件較好，是西藏農業最發達的地區。下游位于林芝一帶。截止國境線，年徑流總量為1.1×1010m3，洪水由強降水形成，持續時間較長。流域東部地區降水量充足，年平均降水量超過900 mm，達到半濕潤地區水平，西部地區降水量少，年均降水量不足100 mm，為干旱地區，整個流域的降水量從東至西呈現遞減的趨勢。
圖1 研究區域位置示意Fig.1 Location of study area
1.2 數據來源
選取流域內9 個氣象站近50 年（1959—2017年）5—9 月逐日降水資料作為研究對象。部分站點開始于1978年和1979年。
2 研究方法
2.1 GPD分布函數和密度函數
GPD 可以直接利用歷年的原始數據，人為設置閾值，在設置好閾值后，以此為標準來抽取每一年超過此閾值的極大或極小值，即“超門限峰值POT”（peaks over threshold），可以提高估算精度［18］。具體為：
GPD的分布函數為：
其相應的密度函數為：
式中：ξ為門限值；a為尺度參數；k為形狀參數，自變量x的取值，取決于k的值，當k< 0時，ξ
表1 雅魯藏布江流域內各站點資料長度Table 1 Basic information of each stations in Yarlung Zangbo River basin
2.1.1 閾值的選取
GPD 模型的核心在于確定閾值，它是正確估計參數的前提。如果閾值選取的過高，會使得超額數據量太少，導致估計出來的參數方差很大；如果閾值選取的過低，則不能保證函數的收斂性，所估計參數有較大偏差。本文主要使用Hill 圖估計、百分位法和年交叉率法來確定各站點日降水量閾值。
Hill 圖法是基于Hill 估計量的一種閾值圖形法［19］，由點（k，1/H（k，n））構成的曲線，通過觀察圖中尾部指數穩定的區域來選擇閾值。其定義為：
百分位法是將該站點的日降水量從小到大排序，并計算相應的累計百分位，某一百分位所對應數據的值就稱為這一百分位的百分位數，文中分別計算了90、95、97、98、99百分位。
年交叉率法是假定超閾值降水極值出現次數服從泊松分布，以一年為時段所量度，極值超過閾值的次數，λ=n/T，其中，λ為年交叉率，T為資料總年數。
2.1.2 參數估計
對GPD 分布進行參數估計的方法有很多，本文主要通過極大似然估計法對參數進行估計。該方法具有很強的靈活性，可以適應不同模型的需求，且統計特性良好，能夠綜合各種有關信息到統計推斷中去。模型估計時，樣本變異可能會導致模型參數的不確定性，由于極大似然方法具有漸進正態性，容易給出估計值及其標準誤差（標準誤差是參數不確定性或變異性的度量之一）。該方法唯一的缺點是計算時迭代繁瑣。具體方法見文獻［20］。
2.1.3 重現期計算
某一指定重現期T時間的降水量分位計算公式：
2.2 獨立性和穩定性檢驗
GPD 擬合需要超閾值數據序列滿足平穩性的條件，因此，擬合之前需對超閾值序列進行穩定性檢驗。文章中使用Mann-Kendall（M-K 檢驗）對序列的變化趨勢和突變點進行檢驗。M-K 檢驗［21］是氣象學、氣候學中經常用來進行突變檢驗的一種非參數檢驗方法，它不要求樣本符合一定的分布。即給定顯著性水平α=0.05，則統計量的臨界值為±1.96。統計量大于0，表示序列呈上升趨勢；反之，表明呈下降趨勢，大于或小于±1.96，表示上升或下降趨勢明顯。該方法能有效區分某一自然過程是處于自然波動還是存在確定的變化趨勢，常用于氣候變化影響下的降水、干旱頻次趨勢檢測。
2.3 擬合度檢驗
通過對超閾值序列進行平穩性檢驗后，即可對GPD 進行擬合，擬合結果需通過Kolmogorov-Simimov（K-S 檢驗）。一般在K-S 檢驗中，先計算需要做比較的兩組觀察數據的累積分布函數，然后求這兩個累積分布函數的差的絕對值中的最大值D。最后通過查表以確定D值是否落在所要求對應的置信區間內。若D值落在了對應的置信區間內，說明被檢測的數據滿足要求。
3 結果與分析
3.1 閾值的選取
通過Hill 圖尾部特征的穩定性來選取閾值（圖2），由于Hill圖具有較強的主觀性，從而會導致選取的閾值不同。以中游地區的拉孜、拉薩和墨竹工卡站為例，拉孜站點的Hill 圖尾部指數趨于穩定大致位于65～70 附近，對應的降水量為7.9～15.4 mm；拉薩和墨竹工卡站的尾部指數均位于160左右時趨于穩定，對應的降水量分別為8.7～18.5 mm 和9.9～20.3 mm。同樣的，下游地區的林芝站在80～120位置時趨于穩定，對應的降水量為14.4～17.6 mm。為更好的確定閾值，在此基礎上，我們結合百分位法和年交叉率法進行閾值確定（表2）。
圖2 各站點Hill圖指數分布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，而對應的是該站點的95 百分位；同樣拉薩站、墨竹工卡站的Hill 圖確定的閾值均位于93 百分位和95 百分位。總體而言，Hill 圖確定的閾值要小于百分位確定的閾值。根據前人研究結果得出，年交叉率為1～2 時對應的閾值可作為GPD 擬合分析時的參考閾值。結合表1 得出，當各站點日降水量達到99 百分位時，年交叉率均穩定在1.5 附近，因此我們確定99 百分位時的閾值為最佳閾值。
3.2 超閾值序列的平穩性檢驗
根據閾值的平穩性和穩定性要求，我們用M-K法在顯著性水平α=0.05 條件下檢測超閾值降水序列的變化趨勢和突變。結果顯示，下游地區的林芝和米林站的統計量略大于1.96，分別為2.3 和2.8，其余大部分站點均位于臨界區域內，通過了顯著性水平檢驗。通過時間序列曲線（UFk曲線）可以看出（圖3），米林、墨竹工卡、南木林站呈顯著增長趨勢，說明三個站點的超閾值序列的日降水量呈逐漸增加趨勢；其余站點呈下降趨勢（此部分只給出拉孜、墨竹工卡、日喀則、米林站的趨勢變化圖，其余圖表省略）。各站點超閾值序列沒有明顯的突變性。
圖3 流域內四個站點的時間序列統計量變化圖Fig.3 Time series statistics of four sites in the basin
3.3 模型擬合效果檢驗
通過對流域內各站點進行GPD 擬合（圖4），并使用K-S檢驗看其是否符合已知理論分布函數。大部分站點的統計值小于0.01顯著性水平，接受原假設。說明雅江流域各站點之間雖存在差異，但由GPD 擬合曲線可知，理論頻數和實測頻數基本相符。另外，雖然西藏地區大部分臺站的觀測記錄起始年代不一致，如拉薩站、拉孜站，資料開始時間分別為1955 年、1977 年，但從圖4 中可看出，數據滿足方程需求，擬合結果表明資料長度并不影響降水極值的統計推斷，且資料年限越長擬合結果越好。
圖4 各站點累積頻率和實際頻率分布曲線對比Fig.4 The distribution of cumulative frequency and empirical frequency over the observation stations
3.4 雅江流域極端降水的分布特征
從各站點的極端降水閾值分布來看，林芝地區的閾值最大，閾值最小的站點為拉孜站。這與我們的觀測事實一致，林芝地區日降水量大，連續降雨日數長，得到的閾值就大。
尺度參數主要是描述極值分布的變率，尺度參數越大，極值波動范圍越大，表明打破極端降水的記錄值也越大。整體而言，雅江流域的尺度參數由下游向上游是逐漸減小的，平均值為5.95。由表3可知，下游地區的尺度參數最大，約為7.00，表明這一區域的極端降水變化幅度很大。從氣候背景來看，該地區位于高原季風區，受印緬槽和西風帶影響，季節性降水較大［22-23］，5—9 月的降水總量可達600 mm，是西藏地區夏季降水量最大的區域，因此可能出現的破極端降水記錄值要高于其他地區。尺度參數最小的區域位于流域中上游地區，受地理位置和大氣環流影響，雨期短，降水量少，且連續降水日數也少，降水極值的范圍比較小，區域打破極端降水的記錄值要比下游地區低。這與前人研究結果一致［24］，我國的干旱地區大部分位于非季風區，降水極值范圍較小，破紀錄的可能性較季風區小。
表3 流域內各站點的GPD模型參數估計及檢驗Table 3 Estimation and validation of parameters in GPD model
形狀參數作為模型的第二個重要參數，不同的形狀有不同的尾部分布特征，它表示該區域極端降水的破紀錄率。由表3 看出，形狀參數正值區主要位于拉孜地區，這些地區發生破紀錄降水事件的可能性比其他地區大。正是由于該區域降水日數少，所以一旦有降水過程，就可能會打破降水極值。而下游地區的形狀參數為負值，說明這些區域的降水發生破紀錄的概率偏小。因為夏季，降雨越頻繁的區域，極值變率大，較均值離散程度大，則破紀錄的概率較小。有降水的地區，不是有較大的形狀參數就是有較大的尺度參數，不可能兩個參數都大。形狀參數較大的地區，稱為“形狀參數主導區”，該區域多持續降水，極端氣候事件較少；尺度參數較大的區域稱為“尺度參數主導區”，該區域降水較少，且多變，極端氣候事件較多［23］。
3.5 雅江流域極值重現期預估
極值模型建立最重要的目的之一就是預測極端事件的重現期或重現水平。如表4 所示，根據重現期公式可以得出，從5 年一遇和10 年一遇的極端降水值來看，雅江流域除拉孜站外，其他地區降雨極值均超過30 mm，日喀則地區的降水極值達50 mm，其中拉薩、澤當、墨竹工卡和南木林站5 年一遇和10年一遇的極端降水量分別在40 mm 左右。相關研究指出［26-27］，當溫度上升為2 ℃時，青藏高原的強降水距平百分率平均增多44.5%～59.5%，大值區出現在山南附近，這與我們的研究結論相符合。各站點在15年一遇的極端降水值之后，極值水平的增長變得非常緩慢，其中林芝地區的增長最緩慢，以0.32 mm·a-1的速率增長；日喀則地區的降雨極值增長率最快，約0.72 mm·a-1。西藏地區極端降水頻率一般為每年4.3 次，強度在20 mm·d-1以上，林芝地區為極端降雨的高值區，且沿雅江一線極端降水的頻次呈增加趨勢［28-29］。
表4 雅江流域各站點日降水量極大值重現水平（單位：mm）Table 4 The maximum daily precipitation of flood season reappeared in Yarlung Zangbo River（unit：mm）
為檢驗各站點不同重現期水平的合理性，將其帶入各站點逐日降水序列中進行驗證。以拉薩和拉孜站為例，拉薩站5 年重現期水平值為39 mm，在1967—2017 年所有逐日降水中，共有10 次大于39 mm的降水過程，平均每5年一次。拉孜站5年重現期水平值為35.1 mm，在1977—2017 年期間，共有9 次大于35.1 mm 的降水過程，平均每4.4 年一次。由GPD 擬合計算出的極端降水重現期水平基本符合實際。值得注意的是這里的“重現期”并不意味著經過T年之后一定會出現的“周期”，它是概率意義上的“統計周期”。
4 結論
雅魯藏布江作為高原河流，由于強降水的時空分布不均而引起洪澇和干旱，并對流域內的農牧業產生重要威脅。因此在氣候日益增暖趨勢下，評估極端降水規律及其發生概率十分必要。本文通過引進GPD 概率分布模型，對西藏地區汛期強降水規律進行模擬。結果表明：
（1）通過Hill 圖法選取的流域內各站點的閾值序列小于百分位法選取的閾值序列，綜合考慮Hill圖法、百分位法及年交叉率法最終確定99百分位時的閾值為最佳閾值。
（2）各站點閾值序列在M-K 顯著性水平檢驗下，無明顯突變。擬合效果通過K-S檢驗，各站點擬合的理論頻數和實測頻數基本相符，且資料長度并不影響降水極值的統計推斷。
（3）通過分析流域內各站點擬合的極端降水特征可知，尺度參數的大值區位于流域下游，即林芝、米林地區，表明該地區的極值波動大；相反地，小值區位于流域中上游的拉孜站附近，表明極值波動小。形狀參數正值區位于流域中上游地區，說明發生破紀錄的降水事件概率較大，擬合結果與實際觀測一致。
（4）從5 年一遇和10 年一遇的極值水平看，雅江流域除拉孜站外，其他地區降雨極值均超過30 mm，日喀則地區的降水極值達50 mm；各站點在15 年一遇的極端降水值之后，極值水平的增長變得非常緩慢，由GPD 擬合計算出的降水極值具有一定的合理性。
