李 晶，崔綠園，閆蕭蕭，楊 震，董金瑋，鄧曉娟

(1. 中國礦業(yè)大學(xué)(北京)地球科學(xué)與測繪工程學(xué)院，北京 100083； 2. 中國科學(xué)院地理科學(xué)與資源研究所，北京 100101)

植被在地球上的能量交換和生態(tài)循環(huán)中起著重要作用[1]，對(duì)于生態(tài)脆弱的草原礦區(qū)而言，植被受氣候、采礦、放牧等人類活動(dòng)的影響較大。利用多光譜遙感數(shù)據(jù)得到的植被覆蓋度(fraction of vegetation coverage,FVC)作為衡量植被狀況、指示植被動(dòng)態(tài)變化的指標(biāo)[2- 4]應(yīng)用廣泛。文獻(xiàn)[5—9]分別對(duì)沙漠化地區(qū)、陰山北麓草原區(qū)、內(nèi)蒙古大青山自然保護(hù)區(qū)、甘南草原區(qū)、淮南礦區(qū)等地的植被覆蓋度動(dòng)態(tài)變化進(jìn)行了研究。從方法上看，利用線性混合理論對(duì)混合光譜進(jìn)行應(yīng)用和解釋是較為普遍的方法，像元二分模型在構(gòu)建時(shí)NDVIveg(NDVI for “pure” vegetation pixel)和NDVIsoil(NDVI for “pure” bare soil pixel)的取值也得到重視[10- 12]。在典型區(qū)域長時(shí)序植被覆蓋度變化研究中傳統(tǒng)的靜態(tài)多時(shí)相比較研究仍有應(yīng)用[13- 14]，但此種方法對(duì)于草原區(qū)FVC變化分析時(shí)易受到氣溫、降水等年際物候差異等的影響，難以科學(xué)揭示生態(tài)是否退化及其驅(qū)動(dòng)因素等根本問題。已有的研究中對(duì)FVC進(jìn)行了分級(jí)[15- 16]，這些基于分級(jí)的分析能夠在一定程度上表現(xiàn)區(qū)域植被演化特征，但難以進(jìn)一步提高定量表達(dá)的精度。此外，多數(shù)長時(shí)序植被覆蓋度研究中只采用一種趨勢分析方法，如回歸分析、Sen趨勢度和Manner- Kendall趨勢檢驗(yàn)相結(jié)合的方法等[17- 19]，多種分析方法的對(duì)比研究較少。

針對(duì)上述問題，本文以內(nèi)蒙古呼倫貝爾市寶日希勒露天礦區(qū)為例，基于1985—2015年年度最大合成NDVI數(shù)據(jù)，采用一元線性回歸法和Sen+Mann- Kendall法，進(jìn)行FVC的長時(shí)間序列分析，并對(duì)比兩種趨勢分析方法的差異。研究結(jié)果有助于獲取科學(xué)評(píng)價(jià)長時(shí)序煤炭開發(fā)活動(dòng)對(duì)地表生態(tài)影響的基礎(chǔ)數(shù)據(jù)，同時(shí)為類似區(qū)域長時(shí)序植被變化監(jiān)測及FVC變化趨勢分析提供方法參考。

1 研究區(qū)概況與數(shù)據(jù)獲取

1.1 研究區(qū)概況

寶日希勒礦區(qū)位于內(nèi)蒙古自治區(qū)呼倫貝爾市陳巴爾虎旗，地理坐標(biāo)為49°10′N—49°40′N，119°10′E—120°30′E，礦區(qū)面積約為255 km2。

考慮到礦區(qū)開采可能影響周邊區(qū)域，本文以自然地形界線和行政界線劃定研究邊界，即以莫爾格勒河和海拉爾河為研究區(qū)北、西、南的界線，東側(cè)以陳巴爾虎旗和海拉爾區(qū)的鄉(xiāng)鎮(zhèn)行政邊界為界線劃定研究區(qū)，總面積約為1792 km2，地勢由東北向西南逐漸降低，地跨森林草原與干旱草原兩個(gè)地帶。

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

本文原始數(shù)據(jù)為1985—2015年的Landsat TM、ETM+、OLI遙感影像數(shù)據(jù)。為消除年際間物候差異等的影響，經(jīng)中科院地理科學(xué)與資源研究所陸地表層格局與模擬實(shí)驗(yàn)室進(jìn)行影像處理后，得到各監(jiān)測年年度最大合成NDVI數(shù)據(jù)。由于1993年和2012年的圖像質(zhì)量不佳，因此研究中使用了29 a的年度最大合成數(shù)據(jù)。

年平均氣溫和年降水?dāng)?shù)據(jù)來源于資源環(huán)境數(shù)據(jù)云平臺(tái)(http:∥www.resdc.cn/)。

2 研究方法

2.1 植被覆蓋度計(jì)算

采用像元二分模型計(jì)算FVC，計(jì)算方法為[15]

(1)

式中，NDVIsoil為裸土或無植被覆蓋區(qū)域的NDVI值；NDVIveg為完全被植被所覆蓋像元的NDVI值。本文中依據(jù)NDVI頻率統(tǒng)計(jì)表，使用0.5%置信度確定NDVI值取值的置信區(qū)間。

2.2 植被覆蓋度時(shí)間變化特征分析

研究區(qū)在1985—2015年第i年各像元FVC均值采用下式計(jì)算

(2)

研究區(qū)FVC的時(shí)間變化特征分析采用轉(zhuǎn)移矩陣法，其中FVC等級(jí)劃分參考現(xiàn)有成果[16]并結(jié)合研究區(qū)實(shí)際情況確定分級(jí)界限，見表1。

表1 研究區(qū)植被覆蓋度分級(jí)等級(jí)

2.3 植被覆蓋度趨勢變化對(duì)比分析

本文采用一元線性回歸法[20]和Sen+Mann- Kendall法[21]進(jìn)行FVC趨勢變化分析。采用上述兩種方法結(jié)合進(jìn)行研究區(qū)的植被覆蓋度趨勢分析可以增強(qiáng)檢驗(yàn)的抗噪性，并在一定程度上提高檢驗(yàn)結(jié)果的準(zhǔn)確性[22]。

2.3.1 一元線性回歸法

一元線性回歸法計(jì)算如下[4]

(3)

式中，n為監(jiān)測時(shí)間。當(dāng)slope>0時(shí)，說明FVC呈增加趨勢，反之降低。

一元線性回歸法的趨勢檢驗(yàn)采用t檢驗(yàn)法[23]，計(jì)算方法如下

(4)

2.3.2Sen+Mann-Kendall法

Sen趨勢度計(jì)算如下[17]

(5)

式中,10表示時(shí)間序列呈上升趨勢，反之下降。

Mann- Kendall趨勢檢驗(yàn)統(tǒng)計(jì)量S的計(jì)算如下[19]

(6)

(7)

根據(jù)時(shí)間序列長度n值大小的不同，顯著性檢驗(yàn)統(tǒng)計(jì)量的選取有所不同。本文時(shí)間序列長度為29，采用檢驗(yàn)統(tǒng)計(jì)量Z來進(jìn)行趨勢檢驗(yàn)，檢驗(yàn)中取顯著水平α=0.05，Z1-α/2=Z0.975=1.96。

Z值計(jì)算方法如下

(8)

3 結(jié)果與分析

3.1 植被覆蓋度總體變化和多時(shí)相比較

研究區(qū)1985—2015年FVC均值的年際變化情況如圖1所示，F(xiàn)VC空間分布如圖2所示。

圖3選取FVC均值接近的1985、1995、2005和2015年4期數(shù)據(jù)進(jìn)行多時(shí)相FVC比較，并對(duì)1985年和2015年FVC作轉(zhuǎn)移矩陣。結(jié)果表明：研究區(qū)FVC年際差異大，1985—2015年高覆蓋的面積和中高覆蓋的面積共減少251.03 km2，其中157.23 km2變?yōu)榈椭当桓采w區(qū)域，52.87 km2變?yōu)闊o植被覆蓋區(qū)域，占總轉(zhuǎn)移面積的83.70%，F(xiàn)VC呈明顯下降趨勢，F(xiàn)VC低值區(qū)域也在不斷變化。通過氣象數(shù)據(jù)可知，1985、1995、2005和2015年研究區(qū)年平均氣溫分別為-1.814、1.312、0.211、0.743℃，年降水量分別為472.138、349.816、410.352、301.259 mm。結(jié)合圖2和氣象數(shù)據(jù)可得，F(xiàn)VC空間分布情況受氣溫和降水影響大，多時(shí)相的靜態(tài)對(duì)比無法準(zhǔn)確揭示草原的長時(shí)序變化趨勢，據(jù)此判斷草原是否退化等容易引起誤判。相比之下，遙感時(shí)序分析方法將時(shí)序影像視為連續(xù)時(shí)態(tài)記錄的完整檔案，可以更為精確地獲取地物在時(shí)序上的變化趨勢和規(guī)律[24- 26]。

3.2 植被覆蓋度空間演化趨勢分析

1985—2015年FVC變化趨勢程度如圖3所示。一元線性回歸法和Sen+Mann- Kendall法趨勢檢驗(yàn)分析結(jié)果見表2。

表2 一元線性回歸法和Sen+Mann- Kendall法趨勢檢驗(yàn)分析結(jié)果

趨勢分析方法植被覆蓋度變化趨勢植被覆蓋度變化程度像元個(gè)數(shù)百分比/(%)一元線性回歸法slope<0,α<0.01極顯著減少668933.35slope<0,0.01<α<0.05顯著減少680313.42α>0.05變化不顯著97221248.83slope>0,0.01<α<0.05顯著增加84504242.45slope>0,α<0.01極顯著增加387271.95Sen+Mann-Kendall法β<0,|Z|>1.96嚴(yán)重退化1294086.50β<0,|Z|<1.96輕微退化2337254.92-1.96≤Z≤1.96變化不顯著67080940.51β≥0,|Z|<1.96輕微改善89087044.75β>0,|Z|>1.96明顯改善661043.32

從圖3可得FVC下降最集中區(qū)域?yàn)檠芯繀^(qū)中部寶日希勒礦區(qū)所在位置，其他大部分區(qū)域的FVC變化斜率在0值左右，兩種趨勢分析方法得到的研究區(qū)FVC變化程度基本一致，但Sen+Mann- Kendall法在植被退化和植被改善區(qū)域面積有所增加。一元線性回歸法分析結(jié)果表明，植被極顯著增加的區(qū)域占1.95%，植被極顯著減少的區(qū)域占3.35%，48.83%的地區(qū)植被變化不顯著；Sen+Mann- Kendall法趨勢檢驗(yàn)分析結(jié)果表明，植被嚴(yán)重退化區(qū)域占6.50%，40.51%的地區(qū)植被變化不明顯，植被明顯改善區(qū)域占3.32%。

兩種趨勢分析法分析結(jié)果的空間分布差異如圖4所示。

由圖4可得，Sen趨勢度法和一元線性回歸法所模擬的植被覆蓋度變化趨勢基本上是一致的，兩者之差在(-0.035 5，0.023 6)，絕大多數(shù)區(qū)域無差異。兩種趨勢分析法得出植被變化趨勢相同的像元占比78.10%，Sen+Mann- Kendall法較一元線性回歸法而言，植被退化面積占比和植被改善面積占比分別多了4.65%、3.67%；負(fù)值和正值反映出Sen趨勢度法在植被退化和植被改善上要比一元線性回歸法更為敏感，克服了一元線性回歸法易受噪聲影響的缺點(diǎn)。Sen+Mann- Kendall方法對(duì)于植被退化和植被改善的選擇更為敏感，植被變化不顯著的區(qū)域減少，結(jié)果更加精確。兩種方法得出的結(jié)論基本一致，均為30年間研究區(qū)的植被覆蓋度的改善面積大于植被覆蓋度的退化面積；在草原礦區(qū)長時(shí)序植被覆蓋演化趨勢的分析中，Sen+Mann- Kendall方法能更準(zhǔn)確地反映其變化趨勢。

4 結(jié) 論

本文以寶日希勒礦區(qū)及周邊為研究區(qū)，應(yīng)用1985—2015年的Landsat TM、ETM+、OLI遙感影像數(shù)據(jù)進(jìn)行FVC的變化趨勢和空間差異分析，主要結(jié)論如下：

(1) 研究區(qū)平均FVC的線性趨勢并不十分顯著，整體呈波動(dòng)下降趨勢。

(2) 采用兩種趨勢方法分析30年間FVC演化趨勢，兩者結(jié)果基本一致，均為研究區(qū)中部寶日希勒礦區(qū)所在位置FVC退化明顯。

(3) 綜合比較兩種FVC趨勢分析方法的結(jié)果，Sen+Mann- Kendall法較一元線性回歸法而言，植被退化面積占比和植被改善面積占比分別多4.65%、3.67%，該方法在FVC改善和退化的反映上更為敏感，植被變化不顯著的區(qū)域減少，能監(jiān)測出輕微的植被變化趨勢，趨勢變化結(jié)果更加精確，更能準(zhǔn)確地反映草原礦區(qū)長時(shí)序的FVC演化趨勢。