基于SWAT模型的氣候變化條件下諾敏河流域徑流變化研究

汪 飛, 洪 林, 馬田遙, 吐尼亞孜, 熊繼東

(1.武漢大學 水利水電學院, 武漢 430072; 2.余姚市水利局,浙江 寧波 315400; 3.寧波興晟新能源科技有限公司, 浙江 寧波 315400)

近年來,氣候變化及其影響引起了越來越多國內外專家學者的關注和研究[1-2]。IPCC第五次研究報告指出[3],地表溫度因受溫室氣體濃度增加而持續升高,降水量、降水強度和降水時空分布也將發生改變,從而導致水資源分布及循環發生相應變化。氣候變化提升了極端水文事件的發生概率[4],也會對區域徑流的時空特性產生不同效應[5]。諾敏河流域位于我國東北部地區,近年來在以氣候變化為主的環境變化因素影響下,流域內水資源趨向不足,水土資源的不協調也逐漸突出。因此,研究氣候變化對諾敏河流域徑流的影響,不僅能深入了解氣候變化條件下徑流的演變規律,也可為流域水資源的綜合管理提供參考依據。

耦合全球氣候模式與水文模型是研究氣候變化對徑流影響的主要方法之一[5-7]。劉春蓁等[8]采用4種GCM情景模式研究了我國7大流域徑流和蒸發的變化可能性,并得出4種模式下松花江流域徑流趨于增加;郭生練等[9]選用BCC-CSM1-1氣候模式及RCP4.5排放情景,預測了長江流域氣溫、降水以及徑流的變化情況;江介倫等[10]結合GWLF水文模型分析了5種大氣環流模式下未來氣候對徑流的影響,發現豐水期水量上升而枯水期水量減小。以上研究雖然對徑流受氣候變化的影響進行了一些探究,著重得到未來氣候變化條件下徑流量的變化范圍,但并未探究徑流的突變、周期等特征。本文選取諾敏河流域作為研究對象,采用SDSM統計降尺度方法,借助大氣環流模式預測流域未來的氣候變化,并通過分布式水文模型SWAT模型模擬在未來氣候情景下徑流的相應特征。

1 研究區概況

諾敏河流域位于大興安嶺東南山麓,在121°45′—124°35′E,48°00′—50°30′N范圍內,流域面積25 098 km2。流域平均海拔1 102 m,最高處可達1 300 m,河口高程最低約為170 m。流域春季干旱多風,夏季熱而濕,秋季較短,冬季寒冷漫長,歷年最高氣溫曾高達40.1℃,最低氣溫達到-35.4℃。諾敏河流域降水年內分配不均,雨量多集中在7—8月。流域總出口水文測站為小二溝站,控制面積占流域總面積的95%以上。

2 研究方法

2.1 趨勢分析

趨勢特征通常使用線性回歸、斯皮爾曼分析法[11]、M-K檢驗法等[12]。本文以線性回歸為主要的研究方式,并將突變分析時使用的M-K法作為補充。

該方法主要反映自變量和因變量的線性相關關系,公式如下:

Y=a+bX

(1)

式中:X為自變量,此處表示年份;Y為因變量,在本文中表示研究對象,如氣溫、降水、徑流;a,b分別為回歸系數,其值通過最小二乘法求得。

在趨勢分析中,需要對研究對象和時間的相關性進行顯著性檢驗,此處采用F檢驗法,統計量F的計算公式為:

(2)

2.2 突變分析

突變分析采用在水文、氣象等長時間序列的突變識別中應用較為廣泛的M-K突變檢驗法[13-16]。

對于樣本長度為n的長序列X,構造的M-K檢驗秩序列如公式(3)—(4)所示:

(3)

(4)

由于時間序列相互獨立,定義新的統計量:

(5)

式中:E(Sk)為累計數Sk的均值;var(Sk)為Sk的方差。

其中,UF1=0,UFk統計量服從標準正態分布。查閱特征值表來得到給定的顯著性水平α(本文采用α=0.05)所對應的臨界值UFα,如果統計量大于臨界值,則表明該序列存在著明顯的增加、減少的趨勢,從而判斷時間序列是否存在趨勢變化特征,作為趨勢分析的補充。

E(Sk),var(Sk)的計算公式如下:

(6)

(7)

按照時間序列逆序編排,用上述的構造過程構造新的統計序列UBk。繪制UBk和UBk兩條曲線圖,若UFk和UBk相交,且交點處的統計值在檢驗的上、下置信度區間內,可認為交點所對應的節點可能就是突變發生的時間點。

2.3 周期分析

周期特征通常采用Morlet小波分析法,Morlet小波可以凸顯長序列存在的周期特征,同時能夠定性地評估序列未來的變化趨勢。

其基本原理是由一組小波函數基來近似地表示所研究的特定函數,通過小波變化計算得到小波系數,再由小波系數的實部、模值等系列數值來分析時間序列中的變化特征。

(1) 小波函數。小波分析中需要構建震蕩變化且最終能快速衰減到0的小波函數Ψ(t)∈L2(R),該函數需滿足:

(8)

式中:Ψ(t)∈L2(R)為基小波函數,該函數經由尺度伸縮及時間軸平移形成子小波函數系:

(9)

式中:Ψa,b(t)為子小波;a為尺度因子,且a≠0,指小波周期的長度;b為平移因子,指時間軸上的平移。

(2) 小波變換。在小波變換中,分別在時間和頻率兩個尺度上分析一維序列的信息。對于已知的能量有限信號f(t)∈L2(R),子小波Ψa,b(t)的連續小波變換為:

(10)

式中:Wf(a,b)為小波變換系數;f(t)為平方可積函數。

由于獲得的數據是離散的,假定函數f(kΔt),(k=1,2,…,N;Δt為取樣間隔),那么離散的小波變換形式為:

(11)

由公式(11)可知,小波分析借由增減伸縮因子來獲得序列的低頻和高頻信息,進而得到不同時間、空間尺度層面的局部變化特征。

(3) 小波方差。小波方差是將小波系數值取平方,并在整個時間軸領域內積分,如公式(12)所示:

(12)

在不同的時間尺度下繪制小波方差變化過程,形成小波方差圖,它能夠體現信號波動過程中的能量匯集情況,從而凸顯出序列中擾動劇烈的信息,分析判斷主周期。

雖然至今未中過一毛錢，我的情緒倒也沒什么大的起伏。畢竟只是動動手指，也沒付出多大成本。只是成功塑造了一個價格敏感型消費者（窮鬼）的形象而已。

3 SWAT模型構建

SWAT模型是由美國農業研究所開發的適用于流域尺度的半分布式水文模型[17],包含氣候、水文、侵蝕、土地覆蓋變化、管理措施、主河道變化、營養物質和水體運動等多個部分[18-21],能夠模擬流域內多種不同的水文循環過程[22]。模型具有很好的適用性,可基于GIS提供的空間數據,模擬長時間段內水文過程,被廣泛用于不同區域變化環境下的流域水文特征研究[23]。

3.1 數據資料

本文選取諾敏河流域的小二溝站作為水文代表站點,選取流域內的8個氣象站點作為研究對象。SWAT模型所需基礎資料及來源見表1,各站點分布如圖1所示。

圖1 諾敏河流域氣象、水文站點分布Fig. 1 Distribution of meteorological and hydrological stations in Nomin River Basin

表1 諾敏河流域基礎數據信息Table 1 Basic data information of Nomin River Basin

3.2 基礎數據庫構建

3.2.1 土壤數據庫 本文使用FAO土壤數據庫直接建立用戶數據庫[16]。模型所需的土壤物理屬性參數按獲取方式分為兩部分:第1部分可直接從FAO土壤數據庫中查取;其余數據通過SPAW軟件和統計分析方法計算得出[24]。

3.2.2 土地利用數據庫 將SWAT模型中定義使用的土地利用分布類型按照植被覆蓋特征進行二級分類。對于不同的土地利用類型,考慮到屬性表中有十幾種類型,且其中有些類型可以合并為一種,因此先將它們進行重新分類,然后再得到能夠直接使用的土地利用數據。

3.2.3 氣象數據庫 按照模型所需的格式要求,對各站點日值氣象數據進行整理,并計算特征值,將結果錄入SWAT模型的基礎數據庫中。

3.2.4 空間離散數據庫 SWAT模型基于流域DEM提取河網水系并進行子流域劃分,再根據不同的土地使用類型、土壤類型和坡度,將具有相同類型的單元進行歸集,將諾敏河流域劃分成273個水文響應單元。

3.3 模型參數率定與驗證

選用決定系數R2以及納什系數Ens來對SWAT模型的擬合效果進行評價,其計算公式為公式(13)—(14)。模型誤差評價參考表見表2。徑流擬合曲線見圖2—3。

圖2 小二溝站率定期(2006-2010年)年日徑流擬合曲線Fig. 2 Fitting curve of annual daily runoff at Xiaoergou station from 2006 to 2010

圖3 小二溝站驗證期(2011-2015年)年日徑流擬合曲線Fig. 3 Fitting curve of annual daily runoff at Xiaoergou station from 2011 to 2015

表2 模型誤差評價Table 2 Model error evaluation

(13)

(14)

從徑流擬合曲線來看,徑流模擬值和實測值趨勢一致,兩條曲線較為貼合。率定期和驗證期的R2分別為0.76,0.73,Ens分別為0.72,0.70,兩者的模擬效果都表現為良好。因此,徑流的模擬能夠滿足模型所需的精度要求,故可以選用調整后的敏感參數值來進行后續的徑流模擬,調整后的敏感參數見表3。

表3 SWAT模型敏感參數率定結果Table 3 Sensitivity parameter calibration results of SWAT model

4 未來情景模式下徑流變化規律分析

對諾敏河流域而言,氣溫和降水是影響徑流變化的最直接、最重要的氣候因素。本文選擇較常用、在我國東北地區適應性強的CanESM 5大氣環流模式,采用SDSM統計降尺度對不同排放情景模式進行統計降尺度,獲得諾敏河流域未來幾十年的氣象要素序列值,并將其輸入率定后的SWAT模型,模擬研究區域徑流對氣候變化的響應,并分析徑流的趨勢、突變及周期變化特征。

4.1 氣候情景設定

為模擬評估未來氣候變化,IPCC第五次評估報告采用CMIP 5中的氣候模式,預設了溫室氣體排放濃度情景(RCPs)。

其中:RCP2.6情景假設全球能源使用發生改變,溫室氣體排放將顯著減少;RCP4.5情景假設可再生能源和碳捕捉系統廣泛使用、森林儲碳量增加,溫室氣體排放量也將顯著減少;RCP8.5情景假設無政策干預,全球人口大幅增加、化石燃料消耗變大,溫室氣體排放和濃度不斷增加。

本文使用CMIP 5中CGCMs模式的數據集CanESM5。該數據集包含了全球的歷史數據集和全球2006—2100年RCP2.6,RCP4.5和RCP8.5三種情景下的基礎數據。

4.2 氣候要素分析

將大氣環流模式中的未來數據統計降尺度,生成3種排放情景模式下的數據。對2021—2100年進行分析,設定2021—2060年為21世紀中期,2061—2100年為21世紀后期。

4.2.1 氣溫變化趨勢分析 繪制3種模式下氣溫的趨勢變化圖并統計其特征值(圖4、表4),3種情景模式下年均氣溫、年平均最高氣溫和年平均最低氣溫都保持增長的趨勢。對比21世紀中、后期3種氣溫變化趨勢可發現,隨著時間的持續,3種排放導致的氣候變暖結果差距將越來越大,尤其對于RCP8.5模式,整個21世紀都保持很高的增長率,增幅幾乎都在0.65℃/10 a以上。對于RCP2.6和RCP4.5模式,氣候變暖得到了有效控制,雖也呈現一定的增幅,但增幅明顯小于RCP8.5模式,增長速率僅為RCP8.5模式下的1/4～1/2。兩者相比較來看,對于21世紀中期,RCP2.6模式下多年平均氣溫稍高于RCP4.5模式,但相差不大;而對于21世紀后期,RCP2.6模式對氣候變暖效果的控制要好于RCP4.5模式。

圖4 諾敏河流域21世紀中期和后期3種排放情景模式氣溫變化趨勢Fig. 4 Temperature change trends of three emission scenario modes in Nomin River Basin in the mid 21st century and late 21st century

4.2.2 降水變化趨勢分析 繪制3種模式下年降水的趨勢圖并統計特征值。由圖5、表5可知,3種情景模式下年降水都保持增長的趨勢。其中,RCP8.5模式下降水的增幅最大,在21世紀后期增幅已經超過20 mm/10 a。對于RCP4.5模式,在21世紀中期,年降水的增幅為6.05 mm/10 a,到了21世紀后期,降水增幅增大至12 mm/10 a。RCP2.6模式下,降水增幅相對較小,且21世紀后期與中期相近,增幅相比另兩種情景模式有所減緩。

圖5 諾敏河流域21世紀中期和后期3種排放情景模式降水變化趨勢Fig. 5 Precipitation change trends of three emission scenario modes in Nomin River Basin in the mid 21st century and late 21st century

表5 諾敏河流域3種排放情景模式年降水特征值Table 5 Annual precipitation characteristic values of three emission scenario modes in Nomin River Basin

4.3 徑流變化特征分析

4.3.1 徑流趨勢特征分析 繪制3種排放情景模式下年徑流的趨勢變化圖并統計特征值(圖6、表6)。結果顯示:3種情景模式下徑流都保持增長趨勢。

圖6 諾敏河流域21世紀中期和后期3種排放情景模式下徑流變化趨勢Fig. 6 Runoff variation trend of three emission scenario modes in Nomin River Basin in the mid 21st century and late 21st century

表6 諾敏河流域3種排放情景模式下年徑流特征值Table 6 Annual runoff characteristic values under three emission scenario modes in Nomin River Basin

其中,RCP8.5模式下徑流增幅最大,在21世紀中期保持1.15億m3/10 a的增速,后期以1.72億m3/10 a增長,因此到21世紀后期,其多年平均年徑流量也最大,達到40.23億m3。RCP4.5模式其次,在21世紀中、后期分別以每10年0.49億m3,0.96億m3的增速增長,到21世紀后期多年平均徑量達到35.34億m3。RCP2.6模式下徑流量增長相對較少,在21世紀中期和后期增長量僅為每10年0.06億m3,0.27億m3。

徑流的增長特征與流域未來的降雨和氣溫關系密切,一方面,降水增加對徑流有著最直接的影響,徑流隨著降水的增加而增加;另一方面可以看到,RCP8.5模式下,21世紀后期降雨增加值非常大,為21世紀中期的3倍,而21世紀后期的徑流增幅僅為中期的1.5倍,同樣地,RCP4.5模式和RCP2.6模式下,降雨增幅分別為2,1.5倍,而徑流增幅則分別為2倍和4倍。這是由于氣溫升高引起蒸散發增加,進而抑制了降水促進徑流增加的作用。

4.3.2 徑流突變特征分析 繪制3種情景模式下2021—2100年諾敏河流域徑流的突變檢驗曲線(圖7—9)。結果顯示:隨著氣候變化,3種情景模式的徑流突變特征產生了很大的差別。

圖7 諾敏河流域RCP2.6模式下徑流M-K法突變識別、Pettitt法突變識別Fig. 7 Runoff M-K and Pettitt method mutation identification under RCP2.6 mode in Nomin River Basin

圖8 諾敏河流域RCP4.5模式下徑流M-K法突變識別、Pettitt法突變識別Fig. 8 Runoff M-K and Pettitt method mutation identification under RCP4.5 mode in Nomin River Basin

圖9 諾敏河流域RCP8.5模式下徑流M-K法突變識別、Pettitt法突變識別Fig. 9 Runoff M-K and Pettitt method mutation identification under RCP8.5 mode in Nomin River Basin

對于RCP2.6模式,正反向序列曲線存在多個交點,且所處的區間位于上下置信區間內,但從年徑流的趨勢變化來看,序列的正向曲線均未能夠突破上下置信度線。結合Pettitt檢驗曲線看,該曲線并不存在明顯的最低點。因此可認為在RCP2.6模式下諾敏河流域年徑流在未來不會突變。

對于RCP4.5模式,正反向序列曲線也存在多個交點,且所處的區間位于上下置信區間內,年徑流序列的正向曲線突破了上置信度線,但并不顯著。結合Pettitt檢驗曲線,該曲線存在一個明顯的最低點。因此可認為在RCP4.5模式下諾敏河流域年徑流在未來存在很大的可能形成突變。

對于RCP8.5模式,正反向序列曲線存在1個交點,且交點位于上下置信區間內,年徑流序列的正向曲線明顯地突破了置信度線,交點位于2067年附近。結合Pettitt檢驗曲線看,該曲線存在一個明顯的最低點,對應的時間為2067年。因此可認為在RCP8.5模式下諾敏河流域年徑流將在2067年附近發生突變,突變后的徑流表現為顯著的增加。徑流的突變特征與趨勢特征也表現出了較好的一致性。

4.3.3 徑流周期特征分析 對3種模式下2021—2100年的徑流進行小波計算,并繪制小波圖(圖10—12)。對于RCP2.6模式,序列存在3個峰值,對應7,19,56 a時間尺度。其中56 a時間尺度為峰值的最高點,可認為是該序列的第一主周期,7 a時間尺度次之,為序列的第二周期,19 a則對應第三周期。結合實部等值線圖可知,雖然30,45 a的時間尺度在方差圖上沒有明顯的峰值特征,但也存在一些周期特征。

圖10 諾敏河流域RCP2.6模式下徑流小波系數實部等值線圖、小波方差Fig. 10 Real part contour map and wavelet variance map of runoff wavelet coefficients under RCP2.6 mode in Nomin River Basin

圖11 諾敏河流域RCP4.5模式下徑流小波系數實部等值線圖、小波方差Fig. 11 Real part contour map and wavelet variance map of runoff wavelet coefficients under RCP4.5 mode in Nomin River Basin

圖12 諾敏河流域RCP8.5模式下徑流小波系數實部等值線圖、小波方差Fig. 12 Real part contour map and wavelet variance map of runoff wavelet coefficients under RCP8.5 mode in Nomin River Basin

對于RCP4.5模式,序列存在4個峰值,對應19,32,44,56 a時間尺度。其中56 a時間尺度為峰值的最高點,可認為是該序列的第一主周期,44 a時間尺度次之,為序列的第二周期,32 a和19 a則對應了第三、第四周期。相比于RCP2.6模式,7 a時間尺度在周期特征變現上產生了波動,而30 a和45 a附近的周期特征得到了更明顯的凸顯。

對于RCP8.5模式,其周期特征和RCP4.5模式相似,峰值對應的時間分別為19,32,44,56 a。其中56 a時間尺度為峰值的最高點,可認為是該序列的第一主周期,44 a時間尺度次之,對應為序列的第二周期,32 a和19 a則對應了第三、第四周期。

5 討論與結論

本文用未來氣候數據驅動SWAT模型進行氣候變化條件下徑流特征分析。從趨勢上看:3種情景模式下徑流都保持增長趨勢。其中,RCP8.5模式下增幅最大,21世紀中期保持1.15億m3/10 a的增速,后期以1.72億m3/10 a增長。RCP4.5模式其次,在21世紀中、后期分別以0.49億m3/10 a,0.96億m3/10 a的增速增長。RCP2.6模式下徑流量增長相對較少,在21世紀中期和后期增速僅為每10年0.06億m3,0.27億m3。

從突變上看:RCP2.6模式下,諾敏河流域年徑流依舊不會突變;RCP4.5模式下,年徑流在未來存在很大的可能形成突變。而RCP8.5模式下,年徑流將在2067年附近發生突變,突變后的徑流表現為顯著的增加。

從周期上看:RCP2.6模式下,年徑流變化周期主要由56 a時間尺度的主周期控制,7 a和19 a時間尺度對應的次周期控制效果一般;RCP4.5模式下,年徑流變化周期主要由56 a和44 a的第一、第二主周期控制,32 a和19 a對應的次周期控制效果不強。RCP8.5模式下,其周期特征和RCP4.5模式相似,相比之下,56 a時間尺度下能量密度更為密集,該周期的控制效果得到進一步增強,而其余幾個周期下的控制效果則進一步減弱。

本文在研究的過程中,水文模型、氣候模式、降尺度模型中都存在不確定性,怎樣量化這些模型在計算中的不確定性,了解計算過程中相互之間存在的分歧,需要更為深入地研究。同時,除氣候變化環境因素外,人類活動如砍伐森林、城市化發展等,都會對徑流產生較大的影響[25-27],如何定量地區分其中氣象要素和人類活動的貢獻率也有待進一步分析。