基于MODIS數(shù)據(jù)的黃河源區(qū)植被覆蓋度反演研究

劉 晶，黨星海,2，陳麗麗，邊雁君，黃立鑫，李志紅，周兆葉*

（1.蘭州理工大學(xué) 土木工程學(xué)院，甘肅 蘭州 730050；2.甘肅省應(yīng)急測(cè)繪工程研究中心，甘肅 蘭州 730050；3.中鐵二十局集團(tuán)市政工程有限公司，甘肅 蘭州 730030）

高寒草地是維持高原地區(qū)生態(tài)平衡的關(guān)鍵物種，高寒草地的退化將導(dǎo)致不可逆轉(zhuǎn)的生態(tài)問(wèn)題，因此對(duì)于高寒草地的監(jiān)測(cè)至關(guān)重要。高寒草地的監(jiān)測(cè)方法主要包括傳統(tǒng)的地面調(diào)查和利用植被指數(shù)進(jìn)行遙感估 測(cè)[1]，地面調(diào)查方法費(fèi)時(shí)費(fèi)力，樣本數(shù)量少，空間代表性差，得出的結(jié)論不足以很好地反映高寒草地覆蓋度的真實(shí)空間分布狀態(tài)[2-3]；遙感估測(cè)的常用手段包括衛(wèi)星遙感、無(wú)人機(jī)遙感等[4]。利用無(wú)人機(jī)能低空飛行、懸停等特點(diǎn)進(jìn)行取樣[5-6]，并配合少量地面取樣，可成為樣方尺度高寒草地覆蓋度、高度、生物量調(diào)查的主要手段，在很大程度上解決傳統(tǒng)觀測(cè)方法的不足[7-8]。 無(wú)人機(jī)操作靈活，可通過(guò)設(shè)置其飛行高度和飛行模式獲取與遙感像元尺度匹配的圖像，能架起樣地實(shí)測(cè)數(shù)據(jù)與遙感影像之間的橋梁。目前基于植被指數(shù)的遙感估測(cè)方法包括回歸方程法、像元二分法和機(jī)器算 法等[1,7]。

本文以黃河源區(qū)的高寒草地為研究對(duì)象，以無(wú)人機(jī)航拍為獲取數(shù)據(jù)的主要手段，結(jié)合地面樣方調(diào)查數(shù)據(jù)和MODIS數(shù)據(jù)提取了植被指數(shù)；再利用回歸模型法和像元二分法反演了源區(qū)高寒草地覆蓋度；最后篩選最優(yōu)指數(shù)和算法進(jìn)行黃河源區(qū)高寒草地覆蓋度估算以及變化趨勢(shì)分析，以期為高寒草地的合理利用與生態(tài)治理提供數(shù)據(jù)支撐。

1 數(shù)據(jù)來(lái)源與研究方法

1.1 數(shù)據(jù)來(lái)源

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

本文所指的黃河源區(qū)為唐乃亥水文站以上區(qū)域，面積為1.2×105km2，主要植被類(lèi)型為高寒草地、高寒灌叢化草地和高寒草甸等[8-9]。該區(qū)域植被生長(zhǎng)季節(jié)性變化明顯，對(duì)氣候變化極其敏感，屬于生態(tài)脆弱帶和氣候變化敏感區(qū)[10-12]。

本文采用大疆精靈Pro四旋翼無(wú)人機(jī)，航高設(shè)置為20 m，樣地大小布設(shè)為250 m×250 m（以便與MODIS數(shù)據(jù)像元大小匹配），樣地內(nèi)以矩形的方式獲取照片，每塊樣地共獲取13張照片，單張照片覆蓋范圍約為 35 m×35 m[3]。2016－2018年7月底8月初在黃河源區(qū)共拍攝326塊樣地，4 238張照片。樣地位置分布如圖1所示。

圖1 研究區(qū)航拍樣地位置分布

1.1.2 MODIS數(shù)據(jù)下載

根據(jù)黃河源區(qū)的經(jīng)緯度范圍，確定行列號(hào)為h25v05、h26v05的兩景數(shù)據(jù)可覆蓋全區(qū)域。該數(shù)據(jù)的時(shí)間分辨率為16 d，空間分辨率為500 m。本文下載數(shù)據(jù)時(shí)間序列為2001－2018年7月11日（以便與野外航拍時(shí)間一致）的MOD11A1和MOD13A1數(shù)據(jù)，共72景影像。

1.2 研究方法

無(wú)人機(jī)樣地植被覆蓋度采用項(xiàng)目團(tuán)隊(duì)自主開(kāi)發(fā)的覆蓋度提取軟件進(jìn)行計(jì)算。首先利用MODIS Tool軟件對(duì)下載的源區(qū)MOD13A1數(shù)據(jù)進(jìn)行預(yù)處理，包括定義投影、裁切等，并提取、計(jì)算得到歸一化植被指數(shù)（NDVI）、增強(qiáng)型植被指數(shù)（EVI）、土壤調(diào)節(jié)植被指數(shù)（SAVI）和修正的土壤調(diào)節(jié)植被指數(shù)（MSAVI）[3]； 然后利用回歸模型法和像元二分模型，分別基于4種植被指數(shù)對(duì)黃河源區(qū)植被覆蓋度進(jìn)行反演，并比較各植被指數(shù)采用上述兩種方法的反演精度；最后比較采用同種反演方法，哪種植被指數(shù)的反演精度更高就采用該植被指數(shù)進(jìn)行逐年反演，并分析黃河源區(qū) 2001－2018年植被覆蓋度的變化情況。

2 黃河源區(qū)植被覆蓋度反演

2.1 黃河源區(qū)樣地植被蓋度提取

本文采用研究團(tuán)隊(duì)自主開(kāi)發(fā)的Pixel Based Manual Classifier V1軟件提取航拍照片的植被覆蓋度和非植被等信息。該軟件能利用閾值法快速、高效、準(zhǔn)確地區(qū)分植被和非植被，并根據(jù)分類(lèi)結(jié)果進(jìn)行植被覆蓋度的計(jì)算。

2.2 植被覆蓋度反演

2.2.1 植被指數(shù)數(shù)據(jù)獲取

首先從經(jīng)過(guò)預(yù)處理的MOD13A1數(shù)據(jù)中提取NDVI、EVI，再提取紅光波段和近紅外波段，用以獲取SAVI和MSAVI[4]。

2.2.2 基于植被指數(shù)的回歸模型法

基于植被指數(shù)的回歸模型法首先建立植被指數(shù)與實(shí)測(cè)植被覆蓋度的關(guān)系，再結(jié)合衛(wèi)星遙感數(shù)據(jù)將該關(guān)系反推到整個(gè)黃河源區(qū)，從而獲得黃河源區(qū)的植被覆蓋度數(shù)據(jù)。

本文首先分別計(jì)算2001－2018年各種植被指數(shù)的均值，再采用“提點(diǎn)法”得到各植被指數(shù)的像元值，并利用其中75%的樣本點(diǎn)數(shù)據(jù)與植被覆蓋度建立擬合關(guān)系，然后通過(guò)擬合關(guān)系進(jìn)行反演，最后利用25%的樣本點(diǎn)進(jìn)行驗(yàn)證，比較其反演精度。

1）利用75%的植被覆蓋度對(duì)各植被指數(shù)均值進(jìn)行采樣，并建立植被指數(shù)與實(shí)測(cè)植被覆蓋度之間的散點(diǎn)圖，得到擬合關(guān)系式，如圖2所示。

圖2 植被指數(shù)與實(shí)測(cè)FVC建立的擬合關(guān)系

2）在ArcGIS中，分別利用各植被指數(shù)的擬合關(guān)系式反演整個(gè)黃河源區(qū)的植被覆蓋度，結(jié)果如圖3所示。

圖3 回歸模型法各植被指數(shù)的反演結(jié)果

3）利用剩余的25%的樣本點(diǎn)數(shù)據(jù)進(jìn)行精度驗(yàn)證，結(jié)果如圖4所示。

圖4 回歸模型法精度驗(yàn)證

2.2.3 像元二分模型法

像元二分模型法是將像元中的地物組成分為植被與非植被兩部分，通過(guò)獲取純植被像元與非純植被像元的光譜信息來(lái)達(dá)到分解混合像元的目的。通過(guò)分解可計(jì)算得到該混合像元中植被因子所占的比例。以NDVI為例，植被覆蓋度的計(jì)算公式為：

式中，NDVIsoil為純土壤像元對(duì)應(yīng)的NDVI值；NDVIveg為純植被像元對(duì)應(yīng)的NDVI值；NDVI為需要計(jì)算植被覆蓋度的混合像元的NDVI值。

具體取值方法為：若有實(shí)測(cè)數(shù)據(jù)，則分別取實(shí)測(cè)植被覆蓋度的最大值和最小值作為FVCmax和FVCmin，遙感影像中的NDVI值作為NDVImax和NDVImin，與實(shí)測(cè)的FVCmax和FVCmin數(shù)據(jù)相對(duì)應(yīng)[6]。

與回歸模型法一樣，采用75%的樣本點(diǎn)進(jìn)行計(jì)算，25%的樣本點(diǎn)進(jìn)行精度驗(yàn)證，擬合結(jié)果如圖5所示。結(jié)合樣地所在像元擬合得到的植被覆蓋度與樣地實(shí)測(cè)植被覆蓋度進(jìn)行精度驗(yàn)證，結(jié)果如圖6所示。

圖5 像元二分模型法各植被指數(shù)反演結(jié)果

圖6 像元二分模型法精度驗(yàn)證

2.3 逐年反演植被覆蓋度

兩種方法精度驗(yàn)證的相關(guān)系數(shù)如表1所示，可以看出，對(duì)于4種植被指數(shù)來(lái)說(shuō)，回歸模型法的模擬精度均優(yōu)于像元二分模型法；對(duì)于回歸模型法來(lái)說(shuō)，4種植被指數(shù)中MSAVI的反演精度最高，R2達(dá)到了0.815 6。 因此，本文采用回歸模型結(jié)合MSAVI建立的關(guān)系：y=146.4x+0.853 3，作為反演2001－2018年黃河源區(qū)植被覆蓋度的依據(jù)。

表1 兩種方法精度驗(yàn)證的相關(guān)系數(shù)R2

3 研究結(jié)果分析

3.1 2001－2018年植被覆蓋度變化分析

本文采用最小二乘線性回歸模型計(jì)算黃河源區(qū)2001－2018年植被覆蓋度的變化情況。該方法常用于分析長(zhǎng)時(shí)間序列的植被變化趨勢(shì)，將時(shí)間變量t作為獨(dú)立變量，與各年的NDVI時(shí)序數(shù)據(jù)進(jìn)行最小二乘回歸分析，從而得到一個(gè)線性方程；再利用該方程，根據(jù)其斜率數(shù)值，獲得源區(qū)2001－2018年植被覆蓋度的變化趨勢(shì)。若斜率小于零，則植被生長(zhǎng)呈下降趨勢(shì)；反之，則呈增長(zhǎng)趨勢(shì)。時(shí)間變化公式為：

式中，y為各年的平均植被覆蓋度；a為擬合參數(shù)；b為線性回歸方程的斜率[6]。

2001－2018年黃河源區(qū)植被覆蓋度變化趨勢(shì)如圖7所示，可以看出，2001－2018年源區(qū)的植被覆蓋度整體呈增長(zhǎng)趨勢(shì)，但增長(zhǎng)不是很明顯，且自2001年 以來(lái)，黃河源區(qū)植被覆蓋度在50%～70%之間。

圖7 2001－2018年黃河源區(qū)植被覆蓋度變化趨勢(shì)分析

3.2 黃河源區(qū)多年植被覆蓋度均值統(tǒng)計(jì)分析

本文根據(jù)植被覆蓋度劃分標(biāo)準(zhǔn)（表2）對(duì)黃河源區(qū)植被覆蓋度多年均值進(jìn)行分級(jí)[13-14]，結(jié)果如圖8所示，并統(tǒng)計(jì)了不同覆蓋度等級(jí)所占的比例（表3）。

表2 植被覆蓋度劃分標(biāo)準(zhǔn)

結(jié)合圖8和表3可知，整個(gè)黃河源區(qū)內(nèi)植被覆蓋度＜10%的裸地，占比為1.41%；植被覆蓋度為10%～30%的低覆蓋區(qū)，占比為11.69%；植被覆蓋度為30%～45%的中低覆蓋區(qū)，占比為17.66%；植被覆蓋度為45%～60%的中覆蓋區(qū)，占比為19.59%；植被覆蓋度＞60%的高覆蓋區(qū)，占比為49.65%；因此黃河源區(qū)約50%的地區(qū)其植被覆蓋度在60%以上。

圖8 黃河源區(qū)植被覆蓋度分級(jí)結(jié)果

表3 黃河源區(qū)植被覆蓋度不同等級(jí)的面積統(tǒng)計(jì)

4 結(jié) 語(yǔ)

1）2001－2018年黃河源區(qū)的植被覆蓋度整體呈增長(zhǎng)趨勢(shì)，但該趨勢(shì)不是很明顯，且自2001年以來(lái)，黃河源區(qū)植被覆蓋度在50%～70%之間。

2）整個(gè)黃河源區(qū)內(nèi)植被覆蓋度＜10%的裸地，占比為1.41%；植被覆蓋度為10%～30%的低覆蓋區(qū)，占比為11.69%；植被覆蓋度為30%～45%的中低覆蓋區(qū)，占比為17.66%；植被覆蓋度為45%～60%的中覆蓋區(qū)，占比為19.59%；植被覆蓋度＞60%的高覆蓋區(qū)，占比為49.65%；因此，黃河源區(qū)約50%的地區(qū)其植被覆蓋度在60%以上。

本文只選用了MODIS數(shù)據(jù)進(jìn)行研究，沒(méi)有結(jié)合其他較長(zhǎng)時(shí)間序列數(shù)據(jù)（如GIMMS）進(jìn)行長(zhǎng)時(shí)間序列的黃河源區(qū)植被覆蓋度的變化趨勢(shì)分析；計(jì)算得到了黃河源區(qū)2001－2018年的植被覆蓋度變化趨勢(shì)，但未進(jìn)行影響因素分析，這也是今后需繼續(xù)研究的內(nèi)容。