基于Sentinel-2影像的廈門市茶園遙感提取

李 艷，張 帆

(河口生態(tài)安全與環(huán)境健康福建省高校重點(diǎn)實(shí)驗(yàn)室 / 廈門大學(xué)嘉庚學(xué)院環(huán)境科學(xué)與工程學(xué)院，福建 漳州 363105)

茶Camelliasinensis是常綠多年生木本植物，也是一種重要的經(jīng)濟(jì)作物，主要分布在我國(guó)熱帶和亞熱帶山區(qū)。福建省是我國(guó)茶葉主要產(chǎn)區(qū)之一，近年來(lái)，隨著茶葉經(jīng)濟(jì)價(jià)值的不斷提升和低碳農(nóng)業(yè)發(fā)展理念的不斷深入，茶園種植面積發(fā)生了明顯的變化，茶園碳匯也成為農(nóng)林業(yè)生態(tài)環(huán)境的重要組成部分[1]。因此，對(duì)茶園范圍的監(jiān)測(cè)具有重要的生態(tài)環(huán)境保護(hù)意義。

由于茶園分布廣泛，且多分布于山區(qū)和丘陵地帶，通過(guò)傳統(tǒng)的實(shí)地調(diào)查獲取茶樹空間分布需花費(fèi)大量時(shí)間和人力，具有較大的局限性。遙感技術(shù)具有宏觀性、動(dòng)態(tài)性、速度快、成本低等特點(diǎn)，是實(shí)現(xiàn)大面積、大尺度茶園種植區(qū)提取的重要手段[2—3]。光學(xué)遙感影像具有較高的空間分辨率和時(shí)間分辨率以及合適的光譜分辨率等優(yōu)勢(shì)，在植被遙感分類中被廣泛使用[4—5]。馬超等[6]基于Landsat8 和時(shí)序的MODIS 影像的光譜特征、植被指數(shù)和地形特征，采用決策樹分類模型對(duì)研究區(qū)的茶園進(jìn)行提取，總體分類精度達(dá)85.71%。趙曉晴等[7]利用時(shí)序Sentinel-2A 影像分析研究區(qū)典型地物的時(shí)序光譜特征與歸一化植被指數(shù)(Normalized Difference Vegetation Index,NDVI)的時(shí)序特征，獲得可用于茶園提取的特征波段，并通過(guò)這些特征波段結(jié)合NDVI 作為特征參數(shù)構(gòu)建決策樹實(shí)現(xiàn)研究區(qū)茶園的提取，但對(duì)于種植分散的茶園，該方法的普適性有待考證。李龍偉等[8]根據(jù)多時(shí)相Sentinel-2 數(shù)據(jù)的紅邊與短波紅外波段構(gòu)建了歸一化茶園指數(shù)(Normalized Difference Tea Plantations Index,NDTI)，并基于NDTI 通過(guò)決策樹模型提取了浙江西北部的茶園，提取總精度達(dá)93.83%。

由于茶樹的光譜反射率值往往與其他一些植物如灌木、果園等重疊，引起提取誤差，而使用基于空間信息的紋理特征可有效解決光譜混合問題[9]。楊艷魁等[10]利用高分二號(hào)(GF-2)影像數(shù)據(jù)構(gòu)建差異歸一化差分植被指數(shù)(Modified Normalized Difference Vegetation Index,DNDVI)和灰度共生紋理(Gray-level Co-occurrence Matrix,GLCT)，然后結(jié)合光譜和紋理特征運(yùn)用支持向量機(jī)(Support Vector Machine,SVM)的方法對(duì)福建省安溪縣茶園進(jìn)行分類。結(jié)果表明，結(jié)合光譜特征和紋理特征提取的茶園精度(89.8%)高于僅用光譜特征提取的茶園精度(79.6%)。柏佳等[11]利用高分1 號(hào)(GF-1)的光譜、紋理特征和 Sentinel-2 時(shí)序特征，采用隨機(jī)森林(Random Forests,RF)的算法提取了茶園，分類結(jié)果的準(zhǔn)確率達(dá)到96.91%。

基于前人的研究，本研究以廈門市同安區(qū)高海拔地區(qū)的蓮花鎮(zhèn)為研究區(qū)，采用多時(shí)相的Sentinel-2遙感影像數(shù)據(jù)，通過(guò)光譜特征選擇植被指數(shù)，并結(jié)合紋理特征，研究茶園提取方法，以期為茶園的管理及政策制定提供參考。

1 研究區(qū)及數(shù)據(jù)

1.1 研究區(qū)概況

廈門市(24°23'～24°54'N、117°53'～118°26'E)位于福建省東北部，總面積1700.61 km2，地形以濱海平原、臺(tái)地和丘陵為主，西北部多中低山，由西北向東南傾斜。研究區(qū)位于同安區(qū)蓮花鎮(zhèn)，地處廈門西北部，最高海拔近千米，屬于山區(qū)，全域植被覆蓋率較高，茶樹為當(dāng)?shù)刂饕?jīng)濟(jì)作物。通過(guò)實(shí)地野外調(diào)查并結(jié)合同時(shí)期Google Earth 影像判讀，研究區(qū)主要地物類型大致可分為茶園、林地、農(nóng)業(yè)用地、建筑用地、裸地以及水體等，地物類型豐富。

1.2 遙感數(shù)據(jù)及預(yù)處理

數(shù)據(jù)源來(lái)自歐洲航天局研發(fā)的Sentinel 系列衛(wèi)星之一的Sentinel-2。Sentinel-2 采用天體平臺(tái)-L(AstroBus-L)，由Sentinel-2A 和Sentinel-2B 兩顆衛(wèi)星聯(lián)合工作，兩顆衛(wèi)星運(yùn)行于同一條軌道上，相位相差180°，雙星座重訪周期可達(dá)5 d。Sentinel-2 攜帶高分辨率多光譜成像儀(MSI)，采用堆掃式成像模式，擁有從可見光到短波紅外的13 個(gè)光譜波段?？臻g分辨率高達(dá)10 m，并在紅邊范圍具有3 個(gè)波段(表1)，為區(qū)域植被生態(tài)環(huán)境特征信息的提取分析提供了全新的解決方案[12]。

表1 Sentinel-2 多光譜影像主要參數(shù)Table 1 Main parameters of Sentinel-2 multispectral images

Sentinel-2 L2A 影像經(jīng)過(guò)輻射定標(biāo)和大氣校正，因此采用歐洲航天局的SNAP 軟件對(duì)原始影像進(jìn)行重采樣與格式轉(zhuǎn)換。

1.3 實(shí)測(cè)數(shù)據(jù)

于2022年5月、10月對(duì)研究區(qū)進(jìn)行野外考察。利用GPS 手持機(jī)和谷歌地球影像，進(jìn)行野外地面樣點(diǎn)布設(shè)及環(huán)境要素的考察與識(shí)別，標(biāo)定764 個(gè)不規(guī)則多邊形地塊，樣本分布均勻，其中80%作為訓(xùn)練樣本，剩余作為驗(yàn)證樣本(圖1)。

圖1 研究區(qū)位置及樣本分布Fig.1 Study area location and the distribution of samples

2 研究方法

基于Sentinel-2 多光譜影像的光譜數(shù)據(jù)和紋理特征，采用隨機(jī)森林(RF)算法進(jìn)行分類，進(jìn)而對(duì)茶園種植區(qū)進(jìn)行提取。首先分析多時(shí)相的Sentinel-2光譜在不同地物的光譜特征差異，通過(guò)敏感光譜建立光譜指數(shù)并分析，然后利用灰度共生矩陣(GLCM)構(gòu)建影像相關(guān)紋理特征，最后應(yīng)用RF 分類方法對(duì)圖像進(jìn)行分類，并進(jìn)行精度檢驗(yàn)和分析。

2.1 分類特征提取

2.1.1 光譜特征

不同地物的光譜曲線具有獨(dú)特的光譜特征，廣泛應(yīng)用于遙感影像分析[13]。通過(guò)Sentinel-2 影像對(duì)研究區(qū)各地物樣本統(tǒng)計(jì)分析得到光譜曲線特征，為茶園提取提供依據(jù)。

綜合考慮研究區(qū)生態(tài)茶園的物候特征以及可獲取的高質(zhì)量遙感數(shù)據(jù)，選擇茶樹春季采摘修剪期(2022年4月8日)、生長(zhǎng)旺盛期(2022年7月22日)、秋季采摘修剪期(2022年10月20日)和生長(zhǎng)減緩期(2022年12月9日)4 個(gè)時(shí)期的 Sentinel-2 數(shù)據(jù)中地表反射率圖像，提取研究區(qū)內(nèi)6 種典型地物樣本(茶園、林地、農(nóng)業(yè)用地、建筑用地、裸地、水體)在不同光譜波段上的光譜反射率，繪制光譜特征曲線圖(圖2)，并進(jìn)行光譜特征分析。

圖2 研究區(qū)典型地物不同時(shí)相光譜曲線Fig.2 The spectral curves of features in several times

2.1.2 光譜指數(shù)

利用光學(xué)遙感的兩個(gè)波段或多個(gè)波段組合建立的植被指數(shù)不僅形式簡(jiǎn)單，而且能夠有效地描述植被狀態(tài)，被廣泛應(yīng)用于植被統(tǒng)計(jì)和提取、植被時(shí)空變化分析等[14]。

根據(jù)光譜在波段上的反映，選擇三種光譜指數(shù)，即歸一化植被指數(shù)(NDVI)、歸一化差值水體指數(shù)(Normalized Difference Water Index,NDWI)和倒紅邊葉綠素指數(shù)(Inverted Red-Edge Chlorophyll Index,IRECI)。

NDVI 是近紅外波段和紅光波段光譜反射率之差與之和的比率，是反映植被生長(zhǎng)狀態(tài)的重要指示因子[15]，表達(dá)式為：

式中，ρNIR和ρred分別為近紅外波段和紅光波段的光譜反射率。

NDWI 是短波紅外波段和近紅外波段的歸一化比值指數(shù)[16]，表達(dá)式為：

式中，ρSWIR為短波紅外波段光譜反射率。根據(jù)典型地物在短波紅外波段的表現(xiàn)，本研究采用的是Sentinel-2 的短波紅外1 波段，即SWIR1。

Sentinel-2 在紅光與近紅外波段的紅邊區(qū)域(670～760 nm)內(nèi)具有3 個(gè)波段，IRECI 包含了Sentinel-2 的3 個(gè)紅邊波段[17]，表達(dá)式為：

式中，ρ783、ρ665、ρ705和ρ740分別對(duì)應(yīng)Sentinel-2 的Band7、Band4、Band5 和Band6 的光譜反射率。

2.1.3 紋理特征

Haralick 等[18]基于統(tǒng)計(jì)的方法提出灰度共生矩陣(GLCM)，可以反映圖像灰度關(guān)于方向、相鄰間隔、變化幅度以及快慢的綜合信息。由于Sentinel-2 包含波段較多，若對(duì)每個(gè)波段都提取紋理特征圖像，容易產(chǎn)生數(shù)據(jù)信息的冗余。因此，本研究對(duì)反射率數(shù)據(jù)采用主成分分析(Principal Components Analysis,PCA)，選取前3 個(gè)主成分。

紋理提取時(shí)窗口尺寸大小會(huì)影響GLCM性能和紋理特征完整性，若窗口較大，會(huì)降低區(qū)分不同地物的靈敏度，而窗口較小會(huì)受像元亮度不均勻性的影響。7×7 像素和9×9 像素的窗口尺寸對(duì)于茶園提取是較合適的[19]。通過(guò)實(shí)驗(yàn)采用7×7 像元窗口，計(jì)算紋理特征值中常用的6 個(gè)紋理：方差(Variance)、同質(zhì)性(Homogeneity)、對(duì)比度(Contrast)、相關(guān)性(Correlation)、角二階矩(Angular Second Moment,ASM)和熵(Entropy)。

由于GLCM 還具有方向性，本研究獲取的紋理值是0°、45°、90°和135°四個(gè)方向的平均值。

綜上所述，研究中將每個(gè)時(shí)期的3 種光譜指數(shù)和主成分變換第一、第二、第三分量的6 種紋理特征作為分類提取特征。

2.2 茶園提取方法

結(jié)合實(shí)地調(diào)研及更高分辨率的 Google Earth圖像，選取茶園(183 個(gè)地塊，共3796 個(gè)像元)與非茶園(423 個(gè)地塊，共28838 個(gè)像元)訓(xùn)練樣本采用隨機(jī)森林(RF)算法進(jìn)行分類，分類工作基于SNAP 實(shí)現(xiàn)。

3 結(jié)果與分析

3.1 典型地物光譜特征

將研究區(qū)6 種主要典型地物樣本分為植被樣本(茶園、林地、農(nóng)業(yè)用地)和非植被樣本(建筑用地、裸地、水體)。從圖2 可以看出，植被樣本和非植被樣本的光譜曲線在不同時(shí)相下，同一地物的變化趨勢(shì)相似，這與研究區(qū)氣候特征有關(guān)，研究區(qū)植被多是常青植被。但在同一時(shí)相下，植被樣本和非植被樣本的光譜曲線有較大差異。

非植被樣本在同一時(shí)相下，不同樣本的光譜曲線趨勢(shì)較為一致，且不同樣本在相同波段的光譜值相差較大，如圖2(b)在近紅外(NIR)波段，裸地的光譜值達(dá)到5300，建筑用地的光譜值接近4000，而水體的光譜值在1100 左右。

植被樣本在同一時(shí)相下，不同樣本的光譜曲線趨勢(shì)較為一致，呈現(xiàn)先上升后下降的趨勢(shì)，在藍(lán)(B)、綠(G)、紅(R)波段的變化比較平緩，且光譜值比較接近；從植被紅邊1(RE1)波段開始，植被樣本的光譜值呈明顯上升趨勢(shì)，茶園在近紅外(NIR)波段處達(dá)最大，隨后光譜值下降。林地和農(nóng)業(yè)用地在不同時(shí)相下，最大光譜值的波段并不一致。在春季采摘修剪期(4月)，茶園與林地在多個(gè)波段的光譜值比較接近；在生長(zhǎng)旺盛期(7月)和秋季采摘修剪期(10月)，除短波紅外2(SWIR2)波段，茶園在其他波段的光譜值均高于其他植被，且有一定的差值；生長(zhǎng)減緩期(12月)，茶園從植被紅邊3(RE3)到水汽(WV)波段的光譜值高于林地和農(nóng)業(yè)用地。總體來(lái)看，茶園與其他植被具有一定光譜可分性的波段主要有 RE2、RE3、NIR 和SWIR1。

3.2 典型地物光譜指數(shù)

研究區(qū)典型地物在多時(shí)相NDVI、NDWI 和IRECI 的特征圖如圖3～圖5所示。

圖3 研究區(qū)典型地物多時(shí)相的NDVI 特征Fig.3 Characteristics of typical land cover in multi-temporal NDVI

在NDVI 特征圖中，非植被類地物(水體、建筑用地)在各個(gè)時(shí)期均呈現(xiàn)暗黑色，植被類地物(茶園、林地、農(nóng)業(yè)用地)各個(gè)時(shí)期均呈現(xiàn)亮白色或灰色(圖3)。建筑用地附近的農(nóng)業(yè)用地在7月和12月特征值較低，呈灰色。大部分林地由于較密集而呈高亮色，部分林地在坡度較大的區(qū)域也呈暗灰色。而茶園在各個(gè)時(shí)期的特征值變化不大，多呈現(xiàn)高亮色。與林地易混淆，無(wú)法實(shí)現(xiàn)很好的區(qū)分。

在NDWI 特征圖中，非植被類與NDVI 特征圖類似，呈暗黑色，但植被類也有部分區(qū)域呈現(xiàn)暗灰色(圖4)。這是由于植被在短波紅外波段的反射率與其葉片含水量相關(guān)，而研究區(qū)水資源較為匱乏，林地、茶園和部分農(nóng)田無(wú)灌溉條件。茶樹與高大林木相比較，根系較淺，在雨水較為匱乏的4月和12月，茶樹冠層含水量低于高大林木，茶園所在區(qū)域灰度值較低。

圖4 研究區(qū)典型地物多時(shí)相的NDWI 特征Fig.4 Characteristics of typical land cover in multi-temporal NDWI

IRECI 對(duì)葉綠素含量變化比較敏感，植被區(qū)域也呈現(xiàn)亮白色和灰色(圖5)。在10月和12月，隨著植被光合作用降低，IRECI 指數(shù)特征值下降，與NDVI 和NDWI 特征圖相似，較難實(shí)現(xiàn)茶園和林木的區(qū)分。

圖5 研究區(qū)典型地物多時(shí)相的IRECI 特征Fig.5 Characteristics of typical land cover in multi-temporal IRECI

3.3 茶園分布

圖6 顯示多時(shí)相下茶園的分類?？傮w而言，4個(gè)時(shí)期提取的茶園空間分布比較一致。研究區(qū)茶園主要分布在北部，中部地區(qū)有零散的茶園。經(jīng)實(shí)地考察，研究區(qū)北部與福建省泉州市的安溪縣交界，茶園種植范圍廣泛。受影像分辨率以及植被光譜特征相似性的限制，在研究區(qū)東南部有部分耕地被錯(cuò)分為茶園。

圖6 研究區(qū)多時(shí)相茶園分布圖Fig.6 Extraction result of tea plantations distribution

3.4 精度評(píng)價(jià)

為了驗(yàn)證分類精度，利用目視判讀結(jié)果隨機(jī)生成921 個(gè)茶園和 8321 個(gè)非茶園驗(yàn)證像元，通過(guò)建立的混淆矩陣，計(jì)算總體精度和總體Kappa 系數(shù)(表2)。研究區(qū)內(nèi)6 種地物在4～12月4 個(gè)時(shí)期的分類總體精度依次為 94.28%、91.05%、95.53%和94.41%；Kappa 系數(shù)依次為0.8889、0.8254、0.9088和0.8881。茶園在4月和10月的用戶精度較好，分別為82.92%和82.25%，這兩個(gè)時(shí)期也是研究區(qū)采茶期；茶園在7月的用戶精度較低，僅51.57%。總體來(lái)說(shuō)，茶園提取的準(zhǔn)確率較高。

表2 茶園提取精度Table 2 Accuracy assessment for extraction of tea plantations

4 結(jié)論與討論

4.1 討論

研究區(qū)植被類型具有相似的光譜特征，因此本研究基于Sentinel-2 豐富的光譜信息，結(jié)合茶園紋理特征，進(jìn)行多時(shí)相的地物分類和茶園提取分析。結(jié)果表明，將光譜特征、植被指數(shù)和紋理特征結(jié)合可有效提取茶園，同時(shí)提取精度受季節(jié)因素影響。趙曉晴等[7]利用Sentinel-2 多時(shí)相影像提取茶園時(shí)，發(fā)現(xiàn)5月和12月是杭州市西湖區(qū)茶園提取特征的最佳組合。Sentinel-2 具有多個(gè)適用于植被監(jiān)測(cè)的紅邊波段，李龍偉等[8]采用Sentinel-2 的紅邊波段構(gòu)建植被指數(shù)(NDTI)作為茶園提取的分類特征。本研究也利用紅邊波段構(gòu)造的植被指數(shù)(IRECI)作為茶園的分類特征之一，由于研究區(qū)多云霧和陰雨時(shí)期較多，本研究未使用Sentinel-2 時(shí)間序列信息進(jìn)行茶園提取，但柏佳等[11]基于 Sentinel-2 與MODIS 數(shù)據(jù)融合獲得完整的 NDVI 時(shí)序數(shù)據(jù)，采用隨機(jī)森林算法提取茶園準(zhǔn)確率高達(dá)96.91%。

盡管本研究結(jié)果可以有效提取茶園范圍，但存在一些不確定性。如林地和農(nóng)業(yè)用地與茶園具有相似的光譜特征，“異物同譜”現(xiàn)象增加了高精度提取茶園的難度；Sentinel-2 的空間分辨率為10 m，對(duì)于小范圍的茶園，無(wú)法準(zhǔn)確辨識(shí)。

4.2 結(jié)論

本文對(duì)廈門市同安區(qū)蓮花鎮(zhèn)茶園進(jìn)行提取，基于影像的光譜特征、植被指數(shù)和紋理特征，構(gòu)造隨機(jī)森林分類模型，研究發(fā)現(xiàn)：

(1)研究區(qū)植被的光譜曲線呈現(xiàn)先上升后下降的趨勢(shì)，而非植被的光譜曲線趨勢(shì)較平穩(wěn)。根據(jù)光譜特征可知，茶園與其他植被具有一定光譜可分性的波段主要有RE2、RE3、NIR 和SWIR1。

(2)在植被類型容易混淆情況下，基于隨機(jī)森林分類方法提取茶園，能夠較好地滿足分類要求。

受多種因素影響，“異物同譜”現(xiàn)象仍然存在，茶園與研究區(qū)其他地物有混淆的情況，如少量林地和農(nóng)地。下一步將結(jié)合高光譜數(shù)據(jù)對(duì)研究區(qū)分類精度開展更加深入的分析和探討。