關中地區氣象干旱時空變化特征分析

薛瑞華, 羅軍剛, 李少軒, 景 鑫

(西安理工大學 西北旱區生態水利國家重點實驗室, 陜西 西安 710048)

干旱是一種在世界各地頻繁發生的慢性自然災害,其發生與緩解都不如其他自然災害那么迅速[1],且干旱的發生通常由地區長期異常降水造成[2],這導致與其他大多數自然災害相比,干旱的影響更大、范圍更廣、持續時間更長。干旱通常分為四類:氣象、農業、水文和社會經濟干旱[3]。氣象干旱主要表現為降雨量不足[4]。當一定時期內的降水量持續少于正常狀態時,河流徑流量會相應減少形成水文干旱,當干旱持續發展時,土壤含水量將繼續減少直到不能滿足生物生長的需水量而形成農業干旱,最終導致作物減產造成經濟損失形成社會經濟干旱[5]。因此,氣象干旱是形成其他類型干旱的原因,監測氣象干旱對減輕干旱危害至關重要[6]。

關中地區是陜西省最典型、最重要的農業區,受氣候地形的影響,關中地區干旱頻發。據統計,近百年來關中地區河流徑流量嚴重減少,水資源開發難度大,但隨著經濟的發展以及人口的增長,對水的需求迅速增加,如今缺水率約為25%且仍在上升[7]。最近一次大型旱災是2014年8月上旬,關中地區遭遇了57年以來最嚴重的大旱,渭河、洛河等河流的徑流比往年同期偏枯一至四成,近20萬人用水緊張,咸陽65萬畝作物受旱,西安97萬畝作物受旱、11萬畝絕收[8]。因此,需密切關注干旱持續時間、嚴重程度等干旱特征,分析其變化趨勢與重現期,并根據干旱特征采取合理的應對措施以減少旱災造成的影響。

目前,關于關中地區的氣象干旱研究已有許多。Zhang等[9]通過標準化降水蒸散指數分析了1955—2015年間關中盆地的干旱變化趨勢;粟曉玲等[10]采用Gumbel Copula函數構建氣象水文綜合干旱指數,探討了該指數的適用性,并使用該指數分析了關中地區的干旱周期與頻率;Chen等[11]探討了1928—1931年陜西關中地區極端干旱事件的社會影響。綜上,已有的關中地區氣象干旱研究,主要是基于干旱指數進行干旱頻率、干旱強度等特征變量的單獨分析,少有對多特征變量進行聯合分析的。而干旱的發生與多個變量相關,其特征以多種方式體現,多特征變量聯合分析可以更好地表征復雜的干旱條件及其影響[12],因此,分析干旱多特征變量的聯合關系對研究干旱有著十分重要的意義。

Copula是最常用的多元分析方法之一,適用于多變量分析,已在極端降水、洪水頻率分析和干旱特征分析等方面得到廣泛應用[13]。劉瑞琳等[14]利用游程理論、單變量和多變量Copula方法,得到了橫斷山地區的干旱特征在不同重現期下的聯合分布特征。李京芳等[15]運用Copula函數計算干旱特征變量之間的聯合累積概率,分析了聯合重現期和同現重現期以及漢江流域水文干旱風險及驅動因素。因此,基于Copula的干旱特征分析在單獨考慮邊際分布或多變量相關性建模方面具有優勢,然而在使用Copula進行分析時,存在參數計算、邊際函數選擇和Copula函數選擇等多個方面的不確定性[16]。

貝葉斯網絡利用可信區間解決了參數的不確定性問題[17],因此本文采用包括基于貝葉斯框架的馬爾可夫鏈蒙特卡羅方法(Markov chain Monte Carlo,MCMC)與局部優化算法在內的多種方法進行參數計算,優選最合適的參數,并通過比較廣泛的邊際函數與Copula函數來有效解決不確定性問題,為關中地區多特征變量聯合分析提供有效方法。

本文基于SPI12時間序列,通過游程理論識別干旱特征,采用局部優化算法與基于貝葉斯框架下的馬爾可夫鏈蒙特卡羅模擬兩種算法進行Copula參數估計,得到最優的Copula函數及其參數;對干旱歷時與干旱烈度進行聯合分析得到關中地區的干旱分布特征,可為關中地區干旱風險評估與水資源管理規劃提供依據。

1 研究區概況與數據

關中地區總面積5.55萬km2,位于陜西省中部,屬于渭河流域,由南向北從秦嶺延伸至北山,由西向東從寶雞延伸至黃河,處于溫暖帶半濕潤與半干旱氣候的過渡地帶,屬大陸性季風氣候[18],西部、北部和南部海拔較高,東部和中部地區海拔較低。地貌類型主要由漫灘、階地和黃土臺地組成,行政區域包括楊凌示范區以及西安、寶雞、咸陽、渭南、銅川5市,是陜西省經濟最發達的地區。

本文收集到關中地區20個氣象站1979—2018年的日降水資料(來源于國家青藏高原科學數據中心http://data.tpdc.ac.cn/zh-hans/),在對數據進行計算前,對缺測的數據采用鏈式方程多重插補法進行插補。研究區域與氣象站的分布如圖1所示。

圖1 關中地區氣象站點分布

2 研究方法

2.1 標準化降水指數

標準化降水指數(standardized precipitation index, SPI)是用來表示某時段降水量出現概率多少的指標,適用于不同地區的干旱檢測與評估[19]。

(1)

(2)

(3)

其中,H(x)=q+(1-q)G(x),q是降雨量為0的概率,G(x)是服從伽馬分布的累計概率;C0=2.515 517,C1=0.802 853,C2=0.010 328,d1=1.432 788,d2=0.189 29,d3=0.001 308。

SPI指數計算結果:SPI>-0.5表示無旱,-1.0

SPI指數用于度量多種時間尺度上的降水不足,關中地區干旱頻發,幾乎年年有旱,因此采用12個月時間尺度的SPI值(記為SPI12)反映長期降雨模式,分析關中地區氣象干旱的年際變化。

2.2 游程理論

游程理論是由Yevjevich[20]提出的研究干旱特征最常用的方法之一[21]?？紤]到實際中會突然出現較大降雨,導致一次嚴重干旱事件的過程分割,從而弱化干旱的實際嚴重程度,故本研究采用王曉峰等[22]提出的三截取游程理論,并在此基礎上進行調整,得到適用于SPI指數的干旱事件識別游程理論。

1) 根據《氣象干旱等級》(GB/T 20481—2017),當SPI值小于或等于-0.5時,有干旱發生。

2) 若只識別出單獨一個月的事件,且這個月的SPI值大于-0.75,就認為這個月沒有發生干旱。

3) 若兩次干旱事件之間僅間隔一個月,且間隔月的SPI值小于0,則將兩次干旱事件合并為一次干旱事件。

通過確定出的干旱事件,基于三截取游程理論可提取出持續時間和烈度等干旱事件的基本特征[23]。其中,干旱歷時表示單次干旱過程的持續時間,干旱烈度表示單次干旱過程持續時間內SPI的累加值。

2.3 自適應Copula函數選擇

自適應Copula函數選擇是依據擬合優度的評價結果選取其中最優的方法作為Copula參數的估計方法。本文采用了兩種算法進行Copula的參數估計,分別是局部優化算法與基于貝葉斯框架下的馬爾可夫鏈蒙特卡羅模擬算法。

2.3.1邊緣分布與相關性檢驗

在Matlab中采用MvCAT工具箱所包含的指數分布、伽馬分布、威布爾分布等在內的16種函數進行單變量邊緣分布函數擬合,并對擬合后的函數進行K-S檢驗,通過赤池信息準則(Akaike information criterion,AIC)進行擬合優度檢驗,根據AIC最小值原則選出合適的邊緣分布函數。

采用Kendall秩相關系數法、Pearson相關系數法和Spearman秩相關系數法對干旱歷時與干旱烈度進行關聯性檢驗。

2.3.2Copula函數

Copula函數通常用于從單變量邊際分布生成多變量分布,該技術的理論基礎是由Sklar[24]定義的。根據Sklar定理,如果FXY是相依隨機變量X和Y的二維聯合分布函數,其邊緣分布函數為FX和FY,則存在一個CopulaC,即:

P(X≤x,Y≤y)=FXY(x,y)=C(FX(x),FY(y))

(4)

上面定義的函數FXY是一個具有邊緣分布函數FX和FY的二維函數。

本文采用了MvCAT工具箱中的Gaussian、t、Clayton、Frank、Gumbel等26種Copula函數進行計算[25]。

2.3.3局部優化算法

本文采用的局部優化算法是一種基于梯度的內點優化算法[26]。此算法通過一個稠密的擬牛頓近似方法(dense quasi-Newton approximation)估算Hessian矩陣,并使用內點法尋找可行解空間。

在求解參數的過程中,每次的參數估計將從隨機的30個點出發尋找最優解,以此減少陷入局部最優的概率。

2.3.4貝葉斯框架

貝葉斯統計推斷對于不確定性評估和參數估計是有效的,其中假定參數為具有概率分布的隨機變量[16]。貝葉斯推理的基本定理是:

(5)

2.3.5馬爾科夫鏈蒙特卡洛模擬

本研究使用一種新混合進化的MCMC方法,該方法采用自適應采樣分布來描述貝葉斯環境中的后驗參數區域。首先通過拉丁超立方抽樣(Latin hypercube sampling, LHS)算法對整個先驗空間進行隨機搜索,通過樣本的隨機分配從中選擇出可能性值最高的樣本作為馬爾科夫鏈的起點,這種智能起點的選擇使不同吸引區域的點作為起始點成為可能,從而避免了退化。在進行后驗參數空間搜索時,采用10%概率的斯諾克方向更新和90%概率的并行方向更新,在并行方向中有NAM條馬爾科夫鏈采用自適應采樣算法(adaptive metropolis, AM) 豐富跳躍方向,增強MCMC早期的搜索能力,其他鏈則采用差分進化算法(differential evolution, DE)加快收斂速度。

2.3.6擬合優度評價

在本研究中采用均方根誤差(root mean square error,RMSE)、納什效率系數(Nash-sutcliffe efficiency coefficient,NSE)以及AIC進行Copula函數擬合優度檢驗。

1) 均方根誤差

(6)

(7)

(8)

2) 納什效率系數

(9)

3) 信息準則

AIC=nln (MSE)+2

(10)

式中:Pei為聯合分布經驗頻率;Pi為聯合分布理論頻率;n為樣本容量;mi為干旱特征小于某一閾值的干旱事件發生的次數。

2.4 重現期計算

單變量重現期通常會高估或低估某一特定事件的風險率。因此,本研究計算了氣象干旱特征的二變量重現期。該方法假設X1、X2代表與干旱事件相關的特征變量值,x1、x2代表與干旱事件相關的特征變量序列對應干旱特征的給定閾值。

二維聯合回歸周期可以表示為:

(11)

二維共現重現期可以表示為:

(12)

式中:N表示總觀測次數;n表示樣本容量。

如圖2所示,本文主要包括干旱指數計算、干旱特征識別、邊緣分布函數篩選、構建Copula聯合分布以及關中地區干旱特征分析這五個部分。

圖2 干旱特征分析流程圖

3 結果分析

3.1 干旱特征

根據各測站的月累計降雨量數據,采用算數平均法計算出關中地區五個市的逐月累計降雨量,并根據逐月累計降雨量得到各市的SPI12時間序列。關中地區及五個市1980—2018年逐月累計降雨量與SPI12時間序列的變化如圖3所示。從圖中可以看出,關中地區整體的SPI12指數呈現下降趨勢,寶雞市與渭南市的SPI12指數有較明顯的下降趨勢,西安市和銅川市的SPI12指數呈現不明顯的下降趨勢,咸陽市的SPI12指數呈現不明顯的上升趨勢,表明在研究時段內關中地區整體的干旱愈發嚴重。

圖3 關中地區月累計降水量與SPI12指數變化圖

從圖3可知,關中地區經常發生輕旱與中旱,但重旱與特旱的情況較少。關中地區在1997年間干旱情況十分嚴重,這與1994—1997年間降雨量偏少有著十分密切的關系,根據《中國氣象災害記錄》記載,關中地區1995年全年降水嚴重偏少,其中上半年偏少達6～8成,9月初至12月底降水偏少達5～8成,1997年4～7月降水比常年偏少3～9成,8、9月降水偏少6～9成[27],與研究結果一致。從圖3中還可知,近年來較為嚴重的干旱發生在2014年,這與報道的多個縣區降水量偏少5成以上,出現大面積未灌或少灌區域,部分地區的獼猴桃出現枯死現象的結果一致。

通過游程理論,對各市SPI12時間序列進行識別,得到寶雞市、咸陽市、西安市、銅川市以及渭南市的干旱次數分別為16次、14次、17次、21次以及18次。圖4為不同地區干旱歷時與干旱烈度的分布情況。關中地區整體的平均干旱歷時為7.55月,其中寶雞、咸陽、西安的平均干旱歷時分別為8.94月、9.86月以及9.24月,高于關中地區的平均干旱歷時;渭南的平均干旱歷時為7.83月,與關中地區的平均干旱歷時相接近;銅川的平均干旱歷時為6.38月,低于關中地區的平均干旱歷時。關中地區的平均干旱烈度為8.10,寶雞、咸陽、西安、渭南的平均干旱烈度分別為9.90、11.06、9.70以及9.24,均高于關中地區的平均干旱烈度;銅川的平均干旱烈度為7.10,低于關中地區的平均值。寶雞、咸陽以及西安的平均干旱歷時與平均干旱烈度均較大,干旱情況較為嚴重;渭南市平均干旱歷時較小,平均干旱烈度較大,且烈度的中位數遠遠小于平均值,表明渭南市多發小干旱事件,但也有極端干旱情況的出現;銅川市的平均干旱歷時與平均干旱烈度均較小,表明銅川市的干旱情況較輕。

3.2 邊緣分布

表1給出了20個測站識別出的干旱歷時與干旱烈度的邊緣分布函數與相關性檢驗。由表1可知,關中地區的干旱歷時主要服從指數分布和廣義帕累托分布,占所有測站的70%,其次是疲勞壽命分布以及逆高斯分布;干旱烈度主要服從指數分布和廣義帕累托分布,占所有測站的80%,其次是疲勞壽命分布以及逆高斯分布。采用Kendall rank、Spearman’s rank-order以及Pearson product-moment三種相關系數進行檢驗,見表1的第4～6列。從表1中可以看出,除了扶風站,所有測站的Kendall值都在0.8以上,Pearson相關系數都在0.9以上,Spearman相關系數除長武站都在0.9以上,且所有站點的擬合結果均通過顯著性為5%的K-S檢驗,表明所有測站的干旱歷時與干旱烈度均有密切的相關關系。

表1 干旱歷時與干旱烈度邊緣分布情況與相關性檢驗

3.3 Copula聯合分布

采用MvCAT工具箱內的局部優化算法與貝葉斯框架下的馬爾科夫鏈蒙特卡洛模擬兩種方法進行參數估計,并結合RMSE、NSE、AIC三種擬合優度評價指數選擇出最優的Copula函數,結果見表2,其中θ1、θ2分別對應Copula函數中待定參數的最優值。從表2中可以看出,有30%的測站采用局部優化算法,有55%的測站采用馬爾科夫鏈蒙特卡洛方法,有15%的測站兩種方法的最優情況相同。適合關中地區干旱歷時與干旱烈度聯合分布的Copula函數主要為Gaussian、Clayton以及BB1這三種函數。其中對BB1函數的參數估計更適合使用局部優化算法,對Gaussian以及Clayton函數的參數估計更適合使用馬爾科夫鏈蒙特卡洛模擬。

3.4 重現期分析

干旱重現期的計算結果結合了干旱發生的概率,是一種可用來評估干旱風險的簡單實用的方法[28]。根據計算結果與克里金插值法繪制干旱歷時與干旱烈度的聯合重現期與同現重現期,如圖5、圖6所示。其中(a)、(b)、(c)、(d)4個圖分別為根據關中地區總體干旱特征分布情況確定的小干旱特征事件(干旱歷時為2個月,干旱烈度為1.5)、平均干旱特征事件(干旱歷時為8個月,干旱烈度為8)、嚴重干旱特征事件(干旱歷時為12個月,干旱烈度為12)以及極端干旱特征事件(干旱歷時為26個月,干旱烈度為20)下的重現期。

圖6 關中地區干旱歷時與干旱烈度同現重現期分布

從圖5中可以看出,對于小干旱事件,關中地區整體的聯合重現期均在一年多;平均干旱事件的短歷時聯合重現期主要分布在寶雞以及西安市東部地區,為1 a多,長歷時聯合重現期主要分布在銅川、咸陽東部、渭南西部以及西安北部地區,最長為4.7 a;嚴重干旱事件的聯合重現期自西向東呈增加趨勢,聯合重現期最短出現在寶雞的扶風,為2.19 a,聯合重現期最長出現在銅川市的洛川測站,為5.3 a;極端干旱事件長歷時聯合重現期主要分布在銅川東部與渭南西北部,最長聯合重現期出現在銅川市的洛川測站,為11 a,短歷時聯合重現期主要分布在寶雞東南部與西安西部地區,最短聯合重現期出現在西安市的武功測站,為3.7 a。

從圖6中可以看出,長同現重現期主要分布在銅川市,短同現重現期主要分布在寶雞市。小干旱事件的短同現重現期主要集中在銅川市以及咸陽市東部地區,最短同現重現期出現在寶雞市扶風站,為1.1 a,長同現重現期主要集中在寶雞市、西安市西部以及渭南市東部地區,最長同現重現期出現在銅川市的宜君站,為1.6 a;平均干旱事件的短同現重現期主要集中在寶雞市以及渭南市東部,最短同現重現期出現在寶雞市扶風站,為2 a,長同現重現期主要分布在銅川市,最長同現重現期出現在銅川的宜君站,為5.5 a;嚴重干旱事件的短同現重現期主要分布在寶雞市、西安市西部、咸陽市西南部以及渭南市東部地區,最短同現重現期出現在寶雞市的扶風站,為3.1 a,長同現重現期主要分布在銅川市以及渭南市西部,最長同現重現期出現在銅川市白水站,為10 a;極端干旱事件的同現重現期自西向東呈增加趨勢,最短同現重現期出現在咸陽市的永壽站,為10.6 a,最長同現重現期出現在銅川市白水站,為138.9 a。

從圖5、圖6中可以看出,聯合重現期與同現重現期有著相似的分布規律,由于同現重現期的出現有著更高的要求,因此同現重現期隨著烈度歷時的增大而快速增加。關中地區最頻發的為小干旱特征事件,幾乎每年都會出現,其他干旱事件自西向東重現期逐漸增加,長重現期主要分布在銅川市以及渭南市西北部地區,短重現期主要分布在寶雞市、西安市西部以及渭南市東部少部分地區。

關中地區干旱特征呈現此種變化的原因與降水分布有關,通過馮晶等[29]對關中地區降雨特征的分析可知,除咸陽市外,其他城市的降雨量均呈下降趨勢,這與本文的研究結果一致,同時馮晶還指出,降雨自西到東呈增加趨勢,因此采用基于降水數據的SPI12指數計算出的重現期也呈現西邊重現期短、東邊重現期長的變化特征。Zhang等[9]通過計算干旱頻率得到關中地區頻發輕旱事件,這也與本文的研究結果一致。

本文通過重現期計算,得到關中西部、西安周邊以及渭南市東部地區的重現期較短,干旱較為嚴重;喬亮[30]通過對關中地區旱災累積程度的分析,認為寶雞、西安、藍田、韓城、潼關(寶雞位于關中西部,西安和藍田位于西安,韓城以及潼關位于渭南市東部地區)的旱災累積程度最為嚴重,這也與本文的研究結果一致。

4 結 論

本文基于SPI12指數,結合游程理論,分析了關中地區的干旱特征變化情況,并利用Copula函數計算出聯合特征重現期,得出以下結論。

1) 適合關中地區干旱歷時與干旱烈度聯合分布的Copula函數主要為Gaussian、Clayton以及BB1這三種函數。其中BB1函數的參數估計適合使用局部優化算法,Gaussian以及Clayton函數的參數估計更適合使用馬爾科夫鏈蒙特卡洛模擬算法。

2) 關中地區整體的SPI12指數值呈現下降趨勢,寶雞市與渭南市的SPI12指數有較為明顯的下降趨勢,西安市和銅川市的SPI12指數呈現不明顯的下降趨勢,咸陽市的SPI12指數呈現不明顯的上升趨勢,說明近年來關中地區的干旱情況有加重的趨勢。

3) 寶雞市、咸陽市以及西安市的干旱烈度與歷時均值較大,干旱情況較為嚴重,渭南市干旱特征平均值大、中位數小,表明有極端干旱的情況發生,但小干旱事件居多,銅川市干旱情況較輕,干旱特征值較小。

4) 基于干旱歷時與干旱烈度的聯合分布計算,關中地區頻發小干旱特征事件,幾乎每年都會出現,其他干旱事件的重現期自西向東逐漸增加,聯合重現期與同現重現期有相同的變化趨勢,長重現期主要分布在銅川市以及渭南市西北部地區,短重現期主要分布在寶雞市、西安市西部以及渭南市東部少部分地區。