北大沙河流域徑流變化歸因分析

何玉芬，楊漢波，唐莉華，雷慧閩，楊大文

(1.清華大學(xué)水利水電工程系，北京 100084； 2.清華大學(xué)水沙科學(xué)與水利水電工程國家重點實驗室，北京 100084)

徑流是人類可利用水的主要來源之一。我國多個地區(qū)觀測到徑流發(fā)生了顯著的變化[1-2]，越來越多的研究關(guān)注徑流變化的驅(qū)動因子，其中，氣候變化和人類活動被認(rèn)為是兩大主要因素[3]。氣候變化通過降水、溫度和日照等氣象因子影響徑流過程，人類活動的影響主要體現(xiàn)在土地利用和覆被變化、水庫運行和地表水、地下水的直接取水[4]。研究氣候變化和人類活動對徑流的影響，分離和量化二者的貢獻(xiàn)，一方面有助于認(rèn)識水文過程與氣候和環(huán)境變化的相互關(guān)系，另一方面可為水資源管理等問題提供科學(xué)依據(jù)。

量化氣候變化和人類活動對徑流的影響主要有3種方法：一是流域試驗法，但局限于小流域的研究，需要大量的人力物力，代價昂貴；另外兩種方法分別基于水文模型和統(tǒng)計分析[5]。常見的水文模型有集總式水文模型和分布式水文模型。Chang等[4]分別用TOPMODEL和可變滲透能力(variable infiltration capacity, VIC)模型分析了涇河流域徑流變化的原因，并比較了水文模型法和氣候彈性法歸因結(jié)果的差異；付曉花等[6]基于SWAT模型分析了氣候變化和人類活動對漳衛(wèi)南運河流域徑流變化的影響，尤其量化了水土保持、城鎮(zhèn)化等不同人類活動因素對徑流變化的貢獻(xiàn)率。基于統(tǒng)計分析的方法常采用基于Budyko框架的彈性系數(shù)法[7]。薛聯(lián)青等[8]基于Budyko框架，對塔里木河流域徑流變化對氣候變化和人類活動的響應(yīng)進(jìn)行了敏感性分析；Wang等[9]通過Budyko框架量化分析了黃河源地區(qū)徑流對氣候變化的響應(yīng)，認(rèn)為該地區(qū)的徑流下降主要受凍土退化的間接影響。

在量化歸因分析時，過去的研究多集中于年尺度，而水資源季節(jié)分布不均直接導(dǎo)致了不同月份用水的差異，因此，探討氣候變化和人類活動對月徑流的影響至關(guān)重要。Thomas[10]提出的ABCD模型由4個參數(shù)(a、b、c、d)控制，明確地包含了主要水文過程，近年來被廣泛用于年尺度和月尺度的水文過程分析中。Li等[11]利用ABCD模型模擬漢江上游土壤蓄水量和實際蒸散發(fā)的變化，定量評估了氣候變化和人類活動對不同時期徑流的影響；Wang等[12]通過考慮地下水蒸發(fā)對ABCD模型進(jìn)行了修正，研究了鄂爾多斯高原干旱到半干旱地區(qū)的年徑流變化。此外，還有不少學(xué)者開展了ABCD模型的應(yīng)用研究[5,13-14]。

目前針對徑流變化歸因分析的研究主要存在兩個問題：一是對月徑流變化的歸因分析較少；二是對人類活動的分析不夠深入，許多研究將總徑流變化量與氣候變化影響量的余項作為人類活動影響量。因此，本文采用ABCD模型建立分離氣候變化和人類活動對徑流影響的分析框架，結(jié)合基于Budyko框架的彈性系數(shù)法，開展北大沙河流域徑流變化原因研究，分析人類活動影響下的徑流季節(jié)性變化。

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

北大沙河流域(36.24°N～36.51°N，116.83°E～117.09°E)位于山東省濟南市長清區(qū)，河流發(fā)源于武家莊鄉(xiāng)摩天嶺西麓，屬于黃河下游右岸支流，全長54.3 km，流域面積584 km2，多年平均年徑流量為0.76億 m3。流域海拔30～1 000 m，地勢東南高、西北低，呈現(xiàn)山地-丘陵-平原的地形變化。該流域為溫帶大陸性季風(fēng)氣候，夏季炎熱多雨，冬季寒冷干燥，多年平均氣溫13.7 ℃，多年平均降水量845 mm，降水量年內(nèi)分布不均，汛期降水占全年的76%，多年平均蒸發(fā)皿(直徑為20 cm)蒸發(fā)量達(dá)2 300 mm。流域土壤類型以棕壤土和褐土類為主，植被類型以溫帶落葉闊葉林為主。作為濟南市重要的地下水補給區(qū)、飲用水源地和生態(tài)走廊，北大沙河在改善當(dāng)?shù)厣鷳B(tài)環(huán)境、保障城市用水安全等方面具有重要作用。

同海河流域、黃河中下游等流域類似，北大沙河流域受人類活動干擾較大。流域上游段(張夏大橋以上)建有中型水庫1座(石店水庫)、小(1)型水庫7座、小(2)型水庫14座、塘壩105座，總攔蓄庫容為3 295萬m3。統(tǒng)計資料顯示，截至2018年，山東省共建成大型灌區(qū)50處、中型灌區(qū)444處、小型灌區(qū)14.77萬處，相比于1980年，有效灌溉面積增加了78萬hm2，人類活動逐漸加強，對水資源的需求逐漸加大。在水資源短缺、供需矛盾加劇的大背景下，水庫的調(diào)蓄作用越發(fā)突顯。

本研究中，徑流數(shù)據(jù)采用崮山水文站(36.48°N，116.87°E)記錄數(shù)據(jù)，來自《中華人民共和國水文年鑒》。崮山水文站1979年建站，1980年起有數(shù)據(jù)記錄，1991—2005年沒有數(shù)據(jù)記錄。因此，提取1980—1990年和2006—2018年兩個時期的數(shù)據(jù)，并以前者為基準(zhǔn)期，后者為擾動期。

降水?dāng)?shù)據(jù)和20 cm小型蒸發(fā)皿觀測數(shù)據(jù)采用濟南站(36.36°N，117.03°E，高程為170 m)和泰山站(36.15°N，117.06°E，高程為1 534 m)的實測數(shù)據(jù)，來源于國家氣象科學(xué)數(shù)據(jù)中心。潛在蒸散發(fā)量為蒸發(fā)皿蒸發(fā)量與蒸發(fā)皿系數(shù)的乘積，參考劉波[15]在全國范圍內(nèi)確定的蒸發(fā)皿系數(shù)分布圖，取山東省的蒸發(fā)皿系數(shù)為0.65。對各站點降水量和潛在蒸散發(fā)量取算術(shù)平均作為流域面平均降水量和潛在蒸散發(fā)量。

2 研究方法

2.1 基于ABCD模型的歸因分析

ABCD模型是一個非線性水量平衡模型，以潛在蒸散發(fā)和降水?dāng)?shù)據(jù)作為輸入，可用于模擬基流、地下水、土壤水以及實際蒸散發(fā)量的變化[16]。模型結(jié)構(gòu)如圖1所示。

圖1 ABCD模型結(jié)構(gòu)

模型將流域分為土壤水和地下水兩層，設(shè)時間步長為t，定義可用水量Wt和可能蒸散發(fā)量Yt兩個狀態(tài)變量，假設(shè)Yt是Wt的非線性函數(shù)：

(1)

式中：a為土壤飽和前形成徑流的概率；b為流域?qū)嶋H總蒸發(fā)量與土壤水分儲量之和的上限，mm。

可用于產(chǎn)生徑流的水量Wt-Yt被分為直接徑流和地下水補給：

Dt=(1-c)(Wt-Yt)

(2)

Rt=c(Wt-Yt)

(3)

式中：Dt為直接徑流，mm；Rt為地下水補給量，mm；c為地下水補給系數(shù)。

假設(shè)地下水層為一個線性水庫，通過參數(shù)d進(jìn)行控制：

Ft=dGt

(4)

式中：Ft為基流，mm；Gt為地下水儲量，mm；d為地下水儲放系數(shù)。

總徑流為直接徑流與基流之和：

Qt=Dt+Ft

(5)

式中Qt為總徑流，mm。

采用遺傳算法率定ABCD模型的參數(shù)，選取年水量平衡誤差作為率定的目標(biāo)函數(shù)，計算公式為

(6)

式中：e為年水量平衡誤差；Qobs,i、Qsim,i分別為第i個觀測和模擬的徑流值，mm；m為數(shù)據(jù)序列長度。通過多次迭代與模型預(yù)熱，確定初始土壤水和地下水儲量。

在ABCD模型的基礎(chǔ)上，提出以下歸因分析的框架：利用基準(zhǔn)期率定的參數(shù)重建擾動期的徑流過程，總徑流變化為兩個時期的觀測徑流之差，氣候變化導(dǎo)致的徑流變化為兩個時期的模擬徑流之差，人類活動影響量為總徑流變化減去氣候變化影響量，即：

ΔQ=Qobs2-Qobs1

(7)

ΔQc=Qsim2-Qsim1

(8)

ΔQh=ΔQ-ΔQc

(9)

式中：ΔQ為總徑流變化，mm；ΔQc、ΔQh分別為氣候變化和人類活動引起的徑流變化，mm；Qobs1、Qobs2分別為基準(zhǔn)期和擾動期的觀測徑流，mm；Qsim1、Qsim2分別為基準(zhǔn)期和擾動期的模擬徑流，mm。

氣候變化和人類活動的貢獻(xiàn)率為

(10)

(11)

式中ηc、ηh分別為氣候變化和人類活動對徑流變化的貢獻(xiàn)率。

2.2 基于Budyko框架的彈性系數(shù)法

與水文模型法相比，基于Budyko框架的彈性系數(shù)法不僅具有清晰的物理意義，且對參數(shù)選擇和模型構(gòu)建的要求相對較低。Budyko框架經(jīng)過多年的發(fā)展，衍生出了許多改進(jìn)的模型，本研究采用Choudhury等[17-18]提出的方程，并結(jié)合水量平衡方程，得到：

(12)

式中：Q為多年平均徑流深，mm；P為多年平均降水量，mm；E、E0分別為多年平均實際蒸散發(fā)量和潛在蒸散發(fā)量，mm；n為反映流域特征的下墊面參數(shù)。

為方便起見，式(12)可以表示為

(13)

根據(jù)Milly 等[19]的研究，氣候變化的影響量可以表示為

(14)

3 結(jié)果與分析

3.1 氣候和徑流變化

圖2為流域1980—1990年和2006—2018年徑流深和降水量變化。基準(zhǔn)期和擾動期的多年平均降水量分別為830 mm和870 mm，均呈現(xiàn)較大的年際波動，但未發(fā)生顯著趨勢變化(顯著性水平0.05)。同基準(zhǔn)期相比，擾動期的降水量呈現(xiàn)相對明顯的增加趨勢。基準(zhǔn)期和擾動期的多年平均蒸散發(fā)量分別為1 218 mm 和1 136 mm。

圖2 基準(zhǔn)期和擾動期年徑流深和降水量變化

1980—1990年北大沙河流域年徑流深在0～270 mm之間，多年平均徑流深為101 mm；2006—2018年年徑流深在5～258 mm之間，多年平均徑流深為103 mm。兩個時期的徑流總量接近，且由Mann-Kendall趨勢檢驗發(fā)現(xiàn)徑流均未發(fā)生顯著增加或減少的趨勢，說明在氣候變化和人類活動的共同作用下，徑流呈現(xiàn)較大的年際波動，但未發(fā)生顯著的趨勢變化。在月尺度上，北大沙河流域徑流具有較大的季節(jié)變異性，集中在汛期7—9月；基準(zhǔn)期汛期徑流占全年徑流的60%以上，擾動期汛期徑流占全年徑流的90%以上。非汛期徑流量較小，2006—2018年出現(xiàn)了斷流現(xiàn)象，來水幾乎全部為汛期徑流。

ABCD模型預(yù)熱得到的初始地下水和土壤水含量分別取190 mm和400 mm。得到的a、b、c、d參數(shù)的取值如表1所示，符合類似研究給出的參數(shù)范圍。計算得到的年水量平衡誤差為4.5%，意味著模擬結(jié)果較為可靠。

表1 ABCD模型參數(shù)取值

3.2 氣候變化和人類活動對徑流的貢獻(xiàn)

利用彈性系數(shù)法評估徑流對降水和潛在蒸散發(fā)的敏感性，由式(12)得到的下墊面參數(shù)n為2.49。徑流對降水和潛在蒸散發(fā)的敏感性分別為0.40和-0.19，表明徑流對降水的變化更為敏感。年尺度上，彈性系數(shù)法得到的氣候變化和人類活動對徑流變化的貢獻(xiàn)量分別為32.0 mm和-30.0 mm。

表2給出了利用ABCD模型分析年尺度、季節(jié)(春季為3—5月,夏季為6—8月，秋季為9—11月，冬季為12—2月)尺度和月尺度上徑流對氣候變化和人類活動響應(yīng)的結(jié)果。無論在年尺度還是季節(jié)尺度上，氣候變化導(dǎo)致徑流增加，而人類活動導(dǎo)致徑流減少。年尺度上，氣候變化和人類活動對徑流變化的貢獻(xiàn)量分別為43.0 mm和-41.0 mm；季節(jié)尺度上，春、冬兩季的徑流變化受人類活動的主導(dǎo)，夏、秋兩季的徑流變化則受氣候變化的主導(dǎo)。

表2 不同時間尺度ABCD模型歸因分析結(jié)果

月尺度上，氣候變化使月徑流增加了1.3～18.1 mm，在1—6月的貢獻(xiàn)率較低，8月最高。擾動期與基準(zhǔn)期相比，8月平均降水量增加了49 mm，平均潛在蒸散發(fā)量減少了21 mm(圖3)，二者的共同作用使氣候變化在這一階段的貢獻(xiàn)尤為突出。人類活動在1—6月的貢獻(xiàn)率較低，7—12月較高，說明水庫的調(diào)蓄主要發(fā)生在7—12月。值得注意的是，除9月外，人類活動均導(dǎo)致徑流的減少(1.9～10.1 mm)，而9月人類活動使得徑流增加了11.3 mm。

在年尺度上，兩種方法的歸因分析結(jié)果較為一致。氣候變化導(dǎo)致年徑流增加，而人類活動導(dǎo)致年徑流減少，增加和減少量相當(dāng)，即徑流對氣候變化和人類活動的響應(yīng)程度相當(dāng)。

4 討 論

4.1 結(jié)果對比

在年尺度的歸因分析中，ABCD模型和彈性系數(shù)法計算得到的氣候變化對徑流變化的貢獻(xiàn)量分別是43.0 mm和32.0 mm，二者相差較大，可能的原因是彈性系數(shù)法未考慮降水等因子的年內(nèi)分配過程(即季節(jié)性)的變化，而不少研究指出氣候季節(jié)性變化對流域水量平衡有顯著影響[22]。圖3反映了兩個時期的降水量和潛在蒸散發(fā)量季節(jié)性變化，從基準(zhǔn)期到擾動期，降水量和潛在蒸散發(fā)量的季節(jié)變異性較大，尤其是8月，平均降水量增加了49 mm，平均潛在蒸散發(fā)量減少了21 mm，二者的共同作用使夏季徑流增加。采用Woods定義的氣候季節(jié)性指數(shù)來量化降水量和潛在蒸散發(fā)量差值的季節(jié)性波動[23]，其值越大說明兩者差值的季節(jié)性波動越大，即水分和能量的匹配性越差。本文計算得到兩個時期的氣候季節(jié)性指數(shù)分別為0.33和0.43，根據(jù)楊漢波等[22]的結(jié)論，氣候季節(jié)性指數(shù)的增加將導(dǎo)致參數(shù)n的減小，而式(13)反映出徑流量是n的減函數(shù)，因此氣候季節(jié)性指數(shù)的增加將導(dǎo)致徑流增加。常用的彈性系數(shù)法，由于沒有考慮氣象因子季節(jié)性變化的影響，可能低估了氣候變化的影響；而ABCD模型以月過程為輸入，包括了氣候季節(jié)性變化的影響，可能得到更為合理的結(jié)果。

圖3 降水量與潛在蒸散發(fā)量年內(nèi)分配過程

本文利用ABCD模型分離氣候變化和人類活動的影響，氣候變化導(dǎo)致的徑流變化量為兩個時期的模擬徑流之差(式(8))，人類活動導(dǎo)致的徑流變化量為總的徑流變化減去氣候變化導(dǎo)致的徑流變化量(式(9))。此前不少研究采用ABCD模型進(jìn)行徑流變化歸因分析時[13,24]，將擾動期的模擬徑流量和觀測徑流量之差作為人類活動的影響，沒有考慮到模型本身可能產(chǎn)生的系統(tǒng)誤差。事實上，擾動期的模擬徑流量與觀測徑流量間的差異包括了人類活動影響和模型誤差，因此，本文提出的框架分析結(jié)果更為合理。

4.2 人類活動對北大沙河流域徑流變化的影響

由表2可知北大沙河流域夏、秋兩季的徑流主要受氣候變化的影響，體現(xiàn)為夏季降水量的增加和潛在蒸散發(fā)量的減少；而春、冬兩季人類活動的影響占主導(dǎo)地位，主要體現(xiàn)在水庫調(diào)蓄上。流域上游段建成的水庫通過調(diào)節(jié)流量的大小和持續(xù)時間來改變自然流動狀態(tài)，使月徑流分配發(fā)生變化，汛期徑流占年徑流總量的比例減少, 枯季時被用于農(nóng)田灌溉，洪峰被削弱，月徑流過程坦化。

人類活動對徑流的影響表現(xiàn)出明顯的動態(tài)和季節(jié)性特征，枯季降水量和蒸發(fā)量均較少，人類的直接取水使其調(diào)蓄作用被突顯，由表2可知，水庫在下半年的調(diào)蓄量更大。汛期(7—9月)石店、武家莊水庫和小崮山水庫發(fā)揮蓄水功能，導(dǎo)致出口段夏季徑流減少。9月人類活動導(dǎo)致徑流增加，可能有兩方面原因：一是人類通過河道或渠道取用水庫水進(jìn)行灌溉，二是上游水庫陸續(xù)向下游河道供水，兩種方式均可能導(dǎo)致河道斷面觀測到的徑流在短期內(nèi)增加。

研究發(fā)現(xiàn)北大沙河在擾動期出現(xiàn)了枯水期斷流現(xiàn)象。1960年以來，多個攔水大壩的興建以及引黃灌溉規(guī)模的擴大, 使得實測黃河年徑流量遠(yuǎn)低于天然徑流量，導(dǎo)致下游山東等地長時間斷流。Liu等[25]重建了黃河中游過去500多年天然徑流量變化過程, 并定量計算出人類活動耗水量，說明在研究徑流變化時有必要考慮人類活動的大量耗水。王樂揚等[26]也指出人類活動已成為我國北方主要河流徑流變化的重要原因。

本文研究結(jié)果表明，除9月外，人類活動導(dǎo)致徑流減少，而觀測徑流仍能維持不變，這得益于降水量的增加和蒸發(fā)能力的下降。王國慶等[27-28]基于未來氣候情景分析得出，未來黃河下游可利用水資源量有減少趨勢，且華東地區(qū)北部表現(xiàn)出水資源對氣候變化的高度脆弱性，如果人類活動的負(fù)面影響持續(xù)，流域水資源供需矛盾可能進(jìn)一步加劇[29]。因此，相關(guān)部門需加強可持續(xù)的水資源管理，進(jìn)行更加合理的規(guī)劃配置。

4.3 研究存在的不足

彈性系數(shù)法中，人類活動影響為觀測到的徑流變化與氣候變化影響量的余項，包括流域蓄水量、植被等因素變化的影響，難以單獨分離出每個因子的變化對徑流產(chǎn)生的影響；氣候季節(jié)性指數(shù)變化的影響被歸結(jié)到人類活動影響中，低估了氣候變化的影響。此外，彈性系數(shù)法忽略了徑流對氣象要素的高階偏導(dǎo)[30]，在氣候年內(nèi)變異性較大的地區(qū)可能會產(chǎn)生較大的誤差。

本研究利用ABCD模型進(jìn)行歸因分析時，將人類活動視為總徑流變化與氣候變化影響量的余項，沒有考慮流域土地利用類型、植被的覆蓋度等因素，認(rèn)為流域下墊面特征在兩個時期保持一致，會帶來一定的誤差[31]；ABCD模型中地下水和土壤水的初始值設(shè)置等因素也會帶來一定的不確定性。此外，由于水庫分布的空間復(fù)雜性和調(diào)度的時間復(fù)雜性，本文沒有研究小型水庫對徑流的時空影響，也沒有考慮詳細(xì)的水庫規(guī)模和運作方式，若要進(jìn)行更加細(xì)致的歸因分析，需分離水庫的影響。

5 結(jié) 論

a.在氣候變化和人類活動的共同作用下，北大沙河流域徑流呈現(xiàn)較大年際波動，但未發(fā)生顯著趨勢變化，年內(nèi)具有較大的季節(jié)變異性，徑流集中在汛期7—9月，2006—2018年出現(xiàn)了枯水期斷流的現(xiàn)象。

b.從基準(zhǔn)期(1980—1990年)到擾動期(2006—2018年)，年尺度上，ABCD模型和彈性系數(shù)法得到的氣候變化導(dǎo)致徑流的增加量分別為 43.0 mm和32.0 mm，彈性系數(shù)法由于未考慮氣象要素年內(nèi)分配變化的影響，可能低估了氣候變化的影響；月尺度上，氣候變化使徑流在汛期增加了4.4～18.1 mm，非汛期增加了1.3～2.2 mm，在夏、秋兩季占主導(dǎo)地位，人類活動使徑流9月增加了11.3 mm，其他月份減少了1.9～10.1 mm，在春、冬兩季占主導(dǎo)地位。

c.人類活動主要體現(xiàn)在水庫對月徑流過程的調(diào)蓄上，在下半年的調(diào)蓄量更大。水庫的修建為解決水資源時空分布不均問題提供了解決方案，但改變了天然水循環(huán)過程，在未來降水有減少趨勢的情形下，應(yīng)對水資源進(jìn)行可持續(xù)規(guī)劃與配置。

d.ABCD模型在小流域徑流變化歸因分析中具有適用性，為年、月尺度的徑流變化研究提供了一種更為合理、靈活的手段。