2005—2030 年無定河流域土地利用變化對徑流的影響

趙雪巖，張 鑫，孫 媛

（西北農(nóng)林科技大學(xué) 水利與建筑工程學(xué)院，陜西 楊凌 712100）

隨著社會發(fā)展，土地利用變化通過調(diào)節(jié)產(chǎn)匯流和降水蒸發(fā)過程來影響水文循環(huán)和水資源量的變化［1－3］。自1960 年以來，無定河流域受人類活動（如退耕還林還草、防風(fēng)治沙、開荒種田等）的影響，其水文循環(huán)發(fā)生了顯著變化［4－7］。前人大多聚焦于無定河流域植被恢復(fù)前后水文通量的大幅度變化進(jìn)行研究，如王計(jì)平等［4］認(rèn)為1995 年前后無定河流域土地流轉(zhuǎn)頻繁；任宗萍等［5］認(rèn)為水土治理措施使無定河支流大理河徑流量在1996 年出現(xiàn)拐點(diǎn)；張守紅等［6］認(rèn)為水保措施是無定河流域徑流量減少的主要原因。然而，水保措施實(shí)施后徑流、生態(tài)的大幅變化趨勢是不可持續(xù)的，國家山水林田湖草沙等自然資源共同發(fā)展演化的趨勢是從大修大補(bǔ)的退耕還林（草）轉(zhuǎn)向小微型、質(zhì)量型發(fā)展［8］。

SWAT（Soil and Water Assessment Tool）模型是一個基于物理機(jī)制的、復(fù)雜的半分布式水文模型［9－16］，該模型能夠模擬流域水文循環(huán)過程和土地利用變化過程等地表自然物理化學(xué)過程，并在國內(nèi)外得到了廣泛應(yīng)用［17－20］。Yang 等［9］認(rèn)為灤河流域退耕還林使夏季流量增大，其他季節(jié)流量減小；傅春等［10］研究表明撫河流域林地增加、耕地減少使徑流量減少；楊倩等［11］研究表明漢江上游退耕還林使年徑流量減少；Ma 等［12］研究表明岷江流域耕地增加、草地減少使產(chǎn)水量減少；李爍陽等［13－16］在渭河流域和汾河上游均開展不同角度的研究。

本文以無定河流域?yàn)檠芯繉ο?，運(yùn)用SWAT 模型，利用五期土地利用數(shù)據(jù)（2005 年、2015 年、2020年、2025 年、2030 年，其中2025 年和2030 年為CAMarkov 模型預(yù)測結(jié)果）驅(qū)動模型，分析4 種土地利用變化及其對徑流的影響。

1 研究區(qū)概況與數(shù)據(jù)來源

1.1 研究區(qū)概況

無定河流域位于北緯37°2′—39°1′、東經(jīng)107°47′—110°34′（見圖1），多年平均徑流量為11.51 億m3。流域處在東部季風(fēng)氣候和西北干旱氣候的過渡帶，多年平均氣溫為7.8～9.6 ℃，多年平均降水量為300～500 mm，降水量由東南向西北遞減，年際變化大，降水多集中在夏秋兩季，年內(nèi)分配不均［4－7］。

圖1 無定河流域高程及水文站和氣象站分布

1.2 數(shù)據(jù)來源與處理

1.2.1 空間數(shù)據(jù)處理

空間數(shù)據(jù)介紹詳細(xì)見表1。其中，2005 年、2015年、2020 年三期土地利用數(shù)據(jù)用于搭建SWAT 模型的土地利用數(shù)據(jù)庫，原始數(shù)據(jù)中的二級分類需要重分類為草地、林地、耕地、未利用地、居民用地和水域6 種類型（見表2，其中土地利用分類屬性前的數(shù)字為代號）；2010 年土地利用數(shù)據(jù)用于構(gòu)建CA－Markov 模型中的土地轉(zhuǎn)移概率矩陣。HWSD 用于搭建土壤屬性數(shù)據(jù)庫，同時還需輸入默認(rèn)參數(shù)值，如上層鹽度默認(rèn)0.1等；依托HWSD，利用SPAW 軟件進(jìn)行計(jì)算后獲取其他參數(shù)。另外，石灰性雛形土和石灰性沖積土根據(jù)屬性的不同進(jìn)一步細(xì)分為1 類和2 類。土壤類型和面積統(tǒng)計(jì)詳見表3 和圖2（空間參考坐標(biāo)系統(tǒng)一設(shè)置為WGS 1984－utm－48n）。

表1 空間數(shù)據(jù)介紹

表2 土地利用數(shù)據(jù)重分類

表3 不同類型土壤的面積統(tǒng)計(jì)

圖2 無定河流域土壤分類

1.2.2 氣象水文數(shù)據(jù)處理

氣象數(shù)據(jù)包括逐日降水量、最高／最低氣溫、太陽輻射、風(fēng)速、相對濕度等，來源于中國氣象數(shù)據(jù)網(wǎng)，氣象站分布見圖1。水文數(shù)據(jù)為逐日平均流量，來源于黃河流域水文年鑒（中游區(qū)上段），水文站分布見圖1。兩者時間段均為2006—2018 年。另外，各類氣象數(shù)據(jù)需要按照站編號、站名、經(jīng)緯度、高程做成索引表（txt文件），缺測數(shù)據(jù)使用天氣發(fā)生器模擬進(jìn)行插補(bǔ)。

2 研究方法

2.1 土地利用數(shù)據(jù)模擬

CA－Markov 模型結(jié)合了CA 模型模擬復(fù)雜系統(tǒng)空間變化的能力和Markov 模型的定量預(yù)測優(yōu)勢，可作為土地利用模擬工具，并在IDRISI 軟件中實(shí)現(xiàn)，其實(shí)質(zhì)是建立土地轉(zhuǎn)移概率矩陣并利用卡帕系數(shù)Kappa評估概率矩陣模擬出的土地利用數(shù)據(jù)和實(shí)測數(shù)據(jù)的差異，從而建立模型用于輸出預(yù)測數(shù)據(jù)［21－23］。

式中：PO為恰當(dāng)被模擬的部分；PC為由偶然因素導(dǎo)致的被希望糾正的部分；PP為觀測和模擬完全匹配的理想比例，取為1；Lij為i類土地利用類型轉(zhuǎn)變?yōu)閖類土地利用類型的概率；n為土地利用類型總數(shù)。

如果Kappa ＝1，說明二者是一致的；如果0.75 ≤Kappa ＜1，說明高度一致；如果0.5 ≤Kappa ＜0.75，說明中度一致；如果Kappa＜0.5，說明不一致［21］。

本文輸出2025 年和2030 年的無定河流域土地利用數(shù)據(jù)的主要步驟［21－22］：①利用Markov 模型建立時間間隔為5 a（2005—2010 年）的過渡概率矩陣；②設(shè)定預(yù)測的起始年份為2010 年，并利用CA 模型設(shè)置循環(huán)迭代次數(shù)為5 次，分別用來模擬2015 年土地利用數(shù)據(jù)。③模擬2015 年土地利用數(shù)據(jù)與實(shí)際2015 年土地利用數(shù)據(jù)進(jìn)行比較，以Kappa來評估模擬效果。

2.2 SWAT 模型

在SWAT 模型中，一個流域被劃分為由河流網(wǎng)絡(luò)相互連接的許多子流域，每個子流域根據(jù)土地利用類型、土壤類型和坡度進(jìn)一步劃分為水文響應(yīng)單元（HRUs），模擬水的流動軌跡從HRUs 到子流域再到流域出口斷面。水文過程基于如下水量平衡方程［18－20］：

式中：SWt為最終土壤含水量，mm；SW為初始土壤含水量，mm；R、Qsurf、ET、P、QR分別為滲漏量、徑流量、蒸散發(fā)量、日降雨量、地下水含量，mm。

常用Nash－Suttcliffe 系數(shù)（Ens）、相關(guān)系數(shù)（r）［18］和百分比偏差（PBIAS）［24］來評估SWAT 模型的模擬效果。

式中：Qoi為i時的觀測徑流量；Qsi為i時的模擬徑流量；為觀測徑流量平均值；為模擬徑流量平均值；n為觀測數(shù)據(jù)總數(shù)［24］。

當(dāng)r ＝1 時，結(jié)果非常吻合；r≥0.5 時，模擬結(jié)果可以接受。當(dāng)Ens＝1 時，模擬結(jié)果和實(shí)測值完全吻合；0.5＜Ens＜1 時，模擬結(jié)果可以接受。當(dāng)PBIAS ＝0 時，模擬結(jié)果非常好；PBIAS ＜0 時，模擬結(jié)果偏大；PBIAS ＞0 時，模擬結(jié)果偏?。?4］。

年內(nèi)分配完全調(diào)節(jié)系數(shù)Cr可表示徑流量年內(nèi)分配不均勻程度，值越大不均勻程度越大［25］。

式中：t為月份；R（t）為t月的徑流量；為多年平均徑流量。

3 結(jié)果與分析

3.1 無定河流域土地利用變化分析

CA－Markov 模型運(yùn)行結(jié)果表明：模擬的2015 年土地利用數(shù)據(jù)和實(shí)測數(shù)據(jù)的Kappa系數(shù)為0.930 7，說明二者高度一致。所以用2005—2010 年土地轉(zhuǎn)移概率矩陣和2010 年土地利用數(shù)據(jù)分別迭代15、20 次來模擬2025 年、2030 年土地利用數(shù)據(jù)。

本文用于土地利用分析的數(shù)據(jù)分辨率為1 km。由表4 可知：①無定河流域的土地利用類型主要是草地、耕地和未利用地，其中草地占流域總面積的40%左右，耕地占30%左右，未利用地占22%左右，其余為林地、水域和居民用地，共占8%左右；②草地面積先增加后減小，耕地先增加后減小而后保持維定，未利用地先減少后增加，居民用地持續(xù)擴(kuò)張且偶有波動；③水域面積變化很小，說明流域25 a 內(nèi)的水庫、河流、灘地基本保持不變。由圖3 可知：①草地主要集中在流域中上游地區(qū)，并在下游的耕地中有少量零星草地；②耕地主要集中在流域中下游；③未利用地大多分布在中上游地區(qū)，呈點(diǎn)狀分布，并且2030 年未利用地在東北部大塊狀出現(xiàn)，應(yīng)引起關(guān)注；④居民用地增加主要集中在河道附近和流域東北部。

表4 不同時期土地類型面積統(tǒng)計(jì)

圖3 各土地利用類型分布情況

本文涉及4 種土地利用變化方式，見表5（記為C1、C2、C3、C4）。為探究C1～C4 的規(guī)律，建立4 個土地轉(zhuǎn)移矩陣（見表6～表9）。綜合圖4、表6～表9 可知：①居民用地雖然在流域內(nèi)占比小，但是增長幅度最大，C1～C4 分別增長110%、220%、120%、210%；②C1～C4的轉(zhuǎn)移矩陣中，互相轉(zhuǎn)化的土地利用類型的面積遠(yuǎn)小于沒有發(fā)生轉(zhuǎn)化的面積，說明2005—2030 年土地流轉(zhuǎn)僅在小范圍內(nèi)進(jìn)行；③C1，未利用地面積的1.2%和草地面積的0.6%共同轉(zhuǎn)向耕地，使耕地增加1.3%；④C2，4.7%未利用地和4.3%耕地轉(zhuǎn)向草地，使草地增加3%，同時1.27%的草地轉(zhuǎn)出為居民用地；⑤C3，未利用地減少10.4%，轉(zhuǎn)出為耕地、居民用地和草地；⑥C4，1.1%草地和2.1%未利用地轉(zhuǎn)向耕地，使耕地增加2.9%，與C1 相比變化幅度更大，同時4.3%的草地轉(zhuǎn)為耕地、未利用地和居民用地。

表5 不同土地利用變化方式

表6 C1 土地利用轉(zhuǎn)移矩陣 km2

表7 C2 土地利用轉(zhuǎn)移矩陣 km2

表8 C3 土地利用轉(zhuǎn)移矩陣 km2

表9 C4 土地利用轉(zhuǎn)移矩陣 km2

3.2 SWAT 模型校準(zhǔn)和驗(yàn)證

根據(jù)本文研究需要，并結(jié)合前人研究經(jīng)驗(yàn)［26－30］，采用1 km 的土壤數(shù)據(jù)集和土地利用數(shù)據(jù)集構(gòu)建模型，共劃分18 個子流域（見圖5）和372 個HRUs。把2006—2013 年設(shè)定為模型率定期，把2014—2018 年設(shè)定為模型驗(yàn)證期。2005 年土地利用數(shù)據(jù)用于模型初始構(gòu)建，利用白家川水文站還原后的逐月平均流量（以下稱為還原值）進(jìn)行校準(zhǔn)。同時使用SWAT－CUP自動校準(zhǔn)和驗(yàn)證參數(shù)，見表10。由圖6 可知：①率定期和驗(yàn)證期的多年月平均流量模擬值和還原值的變化趨勢有較好的一致性，擬合的線率分別為0.913、0.853 3；②率定期r、Ens、PBIAS分別為0.85、0.84、5.4%，驗(yàn)證期r、Ens、PBIAS分別為0.89、0.88、5.0%。綜上，在無定河流域建立的SWAT 模型有良好適用性，可以用來研究C1～C4 對徑流的影響。

圖5 子流域分布

表10 參數(shù)率定結(jié)果

圖6 模型率定期和驗(yàn)證期流量的模擬值和還原值

3.3 分析土地利用變化對徑流影響的試驗(yàn)設(shè)計(jì)

氣象數(shù)據(jù)不變，分析C1～C4 對徑流的影響，具體試驗(yàn)設(shè)計(jì)見表11。其中，定義Q1～Q5 為情景1～5 下得到的徑流模擬值的時間序列，R1～R5 為由流量時間序列得到的徑流深，ΔR1～ΔR4 為C1～C4 引起的徑流深變化。

表11 試驗(yàn)設(shè)計(jì)

3.4 不同土地利用變化方式對徑流的影響分析

由圖7 可知：①ΔR1～ΔR4 的變化趨勢基本一致，且C1～C4 都會導(dǎo)致流域徑流量增加，徑流深多年總增加量從大到小排序?yàn)棣2＞ΔR4＞ΔR3＞ΔR1，說明C1～C4 均對流域徑流起促進(jìn)作用，雖然C1～C4 中沒有進(jìn)行流轉(zhuǎn)的土地占各自土地利用類型面積的88%以上（見圖4、表6～表9），且草地、耕地和未利用地面積各有增減，對徑流有不同程度的促進(jìn)或抑制作用，但是居民用地面積均增加，說明居民用地面積增加對徑流的促進(jìn)作用占主導(dǎo)地位，城市化程度高的流域，不透水地面的面積大，下滲量和下滲強(qiáng)度降低，屬于超滲產(chǎn)流，且居民用地主要分布在河道兩側(cè)，從而使匯流時間變短，徑流量增大，這與楊滿根［31］在淮河流域得出的結(jié)論一致。②ΔR2 多年變化總量達(dá)9.54 mm，共計(jì)2.74 億m3，占流域內(nèi)現(xiàn)有水量的23.81%，說明C2 對徑流的促進(jìn)作用最大。變化量最小的是ΔR1，多年變化總量達(dá)3.97 mm，C1 中居民用地面積增加1.1 倍，說明C1 對徑流的促進(jìn)作用最小。③ΔR3 整體上略大于ΔR1，變化總量為4.40 mm；ΔR4 整體上略小于ΔR2，變化總量為8.40 mm，針對C4 使ΔR4 增加較大的情況，應(yīng)引起關(guān)注。

圖7 土地利用變化引起的徑流深變化

為進(jìn)一步探究ΔR的年內(nèi)變化情況，本文從多年月平均徑流深和年內(nèi)分配完全調(diào)節(jié)系數(shù)（Cr）兩方面入手。對R1～R5 中每年相同月份求平均值，而后分別和情景1 做差比較，得ΔR的年內(nèi)變化，由圖8 可知：①ΔR1～ΔR3 中，多年月平均徑流深的年內(nèi)增加趨勢一致，呈先減小后增大而后減小的趨勢。多年月平均徑流深的增大量排序?yàn)棣2＞ΔR3＞ΔR1＞ΔR4，ΔR2 的年內(nèi)變化總量達(dá)0.92 mm，說明C1～C3 對ΔR的各月流量均起促進(jìn)作用，且C2 對各月流量的促進(jìn)作用最顯著，原因是C2 草地增加最多，由耕地和未利用地轉(zhuǎn)入，使地表粗糙度增大，這與孫鐵軍等［17］的結(jié)論一致。且經(jīng)過多年水保措施改造后土壤含水量高［17］，雨水落在草地上形成超滲產(chǎn)流；居民用地面積增加最多，形成不透水面積最大［31］，兩方面原因共同促使各月徑流深增加量最多。②ΔR4 中，1—6 月徑流深減小總量共計(jì)0.26 mm，7—12 月徑流深增大總量共計(jì)0.68 mm，說明C4 對上半年徑流起抑制作用，對下半年徑流起促進(jìn)作用，且促進(jìn)程度比抑制程度更高，原因可能是C4 未利用地和草地向耕地轉(zhuǎn)移，使耕地面積大量增加，農(nóng)作物在春天需要大量使用河水灌溉［7］，河水從河道進(jìn)入土壤和地下，從而起到抑制徑流的作用；下半年降水量增加，灌溉用水量少，從而起到促進(jìn)徑流的作用，針對草地和未利用地大量向耕地轉(zhuǎn)移，導(dǎo)致上半年各月徑流減少的情況應(yīng)引起關(guān)注。③從年內(nèi)分布的整體來看，ΔR1～ΔR4 的最大月份均在8 月，最枯月份均在3 月；C1～C4 均使7—12 月ΔR遠(yuǎn)大于1—6 月的，降水是本流域水資源的主要補(bǔ)給方式，且多集中在夏季和秋季，河流進(jìn)入汛期，C1～C4 對汛期徑流的影響更大。

圖8 C1～C4 引起的徑流深變化（ΔR）的年內(nèi)分布

依據(jù)Q1～Q5 中每月Q值按式（7）和式（8）計(jì)算每年的Cr值，而后分別和情景1 的Cr值做差得ΔCr，繪制2006—2018 年ΔCr分布圖。由圖9 可知：①2006—2018 年的ΔCr1～ΔCr4 均大于0，變化幅度在0.14 之間，說明一定程度上C1～C4 均使徑流深的年內(nèi)分配不均勻，原因是C1～C4 各種土地類型面積的12%左右發(fā)生流轉(zhuǎn)，造成下墊面在小范圍內(nèi)變化。②2006—2018 年的ΔCr1、ΔCr3 值均小于ΔCr2、ΔCr4，說明C1、C3 使徑流深年內(nèi)分配相對來說更均勻些，雖然一年內(nèi)大部分降水集中在夏秋兩季，使河水上漲，但是C1、C3 使耕地面積不同程度增加，無定河流域地處陜北，多種玉米和馬鈴薯［7］，一年內(nèi)均需消耗不同程度的灌溉用水量，因此使徑流深年內(nèi)分配不均勻程度較小。③ΔCr1～ΔCr4 中，2014—2018 年的差值整體均呈現(xiàn)上升趨勢，ΔCr2 的上升趨勢最為顯著，說明C2 對2014—2018 年的徑流深年內(nèi)分配不均勻程度影響最大。

圖9 C1～C4 引起的年內(nèi)分配不均勻系數(shù)的分布

4 結(jié) 論

為研究無定河流域2005 年后的土地利用的小幅度變化對徑流的影響，采用2005 年、2015 年、2020 年三期土地利用數(shù)據(jù)和基于CA－Markov 模型預(yù)測的2025 年、2030 年兩期數(shù)據(jù)，分析4 種土地利用變化方式，并構(gòu)建SWAT 模型，研究不同土地利用變化方式對徑流的影響。

（1）2005 年后無定河流域草地、耕地和未利用地占總流域面積的90%以上。2005—2030 年耕地、草地在小范圍內(nèi)互相轉(zhuǎn)化；未利用地持續(xù)減少，主要轉(zhuǎn)化為耕地和草地；居民用地占流域面積較少，但擴(kuò)張極其迅速。

（2）基于2005 年土地利用數(shù)據(jù)構(gòu)建的SWAT 模型，經(jīng)過參數(shù)率定后，在率定期和驗(yàn)證期的r和Ens均在0.8 以上，說明對無定河流域水文過程的模擬效果良好，適用于研究本流域土地利用變化對徑流的影響。

（3）4 種土地利用變化均使本流域2006—2018 年徑流量有不同程度增加，且居民用地面積增加在促進(jìn)徑流增加中占主導(dǎo)作用。

（4）4 種土地利用變化方式對7—12 月徑流的促進(jìn)作用均遠(yuǎn)大于1—6 月；草地增加地表粗糙度，居民用地增加不透水面積，二者是使ΔR2 最大的主要原因；第4 種土地利用變化對上半年徑流起抑制作用，對下半年徑流起促進(jìn)作用。

（5）4 種土地利用變化均使2006—2018 年的徑流年內(nèi)分配不均勻，第1 和第3 種土地利用變化使徑流年內(nèi)分配相較于其他兩種變化來說均勻些。

（6）第4 種土地利用變化雖然未利用地總面積減少，卻在東北部大塊出現(xiàn)，并且引起徑流年際增加量大、上半年徑流減少下半年徑流增加、歷年徑流深的年內(nèi)分配極不均勻，該情況應(yīng)引起關(guān)注。

本文利用土地轉(zhuǎn)移概率矩陣進(jìn)行未來土地利用預(yù)測，沒有考慮政策導(dǎo)向和氣候變化的影響，同時所采用土地利用數(shù)據(jù)的分辨率也有待于進(jìn)一步提高，進(jìn)而提高成果的精度，這將是下一步研究的重點(diǎn)。