河套灌區出讓水權對天然植被影響研究

常布輝 李根東 蘇小飛 王會永 王艷明 溫樂

關鍵詞：ＭＩＫＥＳＨＥ；生態環境約束；水權交易；沈烏灌域；河套灌區

０引言

水權交易的本質是對水資源進行時空的二次分配［１－２］，因此可能會引起一系列生態負外部效應，尤其是在降水量少、以地表水為主要水源的西北內陸灌區［３］。在水資源短缺加劇的背景下必然存在著復雜的社會利益沖突［４－５］，因此國內外學者對水權交易的外部生態環境效應進行了大量研究。

張忠［６］結合東陽義烏地區水權交易的典型案例，圍繞水權交易過程中產生的外部效應進行了闡釋，研究結果表明當初始水權已經被不同的水資源利用主體完全分配，針對生態安全等新增的生態用水問題就必須通過水權的轉讓或交易解決。基于生態經濟學的Ｅｍｅｒｇｙ理論和外部理論，Ｌｙｕ等［７］在分析農業水權交易導致出讓水權者產生生態經濟損失后，提出了一種包括社會、經濟和生態環境補償的突發事件量化方法，同時發現中國水權交易市場的生態經濟補償呈持續上升趨勢，２０１３年至２０１７年由９４００萬元上升到了１．１４億元。

２０世紀８０年代國外出現了水權交易與生態環境的關系研究［８］。例如：美國加利福尼亞州的ＳａｎＪｏａｑｕｉｎ河谷將農業水權轉向工業，導致灌溉補水不足，使得砷等有毒物質得不到有效的降解，最終導致Ｋｅｓｔｅｒｓｏｎ國家野生動物園庇護所發生污染并關閉［９］。洛杉磯通過水權交易獲得蒙諾湖支流取水權后，導致該湖的水資源量大幅度下降，生態平衡被打破，出現了大量野禽死亡的現象［１０］。Ｅｔｃｈｅｌｌｓ等［１１］發現水資源在維持河流生態環境方面的價值相當于或大于其消費性利用的價值，水權交易不僅降低了水資源的質量，而且可能引起地下水位下降、土壤鹽堿化等問題，從而降低水資源的再生能力。Ｌａｄｓｏｎ等［１２］從環境角度解釋了改變水權交易取水的時間和地點所帶來的外部環境效應：在空間上，由于不同取水點和運輸渠道的漏損率不同，因此在進行水權交易時對第三方產生的影響程度也不同；在時間上，干旱期取水比非干旱期取水所帶來的外部影響大。Ｈｅａｎｅｙ等［１３］考慮時間因素的影響，把水權交易外部效應分為供給可靠性、取水時效性、儲蓄水及取水費用、水質４類。Ｈａｎａｋ［１４］討論美國加利福尼亞州地下水區域間交易對第三方回流的影響，認為出口禁令可以有效限制地下水輸出，減輕地下水位的下降程度。

基于以上研究可知，雖然目前國內外學者在水權確權方法、水權交易模式和水權交易外部影響等研究領域已經獲得了較為豐富的研究成果，但在水權交易外部影響方面，國內對水權交易的經濟社會效應關注較多，對生態環境效應的研究相對較少，水權交易對生態環境的作用機理、影響程度及應對措施等仍需要深入研究［１５］。

１研究區概況與數據來源

１．１研究區概況

河套灌區位于內蒙古自治區巴彥淖爾市，是全國３個特大型灌區之一，也是我國最大的一首制自流引黃灌區。沈烏灌域位于三盛公水利樞紐西北部，烏蘭布和沙漠穿沙公路以北，磴口縣與杭錦后旗行政邊界以南，河套總干渠及烏拉河干渠以西，狼山沖洪積坡地邊界以東，總面積約１８．６８萬ｈｍ２。沈烏灌域屬于溫帶大陸性干旱氣候區，降水稀少，蒸發強烈，干燥多風，晝夜溫差大，日照時間長。２０１９年沈烏灌域通過實施節水改造向外轉讓水權１．２億ｍ３，現狀引水由原來的４．２６億ｍ３下降到３．０６億ｍ３。

１．２數據來源

１）地下水埋深。課題組在沈烏灌域布設了４７眼地下水埋深監測井，監測獲取了２０１６—２０１９年的五日地下水埋深。

２）歸一化植被指數（ＮＤＶＩ）。選用２０１６年８月２１日Ｌａｎｄｓａｔ＿５和２０１９年９月４日Ｌａｎｄｓａｔ＿７衛星數據，分辨率為３０ｍ，云量少，使用ＥＮＶＩ５．２經過輻射定標、大氣校正、裁剪等預處理后計算得到研究區自然植被分布區的ＮＤＶＩ值。

２技術路線

技術路線見圖１。

３地下水埋深與ＤＮＶＩ之間的關系

３．１數據處理

ＮＤＶＩ數據是３０ｍ×３０ｍ分辨率的柵格數據，而地下水埋深數據來源于分散在研究區的地下水監測井，為分析二者之間的關系，需要對數據進行預處理，以保證兩個變量數據的匹配性。數據預處理的方法有３種。

１）空間網格作為基本單元。本研究分析空間網格設置為１ｋｍ×１ｋｍ，需要利用ＡｒｃＧＩＳ工具將３０ｍ×３０ｍ的ＮＤＶＩ數據重采樣至１ｋｍ×１ｋｍ。同時，利用克里金方法插值到對應空間網格。通過匹配共獲取５１７３６組數據樣本，基于空間網格的地下水埋深與ＮＤＶＩ散點圖（圖略）可以看出，兩者之間的關系沒有明顯規律。

２）泰森多邊形法劃分基本單元。根據離散分布的地下水埋深監測井，得到地下水埋深監測井的泰森多邊形，如圖２所示。利用ＡｒｃＧＩＳ的空間統計工具計算每個泰森多邊形內ＮＤＶＩ平均值，２０１６—２０１９年３９個（部分監測井數據不連續）多邊形區域一共確定１５６組樣本。

利用泰森多邊形獲取的１５６組樣本繪制地下水埋深與ＮＤＶＩ關系圖（見圖３），由圖３可以看出，地下水埋深與ＮＤＶＩ成較明顯的線性關系，隨著地下水埋深的增大，ＮＤＶＩ明顯下降。

３）土質單元劃分基本單元。土壤類型是影響植被生長的重要因子，根據世界土壤數據庫土質類型，沈烏灌域可劃分為６６個土質單元。利用ＡｒｃＧＩＳ的空間統計工具計算每個土質單元ＮＤＶＩ平均值和地下水埋深，２０１６年各土質單元ＮＤＶＩ平均值和地下水埋深如圖４所示。

２０１６—２０１９年６６個土質單元（其中１個土質單元無ＮＤＶＩ數據）一共確定２６０組樣本。圖５為２６０組樣本的地下水埋深與ＮＤＶＩ關系散點圖，由圖５可知，地下水埋深小于３．４ｍ時，地下水埋深與ＮＤＶＩ沒有關系；地下水埋深大于３．４ｍ時，地下水埋深與ＮＤＶＩ成明顯的線性關系。

３．２結果分析

通過數據預處理可知，基于泰森多邊形基本單元和土質基本單元均得到地下水埋深與ＮＤＶＩ的線性關系，利用皮爾遜相關系數法和多元線性回歸法分別對其進行相關關系分析。

１）泰森多邊形基本單元的分析結果。２０１６—２０１９年３９個泰森多邊形基本單元的地下水埋深與ＮＤＶＩ的皮爾遜相關系數為－０．３７５，且通過顯著性水平為１％的顯著性檢驗。地下水埋深與ＮＤＶＩ的相關性屬于負相關，即地下水埋深越小，泰森多邊形基本單元內的ＮＤＶＩ平均值越大。

以泰森多邊形基本單元的ＮＤＶＩ平均值作為被解釋變量，將地下水埋深、土壤含鹽量、土壤沙含量、土壤淤泥含量、土壤黏土含量、土壤有機碳含量、土壤電導率、土壤碎石體積百分比等要素因子作為解釋變量，構建多元線性回歸分析模型，得到以下結論：多元線性回歸模型的擬合優度為０．１４，表明選取的解釋變量對ＮＤＶＩ值的解釋程度為１４％；地下水埋深的檢驗統計量為－３．８２，Ｐ值小于０．０５（Ｐ在ｔ檢驗中用于對各變量系數顯著性進行檢驗，Ｐ＜０．０５表示顯著，Ｐ＞０．０５表示不顯著），表明泰森多邊形單元內的地下水埋深對ＮＤＶＩ平均值具有顯著的影響；地下水埋深的偏回歸系數為－０．００８６，說明在其他解釋變量不變的條件下，地下水埋深每下降１ｍ，泰森多邊形單元內ＮＤＶＩ平均值會降低０．００８６。

２）土質基本單元的分析結果。當地下水埋深小于３．４ｍ時，土質基本單元地下水埋深與ＮＤＶＩ沒有關系；當地下水埋深大于３．４ｍ時，土質基本單元地下水埋深與ＮＤＶＩ成顯著線性關系。因此，將地下水埋深３．４ｍ作為臨界值，將２０１６—２０１９年６６個土質單元的２６０組樣本分別進行分析。在２６０組樣本中，有２２０組樣本地下水埋深小于３．４ｍ，記為Ａ組；４０組樣本地下水埋深大于３．４ｍ，記為Ｂ組。

Ａ組地下水埋深與ＮＤＶＩ的皮爾遜相關系數為０．０２，且未通過顯著性檢驗，表明當地下水埋深小于３．４ｍ時，地下水埋深與ＮＤＶＩ不存在相關性。Ａ組多元線性回歸的結果顯示，地下水埋深的檢驗統計量為０．０５（Ｐ＞０．０５），沒有通過顯著性檢驗；而土壤鹽分和土壤沙含量、土壤黏土含量等變量通過了顯著性檢驗（Ｐ＜０．０１）。結果表明：當地下水埋深小于３．４ｍ時，地下水埋深對ＮＤＶＩ沒有顯著的影響，而土壤鹽分和土壤質地對ＮＤＶＩ具有顯著的影響。

Ｂ組地下水埋深與ＮＤＶＩ的皮爾遜相關系數為－０．５５９，且通過了顯著性水平為１％的顯著性檢驗，表明當地下水埋深大于３．４ｍ時，地下水埋深與ＮＤＶＩ存在顯著負相關性。Ｂ組多元線性回歸模型的擬合優度為０．５８，表明選取的解釋變量對ＮＤＶＩ值的解釋程度為５８％；地下水埋深的檢驗統計量為－２．３５（Ｐ＜０．０５），表明當地下水埋深大于３．４ｍ時，地下水埋深對土質基本單元ＮＤＶＩ平均值具有顯著的影響；地下水埋深的偏回歸系數為－０．０８０７，說明在其他解釋變量不變的條件下，當地下水埋深大于３．４ｍ時，地下水埋深每下降１ｍ，土質基本單元ＮＤＶＩ平均值會降低０．０８０７。

通過對比泰森多邊形基本單元和土質基本單元Ｂ組的多元線性回歸結果可知，土質基本單元Ｂ組的多元線性回歸模型的擬合優度遠高于泰森多邊形基本單元。由土質基本單元Ａ組的多元線性回歸結果可知，土壤質地在一定條件下對ＮＤＶＩ具有顯著影響，因此泰森多邊形基本單元的分析結果準確性會受區域內土壤質地差異性的影響。同時考慮泰森多邊形基本單元多元線性回歸模型的擬合優度相對較低，本研究認為土質基本單元的分析結果具有更高的可信度。

綜上所述，利用皮爾遜相關系數法和多元線性回歸法分析地下水埋深與ＮＤＶＩ的關系，得到以下結論：當地下水埋深小于３．４ｍ時，地下水埋深與ＮＤＶＩ沒有相關性；當地下水埋深大于３．４ｍ時，地下水埋深與ＮＤＶＩ具有顯著負相關性，在其他因素不變的條件下，地下水埋深每下降１ｍ，土質基本單元的ＮＤＶＩ平均值會降低０．０８０７。

４模型構建及驗證

結合沈烏灌域水循環特點，構建基于ＭＩＫＥＳＨＥ模型的沈烏灌域分布式水文模型。鑒于灌域內水分水平運動比較弱，而垂向運動較強的特點，模型主要由氣象模塊、地表模塊（簡化處理）、不飽和帶、飽和帶以及一維河流模塊組成。

４．１模型結構及原理

ＭＩＫＥＳＨＥ模型由數個獨立的且相互聯系的基于物理過程的模塊構成，每個子模塊用于一個主要的水文過程的描述，根據不同的模擬要求，這些模塊可以互相分離也可以綜合起來應用，分離開來可以分別描述水文循環的各個過程，綜合起來可以描述整個流域的水文循環過程。

１）蒸散發計算模塊。在ＭＩＫＥＳＨＥ蒸散發模塊的計算中，需要運用氣象和植被覆蓋數據對降雨蒸發量進行預測。其中包括植物冠層對降雨的截留量、由植物冠層降落到土壤層的水量以及植物冠層、土壤表面和根系的蒸發量等決定因素。

２）坡面流模塊。當流域地表的凈雨強度超過土壤的下滲能力時，就會產生地表徑流，地表徑流沿坡面匯入河網，最后由河網內的流域出口斷面流出。而徑流在河網內的流經路徑及徑流量主要由流域地形、地表的摩阻力以及通過蒸散發和下滲失去的水量決定。在ＭＩＫＥＳＨＥ的坡面流模塊中，主要采用擴散波對圣維南方程組進行近似模擬，運用有限差分的方法對方程組進行計算。

３）不飽和帶模塊。不飽和帶是流域水文系統中一個至關重要的部分，同時也是ＭＩＫＥＳＨＥ的一個核心計算模塊，在模型模擬應用中起著重要作用。受重力作用影響，不飽和帶出流在下滲過程中主要作垂直方向上的運動。ＭＩＫＥＳＨＥ計算模塊包含了反復迭代的計算過程，用于校正土壤含水量以及計算土壤淺層地下水位的動態變化。

ＭＩＫＥＳＨＥ在計算不飽和帶模塊時主要采用３種方法：Ｒｉｃｈａｒｄｓ方程、模擬簡單的重力水出流過程、兩層水量平衡模型。

４）飽和帶計算模塊。在ＭＩＫＥＳＨＥ的水流運動模塊中，飽和帶模塊用于計算流域中的飽和地下水流運動。ＭＩＫＥＳＨＥ可以描述三維水流在異質含水層不定邊界上的運動。模型對于時空變換的獨立參數，采用三維達西方程進行模擬，并用隱式有限差分進行數值求解。

ＭＩＫＥＳＨＥ對上述地下水模擬提供兩種解決辦法，一種是逐次超松弛法，另一種是前承條件共軛梯度法。

５）河流計算模塊。模型中的不同水量最后通過河網進行匯流，最后到達流域出口。坡面流、壤中流和基流通過格網向流域較低的方向傳播，最后到達流域相鄰格網的最低點，也就是河流的位置，然后匯入河流，進行下一步的河流匯流計算。由于ＭＩＫＥＳＨＥ本身不包含河道匯流的計算，因此需要耦合ＤＨＩ公司研發的另外一個水動力演算模型ＭＩＫＥ１１進行這一部分的計算。耦合了ＭＩＫＥ１１水動力模型之后，能夠對地表和地下的水流動力過程進行完整的、動態的耦合模擬。

ＭＩＫＥ１１河流模型采用的是水動力學模型，即明渠不穩定流隱格式有限差分解，其差分格式采用六點中心隱式差分格式，其數值計算采用傳統的“追趕法”。ＭＩＫＥ１１計算參數包括兩類：數值參數，主要是方程組迭代求解時的有關參數，如迭代次數及迭代計算精度；物理參數，主要是河網的阻力參數。

４．２模型構建

１）地表高程。本文使用的地表高程數據來源于地理空間數據云提供的ＳＲＴＭＤＥＭＵＴＭ３０ｍ分辨率數字高程數據產品。

２）氣象數據。搜集沈烏灌域附近氣象站數據，用Ｐｅｎｍａｎ－Ｍｏｎｔｅｉｔｈ公式計算研究區域的參考蒸散發量。

３）土地利用。通過遙感影像解譯，獲取區域及各分干渠范圍內耕地、水域、草地、居民區、道路、裸地等土地利用類型。土地利用中不同下墊面對應的植被葉面積指數、根系深度以及蒸散發參數都是不同的。葉面積指數的獲取在沒有實測數據的前提下，借鑒已有研究中ＮＤＶＩ與葉面積指數之間的經驗關系。

４）農田灌溉。灌溉模塊是對灌域從引水到田間灌水過程的模擬，引水經過渠道逐級損失進入田間，一部分形成下滲補給土壤水和地下水，一部分消耗于蒸散發。由于每年種植結構以及總引水量受人為影響大，總引水在各分干渠及主要支渠間的配水過程非常復雜，因此沒有客觀規律可循。為了盡可能真實反映灌溉過程中的渠系滲漏、田間灌溉和滲漏，在渠道引水測流數據的基礎上，結合渠道滲漏試驗和田間試驗確定的渠系水利用系數、田間水利用系數，分渠段計算得出相關渠段控制范圍內的渠系滲漏水量和田間灌溉水量。

５）農田排水溝。模型中河流模塊用于模擬排水溝的排水過程。排水溝位置確定后，通過實測的溝底高程以及末端測流斷面特征和水位流量關系，結合排水溝比降，每２ｋｍ設置一個計算斷面，以滿足計算精度的需要。最終通過河流模塊與ＭＩＫＥＳＨＥ地表、土壤、地下水模塊之間的耦合完成設置。

６）土壤數據及參數。土壤數據主要來源于ＦＡＯ網站上的ＨＷＳＤ數據庫。數據庫中土壤分為上層（０～３０ｃｍ）和下層（３０～１００ｃｍ）。運用ＳＰＡＷ軟件的ＳＷＣ模塊，結合研究區土壤柵格屬性表，按對應的屬性尋找出變量。根據屬性表中查出的變量，運用ＳＰＡＷ軟件計算得出ＭＩＫＥＳＨＥ模型需要的土壤參數。

４．３模型參數率定及驗證

２０１５年為率定期，２０１６年為驗證期。通過調整地下水模型結構與含水層參數，利用不同時期地下水實測數據及二維流場，對模擬水位和流場進行校驗。通過修正模型結構與參數，控制并再現地下水徑流條件，使地下水達到了較好擬合。由擬合結果（圖略）可以看出，模擬結果與實測數據的擬合情況較好。對已有的４４眼地下水觀測井擬合結果進行統計，見表１。由表１可以看出，平均誤差在０．５ｍ以內的井數為３４眼，占比７７．２７％。綜上所述，模型精度滿足要求，可以用于地下水埋深變化預測。

５地下水埋深預測

通過模型模擬，對沈烏灌域在３．０６億ｍ３總引水量下的地下水埋深進行了預測分析。２０２５年沈烏灌域地下水埋深預測結果見圖６。由圖６可以看出，２０２５年地下水埋深在３．４ｍ以上的區域明顯擴大，由２０２０年的４．８６萬ｈｍ２增加到了９．３４萬ｈｍ２。由圖７可以看出，到２０３０年，沈烏灌域地下水埋深持續增加，大于３．４ｍ的區域持續擴大，僅在東南部、北部和西部局部范圍內存在小于３．４ｍ的區域。大于３．４ｍ的區域由２０２５年的９．３４萬ｈｍ２持續增加到了２０３０年的１４．５６萬ｈｍ２，占整個灌域面積（１８．６８萬ｈｍ２）的７８％。

６結果分析討論

通過上述分析可知，在保持３．０６億ｍ３引水量的情況下，沈烏灌域的地下水位會持續下降。預測結果顯示，到２０２５年將有５０％的區域地下水埋深超過３．４ｍ，到２０３０年將有７８％的區域地下水埋深超過３．４ｍ。在地下水埋深超過３．４ｍ的區域，天然植被狀況將會受到明顯影響，主要體現在植被覆蓋度下降。

７結論

基于ＭＩＫＥＳＨＥ分布式水文模型，采用實測水文數據和基于遙感的ＮＤＶＩ數據，運用數值模擬和統計分析等手段，對沈烏灌域出讓水權后的地下水埋深和ＮＤＶＩ進行了預測分析，得到的結論如下。

當地下水埋深小于３．４ｍ時，ＮＤＶＩ與土壤質地存在一定的相關性；當地下水埋深大于３．４ｍ時，地下水埋深與ＮＤＶＩ具有顯著負相關性，在其他因素不變的條件下，地下水埋深每增加１ｍ，土質基本單元的ＮＤＶＩ平均值會降低０．０８０７。

根據模型預測結果，相較于２０２０年，２０２５年地下水埋深在３．４ｍ以上的區域將增加到９．３４萬ｈｍ２；到２０３０年，沈烏灌域地下水埋深將持續增加，大于３．４ｍ的區域將增加到１４．５６萬ｈｍ２，僅在東南部（引黃總渠首）、北部（地下水排泄區）和西部（山洪補給區）局部范圍內存在小于３．４ｍ的區域。針對天然植被可能存在的狀況，灌域管理單位要加強灌域歷年凌汛期以及灌溉間隙的生態補水調度，針對地下水埋深增加區域的補水水體進行補水級別劃分，以減緩重點區域地下水埋深增加趨勢。

農業作為黃河流域的用水大戶，隨著黃河流域水資源短缺與經濟社會和生態環境發展矛盾的加劇，水權交易作為推動區域水資源優化配置的有力措施，需要在注重經濟社會效益的同時，更加重視生態環境的響應，全面、系統地做好水權交易引起的水資源重新配置后的外部效應將是今后一個重要研究方向。

本文的不足之處在于地下水觀測數據序列不長，對灌域內天然植被種類的變化考慮不足。下一步將在灌域地下水埋深變化頻繁的區域內設置典型區域生態樣方，針對天然植被覆蓋度變化以及植被種類更替規律等進行系統監測研究。