東勝礦區(qū)近50年地表干濕狀況時(shí)程變化分析

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

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

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

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

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

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

1 研究方法

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

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

1.2 干燥指數(shù)

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

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

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

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

1.3 Mann-Kendall趨勢(shì)檢驗(yàn)和突變檢驗(yàn)

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

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

1.4 水文頻率分析

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

1.5 小波分析

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

1.6 主導(dǎo)因子和貢獻(xiàn)值分析

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

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

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

2 結(jié)果與分析

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

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

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

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

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

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

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

2.3 重現(xiàn)期分析

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

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

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

2.4 周期特征分析

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

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

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

圖6 干燥指數(shù)的Morlet小波實(shí)部圖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個(gè)較明顯的峰值:5、7、12和24年,其中第1震蕩主周期為12年,第2震蕩主周期為7年,第3、第4主周期分別為24和5年。這說明上述4個(gè)周期的波動(dòng)控制著該地區(qū)干濕狀況在整個(gè)時(shí)間域內(nèi)的變化特征。

2.5 主導(dǎo)因子分析

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

3 討論

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

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

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

4 結(jié)論

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

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個(gè)循環(huán)。中尺度和小尺度則分別以7～12年和3～5年震蕩周期為主。第1震蕩主周期為12年,第2震蕩主周期為7年。

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

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