基于MODIS數(shù)據(jù)的羌塘國家自然保護(hù)區(qū)荒漠化遙感反演與評價

侯富強(qiáng)，李國明，闞璦珂

（1.成都理工大學(xué) 信息科學(xué)與技術(shù)學(xué)院，成都 610059；2.成都理工大學(xué) 地球物理學(xué)院，成都 610059）

基于MODIS數(shù)據(jù)的羌塘國家自然保護(hù)區(qū)荒漠化遙感反演與評價

侯富強(qiáng)1，李國明2，闞璦珂2

（1.成都理工大學(xué) 信息科學(xué)與技術(shù)學(xué)院，成都 610059；2.成都理工大學(xué) 地球物理學(xué)院，成都 610059）

基于MODIS數(shù)據(jù)采用因子分析法，選取MSAVI、NDVI、LST、TVDI四個指標(biāo)反演羌塘國家自然保護(hù)區(qū)荒漠化動態(tài)變化過程，并進(jìn)行分級評價。結(jié)果表明：①近十年間保護(hù)區(qū)荒漠化面積總體呈增加趨勢，中度和極重度荒漠化面積呈減小趨勢，輕度、重度荒漠化面積呈增加趨勢；②若拉核心區(qū)、緩沖區(qū)東部無明顯荒漠化現(xiàn)象，實驗區(qū)與緩沖區(qū)西部荒漠化較為嚴(yán)重，美馬錯核心區(qū)、瑪依核心區(qū)、緩沖區(qū)中部荒漠化程度較輕。

荒漠化；因子分析法；動態(tài)變化；分級評價

0 前言

荒漠化作為當(dāng)今世界十大環(huán)境問題之首，倍受國際社會的關(guān)注［1］。我國是世界上荒漠化面積大、分布廣、受其危害最嚴(yán)重的國家之一，西藏又是全國荒漠化較為嚴(yán)重的省份。藏北羌塘地區(qū)作為荒漠化監(jiān)測和防治的典型區(qū)域，研究評價其荒漠化程度對高寒生態(tài)脆弱地區(qū)具有重要意義。作者在本文基于MODIS數(shù)據(jù)，采用因子分析法，選取了MSAVI、NDVI、LST、TVDI四個指標(biāo)反演羌塘國家自然保護(hù)區(qū)荒漠化動態(tài)變化過程，并進(jìn)行分級評價。

1 研究區(qū)概況

羌塘國家自然保護(hù)區(qū)（以下簡稱“保護(hù)區(qū)”，見下頁圖1）位于西藏西北部，昆侖山、可可西里山以南，岡底斯山以北，面積為29.8×104km2。范圍涉及西藏那曲地區(qū)西部三縣（安多、尼瑪、班戈）和阿里地區(qū)北部三縣（日土、改則、革吉）。根據(jù)《西藏自治區(qū)第三次荒漠化和沙化監(jiān)測報告》資料可知以上區(qū)域是荒漠化較為嚴(yán)重的地區(qū)之一。

2 遙感反演指標(biāo)選取

2.1 改進(jìn)型土壤調(diào)整植被指數(shù)（MSAVI）

植被指數(shù)是根據(jù)植被反射波段的特性，計算反映地表植被生長狀況、覆蓋情況、生物量和植被種植特征的間接指標(biāo)［2、3］。MSAVI作為改進(jìn)型土壤調(diào)節(jié)植被指數(shù)，一般認(rèn)為將土壤背景的影響減至最低，可增強(qiáng)對植被的敏感性。MSAVI指標(biāo)采用植被指數(shù)數(shù)據(jù)產(chǎn)品MOD13Q1在2001年～2010年數(shù)據(jù)：MOD13Q1空間分辨率為250m×250m，時間分辨率為16d。保護(hù)區(qū)覆蓋h24v05、h25v05共兩個區(qū)域。

2.2 歸一化植被指數(shù)（NDVI）

歸一化植被指數(shù)（NDVI）能夠敏感地反映出植被生物量、覆蓋度，以及葉綠素含量等生物理化性質(zhì)，在一定程度上反映了像元所對應(yīng)區(qū)域的土地覆蓋狀況［4］，因而在土地覆蓋變化研究中，常被用來描述土地覆蓋的基本特征以及開展相關(guān)的定量研究。選用MOD13Q1的2001年～2010年數(shù)據(jù)，地面分辨率250m，每16天一次。

圖1 羌塘國家自然保護(hù)區(qū)位置與區(qū)劃Fig.1 The Qiangtang national nature reserve location and regionalization

2.3 陸地表面溫度（LST）

對于植被茂密的地表，遙感反演所得到的地表溫度，是指植被葉冠的表面溫度。對于植被稀疏的地表，地表溫度是地面、植被葉冠等溫度的混合平均值［5］。采用MODIS數(shù)據(jù)的MOD11A2產(chǎn)品（2001年～2010年），時間分辨率8d，空間分辨率為1km的陸地表面溫度產(chǎn)品，包含白天地表溫度、夜間地表溫度、31波段和32波段通道發(fā)射率等資料。

2.4 溫度植被旱情指數(shù)（TVDI）

溫度植被旱情指數(shù)為地表溫度與NDVI的綜合反映，被定義為基于NDVI－LST特征空間提取出的水份脅迫指標(biāo)，較單獨的地表溫度指示度有所提高［6］。由于植被覆蓋度與光譜植被指數(shù)存在一定關(guān)系，而植被覆蓋度決定了傳感器接收到土壤背景和植被冠層的可見光和熱紅外信息，從而影響了遙感影像獲取的輻射溫度［7］。土壤水份狀況與表面溫度之間不存在直接的關(guān)系，但土壤水份無疑是影響植被冠層溫度的重要因素［8］。因此，一定植被覆蓋條件下的冠層溫度，能夠間接反映土壤含水情況。

3 數(shù)據(jù)處理與分析

3.1 方法流程（見下頁圖2）

（1）針對MODIS數(shù)據(jù)的數(shù)據(jù)結(jié)構(gòu)及通道特點，實現(xiàn)遙感數(shù)據(jù)的預(yù)處理、空間幾何校正。

（2）針對MSAVI、NDVI、LST、TVDI四個指標(biāo)進(jìn)行荒漠化信息提取，對荒漠化過程進(jìn)行反演。

（3）對各指標(biāo)2010年與2001年荒漠化結(jié)果進(jìn)行差值分析，完成基于各指標(biāo)的荒漠化變化分布制圖。

（4）分級評價羌塘自然保護(hù)區(qū)荒漠化空間分布和程度。

3.2 數(shù)據(jù)處理

3.2.1 MSAVI數(shù)據(jù)處理

使用ENVI提供的IDL接口，利用綠色植被在不同波段的光譜反射特征進(jìn)行計算，由用戶自定義函數(shù)計算得到。利用歸納法可推導(dǎo)得MSAV表達(dá)式，如公式（1）所示。

圖2 方法流程Fig.2 Workflows of methods

式中 ρnir和ρred分別代表MODIS第1波段和第2波段的反射率。

反演結(jié)果見下頁圖3（a）。

3.2.2 NDVI數(shù)據(jù)處理

利用MRT（MODIS ReproJection Tool）對MOD13Q1－NDVI數(shù)據(jù)進(jìn)行地理幾何校正與重采樣批處理，提取NDVI數(shù)據(jù)。同時，為了更加有效地消除云遮蔽、大氣影響、觀測中的幾何關(guān)系等不利因素的影響，對保護(hù)區(qū)MOD13Q1－NDVI數(shù)據(jù)進(jìn)行MVC（最大值合成法）處理，如公式（2），以獲取每月NDVI最大值。

式中 NDVIi為第i月最大NDVI值；NDVIij為第i月第j個16－day合成NDVI值。然后，根據(jù)式（3），計算植被生長季累積NDVI。

式中 NDVIi為植被生長季累積NDVI；NDVI（i，k）為植被在第i年第k月最大NDVI值，i＝2001、2002、…、2010；k＝1、2、…、12。

反演結(jié)果下頁圖（b）。

3.2.3 LST數(shù)據(jù)處理

陸地表面溫度是通過建立31通道、32通道亮溫線性組合的分裂窗算法計算來獲取。

MODIS數(shù)據(jù)共有八個熱紅外波段，根據(jù)分裂窗算法原理，MODIS數(shù)據(jù)選取第31、32兩個波段來反演地表溫度。采用覃志豪等［8］提出的，適用于MODIS數(shù)據(jù)的地表溫度反演算法，該算法如公式（4）所示。

式中 LST是地表溫度（K）；T31和T32分別是MODIS第31波段和32波段的亮度溫度。根據(jù)這兩個波段的圖像DN值和普郎克公式來計算；A0A1、A2是分裂窗算法的參數(shù)，分別定義如下：

a31、b31、a32、b32是常量，可?。?/p>

a31＝－64.603 63；

b31＝0.440 817；

a32＝－68.725 75；

b32＝0.473 453。

其它中間參數(shù)分別計算如下：

其中 ti（q）是i（i＝31、32）波段視角為q的大氣透過率。

反演結(jié)果如圖3（c）所示（見下頁）。

3.2.4 TVDI數(shù)據(jù)處理

溫度植被干旱指數(shù)（TVDI）的定義為［6］式（10）。

式中 LST為陸地表面溫度；Tmax＝a1＋b1＊NDVI，為某一NDVI對應(yīng)的最高溫度，即旱邊；Tmin＝a2＋b2＊NDVI，為某一NDVI對應(yīng)的最低溫度，即濕邊；a1、b1和a2、b2是旱邊和濕邊的線形擬合方程系數(shù)，由線性擬合得到；NDVI為植被指數(shù)。

反演結(jié)果如圖3（d）所示（見下頁）。

3.3 結(jié)果討論

（1）MSAVI差值。值減小最為嚴(yán)重的地區(qū)在保護(hù)區(qū)的東南部，即實驗區(qū)的東部、緩沖區(qū)的東北部，具體區(qū)域分布在班戈縣，在保護(hù)區(qū)西南部表現(xiàn)為增加。

圖3 2010年與2001年MSAVI/NDVI/LST/TVDI差值圖Fig.3 The difference map of MSAVI/NDVI/LST/TVDI in 2001and 2010

（2）NDVI差值。負(fù)值越小，荒漠化程度越重正值越大，荒漠化程度越小。從下頁圖3中可以看出，其值減小最為嚴(yán)重的地區(qū)在保護(hù)區(qū)的東部，即實驗區(qū)、緩沖區(qū)的東部、若拉核心區(qū)，主要區(qū)域分布在班戈縣。

（3）LST差值。為便于研究，可將差值分為三類：①小于－1.00為溫度下降區(qū)；②大于1.00為溫度升高區(qū)；③之間的區(qū)域為微變區(qū)，其中升溫集中的區(qū)域為尼瑪縣。

（4）TVDI差值。由于小于“0”的下降點非常少，所以不予考慮。在整個保護(hù)區(qū)范圍內(nèi)，TVDI指數(shù)都有不同程度上升，這表明旱情的增強(qiáng)。旱情明顯增大的地區(qū)包括保護(hù)區(qū)的東北部、東南部和中西部地區(qū)，其它地區(qū)也有零星分布，其中數(shù)值較高的區(qū)域分布在班戈縣境內(nèi)。

4 分級評價

通過荒漠化程度樣本圖（國家林業(yè)局荒漠化檢測中心繪制的2009年與2000年中國土地荒漠化圖），作者與選取四個指標(biāo)的MODIS數(shù)據(jù)進(jìn)行GIS疊加分析，統(tǒng)計不同荒漠化程度的各指標(biāo)值。通過對像元數(shù)值進(jìn)行分析，確定保護(hù)區(qū)荒漠化分級評價體系（見下頁表1）。

根據(jù)荒漠化分級評價體系計算各個指標(biāo)的荒漠化分級面積，最后通過加權(quán)求平均值法求得2001年～2010年間保護(hù)區(qū)荒漠化程度分級統(tǒng)計面積及其變化趨勢，如圖4所示。

圖4 羌塘自然保護(hù)區(qū)2001年～2010年荒漠化面積分級變化趨勢Fig.4 The desertification area grading change trends o the reserve since 2001to 2010

保護(hù)區(qū)非荒漠化面積、中度荒漠化面積、極重度荒漠化面積呈減少趨勢，輕度荒漠化面積、重度荒漠化面積呈增加趨勢，且非荒漠化面積、極重度荒漠化面積減小較為明顯，重度荒漠化面積增加較為明顯。通過上述四種指標(biāo)反演結(jié)果評價，保護(hù)區(qū)東部地區(qū)，即班戈縣的大部份地區(qū)荒漠化程度較輕。該區(qū)域主要分布在濕潤、亞濕潤氣候類型區(qū)降雨量較為豐富，但荒漠化程度卻呈現(xiàn)較為嚴(yán)重趨勢，同期該地區(qū)的年降水量也呈逐漸減少的趨勢保護(hù)區(qū)西南部，即改則縣中部、革吉縣北部、日土縣東部荒漠化情況較為嚴(yán)重，地處干旱、半干旱氣候類型區(qū)，降雨較少，植被以草地為主，但荒漠化情況有所好轉(zhuǎn)，年降水量增加。改則縣北部，尼瑪縣的大部，班戈縣的少部份變化較小，其氣候比較濕潤，降雨量均勻。影響荒漠化過程的除生態(tài)和氣候等因素外，還來自于人類開發(fā)活動的干擾，過墾、過牧等，都是影響保護(hù)區(qū)荒漠化的可能因素。

表1 保護(hù)區(qū)荒漠化分級評價體系Tab.1 Desertification grading evaluation system of the reserve

5 結(jié)論

作者在本文中采用MODIS數(shù)據(jù)作為主要數(shù)據(jù)源，選取MSAVI、NDVI、LST、TVDI四項指標(biāo)，對羌塘國家自然保護(hù)區(qū)近十年來荒漠化動態(tài)變化進(jìn)行遙感反演，為高寒地區(qū)荒漠化遙感監(jiān)測提供了可循的依據(jù)和方法。結(jié)果表明：

（1）非荒漠化面積在2001年～2005年間，以每年0.17×104km2的速度增加，在2005年～2010年間，以每年0.12×104km2的速度減少，總體呈減少趨勢；中度和極重度荒漠化面積分別以每年0.03×104km2、0.05×104km2的速度減??；輕度、重度荒漠化面積分別以每年0.05×104km2、0.07 ×104km2的速度增加。

（2）若拉核心區(qū)、緩沖區(qū)的東部降雨量較為豐富，溫度適中，沒有呈現(xiàn)明顯荒漠化現(xiàn)象；實驗區(qū)的西部、緩沖區(qū)的西部降雨較少，溫度較高，荒漠化較為嚴(yán)重；美馬錯核心區(qū)、瑪依核心區(qū)、緩沖區(qū)中部荒漠化程度較輕。

［1］ 王莉雯，衛(wèi)亞星，牛錚，等.基于GIS的沙漠化土地面積遙感分析監(jiān)測［J］.地球信息科學(xué)，2007，9（3）：132

［2］ 王濤.干旱區(qū)綠洲化、荒漠化研究的進(jìn)展與趨勢［J］中國沙漠，2009，29（1）：1.

［3］ 張紅，蘇香玲.克拉瑪依市土地荒漠化評價及動態(tài)監(jiān)測［J］.安徽農(nóng)學(xué)通報，2010，16（13）：235.

［4］ GU J，LI X，HUANG C，et al.A simplified data assimilation method for reconstructing time－serie modis ndvi data［J］.Advances in Space Research 2009，44（7）：501.

［5］ SHENG J F，WILSON J P，LEE S.Comparison o land surface temperature（LST）modeled with a spatially－distributed solar radiation model（SRAD）and remote sensing data［J］.Environmental Modelling ＆Software，2009，42（3）：436.

［6］ 李正國，王仰麟，周清波，等.基于TVDI方法的黃土高原地表水分時空分異研究［J］.中國農(nóng)業(yè)資源與區(qū)劃，2008，29（2）：13.

［7］ 周金星，張華，周澤福，等.大興縣荒漠化生態(tài)安全評價研究［J］.干旱地區(qū)農(nóng)業(yè)研究，2009，27（1）：200.

［8］ 覃志豪，高懋芳，秦曉敏，等.農(nóng)業(yè)旱災(zāi)監(jiān)測中的地表溫度遙感反演方法——以MODIS數(shù)據(jù)為例［J］.自然災(zāi)害學(xué)報，2005，14（4）：64.

book=45,ebook=45

1001—1749（2012）03—0353—05

X 144

A

10.3969/j.issn.1001－1749.2012.03.21

侯富強(qiáng)（1953－），男，工程師，從事計算機(jī)應(yīng)用研究。

國家自然科學(xué)基金項目（40761018）

2012－04－17 改回日期：2012－04－19