王藝璇,劉 夏,沈彥軍**

(1.中國科學院遺傳與發(fā)育生物學研究所農業(yè)資源研究中心/中國科學院農業(yè)水資源重點實驗室 石家莊 050022;2.中國科學院大學 北京 100049)

徑流是地表水資源量的重要表現(xiàn)形式,在水循環(huán)過程中發(fā)揮著重要作用。受氣候變化和人類活動的影響,全球和區(qū)域尺度的徑流發(fā)生了巨大變化,徑流過程模擬及變化歸因分析研究也成為國際水科學領域的研究熱點。永定河是北京的“母親河”,其上游流域不僅是重要的水源地,更是京津冀協(xié)同發(fā)展的生態(tài)大動脈,永定河流域的綜合治理和生態(tài)修復作為國家戰(zhàn)略被提上議程。過去幾十年來,受氣候變化和人類活動等因素的影響,永定河上游流域徑流量呈顯著減少趨勢,水資源短缺成為限制流域社會經濟發(fā)展和生態(tài)文明建設的重要因素。因此,定量評估氣候變化和人類活動對徑流減少的貢獻,對深入理解永定河上游流域水文過程和水資源的科學管理具有重要意義。

目前,徑流變化歸因分析的研究已經從定性分析逐漸發(fā)展到了定量識別階段,研究方法主要包括基于數(shù)理統(tǒng)計的分析法、基于物理機制的水文模型法以及基于Budyko 假設和氣候彈性系數(shù)法。統(tǒng)計分析法計算簡便,但整體上缺乏物理機制; 水文模型法物理機制較強,但水文模型需要大量的數(shù)據(jù)輸入、參數(shù)率定過程復雜、易出現(xiàn)過參數(shù)化或異參同效等問題,在缺資料地區(qū)的應用會有較大的局限性; 彈性系數(shù)法需要長期的水文氣象資料,能考慮簡單的物理現(xiàn)實背景,但實際蒸散發(fā)量難以測量,計算實際蒸散發(fā)的經驗方法也具有一定的不確定性。近年來,受全球大數(shù)據(jù)技術迅猛發(fā)展的影響,機器學習應用成為研究熱點。與傳統(tǒng)的方法相比,機器學習方法能夠挖掘大數(shù)據(jù)的潛在價值及其內在聯(lián)系,在流域徑流預測方面具有一定的實用價值,如集成學習方法、人工神經網(wǎng)絡方法、支持向量機等方法在徑流預測方面也得到了廣泛應用。隨機森林(Random Forest)模型作為一種數(shù)據(jù)挖掘的機器學習方法,自2005年Parkhurst 等首次將其應用于水資源科學解決水質問題以來,隨機森林已被廣泛用于各種水文應用,主要包括: 徑流模擬及預報、水位預測、土壤水勢測定及參考蒸散估算等,然而隨機森林模型是否能夠考慮各影響因素與徑流之間的非線性關系來實現(xiàn)流域徑流模擬,且其是否可以應用于徑流變化歸因分析,相關研究還相對較少。關于永定河上游徑流變化的歸因分析當前已有部分學者開展了研究,結果表明氣候變化和人類活動對徑流減少的貢獻率比例分別為10.5%～39.1%和60.3%～89.5%。然而該流域水環(huán)境復雜,且受到了高強度的人類活動影響,研究也多是基于統(tǒng)計學和彈性系數(shù)法,所得氣候變化與人類活動對徑流變化貢獻率的定量分析結果變幅較大,相關結論仍具有較大的不確定性。

本研究以永定河上游典型流域為研究對象,分別以隨機森林模型和SWAT 模型為基礎對流域徑流變化歸因進行定量分析,并將兩種方法的結果進行對比驗證,旨在探索隨機森林模型在徑流變化歸因分析中的適用性,進而為更可靠地評估變化環(huán)境下的流域徑流變化歸因及水資源合理開發(fā)利用提供科學參考。

1 研究區(qū)概況、數(shù)據(jù)來源及方法

1.1 研究區(qū)概況

永定河屬海河水系,發(fā)源于干旱半干旱區(qū),河流全長約747 km,流域面積約為4.7 萬km,流經山西、內蒙古、河北、北京和天津5 省市(112°～117°45′E,39°～41°20′N)。以官廳水庫和三家店水庫為界,分為上、中、下游3 段。永定河上游流域主要包括洋河和桑干河兩大支流(圖1),流域土地利用整體上以耕地、草地和林地為主,三者占流域總面積的90%以上; 流域屬典型的溫帶大陸性季風氣候,年平均氣溫4.6 ℃,年降水量約400 mm,80%徑流量來源于夏季降水補給。

圖1 洋河流域示意圖Fig.1 Schematic diagram of the Yanghe River Basin

本文主要關注永定河上游受人類活動影響較大的洋河流域。洋河地處張家口壩下南部,是永定河上游的第一大支流,全長262 km,流域面積為15 078 km2,山西省的陽高縣和內蒙古的興和縣是河流的發(fā)源地,流經河北省萬全縣、張家口市、宣化縣、下花園、懷來縣等。上源有東洋河、南洋河、西洋河等支流,并于河北省懷安縣柴溝堡附近匯合稱為洋河。

1.2 數(shù)據(jù)來源

本文選取研究區(qū)及附近6 個氣象站(圖1,表1)的逐日平均氣溫、日最高氣溫、日最低氣溫、風速、日照時數(shù)、相對濕度等作為基礎數(shù)據(jù),數(shù)據(jù)來源于中國氣象數(shù)據(jù)共享服務網(wǎng)(http://data.cma.cn)。部分站點資料略有缺失,缺失數(shù)據(jù)基于相鄰站點采用線性插值法補齊。隨機森林模型中輸入的參考作物蒸散發(fā)由Penman-Monteith 公式計算所得。為定量探究氣候變化和人類活動對永定河上游典型流域徑流變化的影響貢獻,基于洋河流域(柴溝堡東站、柴溝堡南站、響水堡站)的實測徑流數(shù)據(jù)(圖1,表1)進行分析,數(shù)據(jù)來源于《中華人民共和國水文年鑒》。研究所用到的2000年的土地利用數(shù)據(jù)來源于中國科學院地理科學與資源研究所,分辨率為30 m。數(shù)字高程模型(DEM)來自于地理空間數(shù)據(jù)云(http://www.gscloud.cn)。土壤類型和質地數(shù)據(jù)來自于中國1∶100萬土壤數(shù)據(jù)庫及野外實地考察。

表1 洋河流域水文站、氣象站基本信息Table 1 Basic information of hydrological stations and weather stations in the Yanghe River Basin

1.3 研究方法

本研究對永定河上游典型流域的徑流進行模擬及變化歸因分析,主要包括以下3 個步驟: 1) MK 突變分析法識別徑流的突變點; 2)運用SWAT 模型和隨機森林模型對突變點前的徑流量進行模擬校準,利用校準好的模型模擬突變后的徑流量; 3)借助徑流影響評價模型進行徑流變化歸因分析計算。

1.3.1 Mann-Kendall 突變檢驗法

基于Mann-Kendall 非參數(shù)方法對水文序列進行突變檢驗。通過構建序列來計算的不同檢驗統(tǒng)計量S,前提是時間序列是無序的:

式中:表示為小于等于但不等于的變量;a代表統(tǒng)計量;1 時UF=0;(S)和var(S)分別代表樣本的均值和方差; UF符合準正態(tài)分布,給定一個置信水平,若UF＞則序列有明顯的上升或者下降趨勢,UB則是UF的相反數(shù)。UF線與UB線在置信區(qū)間的交點為突變點,其對應的時間則為突變時間。

SWAT(Soil and Water Assessment Tool)模型是20世紀90年代由美國農業(yè)部(USDA)開發(fā)的分布式流域水文模型,用于模擬和預測流域長期的徑流、泥沙和面源污染等物理過程,是當前應用較為廣泛的水文模型。目前,已有大量研究證明SWAT 可以較好地應用于中等尺度流域徑流模擬的研究。SWAT模型構建所需數(shù)據(jù)包括氣象、土地利用屬性、土壤特性屬性數(shù)據(jù)和地形、土壤類型、土地利用等空間數(shù)據(jù)。首先基于DEM 和流域水系分布特點,以徑流觀測站為出口,將流域劃分成若干個子流域且包含1 個或多個水文響應單元(HRU); 再次,輸入模型需要的土壤、土地利用屬性數(shù)據(jù),最后輸入氣象數(shù)據(jù)來驅動模型。SWAT 基本原理是在每個水文響應單元或子流域上計算凈雨和產流,再根據(jù)流域特性,基于數(shù)學方程實現(xiàn)匯流計算和徑流的時間變化模擬。模型中的水量平衡表達式為:

式中: SW表示土壤最終的含水量(mm),SW表示土壤前期的含水量(mm),表示在第天的降水量(mm),表示在第天的地表徑流量(mm),則表示為第天的蒸發(fā)量(mm),表示第天時土壤剖面底層的滲透量和側向流量(mm),表示第天地下水流出量(mm)。

本文利用SWAT-CUP 軟件進行模型參數(shù)的率定和驗證,并選用決定系數(shù)()和納什效率系數(shù)(NSE)來評價模型效果。和NSE 越接近1,表示模型效果越好。

隨機森林模型是一種現(xiàn)代分類與回歸技術,也是一種組合式學習技術,其自助樣本集是由變量和屬性隨機化選擇來構建的,分類回歸決策樹是個體學習機的集成學習模型。存在噪聲和缺失值的實測數(shù)據(jù)用隨機森林模型處理效果好且學習速度快,隨機森林模型也很適合處理以小樣本和高維度為特征數(shù)據(jù)的分類和回歸問題,近年來在解決特征分類、時間序列預測、特征提取以及異常點檢測等問題上得到了廣泛應用。

本研究使用的是R 中的“Party”包開發(fā)的隨機森林模型,模型的輸入因子包括降水、氣溫、蒸發(fā)量、日照時數(shù)、風速、相對濕度、月徑流數(shù)據(jù)及土地利用數(shù)據(jù)。主要包括以下3 個步驟: 1)突變分析; 2)徑流模擬; 3)借助徑流影響評價模型進行貢獻率分析。與線性回歸或雙累積曲線法類似,利用隨機森林模型建立了突變點前氣候或人為因素與徑流之間的關系。

突變點前的氣候或人為因素與徑流的關系可以表示為:

其次,氣候驅動的徑流量可由突變前建立的關系確定:

式中:、、、、及分別代表突變前的降水(mm)、溫度(℃)、蒸散發(fā)(mm)、風速(m·s)、日照時數(shù)(h)及相對濕度(%),代表突變前的徑流量(m·s),、、、、及分別代表突變后的降水(mm)、溫度(℃)、蒸散發(fā)(mm)、風速(m·s)、日照時數(shù)(h)及相對濕度(%),表示突變后氣候驅動的徑流(m·s)。

再次,人類活動驅動的徑流變化由突變前的關系確定。

氣候變化和人類活動對徑流變化的貢獻率計算參考式(9)和(10):

式中:表示突變前的多年平均實測徑流(m3s),表示突變后的多年平均實測徑流(m3s);表示氣候變化對徑流變化的貢獻率,表示人類活動對徑流變化的貢獻率。

通過流域水文模型模擬分析變化環(huán)境對徑流量影響的前提是需要假定氣候變化和人類活動這兩者為相互獨立的因子。此種方法主要包括以下幾個步驟: 1)對長時間序列的實測徑流數(shù)據(jù)進行突變年的識別,以突變點為分界將水文序列劃分為天然時期和人類影響時期,并將天然時期認作基準期。2)以突變點作為分界點,利用SWAT 模型或隨機森林模型對突變后的徑流量進行模擬,將模擬出來的徑流量假定為突變后無人類活動影響的天然徑流量。3)通過天然時期的徑流量、人類活動時期的實測徑流量與天然徑流量計算氣候變化和人類活動對徑流變化的貢獻率。氣候變化和人類活動徑流貢獻率的計算公式如下:

式中:是突變前實測徑流量的平均徑流量(10m),是突變后實測徑流量的平均值(10m),是突變后模擬徑流量的平均值(10m),Δ是氣候變化對徑流的影響值(10m),Δ是人類活動對徑流的影響值(10m),Δ為徑流變化總量(10m),、分別代表氣候變化和人類活動對徑流變化影響的百分比。

2 結果與分析

2.1 基于SWAT 和隨機森林模型的洋河流域徑流過程模擬

利用SWAT 模型對永定河上游受人類活動影響較大的洋河流域進行模擬,研究時段為1984-2015年。選取1984-1985年作為模型的預熱年,1986-1990年為率定期,1991年至突變年的前一年作為驗證期。東洋河(柴溝堡東站)、南洋河(柴溝堡南站)、洋河(響水堡站)徑流的突變年分別發(fā)生在2000年、1998年和1997年,因此,東洋河(柴溝堡東站)的驗證期為1991-1999年,南洋河(柴溝堡南站)的驗證期為1991-1997年,洋河(響水堡站)的驗證期為1991-1996年。根據(jù)從上游至下游,從支流到干流的原則,首先通過柴溝堡東和柴溝堡南兩個水文站的月徑流數(shù)據(jù)對東洋河和南洋河兩大支流進行率定和驗證,再將率定好的參數(shù)進行參數(shù)傳遞,輸入SWAT 模型中對響水堡水文站的月徑流數(shù)據(jù)進行率定和驗證。

經參數(shù)率定后的洋河流域子流域各水文站實測徑流值與模型模擬徑流值如圖2 所示,東洋河、南洋河以及洋河流域的模擬與實測月徑流過程較吻合,3 個流域率定期和驗證期的均在0.65 以上,NSE也大多在0.65 以上,只有洋河流域的NSE 在率定期為0.57。SWAT 模型對部分年份枯水期徑流的模擬值偏低,且由于洋河流域位于干旱半干旱氣候區(qū),降雨多為短歷時的暴雨,產流方式以超滲產流為主,在月尺度上模擬的降雨徑流過程較難準確模擬洪峰。但總體來看,SWAT 模型對研究區(qū)域徑流過程模擬精度較好,說明可以采用SWAT 模型模擬洋河流域的天然徑流量變化過程。

圖2 基于SWAT 模型的東洋河、南洋河、洋河流域突變前月徑流實測值和模擬值對比Fig.2 Observed and simulated monthly runoff before abrupt change in the Dongyanghe River Basin,Nanyanghe River Basin,Yanghe River Basin based on SWAT

基于隨機森林模型的洋河流域子流域各站徑流模擬結果如圖3 所示。隨機森林模型在洋河流域徑流模擬結果上與實測值相對接近?；贜SE 和評價結果來看,東洋河徑流模擬的NSE 為0.77,為0.84; 南洋河徑流模擬的NSE 為0.81,為0.86; 洋河流域徑流模擬的NSE 為0.84,為0.91。模型結果表現(xiàn)相對優(yōu)異,說明隨機森林模型建立的突變前氣候和人類活動與徑流變化之間的關系適用于本研究區(qū)。

圖3 基于隨機森林模型的東洋河、南洋河、洋河流域突變前月徑流實測值和模擬值對比Fig.3 Observed and simulated monthly runoff before abrupt change in the Dongyanghe River Basin,Nanyanghe River Basin,Yanghe River Basin based on random forest model

2.2 基于徑流影響評價模型的洋河流域徑流變化歸因分析

基于洋河流域的每個子流域模擬的天然徑流量和實測徑流量,以突變前的實測徑流值為基準值,借助徑流影響評價模型定量評估了氣候變化和人類活動對東洋河、南洋河、洋河3 個流域徑流變化的貢獻率(圖4、表2)?？傮w來看,1986-2015年東洋河、南洋河、洋河流域徑流突變后的實測徑流量相較于SWAT 模型模擬的天然徑流量明顯減少,模擬徑流與突變前實測徑流量相比有一定程度的下降,但差異不大,人類活動是導致洋河流域徑流減少的主要原因。東洋河流域在突變后較突變前的實測年均徑流值減少0.5 億m3左右,人類活動使年均徑流減少約0.43 億m3,其對徑流下降的貢獻率達84.3%;南洋河流域徑流變化受人類活動的影響最大,突變后較突變前的年均實測徑流值減少了近0.3 億m3,人類活動使年均徑流減少約0.27 億m3,氣候變化和人類活動對徑流減少的貢獻率分別為5.0%和95.0%; 洋河流域徑流突變后較突變前的年均實測徑流值減少了1.2 億m3,人類活動使年均徑流減少1.1億m3,氣候變化和人類活動對徑流減少的貢獻率分別為8.9%和91.1%。

表2 基于SWAT 模型的氣候變化及人類活動對東洋河、南洋河、洋河流域徑流變化影響的貢獻率(1986—2015年)Table 2 Contribution rates of climate change and human activities to runoff changes in the Dongyanghe River Basin,Nanyanghe River Basin,Yanghe River Basin based on SWAT model(1986 to 2015)

圖4 基于SWAT 模型(a,b,c)和隨機森林模型(d,e,f)的東洋河、南洋河、洋河流域徑流影響評價模型結果Fig.4 Observed and simulated monthly runoff after calibration and validation periods based on the SWAT and random forest model in the Dongyanghe River Basin,Nanyanghe River Basin,Yanghe River Basin

同上,在隨機森林模型模擬徑流的基礎上,借助徑流影響評價模型定量評估氣候變化和人類活動對東洋河、南洋河、洋河3 個流域徑流變化的貢獻率,結果如圖4、表3所示。東洋河流域徑流突變后較突變前的實測年均徑流減少約0.51 億m,其中,人類活動使徑流減少約0.45 億m,氣候變化和人類活動貢獻率分別為11.5%和88.5%; 南洋河流域突變后較突變前的實測年均徑流減少約0.28 億m,人類活動使徑流減少約0.27 億m,其對徑流下降的貢獻率約為98%; 洋河流域突變后較突變前的實測年均徑流減少近1.2 億m,人類活動使徑流減少約1.15 億m,氣候變化和人類活動對徑流下降的貢獻率分別為3.4%和96.6%。

表3 基于隨機森林模型的氣候變化和人類活動對東洋河、南洋河、洋河流域徑流變化的貢獻率(1986—2015年)Table 3 Contribution rate of climate change and human activities to the runoff changes in the Dongyanghe River Basin,Nanyanghe River Basin,Yanghe River Basin based on the random forest model(1986-2015)

通過與上文中基于SWAT 模型計算得到的氣候變化和人類活動對東洋河、南洋河、洋河流域徑流量影響的貢獻率進行對比,基于隨機森林模型所得的結果在趨勢和數(shù)值上與SWAT 模型接近,表明隨機森林模型方法在洋河流域徑流變化歸因分析研究中具有一定的適用性。

3 討論

本研究嘗試利用隨機森林模型模擬永定河上游洋河流域徑流變化并量化氣候變化和人類活動對其影響?；谒哪Ｐ秃碗S機森林模型的研究結果相似,都表明人類活動是導致該流域徑流變化的主要驅動因素,對徑流減少的貢獻率為84.3%～97.6%。這一結果與李秀通過Budyko 假設的互補關系法得到的結果基本一致(92.0%)。除此之外,田菲等采用趨勢分析法和降水-徑流關系模型對1966-1999年洋河流域徑流變化歸因進行了分析,結果表明人類活動對洋河流域徑流下降的貢獻為83.2%;Xia 等和侯蕾基于Budyko 的方法對氣候變化及人類活動對永定河上游徑流影響進行了定量評價,研究結果表明人類活動是徑流量減少的主要原因,貢獻率為80%～89.5%。本文中基于SWAT 和隨機森林模型所得的1986-2015年人類活動對洋河流域徑流變化貢獻率高于李秀和田菲等的研究,這可能與所選研究時段差異有關。下墊面變化是人類活動對水文過程影響的一個主要方面,隨著“退耕還林”“退耕還牧”等水土流失治理工程的實施及城市化的發(fā)展,2000年后永定河上游流域土地利用的變化較之前更為顯著,建設用地顯著增加,未利用地向耕地和草地轉化,水域面積衰減,人類活動對該流域徑流量變化的影響更加劇烈,因此不同時段人類活動對該流域徑流變化的影響具有一定的差異。另外,本研究中利用兩種模型對人類活動造成東洋河流域徑流減少貢獻率的計算結果為84.3%～88.5%,對南洋河流域貢獻率都達到95%以上。這可能與東洋河流域和南洋河流域土地利用類型差異有關。通過對比,南洋河流域耕地面積占流域總面積的46%左右,且人口多沿河道分布,東洋河流域耕地面積較少,林地、草地面積占流域面積近70%,這也很好地解釋了南洋河流域徑流受人類活動影響顯著高于東洋河的原因?？傮w來看,利用隨機森林模型和SWAT 模型在永定河上游流域進行徑流貢獻變化歸因分析結果與已有研究大體相同,且具有一定的合理性。

不同的模型和方法在進行徑流預測時會得到不同的結果,也能挖掘到不同的信息,每一種方法都有其局限性和應用價值,根據(jù)研究目的選擇適合的方法是有效的研究手段。SWAT 模型具有明確的物理過程,但模型構建對輸入數(shù)據(jù)的精度和種類要求較高且耗時較長,受其結構復雜性的影響,模型在校準和驗證的過程具有較高不確定性。流域徑流過程具有非線性和高度復雜性,隨機森林模型在進行徑流預測時能夠將年徑流量與氣候變量和流域屬性等一系列預測因子聯(lián)系起來,與傳統(tǒng)方法相比,該模型建模過程簡單,且具備捕捉水文系統(tǒng)的非線性特征,其靈活性比線性回歸具有更好的預測性能。已有研究結果表明相較于人工神經網(wǎng)絡、支持向量回歸和線性模型,隨機森林模型在預測水位日變化時預測結果更準確。本研究在永定河上游的洋河流域進行了模擬和適用性分析,但關于隨機森林模型在徑流變化歸因分析上還需要更加細致和完善的研究,多方法結合可能是降低不確定性的重要途徑。此外,雖然隨機森林模型具有學習速度快和應用方便的特點,且可以在一定程度上避免水文模型調參導致的不確定性問題,但仍具有黑箱模型屬性,其模型輸入和結果的可解釋性還需要在未來研究中深入探討。

4 結論

本文以永定河上游的洋河流域為研究區(qū),基于隨機森林模型對流域徑流變化進行了定量歸因分析的嘗試,同時,基于SWAT 模型對模型結果進行了對比分析,以探討隨機森林模型在徑流變化歸因分析中的適用性,得到如下結論:

1)隨機森林模型在永定河上游流域具有較好的模擬精度。從與SWAT 對比結果來看,其對徑流的模擬精度相對較高,NSE 和均高于SWAT 模型。

2)氣候變化和人類活動對永定河上游流域徑流造成了顯著影響,突變后,年均徑流量下降顯著。氣候變化對徑流量減少的貢獻率為2.4%～15.7%; 人類活動對徑流減少的平均貢獻率為84.3%～97.6%。人類活動是造成徑流減少的主要原因。

3)隨機森林模型的定量歸因結果與SWAT 模型結果及前人研究結果相對接近,其在永定河上游流域徑流變化歸因分析研究中具有一定的適用性。