自然降水條件下關中平原農業水資源短缺風險分析

盧澤華,王君勤

(1.四川省農田水利局,成都 610015;2.四川省水利科學院研究院,成都 610072;3.四川大學 水利水電學院,成都 610041)

隨著全球氣候變化和社會經濟發展，水資源短缺問題日益嚴峻，已逐漸成為人類生存與發展的巨額的“生態資源赤字”[1]。2016年全國用水總量6 040.20億m3，其中農業用水(3 768.00億m3)占比62.38%，灌溉用水占農業用水的90%～95%，農業是我國的用水大戶[2]。我國農業水資源和農業水管理面臨巨大挑戰，對水資源系統進行有效的風險管理已成為水資源科學發展的必然趨勢[3]。農業水資源短缺風險評估作為農業水資源風險管理的基礎，可為區域水資源優化配置提供基礎依據，對區域種植結構調整，干旱風險應對具有重要意義[4]。

水資源短缺風險是指在特定的時空環境下，由于來水和用水存在的不確定性，使區域水資源系統發生供水短缺的概率及由此產生的損失[5]。水資源系統是一個復雜的大系統，廣泛存在隨機性和模糊性[6]。目前，對水資源短缺風險的研究大多從隨機模型或模糊模型的角度探討水資源短缺風險問題。羅軍剛等[3]構建了基于熵權的水資源短缺風險評價模型；王紅瑞等[6]構建了基于模糊概率的水資源短缺風險評價模型；Feng等[7]構建了基于信息擴散理論的水資源短缺風險評價模型；Hsieh等[8]構建了基于多站點徑流時空隨機模擬的灌溉用水短缺風險評價模型；Moursi等[9]構建了基于氣候響應函數和大氣環流模式的未來氣候變化下農業水資源短缺風險評價模型。以上模型具有較好的數理基礎，評價指標意義明確，但沒有充分考慮水資源系統組分間的非線性關系，無法對區域農業供水量和作物需水量之間復雜的非線性關系定量描述。

基于變量間非線性相關關系構建的Copula函數能獨立于隨機變量的邊緣分布反映變量間的相關性，可構造任意邊緣分布的聯合分布函數，且在轉換過程中不會改變原始隨機變量，在干旱風險分析、洪水遭遇風險分析方面已得到廣泛應用[10-11]，也有學者利用Copula函數進行水資源短缺風險分析研究[12-13]，如丁志宏等[12]構建了基于Copula函數的寧夏衛寧灌區降水量(P)和參考作物騰發量(ET0)年際聯合分布模型；Zhang等[13]構建了基于Copula函數的河南省陸渾灌區自然降水條件下的水資源短缺風險模型。然而有效降水量(Pe)和作物需水量(ETc)是描述灌溉系統的基本變量[14]，二者的演變特征與匹配程度直接影響農業水資源管理策略，比P和ET0更直觀反映自然降水條件下農業水資源的供需平衡。關中平原地處內陸、遠離海洋，其年降水量僅500～750 mm，對水資源短缺的響應尤為敏感，然而目前對關中平原Pe和ETc相關關系的研究較少，且大多是闡釋其年內、年際變化和簡單的相關關系[15-16]，對其非線性耦合關系還未報道。本文主要在現有研究基礎上，首次運用Copula函數構建關中平原Pe和ETc的二維聯合概率分布模型，應用該模型分析Pe和ETc的條件聯合概率和條件回歸周期，為關中平原灌溉農業發展規劃制定、作物種植制度優化及區域農業水資源管理等提供理論依據。

1 研究區概況與數據來源

關中平原(106°48′—110°36′E，33°35′—35°51′N)地處陜西省中部的渭河流域，介于秦嶺和渭北北山之間，西起寶雞，東至潼關，長約350 km，平均海拔約500 m，面積約3.6萬km2，屬暖溫帶半干旱氣候，是陜西省糧棉油主要產區[17-18]。根據降水量與氣候特征，關中平原可分為關中東部和關中西部[19-20]，其中關中東部年降水量500～650 mm，關中西部年降水量550～700 mm。關中平原糧食生產以一年兩熟的冬小麥和夏玉米輪作為主，冬小麥一般于10月初播種，翌年6月上旬收獲，夏玉米在冬小麥收獲后即時播種，當年9月底收獲[21]。本文以冬小麥播種—夏玉米收獲(10月—翌年9月底)作為一個單位年，計算單位年有效降水量(Pey)和作物需水量(ETcy)。

結合資料系列的完整性，選取關中平原6個代表性站點1962—2016年56 a逐日氣象數據作為基本資料，各站點氣象概況見表1。氣象數據來源于中國氣象科學數據共享服務網(http:∥data.cma.cn/)，包括日最高、最低和平均氣溫(Tmax，Tmin，Tmean)、相對濕度(RH)、日照時數(n)和距地面2 m高處的風速u2(計算時采用FAO推薦方法由10 m風速換算得出風速u2)，以及站點經緯度和海拔高度。少數的缺失數據(占全部數據的0.56%)采用線性內插法和多年平均值法補全。

2 研究方法

2.1 反距離權重法

反距離權重法IDW(inverse distance weighted)是一種充分考慮各因素之間地域性聯系的空間展布方法，因其原理簡便，結果精確，已被廣泛應用[22]。關中平原地形地貌復雜，不同站點海拔差異明顯，所以直接將區域內各站點Pey和ETcy取平均或選擇典型站點來代表區域Pey和ETcy具有明顯局限性。為增強選擇站點的代表性，本文將關中東部和關中西部各站點經度、緯度、海拔做平均，得到2個虛擬站點，根據反距離權重法分別計算各站點對應虛擬站點的權重，將各站點氣象因子乘以權重后相加，得到2個虛擬站點56 a逐日氣象資料，再計算2個虛擬站點的Pey和ETcy，分別代表關中東部和關中西部的Pey和ETcy，權重計算公式如下：

(1)

(2)

式中:xi，yi，zi分別為各站點的經度(°)、緯度(°)和海拔(m);xm,ym,zm分別為虛擬站點的經度、緯度、海拔;ri為第i個站點到虛擬站點的距離;n為區域站點個數;wi為第i個站點的權重。

2.2 有效降水量(Pey)和作物需水量(ETcy)

關中平原Pey的估算參考段愛旺[23]的研究成果，ETcy采用單作物系數法計算[17]。其中，關中平原冬小麥和夏玉米生育期的劃分及Kc的選取參考康紹忠[17]和時學雙等[19]的研究成果。

2.3 構建Pey和ETcy的邊緣分布函數

采用Copula函數進行Pey和ETcy2個特征變量的聯合之前，首先要確定其邊緣分布函數，同時要考慮變量間的相關性。設X和Y分別代表Pey和ETcy，其邊緣分布函數分別為FX(x)和FY(y)，選用常見的5種單變量分布函數(正態分布、伽馬分布、對數正態分布、廣義極值分布及韋爾布分布)[24-25]分別對其進行擬合，利用K-S檢驗確立最優邊緣分布函數。同時采用Kendall秩相關系數(τ)、Spearman秩相關系數(ρ)和Pearson相關系數(γ)度量關中東部和關中西部1962—2016年Pey和ETcy的相關性。本文中單變量分布函數及聯合分布函數中所含參數均采用極大似然法進行估計[24]。

2.4 構建Pey和ETcy的Copula聯合概率分布

Copula函數主要基于變量間相關性進行邊緣變量聯合，根據其定義可知Pey和ETcy的聯合分布函數F(x，y)為：

F(x,y)=P(X≤x,Y≤y)=C[FX(x),FY(y)]

(3)

常見的Copula函數有Archimedean Copula，橢圓Copula，Plackett Copula以及經驗Copula函數。其中，Archimedean Copula函數結構簡單，計算簡便，可以構造出形式多樣、適應性強的多變量聯合分布函數，能滿足多領域應用要求，在實際應用中占有重要地位。目前常用的Archimedean Copula函數有Frank，Ali-Mikhail-Haq，Clayton和Gumbel-Hougaard Copula函數，常見的相關性測度τ與Archimedean Copula函數的參數θ之間保持著一種對應關系，具體形式見參考文獻[11]。其中τ是可以描述變量之間非線性相關性的Kendall相關系數。由參考文獻中Copula函數參數θ和Kendall秩相關系數τ間關系式可知：Clayton Copula和Gumbel-Hougaard Copula函數秩相關系數τ的范圍為τ>0，Ali-Mikhail-Haq Copula函數秩相關系數τ的范圍為-0.1817≤τ≤0.3333，Frank Copula函數秩相關系數τ的范圍為τ∈R。本文基于Pey和ETcy兩個特征變量，分別采用參考文獻中Copula函數擬合這兩個變量的聯合分布函數，并采用K-S檢驗對Copula函數進行擬合檢驗[10]。

2.5 基于Pey和ETcy的條件概率和回歸周期

本文主要考慮在Pey分別處于高(p≤37.5%)、中(37.5%

(4)

(5)

同理，可以計算ETcy分別處于高(p≤37.5%)、中(37.5%

3 結果與分析

3.1 基于Copula函數的農業水資源短缺風險模型

3.1.1Pey和ETcy的邊緣分布函數確定 選用5種單變量分布函數分別對關中東部、關中西部Pey和ETcy進行擬合，利用K-S檢驗確立最優邊緣分布函數，K-S檢驗結果見表2，可以看出5種單變量分布函數對關中東部、關中西部Pey和ETcy的擬合優度均達顯著水平(α<0.05)。其中，關中東部擬合最優的分布函數分別是對數正態分布和廣義極值分布，關中西部擬合最優的分布函數分別是伽馬分布和廣義極值分布。

表2 5種分布函數對關中平原Pey和ETcy擬合的K-S檢驗結果

關中東部、關中西部Pey和ETcy邊緣分布函數的具體參數值見表3，可以看出關中東部Pey多年平均值小于關中西部，而ETcy多年平均值大于關中東部；關中東部、關中西部均表現為Pey的變異系數大于ETcy，即Pey相比ETcy變化更不穩定；關中東部、關中西部Pey和ETcy的經驗累積概率與理論累積概率擬合較好，均呈極顯著相關(α<0.01)，進一步表明邊緣分布及參數的選擇較為合理。劉俊民等[20]研究發現關中平原P的多年變幅較大，薛璐[26]研究發現關中平原ET0年際變化具有較好的穩定性，均與本文結果一致。

表3 關中平原Pey和ETcy邊緣分布的參數值、經驗概率與理論概率的決定系數

3.1.2Pey和ETcy的二維聯合概率分布 采用Kendall秩相關系數τ分析度量關中東部、關中西部Pey和ETcy間的相關性，根據τ與θ的關系，計算Copula函數的參數θ，并對Copula函數進行擬合優度評價，結果見表5。可以看出關中東部、關中西部Pey和ETcy均呈負相關，τ分別為-0.403和-0.228，因此僅Frank Copula函數適合描述Pey和ETcy間的相關關系。對關中東部、關中西部Pey和ETcy聯合概率分布的Frank Copula函數進行K-S檢驗，取顯著性水平α=0.05，n=55時，查柯爾莫格洛夫檢驗分位數表得出對應分位點D0為0.179 8，表5中檢驗統計量D小于分位點0.179 8，表明擬合優度均達顯著水平(α<0.05)，Pey和ETcy聯合分布的經驗累積概率與理論累積概率擬合較好，R2分別為0.976，0.987，且均呈極顯著相關(α<0.01)，因此關中東部、關中西部Pey和ETcy的聯合分布擬合效果較好，即選用Frank Copula函數是合理的。

3.2 農業水資源短缺風險模型應用

關中平原不同頻率的Pey和ETcy閾值見表5，可以看出關中東部的Pey小于關中西部的Pey，而關中東部的ETcy高于關中西部ETcy，因此關中東部的農業水資源短缺程度明顯高于關中西部。

表4 關中平原Pey和ETcy間的相關性度量及Copula函數的參數和擬合優度評價

表5 關中平原不同頻率的Pey和ETcy閾值

圖1給出了ETcy分別處于高(p≤37.5%)、中(37.5%

圖1 關中平原年作物需水量(ETcy)處于不同水平時年有效降水量(Pey)小于等于某特定值的條件概率和條件回歸周期

圖2 關中平原年有效降水量(Pey)處于不同水平時年作物需水量(ETcy)大于等于某特定值的條件概率和條件回歸周期

Pey分別為頻率37.5%，62.5%，87.5%的閾值時，由圖1可知，關中東部的ETcy在處于高(p≤37.5%)、中(37.5%

ETcy分別為頻率62.5%，37.5%，12.5%的閾值時，由圖2可知，關中東部Pey分別處于高(p≤37.5%)、中(37.5%

因此，關中東部、關中西部均表現為當ETcy處于高水平(p≤37.5%)時，Pey不超過某一特定值的條件概率最大，重現期最短；當Pey處于低水平(p≥62.5%)時，ETcy超過某一特定值的條件概率最大，重現期最短，即當ETcy處于高水平(p≤37.5%)或Pey處于低水平(p≥62.5%)時，關中平原自然降水和作物需水處于不協調狀況的可能性較高，供水不能滿足需水要求的概率較大，重現期較短(1～4 a)，農業水資源短缺風險較高。

因此關中東部和關中西部根據該區域的自然降水情況，有效地調節作物種植結構，充分利用區域自然水資源，以水資源總量控制作物布局，以水布局產業發展，建立與水資源承載能力相適應、與節水增收目標相配套的種植業結構。此外應該在考慮關中東部和西部灌區的特點基礎之上，推進農村特色產業區域化布局、規模化發展、產業化經營，形成產業聚集效應和發展的比較優勢，促進高效節水產業的快速發展。

本文對Pey和ETcy兩個特征變量進行分析只適用于描述自然降水條件下的農業水資源短缺風險，增加灌溉用水量、作物產量等特征變量，可進一步建立針對作物不同生育期不同環境下的農業水資源短缺風險，自然降水和灌溉水量組合條件下的農業水資源短缺風險等也是今后研究的主要內容。

4 結 論

(1)關中東部和關中西部Pey和ETcy擬合最優的邊緣分布函數為正態分布和廣義極值分布，關中西部Pey和ETcy擬合最優的邊緣分布函數為伽馬分布和廣義極值分布，Frank Copula函數可較好反映Pey和ETcy間的聯合分布特性。

(2)關中東部的農業水資源短缺程度明顯高于關中西部，當Pe處于低水平(p≥62.5%)或ETc處于高水平(p≤37.5%)條件下，關中平原自然降水不能滿足作物需水的概率較大，重現期較短(1～4 a)，農業水資源短缺風險較高。

(3)本文基于Pey和ETcy兩個變量來進行關中平原農業水資源短缺風險分析，因為歷史數據資料局限性，且隨機變量Pey和ETcy本身可能具有自相關性，本文所建立的聯合分布模型初步表征了Pey和ETcy在普遍意義下的統計特征，但二者同期空間相依關系仍需進一步探究。對隨機變量Pey和ETcy的自相關性進行研究，并基于此考慮單變量的短期時序相依關系，將是下一階段的主要研究內容。