蔣永泉,史正濤,高書鵬
(云南師范大學 地理學部,昆明 650500)
橡膠因其天然獨特的物化特性而成為我國重要的戰(zhàn)略基礎(chǔ)物資之一[1]。橡膠樹于20世紀初引入云南省種植[2],目前,已成為我國最大的天然橡膠生產(chǎn)基地。20世紀60年代云南省瑞麗市開始大面積種植橡膠林,截至2017年,全市共種植橡膠7 053.55hm2,總產(chǎn)值高達6 000萬元余[3]。長久以來,橡膠林種植多以犧牲天然林為主,導致了當?shù)匚锓N多樣性減少、景觀連通性降低、水土流失等[4]諸多問題。因此,準確監(jiān)測橡膠林的種植狀況對于平衡經(jīng)濟發(fā)展和生態(tài)環(huán)境保護具有重要意義。目前,遙感技術(shù)在橡膠林識別中取得了長足的進展,MODIS數(shù)據(jù)是識別橡膠林使用最廣的遙感數(shù)據(jù)源之一[5-6],但由于MODIS空間分辨率低,識別結(jié)果較差,很難將斑塊較小的橡膠林識別出來。隨著遙感技術(shù)的發(fā)展,多源遙感數(shù)據(jù)識別地物已成為遙感應用的趨勢和熱點。與單一數(shù)據(jù)源相比,多源數(shù)據(jù)可充分利用各種數(shù)據(jù)的高光譜、高時空分辨率構(gòu)建各種特征,提高識別精度以滿足應用需求。高書鵬等[7]基于時空數(shù)據(jù)融合算法獲取的ETM+,OLI,Sentinel-2A時序數(shù)據(jù)提取西雙版納的橡膠林,Chen等[8]和Kou等[9]都基于雷達數(shù)據(jù)和Landsat數(shù)據(jù)分別提取海南省和西雙版納的橡膠林,結(jié)果顯示,同時運用2種數(shù)據(jù)源均比單一數(shù)據(jù)源的提取精度高。因此,多源遙感數(shù)據(jù)可以彌補單一數(shù)據(jù)存在的缺陷,實現(xiàn)不同數(shù)據(jù)源之間的優(yōu)勢互補,將多源遙感數(shù)據(jù)進行融合,在一定程度上優(yōu)化了分類結(jié)果,提高了識別精度[10],在橡膠林監(jiān)測管理中具有較好的應用前景。
原種植于亞馬遜河流域的橡膠樹,要求高溫、多雨、靜風、光照充足的氣候條件,而我國屬于非傳統(tǒng)植膠區(qū),適宜種植區(qū)非常有限[11]。瑞麗市作為我國北緣橡膠分布區(qū),最冷月平均氣溫僅有12.5℃,是世界上種植橡膠氣溫最低的地區(qū),經(jīng)常受到寒潮和寒害的影響[12]。1970—1971年、1975—1976年以及1982—1983年瑞麗市橡膠林都分別遭受了不同程度的寒害,使當?shù)亟?jīng)濟受到重創(chuàng)[13-14]。隨著市場需求量擴大和膠價上漲,2004年后,膠農(nóng)開始大面積砍伐森林種植橡膠林并向不適宜種植區(qū)擴張,導致了橡膠生長速度慢、產(chǎn)膠量低[15-16]。因此,掌握橡膠林的空間分布特征及其種植適宜區(qū)對規(guī)劃管理具有重要的意義。本文綜合利用Sentinel-2和MODIS-NDVI多源遙感數(shù)據(jù)識別橡膠林,該方法數(shù)據(jù)易于獲取,可實現(xiàn)大范圍的監(jiān)測,基于識別結(jié)果,對瑞麗市橡膠林空間分布特征及其種植適宜區(qū)分析,可為相關(guān)部門決策和生態(tài)恢復提供參考。
研究區(qū)位于云南省西部瑞麗市,地處橫斷山脈高黎貢山余脈的向南延伸部分,地勢西北高東南低,位于北緯23.38′~24.14′,東經(jīng)97.31′~98.02′之間,其西北、西南、東南三面與緬甸接壤[17]。該區(qū)屬于南亞熱帶季風性氣候區(qū),冬無嚴寒,夏無酷暑,全年分旱雨兩季,年平均氣溫21℃,年降水量1 394.8mm,年平均日照2 330h[18]。瑞麗市總面積為1 020km2,其中森林面積占總面積的46.67%,主要森林類型有熱帶山地雨林,亞熱帶季風常綠闊葉林,亞熱帶山地落葉闊葉林,針葉林和竹林。高等植物類型豐富,共有高等植物284科,943屬,3 733種,其中有國家級保護植物45種,省級保護植物24種[19]。全市植被覆蓋度較高,其中橡膠林是當?shù)胤N植面積最大的人工林[20]。
1) 選用的Sentinel-2遙感影像數(shù)據(jù)可通過歐空局(European Space Agency,ESA)網(wǎng)站(https://scihub.copernicus.eu/dhus/#/home)免費下載獲取。大氣校正通過Sen2Cor[21]小插件完成,其核心算法是大氣輻射傳輸模型libRadtran[22],目的是將大氣上層表觀反射率轉(zhuǎn)化為大氣下層地表反射率。根據(jù)研究需要提取出10m的紅光波段、綠光波段、藍光波段和近紅外波段,以及4個20m的植被紅邊波段并將其重采樣為10m,然后裁剪出研究區(qū)大小的區(qū)域以備后續(xù)使用。
2) MOD13Q1是由美國國家航空航天局(NASA)提供的16d合成的250m空間分辨率NDVI產(chǎn)品,可在NASA網(wǎng)站(https://search.earthdata.nasa.gov/search)免費獲取。獲取的原始MODIS數(shù)據(jù)格式為HDF,投影類型與Sentinel-2數(shù)據(jù)不一致,因此本文采用MRT工具做投影轉(zhuǎn)換,將其地理坐標和投影類型分別轉(zhuǎn)為WGS84和UTM,然后裁剪出研究區(qū)大小的區(qū)域。為了得到10m×10m像元大小的NDVI數(shù)據(jù),將23個波段的NDVI數(shù)據(jù)用立方卷積重采樣處理,使其與Sentinel-2多光譜數(shù)據(jù)具有相同的像元大小。
3) 樣本數(shù)據(jù)包括訓練樣本和驗證樣本,本文基于Google Earth高清影像,對研究區(qū)進行目視解譯,將土地利用類型分為5類,分別為橡膠林、自然林、耕地、建設用地和水體。根據(jù)各地類的空間分布情況、復雜程度以及比例,一共隨機選取了2 563個點,其中橡膠林527個、自然林764個、耕地681個、建設用地337、水體254個。并將各地類按1∶1的比例分開分別作為訓練樣本和驗證樣本。
為分析瑞麗市橡膠林與其他地類在不同時段的生長變化差異情況,選用2019年的MOD13Q1數(shù)據(jù),分析5個地類在23個時相的NDVI響應關(guān)系(圖1)。從總體上看,不同地類對NDVI的響應不同,橡膠林在整個時序變化過程中與其他地類存在明顯的差異。第49d(2月)后,橡膠林開始萌芽進入返青期,NDVI值上升;第129 d(5月)—177d(6月)橡膠林的NDVI值急劇降低,這是由于夏季研究區(qū)內(nèi)多云雨天氣,云量較大NDVI值偏低引起的;在第241—289 d(9—10月)NDVI達到最大值;第305 d(11月)后,橡膠林的NDVI值開始緩慢下降,此時橡膠林開始落葉。
不同地類對不同波段的吸收和反射光譜存在一定的差異,圖2為各地類光譜反射率均值曲線圖。在藍光波段、綠光波段和紅光波段橡膠林易與水體混淆;在植被紅光波段1,橡膠林的光譜值快速升高,且可與其他地類較好的分離;在植被紅邊2和3波段,橡膠林易與建設用地混淆;在植被紅邊波段4,橡膠林的光譜值達到最高,與其他有較好的可分性。因此,本文選用植被紅邊波段1和4用于區(qū)分橡膠林和其他地類。
圖1 地類NDVI時間序列曲線圖
注:“B2”藍光波段;“B3”綠光波段;“B4”紅光波段;“B5”植被紅邊波段1;“B6”植被紅邊波段2;“B7”植被紅邊波段3;“B8A”植被紅邊波段4
基于Sentinel-2多光譜遙感影像提取分類特征,包括植被指數(shù)、紋理特征和光譜特征共23個特征,以及MOD13Q1時間序列數(shù)據(jù)的23個NDVI值(每16d一個),如表1所示。 歸一化植被指數(shù)(NDVI)是反映植被生長狀態(tài)的重要因子[23],與植被的覆蓋度有密切的相關(guān)性。NDVI的計算公式如下[24]:
(1)
式中:ρNIR和ρred分別表示近紅外波段和紅波段的表光反射率。
灰度共生矩陣(GLCM)是20世紀70年代初由Haralick等[25]提出的,是通過灰度的空間特性來描述紋理的方法,能用來反映地物類別在空間特征的差異[26]?;叶裙采仃嚬灿?4種特征指標來定量描述紋理,本文選用了最常用的5種,分別為均值(Mean)、同質(zhì)性(Homogeneity)、相關(guān)性(Correlation)、熵(Entropy)和對比度(Contrast)。
(2)
(3)
(4)
(5)
(6)
式中:ux,uy,σx,σy,分別為行和列的均值和標準差;n為共生矩陣的階數(shù);i,j為共生矩陣的坐標;P為(i,j)處的共生矩陣數(shù)值。
表1 地物分類特征表
本文采用支持向量機分類器進行分類,分別對Sentinel-2多光譜遙感影像、MOD13Q1-NDVI時序數(shù)據(jù)以及多源遙感數(shù)據(jù)進行分類,結(jié)果如圖3所示。從總體上看,3種不同數(shù)據(jù)源得到的瑞麗市橡膠林分布在空間上趨于一致。僅用單一的MOD13Q1-NDVI分類(圖3(a)),不能較準確地識別出一些面積較小且分布零散的橡膠林。僅用單一Sentinel-2數(shù)據(jù)分類(圖3(b)),雖然能將面積較小的橡膠林識別出來,但得到的結(jié)果比較破碎。圖3(c)為綜合運用多源數(shù)據(jù)的分類結(jié)果,由于綜合了兩種數(shù)據(jù)源在空間和時間上的優(yōu)勢,總的分類結(jié)果得到了很大的改善,能在準確識別出一些破碎地塊的同時,又改善了Sentinel-2數(shù)據(jù)提取橡膠林結(jié)果的破碎程度,有效地避免了單一數(shù)據(jù)源的分類錯誤。
圖3 橡膠林識別結(jié)果
本文基于混淆矩陣分別計算分類的制圖精度(Producer Accuracy,PA)、用戶精度(User Accuracy,UA)、總體精度(Overall Accuracy,OA)以及Kappa系數(shù)來評價結(jié)果的可行性。從表2可以看出,綜合應用多源遙感數(shù)據(jù)能夠?qū)崿F(xiàn)瑞麗市橡膠林的高精度識別,優(yōu)于僅用單一遙感數(shù)據(jù)源分類方法,聯(lián)合多源遙感數(shù)據(jù)能夠有效提高分類精度。綜合運用多源數(shù)據(jù),橡膠林的制圖精度PA為89.73%、用戶精度UA為73.98%、總體精度OA為91.41%、Kappa系數(shù)為0.76;與僅用MODIS數(shù)據(jù)相比,PA,UA,OA和Kappa系數(shù)分別提高了18.25%,11.52%,6.10%和0.19;與僅用Sentinel-2數(shù)據(jù)相比,PA,UA,OA和Kappa系數(shù)分別提高了16.73%,3.39%,3.21%和0.12。
表2 橡膠林識別精度評價對比表
從分類結(jié)果可以看出,瑞麗市橡膠林主要分布在的西南和東北方向上,對結(jié)果進行聚類分析處理去除破碎的錯分對象,基于30m的DEM數(shù)據(jù)分析橡膠林空間分布特征。瑞麗市DEM數(shù)據(jù)的海拔高度范圍為682~1 997m,將其劃分為6個等級:682~900m,900~1 100m,1 100~1 200m,1 200~1 300m,1 300~1 400m,≥1 400m;根據(jù)我國水利部頒發(fā)《土壤侵蝕分類分級標準》[27]坡度范圍,將瑞麗市坡度劃分為6個等級:0~5°,5~8°,8~15°,15~25°,25~35°,≥35°。得到表3、表4的統(tǒng)計結(jié)果。
表3 瑞麗市橡膠林在海拔高度上的分布特征
表4 瑞麗市橡膠林在坡度上的分布特征
截至2019年底,瑞麗市橡膠總種植面積為7 693.83hm2,占土地總面積的7.54%,其中有6 433.11hm2種植在海拔為682~900m高度范圍內(nèi),其次是種植在900~1 000m范圍內(nèi),分別占橡膠林總種植面積的83.61%和15.26%;主要種植的坡度為5~15°,占橡膠林總種植面積的83.61%。
中華人民共和國農(nóng)業(yè)部發(fā)布的《橡膠樹種植基地建設標準》(NY/T221-2016)[28]根據(jù)橡膠樹對溫度、水分、風速和光照等膠樹生長要求以及產(chǎn)膠潛力,將云南省西部海拔≥900m高度以上的地區(qū)劃分為不宜作為常規(guī)膠園宜林地。依據(jù)此標準的海拔高度變化將瑞麗市劃分為橡膠適宜種植和不適宜種植2個區(qū)域,然后對2019年瑞麗市植橡膠林適宜區(qū)進行統(tǒng)計分析(圖4),瑞麗市有6 433.11hm2橡膠種植在適宜種植區(qū),占橡膠林總種植面積的83.61%,已有1 260.72hm2橡膠林向不適宜種植區(qū)擴張,占橡膠林總種植面積的16.39%。
圖4 瑞麗市橡膠林適宜種植區(qū)
1) 綜合Sentinel-2多光譜影像和MOD13Q1 NDVI時間序列數(shù)據(jù)的分類方法,可充分發(fā)揮2種遙感數(shù)據(jù)源的不同特征優(yōu)勢,與僅用單一數(shù)據(jù)的結(jié)果相比,有效提高了橡膠林的識別精度,實現(xiàn)了瑞麗市橡膠林高精度識別,綜合運用多源遙數(shù)據(jù)橡膠林識別結(jié)果的制圖精度PA為89.73%、用戶精度UA為73.98%、總體精度OA為91.41%,Kappa系數(shù)為0.76。
2) 截至2019年底,瑞麗市橡膠種植面積為7 693.83hm2,有83.61%的橡膠林種植在海拔為684~900m高度范圍內(nèi),占比最大;其次是種植在900~1 100m高度范圍內(nèi),占橡膠林總種植面積大的15.26%。主要種植的坡度為5~15°,占橡膠林總種植面積的83.61%。
3) 根據(jù)農(nóng)業(yè)部發(fā)布的《橡膠樹種植基地建設標準》將瑞麗市的土地按海拔高度劃分為適宜種植區(qū)和不適宜種植區(qū),截至2019年底,瑞麗市有6 433.11hm2橡膠種植在適宜區(qū),占橡膠林總面積的83.61%,已有1 260.72hm2橡膠林向不適宜種植區(qū)擴張,占橡膠林總面積的16.39%。