寧婷婷,劉文兆,韓曉陽,林文

(中國科學院水利部水土保持研究所,黃土高原土壤侵蝕與旱地農(nóng)業(yè)國家重點實驗室,712100,陜西楊凌;中國科學院大學,100049,北京)

干燥指數(shù)是表征一個地區(qū)地表干濕狀況的指標,長期應(yīng)用在地理學、生態(tài)學等學科研究中,近年來成為全球變化研究中經(jīng)常涉及到的氣候指標之一,一般以潛在蒸散量與降水量的比值表示[1]。在降水量小于潛在蒸散量的地區(qū),其干燥指數(shù)大于1,氣候比較干燥,相反,降水量大于潛在蒸散量的地區(qū)則較為濕潤[2]。潛在蒸散量是指充分供水條件下的區(qū)域蒸散發(fā)能力[3]。干燥指數(shù)可以看作是區(qū)域的氣候指示器,其變化會影響到區(qū)域水文循環(huán)、生態(tài)系統(tǒng)和水資源的管理[4]。

以氣候變暖為主要特征的全球變化已成為一個不爭的事實[5],氣候變化背景下全球和區(qū)域水分循環(huán)以及未來氣候可能出現(xiàn)的變化成為全球氣候變化研究的熱點[6]。一般認為,全球變暖應(yīng)該使地表空氣更加干燥從而導致水體和陸地的潛在蒸散量增加,但事實卻與此相反。已有大量研究[7-9]表明,世界上許多地區(qū)潛在蒸散量在過去的幾十年呈下降趨勢[10]。盡管潛在蒸散量有減小的趨勢;但是全球降水量的趨勢卻因地區(qū)而異,如近幾十年來,北半球高緯地區(qū)和赤道地區(qū)降水量有增多趨勢,而中國、澳大利亞和一些太平洋島國部分地區(qū)降水量卻是減少的[11]。在全球潛在蒸散量和降水量都發(fā)生變化的同時,地表的干濕狀況又發(fā)生了怎樣的變化,這同樣也引起了廣泛的研究。如H.Tabari等[12]發(fā)現(xiàn),伊朗1966—2005年的干燥指數(shù)有減小的趨勢,Zhang Qiang等[13]對珠江流域的干燥指數(shù)研究發(fā)現(xiàn),1960—2005年流域西部有變干的趨勢,而東南部則變得更濕潤,Huo Zailin等[14]的研究表明,1955—2008年中國西北干旱地區(qū)的干燥指數(shù)有減小的趨勢。由此可見,受潛在蒸散量和降水量地區(qū)間不同變化趨勢的影響,干燥指數(shù)的變化也因地區(qū)而異,因此很有必要展開區(qū)域研究。

內(nèi)蒙古東勝礦區(qū)地處毛烏素沙地與黃土丘陵溝壑區(qū)風、水復合侵蝕的中心地帶(圖1),多年平均降水量約為385 mm,潛在蒸散量為1 128.1 mm。行政上屬伊盟東勝等旗市管轄(E 108°40′～110°50′,N 39°20′～ 40°15′)。 礦區(qū)南北長 100 km,東西最寬110 km,面積8 794 km2。礦區(qū)含煤面積1 030 km2,探明儲量131 Gt,盛產(chǎn)低灰、低硫、低磷的優(yōu)質(zhì)動力煤和氣化用煤[15]。煤礦開采是本區(qū)的支柱產(chǎn)業(yè),但同時也帶來了一系列的生態(tài)環(huán)境問題,如土地沙化及水土流失加重、水資源污染及水資源短缺擴大化等,人類活動的影響表現(xiàn)突出。選擇這一區(qū)域,就長時段地表干濕狀況的變化過程進行分析,查明主導氣象因子,并說明其中人類活動的可能影響,對于東勝礦區(qū)乃至整個晉陜蒙接壤區(qū)的礦業(yè)發(fā)展與環(huán)境治理具有重要意義。

圖1 研究區(qū)位置Fig.1 Location of study area

1 研究方法

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

研究區(qū)為內(nèi)蒙古東勝礦區(qū),選擇站點為國家級臺站——東勝站(E 109.98,N 39.83),研究時段為1957—2012年。氣象資料來源于中國氣象數(shù)據(jù)共享服務(wù)網(wǎng)。所用氣象數(shù)據(jù)包括逐日最高、最低和平均氣溫,日照時間,風速,相對濕度和降水量等。

1.2 干燥指數(shù)

自1990年以來,中外學者陸續(xù)提出了20多種干燥指數(shù)的計算方法[1],其中用潛在蒸散量和降水量之比來表達的較多,計算公式為

式中:φ為干燥指數(shù);E0為潛在蒸散量,mm;P為降水量,mm。

潛在蒸散量的計算采用FAO 56推薦的Penman-Monteith模型,該方法較全面地考慮了影響蒸散的各種因素,而且在氣候條件差異較大地區(qū)的應(yīng)用中都取得了較好的結(jié)果,計算精度優(yōu)于其他方法,目前被認為是計算潛在蒸散量精度最高、應(yīng)用最廣泛的一種方法[16-17]。其計算公式[3]為

式中:E0為潛在蒸散量,mm/d;Rn為作物表面凈輻射量,MJ/(m2·d);G 為土壤熱通量,MJ/(m2·d);γ 為濕度計常數(shù),kPa/℃;Δ為飽和水汽壓與溫度關(guān)系曲線的斜率,kPa/℃;θ為平均氣溫,℃;v為在地面以上2 m處的風速,m/s;es為飽和水汽壓,kPa;ea為實際水汽壓,kPa。

1.3 Mann-Kendall趨勢檢驗和突變檢驗

采用非參數(shù)Mann-Kendall(以下簡稱M-K法)趨勢分析及突變檢驗法來檢測干燥指數(shù)的突變。M-K法是用來評估氣候要素時間序列的檢驗方法,以適用范圍廣、人為性少、定量化程度高而著稱[18],其檢驗統(tǒng)計量公式為

式中sgn()為符號函數(shù)。當xi-xj小于、等于或大于零時,sgn(xi-xj)分別為-1、0或1,而M-K趨勢檢驗的統(tǒng)計值文中用Zc表示。

1.4 水文頻率分析

一般認為,皮爾遜Ⅲ型曲線[19]比較符合水文隨機變量的分布,因此,采用該曲線來擬合干燥指數(shù)的多年變化,從而掌握其頻率變化特征以及不同水文頻率下的重現(xiàn)期。

1.5 小波分析

小波分析在氣候變化研究中顯現(xiàn)了年與年之間、波長與波長之間的波幅變化形態(tài)與特征,已廣泛地應(yīng)用于氣象領(lǐng)域,并成為研究氣象要素長期變化的重要工具[20]。通過Morlet小波分析方法研究內(nèi)蒙古東勝礦區(qū)近56年干燥指數(shù)變化規(guī)律,分析方法參照文獻[21]。

1.6 主導因子和貢獻值分析

由干燥指數(shù)的計算公式可知,干燥指數(shù)是降水量 P、日照時間 th、實際水汽壓 ea、最高氣溫 θmax、最低氣溫θmin和風速v的函數(shù),即φ=f(P,th,ea,θmax,θmin,v)。當自變量只有時間t時,利用全導數(shù)公式可計算研究時段內(nèi)各氣象要素隨時間發(fā)生變化后對干燥指數(shù)變化的貢獻值:

式中,等號左邊代表干燥指數(shù)在研究時段內(nèi)的變化量,等號右邊的6項是干燥指數(shù)對各變量的偏導與各變量的斜率乘積,代表研究時段內(nèi)降水量、日照時間、實際水汽壓、最高氣溫、最低氣溫和風速發(fā)生變化后對干燥指數(shù)年變化量的實際貢獻值,并規(guī)定引起干燥指數(shù)上升為正貢獻,引起干燥指數(shù)下降為負貢獻。

以年為單位,通過多元線性回歸方法,先求出式(4)中所對應(yīng)的各變量偏導部分,作為干燥指數(shù)隨各氣象要素的平均變率,再乘以各要素隨時間的平均變率,即可得到式(5)中各氣象要素對干燥指數(shù)變化的貢獻值。此外,各氣象因子對干燥指數(shù)的相對貢獻值可由 C(P)、C(θmax)、C(θmin)、C(v)、C(ea)和C(th)6項氣象因子的相互比值來表示。

2 結(jié)果與分析

2.1 干燥指數(shù)的時程變化

圖2顯示的是近56年東勝礦區(qū)干燥指數(shù)的多年變化趨勢。可知,該地區(qū)干燥指數(shù)的多年平均值為3.2,干燥指數(shù)的波動范圍為1.5～6.3,波動范圍較大,最干燥的一年是2000年,最濕潤的一年是1961年,極值的出現(xiàn)很大程度上都與降水量有關(guān),1961年是56年降水量最多的年份(709.7 mm),而2000年卻是降水量最少的一年(181 mm)。通過線性回歸分析可知,近56年干燥指數(shù)上升幅度僅為0.004/a。

對該區(qū)56年的干燥指數(shù)月均值作M-K檢驗(圖3),可以看出:1、2、4、7、8 和 10 月的 Zc為正值,表示近56年這些月份地表有逐漸變干的趨勢,而其他月份則有逐漸變濕潤的趨勢,但其中只有12月份達到了顯著水平(P＜0.05);所以,在月尺度上,地表干濕狀況值雖有波動變化,但是趨勢不明顯。根據(jù)鄭景云等[22]對干濕區(qū)劃分的標準,該區(qū)春季(3—5月)和冬季(12月—翌年2月)屬于干旱氣候類型,而到了夏(6—8月)、秋(9—11月)季節(jié),又向半干旱過渡。

圖2 近56年干燥指數(shù)趨勢變化Fig.2 Trend variation of aridity index in recent 56 years

圖3 干燥指數(shù)多年月均值M-K檢驗Fig.3 M-K test results of monthly aridity index

2.2 干燥指數(shù)的突變分析

圖4 是近56年東勝礦區(qū)各季節(jié)和全年干燥指數(shù)M-K突變檢驗結(jié)果。圖中實線表示順序變化曲線UF,虛線表示逆序變化曲線UB,2條臨界線的值為±1.96(為顯著水平0.05下的t檢驗值)。若UF大于0則表明序列呈上升趨勢,反之則為下降趨勢;當超過臨界線時,表明上升、下降趨勢明顯。如果UF線和UB線超過y=±1.96信度線,并且2條曲線交于信度線之間,那么交點對應(yīng)的時間即是突變開始時刻[21,23]。圖中,雖然2條線在信度線之間有交點,但是都未突破信度線。說明近56年干燥指數(shù)并不存在統(tǒng)計上的突變點,但也不能忽視干燥指數(shù)發(fā)生的明顯的階段性變化。結(jié)合干燥指數(shù)多年變化趨勢可知,幾個交點中,1998年前后2個階段干燥指數(shù)均值發(fā)生的變化最明顯,1998年前后干燥指數(shù)分別為2.97和3.62,該年以前地表是逐漸變濕潤的,而之后則逐漸變干。

2.3 重現(xiàn)期分析

采用皮爾遜Ⅲ型頻率分布對1957—2012年東勝礦區(qū)干燥指數(shù)進行重現(xiàn)期分析,確定各干燥指數(shù)的概率分布和各種重現(xiàn)期下的干燥指數(shù)的極值。首先通過56年的干燥指數(shù)計算實際頻率值,點繪在機率格紙上,然后通過“適線法”確定皮爾遜Ⅲ型分布曲線3參數(shù)(均值離勢系數(shù)Cv、偏態(tài)系數(shù)Cs),從而得到皮爾遜Ⅲ型的理論頻率曲線。從圖5可以看出,實際頻率曲線和理論頻率曲線擬合較好,說明用皮爾遜Ⅲ型分布來推算該區(qū)不同干燥指數(shù)的重現(xiàn)期及不同重現(xiàn)期下干燥指數(shù)的極值估計結(jié)果是可信的。

圖4 多年干燥指數(shù)的M-K突變檢驗Fig.4 M-K abrupt change test of annual aridity index

用適線法確定皮爾遜Ⅲ型分布曲線后,即可算出各種概率P對應(yīng)的皮爾遜Ⅲ型分布的干燥指數(shù)理論值,50年一遇的干燥指數(shù)極大值和極小值分別為6.32和1.51,20年一遇干燥指數(shù)值分別為5.5和1.70,10年一遇干燥指數(shù)值分別為4.83和1.91。

2.4 周期特征分析

對干燥指數(shù)的時間序列進行小波分析,結(jié)果見圖6。Morlet小波系數(shù)實部圖表示東勝礦區(qū)干燥指數(shù)的時間尺度特征。圖中橫坐標為時間(年份),縱坐標為時間尺度(a),圖中的等值曲線為小波系數(shù)實部值。當小波系數(shù)實部值為正時,代表干燥指數(shù)高值期,在圖中用實線繪出,H表示正值中心,負值則表示干燥指數(shù)低值期,用虛線繪出,L表示負值中心。信號振蕩的強弱是通過等值線的密集程度來表示,等值線越密集表示比常年越干燥。

由圖6可以看出,干燥指數(shù)系列存在著時間尺度上的復雜嵌套結(jié)構(gòu),小尺度變化嵌套在較大尺度的變化下。較大尺度18～25年上的周期震蕩非常明顯,期間經(jīng)歷了逐漸變干—逐漸變濕的4個循環(huán)。其中,1957—1960年、1967—1974年、1982—1988年、1996—2002年、2009年以后地表有逐漸變干的趨勢,其余年份則有逐漸變濕潤的趨勢。值得注意的是,在18～25年較大尺度上,2009年以后的等值線圖未閉合,據(jù)此可推斷未來一段時間里還將繼續(xù)變干。中尺度變化以1985年為界,可分為2個階段:1985年以前以10～15年的尺度為主,存在3個振蕩循環(huán);1985年以后則以7～12年的尺度為主,存在4個振蕩循環(huán)。小尺度變化則以3～5年的尺度為主,震蕩周期比較短。

圖5 東勝礦區(qū)干燥指數(shù)的皮爾遜Ⅲ型分布曲線Fig.5 PearsonⅢdistribution curve of aridity index in Dongsheng mining area

圖6 干燥指數(shù)的Morlet小波實部圖Fig.6 Real part of Morlet wavelet of aridity index

圖7 干燥指數(shù)小波方差圖Fig.7 Wavelet variogram of aridity index

圖7 是東勝礦區(qū)近56年干燥指數(shù)的小波方差圖?？梢钥闯?小波方差圖中出現(xiàn)了4個較明顯的峰值:5、7、12和24年,其中第1震蕩主周期為12年,第2震蕩主周期為7年,第3、第4主周期分別為24和5年。這說明上述4個周期的波動控制著該地區(qū)干濕狀況在整個時間域內(nèi)的變化特征。

2.5 主導因子分析

近56年,降水量、日照時間、實際水汽壓與風速呈下降趨勢,但經(jīng)M-K趨勢檢驗,除風速外均未達到顯著水平;而最高氣溫和最低氣溫呈上升趨勢,二者均達到了顯著水平(P＜0.01)。近56年干燥指數(shù)的變化趨勢為0.004/a,通過全導數(shù)的方法求得各氣象因子對干燥指數(shù)年變化量的貢獻值大小,其中,降水量對干燥指數(shù)年變化的貢獻值為0.004 7,這可以解釋為降水量與干燥指數(shù)呈負相關(guān),而近56年降水量是呈減少趨勢,所以其對干燥指數(shù)的貢獻是正值。日照時間與干燥指數(shù)呈正相關(guān),而56年其呈減小趨勢,所以其對干燥指數(shù)的貢獻是負值,為-0.000 8。同理,可以得到最高氣溫、最低氣溫和風速對干燥指數(shù)年變化量的貢獻值分別為0.005 6、-0.002 7和-0.003 5。由表1可知,6個氣象因子中,水汽壓貢獻值的絕對值最小,僅為0.000 1,令其為1,計算各氣象因子對干燥指數(shù)變化的相對貢獻大小,從而得到影響近56年干燥指數(shù)變化的主導氣象因子主要為最高氣溫,其次為降水量和風速。

3 討論

始于20世紀80年代中期開發(fā)的神府東勝露天礦區(qū),到90年代逐漸進入開采高峰期。90年代初生產(chǎn)規(guī)模只有3 000萬t/a,在2002年突破5 000萬t之后,2003年飆升到7 384萬t,2005年一舉突破1億t大關(guān)達到1億241萬 t,2011年突破2億萬t[24-25]。90年代以來的這種大規(guī)模開采應(yīng)當對地表生態(tài)環(huán)境和區(qū)域氣候帶來影響,而突變檢驗的結(jié)果也很好地驗證了這點。1998年雖然不是統(tǒng)計意義上的突變年,但是,該年前后干燥指數(shù)的均值差異明顯,影響干燥指數(shù)的主導氣象因子發(fā)生了明顯變化:1998—2012年最高氣溫、降水量和風速這3個要素與1984—1997年相比都出現(xiàn)了較大差異,風速和降水量有所減小(前后2階段風速均值分別為2.97和 2.63 m/s,降水量均值分別為 391.6和367.82 mm),而日最高氣溫卻有所增大(分別為11.90和12.85℃)。

表1 各氣象因子的多年變化率及其對干燥指數(shù)的貢獻值Tab.1 Change rate of climate variables and contributions of climate variables to aridity index

煤礦開采特別是露天煤礦開采對下墊面影響很大,大面積開挖擾亂了地表土層和植被,嚴重破壞了原生地貌,而下墊面的大規(guī)模擾動可能通過破壞地表及地下水的循環(huán)系統(tǒng)而減少地表徑流量。相關(guān)研究發(fā)現(xiàn),1997年是東勝礦區(qū)所在窟野河徑流量減少的突變點[26-27],同時,窟野河流域的基流量占總徑流量的比例約為40%,但其基流分別在1980和1996年發(fā)生了2次減少的突變[28]。這些徑流的變化標示著區(qū)域蒸散量以及潛在蒸散量發(fā)生重要的變化,進而導致干燥指數(shù)發(fā)生改變。地表徑流變化與干燥指數(shù)變化的關(guān)系可能難以簡單量化,從突變年份比較相近來看,其中的關(guān)聯(lián)性是存在的,需要就降水量、徑流量、蒸散量及潛在蒸散量的關(guān)系做進一步的研究。此外,在較小空間尺度上,地形、植被等因素對確定干燥指數(shù)的若干氣象因子的影響會有所加強,這在地表有較大擾動情況下特別重要,其間的反饋與關(guān)聯(lián)機制,也還有待于進一步探討。

4 結(jié)論

1)內(nèi)蒙古東勝礦區(qū)干燥指數(shù)的多年平均值為3.2,上升幅度僅為0.004/a。除夏季有逐漸變干的趨勢外,其他季節(jié)都呈現(xiàn)出變濕的趨勢,但變化都不顯著。1998年可認為是一個突變年,該年以前地表是逐漸變濕潤的,而之后則逐漸變干。

2)50年一遇的干燥指數(shù)極大值和極小值分別為6.32和1.51,20年一遇干燥指數(shù)值為5.5和1.70,10年一遇干燥指數(shù)值分別為4.83和1.91。

3)干燥指數(shù)年值在較大尺度18～25年上的周期震蕩最明顯,期間經(jīng)歷了逐漸變干—逐漸變濕的4個循環(huán)。中尺度和小尺度則分別以7～12年和3～5年震蕩周期為主。第1震蕩主周期為12年,第2震蕩主周期為7年。

4)內(nèi)蒙古東勝礦區(qū)干燥指數(shù)近56年波動較小,干燥指數(shù)最大值和最小值的重現(xiàn)期都超過了50年。

5)影響近56年干燥指數(shù)變化的主導氣象為最高氣溫,其次為降水量和風速。