中國東部草原植被綠度時空變化分析及其對煤電基地建設的響應

包 巖,田 野,柳彩霞,范文義,付 曉

1 東北林業大學林學院,哈爾濱 150040 2 中國科學院遙感與數字地球研究所遙感科學國家重點實驗室,北京 100101 3 中國科學院生態環境研究中心城市與區域生態國家重點實驗室,北京 100085

草原植被是陸地生態系統的重要組成部分[1-2],具有明顯的季節變化特點,植被的動態變化可以從植被綠度變化的趨勢上得以反映[3],地表植被呈現的綠度變化與植被覆蓋率緊密相關[4]。應用遙感技術可以計算植被指數(Vegetation Index,VI),從而了解植被覆蓋情況。內蒙古自治區東部草原,常年降水量少,植被覆蓋率較低,因此運用歸一化植被指數(Normalized Difference Vegetation Index,NDVI)不會存在飽和現象[5],即隨著植被綠度的增大NDVI值也隨之增大。

AVHRR(Advanced Very High Resolution Radiometer)NDVI數據集是目前覆蓋時段最長的連續數據集[6],尤其是GIMMS(Global Inventory Modeling and Mapping Studies)NDVI 3g數據集,具有時間序列長、覆蓋范圍廣等優點[7],是目前應用最廣泛的數據集[8]。鄭藝等[9]利用GIMMS NDVI數據,分析了全球干旱區植被時空變化及其對氣候和人文等非氣候因子的響應,結果表明,1982—2012年干旱區NDVI總體呈顯著增長趨勢；賀振等[10]利用1983—2013年的GIMMS NDVI時序數據,探究了黃河流域植被覆蓋時空格局和演化趨勢,研究顯示,黃河流域32年來植被覆蓋出現了遲緩的增加態勢,以0.018/(10 a)的速度增長；Shen等[11]在GIMMS NDVI數據的基礎上應用多元逐步回歸分析中國溫帶草原1982—2005年生長季期間的NDVI,表明GIMMS NDVI能準確估計溫度對中國溫帶草原地區植被的影響。Anyamba等[12]用NDVI 3g數據集更新了Sahel地區的植被動態和趨勢分析,指出該數據可以用來分析和解釋半干旱地區數10年尺度上的地表植被變化趨勢。可以看出,GIMMS NDVI數據集在不同地區都有廣泛應用,用以揭示長時間草原植被的變化規律和驅動因素[8]。

內蒙古東部草原區坐落在我國生態安全“兩屏三帶”的北部防沙地區,受氣候條件制約和地區產業發展影響,具有酷寒、半干旱、土壤貧瘠等生態脆弱性特征。近年來,也有部分學者對這一地區展開研究。如李云鵬等[13]應用1988年以來的遙感數據,獲得了各時間段內的NDVI值,并且利用遙感數據和光譜特征獲得了土地退化植被,進而探究了內蒙古自治區植被覆蓋度的空間變化和及其荒漠化的狀況。張超等[14]根據GIMMS NDVI數據集,對內蒙古1982—2006年的農地、叢林和草原三種植物類型的NDVI分別進行了探究。研究了不同植被類型下NDVI的演變差異,進而探討了植被覆蓋的空間變化規律。張圣微等[15]應用MODIS數據提取了NDVI值,探究了2003—2012年生長季(4—10月份)錫林郭勒草原不同類型的草地植被的變化趨勢,研究表明,草原植被變化由東北向西南呈現出下降的趨勢。葉永昌等[16]基于最大熵模型(MaxEnt)探究了氣象因子與植被散布的關系,進而分析了1961—2010年內蒙古草原植被的空間散布狀況,并且利用綜合模型模擬了凈第一生產力的變化規律。研究結果顯示,氣候因子中的濕潤指數(MI)、年降水量(P)、最暖月的平均溫度(Tw)和最冷月的平均溫度(Tc)是草原植被遍布的主要制約因素。目前,還沒有相關研究報道利用公里空間分辨率NDVI數據從宏觀尺度分析大型煤電基地建設對草原植被的影響,GIMMS NDVI 3g數據集是否可用于分析區域煤電基地對草原生態系統的影響尚還不清楚。

最近幾年,隨著煤炭礦區建設范圍的持續擴大,地區生態環境壓力持續增長,植被綠度在不停地惡化。東部草原區是我國以露天開采為主的重要大型煤電基地,煤炭產能超過4億噸,為我國東北部能源供應提供了保障,但同時,也造成了地下水下降、水土破壞、區域植被退化等生態問題,煤電開發與生態保護矛盾日漸突出,對我國生態環境造成了不小的威脅[17-18]。研究數據顯示,中國露天煤礦每開采一萬噸煤,破壞的土地面積為0.22 hm2,年均損毀土地面積多達1×104hm2[19]。煤炭的露天和地下開采,都會引起煤電基地含水層水位下降[20-21]。認識煤炭礦區環境的狀態,并科學合理地分析引起煤炭礦區環境變化的首要因素,基于此,可以提出加強煤炭礦區的生態環境改善的有效方法,找出緩解生態壓力、改善生態環境問題的合理途徑。本文應用GIMMS NDVI 3g數據集對東部草原區植被綠度進行長時間序列分析,從公里空間尺度上探究我國東部草原區植被綠度變化與煤電基地開發的響應關系,為東部草原煤電一體化開發過程的草原生態修復提供宏觀信息。

1 研究區概況

內蒙古呼倫貝爾市和錫林郭勒盟,地處溫帶北部,屬于大陸性季風氣候,草原資源豐富,草場種類繁多。研究區總面積為54.1554×104km2,位于115°13′—126°04′ E、43°02′—53°20′ N之間。其中呼盟坐落在內蒙古高原的東北邊,由一些高原、山地和平原低地組成了它的整體部分,其地貌的主體主要由大興安嶺山地組成。錫盟的地勢特點主要是海拔較高的平原,并且伴有不同地貌類型的區域,地勢的走向是南部較高,北部較低,東部和南部則是以低山丘陵為主,中間散落著零星的盆地,西部和北部的地勢較平緩(圖1)。研究區不僅擁有豐富的草原資源,而且礦產資源豐富,據統計,2016年全國原煤累計產量為3.364×1010t,其中內蒙古原煤產量就達到了8.38×109t,占全國的24.9%,內蒙古已超越山西,成為中國煤炭產量的第一大省,而呼盟和錫盟的煤炭產能在省內分別排名第二和第三位。

圖1 研究區位置意圖Fig.1 Location map of the study area and coal mine locations

2 數據和方法

2.1 數據來源與預處理

本論文所使用的遙感數據是1981年7月到2010年12月的GIMMS NDVI 3g數據,該數據來源于美國NASA發布的基于NOAA氣象衛星數據全球數據集,采用經過輻射校正和幾何粗校正的NOAA-AVHRR數據,再進一步對每天、每軌影像進行幾何靜校正、除壞線、除云等預處理[22],進而進行NDVI計算及合成。計算公式為NDVI=1000×(b2-b1)/(b2+b1),其中b1、b2為AVHRR的第1、2通道。用該方法可以得到30 a內720幅NDVI半月合成的影像。該NDVI數據集的時間分辨率為16 d,空間分辨率為8 km。NDVI的取值范圍為-1—1,一般認為大于0的為有植被覆蓋區域[23]。

利用呼盟和錫盟的矢量邊界數據從NDVI序列中提取研究區1981—2010年的年NDVI時間序列。在矢量邊界數據的基礎上,結合Google Earth高清影像,在研究區內確定煤礦24個(圖1),在ArcGIS 10.1中對礦區數據進行轉化與提取,并分別建立半徑10、20 km和50 km的緩沖區,得到緩沖區的矢量圖,進而提取緩沖區1981—2010年的年NDVI時間序列。

2.2 研究方法

2.2.1 生長季最大值提取

采用最大值合成法(Maximum Value Composites,MVC)獲得月NDVI最大值,可以進一步消除云和太陽高度角的影響。在此基礎上,通過月NDVI最大值,獲取研究區內1981—2010年間每個像元年生長季(5月—8月)最大值NDVI,來代表該像元點生長季NDVI的最大值,以探究NDVI年際變化特征。通過分析NDVI的生長季最大值,了解東部草原區1981—2010年30年間的NDVI空間分布整體格局,得到東部草原區植被綠度的時空變化趨勢。計算公式如下：

MNDVIij=MAX(NDVIij1,NDVIij2)

(1)

在此基礎上,求出生長季的NDVI最大值。公式為：

GNDVIi=MAX(MNDVIi5,MNDVIi6,MNDVIi7,MNDVIi8)

(2)

式中i為年序號,取值范圍是1981—2010；j為月序號,取值范圍是1—12；MNDVIij為第i年第j月的NDVI最大值；GNDVIi為第i年的生長季NDVI的最大值；NDVIij1,NDVIij2分別為第i年第j月上半月和下半月的NDVI值。

2.2.2 趨勢線分析

趨勢線分析法可以較好地分析柵格中每個柵格像元的變化規律,可全面系統的體現研究區內植被覆蓋時空的格局變化特點。本論文應用此方法來分析1981—2010年東部草原區生長季NDVI最大值的變化趨勢。

為了探究NDVI隨年份的變化情況,本文利用一元線性回歸分析來分析逐個柵格像元的變化規律。擬合NDVI相對于年份的直線方程,得到每個柵格NDVI與年份的回歸模型,從而獲得一幅30年間變化斜率的影像。對于每一個像元點,采用最小二乘(Ordinary Least Squares,OLS)估計,建立植被指數的生長季最大值(GNDVI)與年份(YEAR)的線性關系,擬合線性方程：

GNDVI=SLOPE×YEAR+B

(3)

直線方程的斜率SLOPE代表了該像元點植被綠度的變化規律,進而用來作為輸出影像中該像柵格的值。若1981—2010年植被NDVI是增加的態勢,則擬合斜率大于0,表示植被綠度往好的方向發展；若1981—2010年植被NDVI是減少的態勢,則擬合斜率小于0,即植被綠度往不好的方向發展。

2.2.3 相關性分析

本文擬將選出的礦區與緩沖區年際GNDVI值進行相關分析,為了檢測礦區對周邊植被綠度變化的影響,設置了3個級別的緩沖區,并獲得了30年間緩沖區的GNDVI值及其變化斜率。采用Pearson相關系數分析礦區與緩沖區GNDVI的相關性,公式如下：

(4)

式中,Pearson相關系數是用協方差與兩個變量的標準差得到的,能夠更好的度量兩個變量之間的相關程度,取值范圍為-1—1,當兩個變量的線性關系越大時,相關系數的絕對值越趨近1。

3 結果

3.1 東部草原區植被綠度變化趨勢

1981—2010年呼盟和錫盟GNDVI線性擬合斜率的灰度影像見圖2,將獲得的影像進行分級統計,以0為界分為兩級。其中呼盟的GNDVI線性擬合斜率的取值范圍為-0.007—0.034,斜率大于0的區域占總面積的40.84%,斜率小于0的占59.16%,綠度變化總體上為減少趨勢,小于0的區域相對比較集中,主要分布在呼盟的西部和東南部地區,即呼倫貝爾草原地區；錫盟擬合斜率的取值范圍為-0.006—0.026,斜率大于0的區域占總面積的26.87%,小于0的區域占總面積的73.13%,植被綠度減少非常明顯,退化嚴重,小于0的區域主要分布錫盟北部和東南部?？梢钥闯?無論是呼盟還是錫盟,植被綠度整體上均呈下降趨勢,錫盟植被綠度下降尤為嚴重,東部草原區植被退化明顯,這與葉永昌等[16]在2016年所研究的結果相似,他的研究表明1961—2010年內蒙古草甸草原、典型草原和荒漠草原分布的面積分別降低了5%、1%和62%,內蒙古草原面積整體下降了11%,說明內蒙古草原植被退化明顯,草原向著荒漠化的趨勢發展。

圖2 1981—2010年東部草原區綠度變化Fig.2 Greenness variation in the eastern grassland from 1981 to 2010

3.2 礦區及緩沖區植被綠度變化趨勢分析

統計分析了1981—2010年30年間24個礦區的植被生長季最大NDVI值及其與年份的變化斜率,如圖3所示。從圖中可以看出,除東明礦區和扎賚諾爾靈泉礦區外,其余礦區的斜率都呈明顯的下降趨勢。大雁1、大雁3、呼盛、天順、白音華、查干淖爾、賀斯格烏拉南、勝利西和勝利錫凌礦區的回歸斜率都小于其緩沖區的回歸斜率,其余礦區的回歸斜率都大于其緩沖區的回歸斜率。本論文在礦區周圍分別設置了半徑為10、20 km及50 km的緩沖區,并分析了緩沖區GNDVI的年際變化,發現緩沖區的GNDVI值也都呈現動態下降的趨勢,斜率均是小于0的。隨著緩沖區范圍的增大,GNDVI值也表現出增加的趨勢,這也說明了礦區對周邊地區植被綠度的影響。

將礦區與相應緩沖區的GNDVI值進行相關性分析(圖4),可以看出,除大雁3礦區外,各組數據的相關系數總體較大,并且在0.05置信水平上顯著相關,說明緩沖區GNDVI值與礦區的GNDVI值整體變化趨勢相同(圖3),并且隨著緩沖區范圍的擴大,相關系數的值呈現明顯降低的趨勢(大雁3礦除外),說明礦區的GNDVI值在一定程度上影響著緩沖區的GNDVI值。圖3中可以看出,呼盟礦區及緩沖區GNDVI總體大于錫盟的GNDVI,說明錫盟的植被綠度比呼盟的差。

圖3 礦區及緩沖區GNDVI變化趨勢Fig.3 GNDVI variation trend of the mining area and buffer zone

圖4 礦區GNDVI均值與緩沖區GNDIV均值的相關性 Fig.4 Correlations of mean GNDVI between coal mines and their buffers 陰影突出區域為井工礦

3.3 東部草原區植被綠度分級

通過對礦區及其緩沖區的分析可知,除東明礦區和扎賚諾爾靈泉礦區外,其余礦區及其相應緩沖區的斜率都呈明顯的下降趨勢。說明礦區周圍植被的綠度在減少,植被呈現退化的趨勢。我們將斜率劃分為5個閾值(表1—2),將劃分的閾值做成分布圖(圖5),可以更直觀地看出東部草原區30年來植被綠度變化情況。

從圖5、表1和表2中可以看出,東部草原區植被綠度總體呈現減少的趨勢。其中,呼盟在1981—2010年30年間,植被綠度減少(顯著減少和輕度減少)的像元所占的百分比分別為0.765%和41.83%,減少的像元散落在呼盟的各個方位,植被綠度增加(輕度增加和顯著增加)的像元所占的百分比分別為26.54%和0.415%；錫盟30年間植被綠度減少(顯著減少和輕度減少)的像元所占的百分比分別為4.31%和51%,主要分布在錫盟的東部和北部,植被綠度增加(輕度增加和顯著增加)的像元所占的百分比分別為21.7%和0.09%,增加的區域主要分布在錫盟的西北部。錫盟減少的像元所占的百分比(55.31%)多于呼盟減少的像元所占的百分比(42.595%),增加的像元所占的百分比(21.79%)少于呼盟增加的像元所占的百分比(26.955%),說明錫盟的植被綠度情況要比呼盟的植被綠度情況差。

表1 呼倫貝爾市GNDVI斜率分級統計

表2 錫林郭勒盟GNDVI斜率分級統計

圖5 東部草原區植被GNDVI時空變化分布Fig.5 Spatial and temporal distribution of GNDVI in the eastern grassland

3.4 GNDVI殘差分析

為了分析礦區的采礦活動對草原植被綠度的影響,本文進行了GNDVI的殘差分析[24]。首先將GNDVI與氣候因子進行回歸分析,在每個柵格上建立GNDVI與氣候因子的回歸模型,然后利用這個回歸模型,可以估計每個柵格上每年的GNDVI值,用實測的GNDVI值減去估計的GNDVI值,這樣就獲得了1981—2010年每年的GNDVI殘差,最后對GNDVI殘差與其對應的年份進行一元回歸方程的計算,獲得的斜率作為輸出影像的像元值,即可以獲得一幅30年間GNDVI殘差的年際變化分布圖,如圖6所示。

在沒有采礦等人為活動擾動的狀況下,隨著時間的變化,殘差值應在0值上下呈現波動的變化趨勢。如果殘差的變化斜率呈現出減少的趨向,即擬合斜率小于0,則說明采礦活動使草原的植被綠度下降,加劇了草原的荒漠化；如果殘差的變化斜率是上升的趨向,則說明礦區的開采改善了草原的生態環境,提高了草原的植被綠度[25]。

從圖6中可以看出,1981—2010年東部草原區GNDVI的殘差斜率小于0的區域占得比重較多,其中呼盟GNDVI殘差斜率小于0的區域占總面積的60.99%,分布在呼盟的各個地域,錫盟GNDVI殘差斜率小于0的區域占總面積的55.02%,主要分布在錫盟的東北部,這與植被綠度的空間變化趨勢相一致(圖5)?？梢钥闯?30年間,東部草原區的植被破壞多于植被恢復,植被綠度在下降。結合礦區的位置,可以發現,寶日希勒礦區、東明礦區和勝利礦區的GNDVI殘差斜率均是大于0的,說明這些礦區的煤礦開采并沒有加劇植被綠度的退化,GNDVI在這些地區主要受氣候因子的影響。而在其他礦區,GNDVI殘差的斜率都是小于0的,說明礦區的開采加劇了植被綠度的退化,這些地區的植被綠度不僅受氣候因子(溫度、降水)的影響,還與礦區開采等人類活動因素有關。

圖6 1981—2010年GNDVI殘差變化時空分布圖Fig.6 Spatial and temporal distribution of GNDVI′s residual in the eastern grassland

4 討論

分析結果顯示,所選取的24個礦區(呼盟15個,錫盟9個)除東明和扎賚諾爾靈泉礦區外,其余礦區GNDVI的擬合斜率均是小于0的,植被的綠度在降低,說明礦區對植被綠度影響為負。礦區的開采在一定程度上影響了植被的綠度變化趨勢[26-27]。植被綠度的增加可能主要和氣候因子有關,氣候適宜、降水增加都可以引起植被綠度的改善。另外,煤礦開采區土地復墾和植被修復對植被綠度的增加也有貢獻,如圖5,寶日希勒礦區和勝利礦區的植被綠度有增加趨勢。而植被綠度的減少主要是自然因素與人為因素的影響結果,自然因素方面,東部草原區常年干旱,降水量減少。造成土地荒漠化,植被覆蓋度降低,從而造成植被綠度的下降；人為因素方面,人口增長過快,過度放牧,不合理地開采煤礦資源等,這些都能引起生態環境的退化,使植被綠度快速下降[28]。

觀察圖3可以看出,煤炭的露天和井工開采,都對草原生態系統植被綠度有負面影響。83.33%的礦區在2005—2007期間GNDVI的值都呈明顯下降趨勢,這與它們的開采年份有關。其中,五牧場和白音華礦區是于2005年開采,在之后的5年間,它們的GNDVI值呈現下降的趨勢；呼盛、蒙西、天順、賀斯格烏拉南礦區是于2006年開始開采,在2007年它們的GNDVI值呈跳躍式下降；除扎賚諾爾靈泉、勝利東2、敏東和西一礦區,其余礦區在開采后的GNDVI值都呈現下降的趨勢,這說明煤礦的露天和井工開采都對草原植被綠度有著負影響。扎賚諾爾靈泉露天礦區30年來GNDVI值呈增長的趨勢,這可能與它開采時間較長有關,靈泉露天礦于1936年的偽滿時期進行開采,與本文選取的時間段相差較大,造成了結果的不確定性；而勝利東2、敏東和西一礦區的開采時間在2008年前后,與本文的時間段重合較小,其結果也不能反映礦區GNDVI值的趨勢變化情況。

從圖3和圖4中可以看出,呼盟礦區GNDVI與不同緩沖區的相關性差異相比錫盟更明顯,可能的原因是呼盟的生長季植被綠度水平高于錫盟,一旦煤礦開采,NDVI下降程度更劇烈,因此對草原的植被影響更敏感。實驗表明(圖4),無論是煤礦的露天開采還是井工開采,都會對礦區周邊的植被綠度產生影響,即離礦區越遠的區域(緩沖區半徑大的區域),植被綠度的減少越少。大雁3礦區情況稍有反常,可能的原因是大雁礦區已經停產,采礦設施均已拆除,在原地及周邊建設了礦山公園。

東部草原區植被綠度的退化可能還與氣候因素有關[9, 29],因此,本研究分析東部草原區15個氣象站點的氣候數據,數據來源于中國氣象數據網(http://data.cma.cn/)。結合高程數據,對15個氣象站點生長季的降水和溫度數據進行克里金插值[30],從而得到東部草原區1981—2010年生長季氣候數據。將呼盟和錫盟各個像元年GNDVI做平均值,得到年平均GNDVI,與氣候數據進行對比分析(圖7),可以看出,整個東部草原區在30年間,生長季的降水量呈下降趨勢,溫度有緩慢升高趨勢,錫盟的降水量明顯低于呼盟的降水量,這也是呼盟礦區的GNDVI整體都比錫盟的GNDVI高的原因之一。

圖7 東部草原區NDVI值與氣候因子關系圖Fig.7 Relationship between NDVI and climatic factors in the eastern grassland

5 結論

本文將東部草原區的呼盟與錫盟分開進行統計分析,從宏觀尺度上更能直觀地反映兩個區域內植被綠度的時空變化趨勢,為東部草原煤電基地開發區植被生態修復提供區域差異信息,做到礦區生態系統修復的因地制宜[31]。綜合比較可以發現,在植被綠度的覆蓋情況、GNDVI的年際變化斜率、生長季降水量各個因素中,呼盟的表現值都要好于或高于錫盟的表現值。這是由于呼盟有充沛的水資源,且降水量較高,溫度適宜,有利于植被的生長,所以植被綠度明顯好于錫盟,因此在呼盟,煤礦開采對草原植被綠度的影響沒有錫盟明顯。但是,不論露天還是井工開采,都對植被綠度有負面影響。本研究公里空間分辨率的遙感數據無法區分開采方式對植被綠度影響的差異,今后還需采用中高分辨率遙感數據深入評估。