基于MODIS植被豐度的沙漠化界線提取方法探究
——以陜甘寧地區(qū)為例

馬靜宇，秦 怡，任宇鵬，王 鑫，谷 娟，郭仲皓，鄒樂君，沈曉華

(浙江大學地球科學學院，浙江 杭州 310000)

沙漠化是指包括氣候變化和人類活動等各種因素造成的干旱、半干旱和亞濕潤干旱地區(qū)的土地退化過程[1]，主要表現為植被面積的減少和沙漠面積的增加[2-3]。植被覆蓋度作為刻畫陸地表面植被數量的重要參數，是沙漠化過程土地自然屬性和生態(tài)變化最直觀的反映，可以被用來大致劃分沙漠區(qū)和非沙漠區(qū)，而且遙感手段已經能夠定量提取其信息，已有不少學者在沙漠化研究中用植被覆蓋度對沙漠化進行監(jiān)測與評價[4-5]。

沙漠化實質是一種退化過程，如果能夠刻畫出沙漠與非沙漠地區(qū)的界線，可以更直觀、準確地反映沙漠化過程，但目前相關研究還非常少。Lamprey[6]最早用100 mm等降水量線作為植被與沙漠的分界線，對撒哈拉沙漠進行監(jiān)測。Tucker等[7]利用NOAA/AVHRR NDVI數據推算降水量，以200 mm等降水量線作為撒哈拉沙漠的邊界。國內學者王澄海[8]通過對我國1981—2001年平均NDVI的0.12等值線和200 mm年均降水等值線位置的定性判斷，用NDVI的0.12等值線作為沙漠化邊界。這種基于降水量數據的沙漠化劃界方法實現了對沙漠化范圍的定性表達，但缺乏對沙漠化形成機制的科學討論，缺少理論依據。Yu等[9]、劉愛霞[10]提出以區(qū)域內是否出現植被萌芽期事件作為沙漠與非沙漠區(qū)的界線，對蒙古戈壁沙漠和我國西部干旱半干旱區(qū)進行沙漠化劃界。由于NDVI是植被長勢的反映，在干旱半干旱區(qū)植被生長易受降水等影響，植被的高波動使NDVI的取值有較大不確定性，導致基于NDVI確定的沙漠化界線不斷波動，誤導人們對沙漠化的真實判斷[11]。

MODIS是搭載在TERRA和AQUA衛(wèi)星上的一個重要傳感器，其空間分辨率為250～1000 m，具有1 d的高時間分辨率，為大面積、長時間的沙漠化實時動態(tài)監(jiān)測研究提供了新的遙感數據源。本研究基于對MODIS進行混合像元分解獲取的像元植被豐度，提出了累積曲線梯度變化法，通過分析每一年的累積曲線的梯度變化，找到反映植被生長狀態(tài)改變的梯度突變點，從而確定能夠表達沙漠化界線的相對指標，進而以陜甘寧地區(qū)為例，開展沙漠化界線的識別與分析，以期對沙漠化的動態(tài)監(jiān)測提供新的途徑。

1 研究區(qū)與試驗數據

研究區(qū)位于北緯35°—40°，東經105°—112°(如圖1所示)，包括陜西的榆林、延安兩市，甘肅的慶陽市，以及寧夏的石嘴山市、銀川市、吳忠市、中衛(wèi)市和固原市。研究區(qū)深處內陸，大部分地區(qū)處于干旱、半干旱氣候區(qū)，降水量少且分布不均勻，地形、地貌復雜。由于生態(tài)環(huán)境脆弱，以及人類不合理的土地利用，使得該地區(qū)土地沙漠化問題極為突出。

試驗數據包括：①MODIS L1B數據，空間分辨率為1 km，時間范圍為2001年1月1日—2015年12月31日，完成輻射定標、大氣校正、幾何校正以及云掩膜。②Landsat 8數據，空間分辨率為30 m，行列號為121/037，獲取日期為2015年10月11日，并進行輻射定標、大氣校正和幾何校正。

圖1 研究區(qū)地理位置

2 研究方法

2.1 提取像元植被豐度

MODIS遙感影像滿足了對沙漠化監(jiān)測的高時間分辨率要求，但由于空間分辨率低、混合像元現象嚴重，需要通過混合像元分解將每個像元分解為它的覆蓋類型組分占該像元的百分含量(即豐度)[12]。

本文采用郭峰博士提出的基于地物波譜特征的混合像元分解模型[13]來獲取植被豐度，即利用三次Bernstein基函數構建水體、植被、裸地三者地物的理論特征波譜曲線，用于表征混合光譜中水體、植被、裸地的相對含量；再以Landsat 8的水體、植被、裸地的分類結果作為參照豐度，建立像元光譜分解指數與實際豐度之間的定標模型，消除基函數曲線特征與地物光譜特征之間的差異，獲取地物的真實豐度。本文研究建立的定標模型通過F檢驗，得到定標參數和檢驗參數的顯著性概率均大于99%，表明可以很好地用于植被豐度提取，選取每月一景無云或少云的MODIS影像數據共180景數據，對每年12個月的植被豐度分別進行疊加分析，剔除5%異常值，以有效去除云、雨等天氣影響，并用年最大植被豐度表示該年的植被覆蓋情況。

2.2 生成不同植被豐度的像元累積頻率年變化曲線

植被豐度反映的是單個像元的植被覆蓋情況，以單個像元為變化單位的影像，變化頻率太高，而且在遙感圖像中，相同植被豐度的像元分布零散，不利于刻畫不同生長環(huán)境下植被分布與生長漸變狀況，因此，需先對相同植被豐度的像元進行累加，獲得不同植被豐度的像元累積頻率，即

(1)

式中，Pi表示植被豐度為i%時的像元累積頻率；ni表示不大于植被豐度值i%的累積像元數；N表示一幅遙感影像的總像元數。

以年為橫坐標，像元累積頻率為縱坐標，可以繪制不同植被豐度的像元累積頻率年變化曲線，簡稱累積曲線(如圖2所示)。

在累積曲線中，同一年內，相鄰植被豐度的像元累積頻率變化量能夠在該環(huán)境(主要是水環(huán)境)影響下的植被生長狀態(tài)：在低植被豐度的沙漠地區(qū)，由于環(huán)境惡劣，植被分布稀疏、成活率低、不易生長，因此單位植被豐度內的像元累積頻率的變化值較大(累積曲線分布較為稀疏)；而在非沙漠地區(qū)，植被分布密集且易于生長，像元累積頻率的單位變化值較小(累積曲線分布較為密集)。通過總結分析多年植被生長規(guī)律，可以初步確定沙漠與非沙漠的界線范圍，有助于下一步對沙漠化界線的判斷。

圖2 陜甘寧地區(qū)不同植被豐度的像元累積頻率年變化曲線

2.3 確定沙漠化界線指標

從沙漠到非沙漠地區(qū)，本質上是水資源空間再分配的過程[14]：在沙漠地區(qū)，植被生長所需要的水分主要來源于天然降水；在非沙漠地區(qū)，受農田灌溉等影響，人為干預占主導地位。這兩種不同生長背景下的植被由于受到的生長動力因素較為單一，因而分別呈現較為穩(wěn)定且相同的生長狀態(tài)和變化趨勢，那么在較難生長和易于生長之間，就存在植被的一個生長過渡狀態(tài)。研究表明，在沙漠與非沙漠過渡帶內，以水為主導因子的環(huán)境梯度產生植被梯度現象[15-16]，打破了植被原有較為穩(wěn)定的生長狀態(tài)，加之過渡區(qū)的植被既不像干旱區(qū)植物有很強抗旱能力，又不像濕潤區(qū)植物有很強繁衍能力，因而植被變化比相鄰區(qū)域更迅速，在數學中往往用梯度來反映這一規(guī)律。梯度是某函數f(x,y)的一階導數，定義為

(2)

根據梯度的性質，在同種物質內部，數值基本相同且變化不大，而在不同物質交界處，數值變化十分迅速，梯度的局部變化率非常大[17]，反映在累積曲線上即在某一范圍植被豐度內形成突變梯度，梯度變化最大處所對應的位置即可近似為沙漠與非沙漠的分區(qū)界線。在實際操作中，用兩條擬合直線的交點近似作為該梯度變化點，取滿足條件的多年x值的平均值，作為研究區(qū)的沙漠化界線指標。

3 結果與討論

3.1 陜甘寧研究區(qū)的沙漠化界線指標

分析圖2可以看出，低植被豐度的像元累積頻率變化基本穩(wěn)定，處于沙漠與非沙漠地區(qū)過渡帶時，由于兩種性質不同的生態(tài)系統(tǒng)的交匯與重疊，植被對氣候的適應和環(huán)境的調控能力弱，加之人類活動也會對該區(qū)域植被產生明顯的影響，如2003年與前一年相比，像元累積頻率的年變化趨勢出現了反向特征，這很可能是由于人們采取的植樹等水土保持措施使該范圍內的像元數量增多。隨著氣候、土壤等條件逐漸變好，灌溉設施更加完善，人為因素逐漸占主導，在該環(huán)境下植被生長又基本保持穩(wěn)定。因此確定全區(qū)沙漠化界線的大致范圍在植被豐度為10%～20%處。

以植被豐度1%為間隔，那么相鄰植被豐度的像元累積頻率之差就是不同植被豐度對應的像元頻率，分析2001—2015年每一年累積曲線的梯度變化，以2001年為例(如圖3所示)，在植被豐度10%～20%范圍內，利用兩條擬合直線的交點估計其梯度變化點在植被豐度為15%處，通過統(tǒng)計各年梯度變化點，取平均值15%作為陜甘寧地區(qū)的沙漠化界線指示指標。

圖3 梯度變化點的確定

3.2 沙漠化界線的野外驗證

在本文研究期間，筆者曾于2016年8月10日—8月11日，對寧夏南部地區(qū)(105°E—107°E，37°N—40°N) 71個采樣點進行實地考察(如圖4所

示)。獲取2016年的MODIS影像數據，繪制2016年寧夏南部地區(qū)植被豐度等值線分布圖，結合采樣點的空間分布，重點分析植被豐度為5%、10%、15%和20%附近采樣點的地表景觀特征。圖4左圖為2016年10月2日MODIS 250 m分辨率的假彩色合成影像。

植被豐度小于10%的地區(qū)(如圖4中采樣點88、75、131所示)，植被種類較為單一且結構簡單，植被覆蓋稀疏；植被豐度介于10%至20%之間的區(qū)域(如圖4中采樣點72、107、108所示)，天然草地和農田并存，天然降水和農田灌溉成為植被生長所需水分來源，植被生長受到自然和人類的雙重作用，其中72號采樣點中整齊分布的樹苗可以反映出為防治風沙運動、減輕沙漠化而采取種樹等人工措施；植被豐度大于20%的地區(qū)(如圖4中采樣點122、125、103所示)，農田廣泛分布，植被覆蓋度高且植被生長旺盛，肥沃的土壤、濕潤的氣候條件等給植被的生長創(chuàng)造了優(yōu)越的環(huán)境。

圖4 寧夏北部野外采樣點分布圖及其地表狀況

通過野外調查，進一步證實了由于過渡區(qū)域水源的變化，人為與自然營力的共同作用，造成了該區(qū)域植被生長狀態(tài)和植被覆蓋情況的明顯變化，確定了以植被豐度15%作為沙漠化界線的合理性，間接驗證了累積曲線梯度變化法對沙漠化評價結果的有效性。

3.3 陜甘寧研究區(qū)沙漠化動態(tài)監(jiān)測

圖5為2001—2015年植被豐度等于15%的等值線平面圖，即沙漠化界線空間分布圖。通過對15%線包圍面積的定性比較可以看出，2001—2015年沙漠化面積整體減少，但是在15年間沙漠化界線不斷波動，且變化區(qū)間主要在內蒙古自治區(qū)與陜西榆林地區(qū)的交界地帶，以及寧夏南部地區(qū)。這些地區(qū)歷來就是農牧交錯演化頻繁的地區(qū)，由于不同景觀要素相互作用，加之人為活動引起的植被變化較大，使該地區(qū)具有高敏感和生態(tài)脆弱性。

圖5 2001—2015年陜甘寧研究區(qū)沙漠化界線空間分布

繪制植被豐度為15%的像元累積頻率年變化曲線(如圖6所示)，其大小反映了所包含沙漠化土地的面積(像元累積頻率與像元面積的乘積)，由于降水對陜甘寧地區(qū)農業(yè)發(fā)展和生態(tài)系統(tǒng)的恢復有著極為重要的影響，因此，結合由美、日聯(lián)合研制和發(fā)射的熱帶測雨衛(wèi)星(tropical rainfall measuring mission，TRMM)反演的TRMM 3B42降水數據繪制研究區(qū)年均降水量曲線，綜合評價陜甘寧地區(qū)的沙漠化過程。

圖6 陜甘寧地區(qū)年降水量和植被豐度15%的像元累積頻率變化

總體上看，2001—2015年間，陜甘寧地區(qū)的沙漠化在波動中總體呈現不斷改善的趨勢。2001—2004年，沙漠化土地面積明顯減少，這段時間恰好是1999年退耕還林生態(tài)工程建設的開始實施階段，生態(tài)治理工程初見成效；2004—2005年間沙漠化土地面積增加，這可能與該年段陜甘寧地區(qū)整體降水條件有關，從2003年開始的低降水量對植被長勢造成深刻影響[18]，植被生長受到限制，部分植被因缺水干枯死亡，抗御風沙的能力被削弱。盡管沙漠化受氣候波動影響，但人類活動的影響也逐漸上升到主導地位[19]，當沙漠化程度加重時，人們一方面通過灌溉等保證農作物和植被的生長，另一方面通過植樹造林進行固沙，可以看到在2005年后沙漠化得到及時、明顯的改善；2006—2011年基本處于穩(wěn)定狀態(tài)，到2012年沙漠化面積明顯減少，由于在該時期降水沒有明顯的變化，沙漠化的逆轉趨勢很可能受人為活動影響所致，但由于缺乏長效保護機制，在2013、2014年沙漠化程度明顯加劇，到2015年才得以進一步的控制。

沙漠與非沙漠地區(qū)的過渡區(qū)域，生態(tài)環(huán)境的脆弱性和過渡特征會使沙漠化程度加劇，而降水增多等有利的自然條件也會使沙漠化有自然逆轉的可能；人工干擾和破壞會導致沙漠化的加劇，但同時人為治理也會促進沙漠化土地的逆轉。2001—2015年陜甘寧地區(qū)沙漠化逆轉趨勢較為明顯，但是植被退化的沙漠化現象依舊存在，人為保護措施有待進一步恢復、建設和加強，從而創(chuàng)造植被保持穩(wěn)定生長的有利條件，鞏固治沙成果。

4 結 語

采用全覆蓋、高時間分辨率、低空間分辨率的MODIS遙感影像，可獲得具有時間序列的植被豐度數據，為實現大范圍沙漠化動態(tài)過程監(jiān)測提供有效途徑。其中的關鍵在于如何定量劃定沙漠與非沙漠地區(qū)的界線。本文提出的基于像元植被豐度的累積曲線梯度變化法為沙漠化劃界提供了一種客觀有效的方法，借助該方法，可利用MODIS影像數據實現對沙漠化時空分布與動態(tài)變化的快速監(jiān)測與定量評價。通過對陜甘寧地區(qū)沙漠化界線的變化過程分析，發(fā)現2001—2015年該地區(qū)的沙漠化總體呈現在波動中不斷改善的趨勢，盡管降水會通過影響植被生長影響沙漠化的發(fā)生與逆轉，人類活動的干預也起到了明顯的作用。

參考文獻：

[1] 朱震達，劉恕.關于沙漠化的概念及其發(fā)展程度的判斷[J].中國沙漠，1984，4(3)：2-8.

[2] 曾永年，向南平，馮兆東，等.Albedo-NDVI特征空間及沙漠化遙感監(jiān)測指數研究[J].地理科學， 2006，26(1)：75-81.

[3] GOUVEIA C M，PASCOA P，RUSSO A，et al.Land Degradation Trend Assessment over Iberia during 1982—2012[J].Cuadernos De Investigacion Geografica，2016，42(1)：89-112.

[4] 劉廣峰，吳波，范文義，等.基于像元二分模型的沙漠化地區(qū)植被覆蓋度提取——以毛烏素沙地為例[J].水土保持研究，2007，14(2)：268-271.

[5] DAWELBAIT M，MORARI F.Monitoring Desertification in a Savannah Region in Sudan Using Landsat Images and Spectral Mixture Analysis[J].Journal of Arid Environments，2012，80(80)：45-55.

[6] LAMPERY H F.Report on the Desert Encroachment Reconnaissance in Northern Sudan[R].Khartoum，Sudan：National Council for Research，Ministry of Agriculture，Food and Resources，1975.

[7] TUCKER C J，DREGNE H E，NEWCOME W W.Expansion and Contraction of the Sahara Desert from 1980 to 1990[J].Science，1991，253(5017)：299-301.

[8] 王澄海，惠小英.以植被指數0.12為指標看我國的荒漠化與草原界限的變化[J].中國沙漠，2005，25(1)：88-92.

[9] YU F，PRICE K P，ELLIS J，et al.Interannual Variations of the Grassland Boundaries Bordering the Eastern Edges of the Gobi Desert in Central Asia[J].International Journal of Remote Sensing，2004，25(2)：327-346.

[10] 劉愛霞.中國及中亞地區(qū)荒漠化遙感監(jiān)測研究[D].北京：中國科學院研究生院(遙感應用研究所)， 2004.

[11] 韋振鋒，任志遠，張翀，等.1999—2010年陜甘寧黃土高原區(qū)荒漠化空間特征[J].中國沙漠，2014，34(5)：1230-1236.

[12] KESHAVA N，MUSTRAD J F.Spectral Unmixing[J].IEEE Signal Processing Magazine，2002，19(1)：44-57.

[13] 郭峰.基于Bernstein基函數的遙感影像混合像元分解研究[D].杭州：浙江大學，2015.

[14] 趙成義，王玉朝，李國振.荒漠-綠洲邊緣區(qū)研究[J].水土保持學報，2001，15(3)：93-97.

[15] 馬瑞，王繼和，屈建軍，等.民勤綠洲-荒漠過渡帶植被梯度特征及其土壤水環(huán)境[J].生態(tài)學雜志， 2010，29(6)：1075-1080.

[16] 何小武，劉廣全，郭孟華.黃土高原植被建設的水資源環(huán)境及對策[J].水利學報，2008，39(7)：843-847.

[17] KINDLMANN G，DURKIN J W.Semi-automatic Gener-ation of Transfer Functions for Direct Volume Rendering[M].Los Alamitos：IEEE Computer Society Press，1998: 79-86.

[18] 盧曉杰，張克斌，李瑞.我國北方農牧交錯帶氣候對植被的影響——以寧夏鹽池為例[J].水土保持研究，2007，14(6): 193-197.

[19] 王濤，朱震達.我國沙漠化研究的若干問題——沙漠化的概念及其內涵[J].中國沙漠，2003，23(3)：209-214.