房家瑋,胡念釗,王懷玉,劉詠梅, 2
(1.西北大學(xué)城市與環(huán)境學(xué)院, 陜西 西安 710127;2.陜西省地表系統(tǒng)與環(huán)境承載力重點實驗室, 陜西 西安 710127)
高寒草地具有氣候調(diào)節(jié)、水源涵養(yǎng)、水土保持、生物多樣性保護(hù)等多種生態(tài)功能,是青藏高原高寒生態(tài)系統(tǒng)價值的重要維系,也是當(dāng)?shù)匦竽翗I(yè)可持續(xù)發(fā)展的重要保障[1]。近年來,高寒草地呈現(xiàn)不同程度的退化趨勢,狼毒(Stellera chamaejasme)等毒害草在退化草地的擴(kuò)散日益加劇,對高寒草地生態(tài)系統(tǒng)的影響日益顯現(xiàn)[2-3]。狼毒是瑞香科狼毒屬多年生草本植物,其為頂生頭狀花序,花白色、淡紅或黃色,花期為6 月下旬至7 月下旬,廣泛分布于海拔2 300~4 200 m 干燥而向陽的高山山坡或河灘臺地,是青藏高原危害最嚴(yán)重的有毒入侵植物之一[4-6]。狼毒根系發(fā)達(dá),對土壤水分的利用效率高且對土壤養(yǎng)分轉(zhuǎn)化速度快,在群落物種競爭中具有明顯優(yōu)勢,已成為中、重度退化草地的優(yōu)勢種或建群種[7-8]。開展青藏高原高寒草地狼毒入侵調(diào)查,監(jiān)測其發(fā)生發(fā)展區(qū)域,對草地資源管理和退化草甸修復(fù)至關(guān)重要。
近年來,隨著影像分辨率的提高以及分類技術(shù)的不斷完善,遙感技術(shù)在草地入侵植物分布識別中的應(yīng)用日漸深入。依據(jù)入侵物種的物候特征并結(jié)合關(guān)鍵生長期的遙感影像,Li 等[9]利用IKONOS 4 m 多光譜影像對狼毒分布區(qū)域進(jìn)行識別;Jana等[10]采用Rapid Eye 5 m 多光譜影像對入侵植物巨型豬草(Heracleum mantegazzianum)的生長進(jìn)行監(jiān)測;楊俊芳等[11]基于GF-1 和GF-2 全色及多光譜波段對黃河三角洲互花米草(Spartina alterniflora)的分布范圍進(jìn)行提取。上述研究使用單時相高分辨率影像取得了較理想的識別效果,但高分辨率影像在大區(qū)域尺度的廣泛應(yīng)用受到一定程度的限制[12]。同時,由于植物群落構(gòu)成復(fù)雜多樣,單時相中低分辨率遙感影像對物種的識別精度相對較低,多時相遙感影像能夠反映植物的物候差異,為其精準(zhǔn)識別提供有效數(shù)據(jù)源。盧元兵等[13]基于Landsat8 OLI的8 期影像構(gòu)建混合三維和二維卷積的神經(jīng)網(wǎng)絡(luò)識別模型,提出了番茄(Solanum lycopersicum)等農(nóng)作物的有效分類方案;史飛飛等[14]使用Landsat 8 OLI和GF-1WFV 4 月-11 月的作物生長期時序數(shù)據(jù)對枸杞(Lycium barbarum)種植區(qū)進(jìn)行分類,取得了較好的識別精度。現(xiàn)有研究側(cè)重于選取關(guān)鍵生長期的中高分辨率遙感影像,以大面積規(guī)律分布的農(nóng)作物作為識別目標(biāo),而將多時相遙感影像與入侵植物獨特的物候信息相結(jié)合,應(yīng)用于毒雜草入侵調(diào)查監(jiān)測,這方面的研究有待于進(jìn)一步開展[15]。
本研究在青海省海北州祁連縣選取狼毒廣泛分布的退化草甸作為研究區(qū),采用狼毒花期前和盛花期的Sentinel-2 多光譜影像,基于Google Earth Engine平臺開展Sentinel-2 影像時間序列去云處理,采用多時相信息篩選狼毒識別敏感指數(shù)并構(gòu)建最優(yōu)特征組合方案,結(jié)合隨機(jī)森林算法提取狼毒分布區(qū)域,評估多時相Sentinel-2 影像在區(qū)域尺度狼毒遙感識別中的適用性,為退化高寒草地毒雜草的綜合防治提供技術(shù)支持。
研究區(qū)位于狼毒入侵廣泛、草地退化嚴(yán)重的青海省海北藏族自治州祁連縣(100°11′15.69″~101°10′34.21″ E,37°45′41.84″~38°16′30.49″ N),面積約4 000 km2,海拔2 596~4 975 m。該區(qū)域為典型的高原大陸性氣候,光照充足、氣溫日較差大、降水量較少、雨熱同期。研究區(qū)草地類型主要為高寒草甸,群落優(yōu)勢種有狼毒、高山嵩草(Kobresia pygmaea)、披針葉黃華(Thermopsis lanceolata)、早熟禾(Poa annua)、異葉米口袋(Gueldenstaedtia diversifolia)、委陵菜(Potentilla chinensis)等。群落平均蓋度約為48.9%,狼毒最大蓋度達(dá)到57.7%,空間上呈現(xiàn)斑塊狀聚集分布的特征。研究區(qū)內(nèi)與狼毒同期的植物群落主要呈現(xiàn)綠色,在盛花期大面積密集分布的狼毒花為粉白色,其他植物群落間或出現(xiàn)零散的黃色或紫色花,對狼毒群落光譜識別影響微弱,狼毒盛花期獨特的物候特征為其遙感識別提供了基礎(chǔ)。
本研究主要采用Sentinel-2A 衛(wèi)星多光譜影像,下載自歐洲航天局(European Space Agency,ESA)數(shù)據(jù)中心(https://scihub.copernicus.eu/dhus),產(chǎn)品級別為L2A 級,選取空間分辨率為10 m 的藍(lán)波段B2(458~523 nm)、綠波段B3 (453~578 nm)、紅波段B4 (650~680 nm)和近紅外波段B8 (785~900 nm)。根據(jù)該區(qū)域狼毒生長的具體物候期[8],選用2019年、2020 年和2021 年狼毒花期前(6 月20 日至30日)與盛花期(7 月5 日至20 日)內(nèi)云量少于20%的影像作為狼毒識別的基礎(chǔ)數(shù)據(jù)。
研究區(qū)在7 月-8 月多云雨天氣存在大量的云霧遮擋區(qū)域,導(dǎo)致影像質(zhì)量下降。本研究利用Google Earth Engine 平臺對所用影像進(jìn)行去云處理,獲取每景影像的可用像元。首先利用Sentinel-2A 質(zhì)量評估(quality assessment,QA) 波段中的云掩膜波段QA60 產(chǎn)品對影像進(jìn)行云檢測并掩膜去云;其次對去云后產(chǎn)生的空值區(qū)采用無云影像取中值進(jìn)行補充,最終得到花期前和盛花期兩個時相的影像。
本研究于2020 年7 月上中旬-2021 年7 月上中旬在研究區(qū)開展了植物群落調(diào)查。根據(jù)區(qū)域內(nèi)植物群落的典型性,沿青羊溝和八寶河谷地布設(shè)兩條樣帶,共布設(shè)37 個樣地。為了和Sentinel-2A 多光譜影像的空間分辨率相匹配,每個樣地選在群落分布均勻且面積大于10 m × 10 m 的區(qū)域,并于區(qū)域中心布設(shè)1 m × 1 m 樣方 2~3 個,共布設(shè)93 個樣方(圖1)。其中,草地樣方(無狼毒生長樣方) 39 個,狼毒生長樣方54 個。樣方調(diào)查內(nèi)容主要包括植物種類、每個物種的數(shù)量和高度,利用數(shù)碼相機(jī)獲取樣方蓋度照片,并采用GNSS RTK 測定樣方坐標(biāo)、坡度和高程。對樣方照片進(jìn)行解譯后得到樣方群落蓋度和物種分蓋度,樣方群落蓋度為29.80%~70.00%。研究中依據(jù)狼毒群落蓋度變化梯度設(shè)置狼毒樣方,狼毒群落蓋度為1.88%~57.68%,不同覆蓋度狼毒群落樣方如圖2 所示。因所選Sentinel-2A 影像波段的分辨率為10 m,當(dāng)樣方的狼毒蓋度小于10%時識別效果不佳,故進(jìn)行剔除,最終研究中采用46 個狼毒樣方和39 個草地樣方作為分類樣本數(shù)據(jù),同期使用GNSS RTK 在研究區(qū)內(nèi)采集139 個狼毒驗證點,147個非狼毒驗證點(包含其他草地群落、道路、河流、居民點及農(nóng)田)作為分類結(jié)果驗證數(shù)據(jù)。地面驗證點主要記錄經(jīng)緯度坐標(biāo)及類別信息。
圖1 研究區(qū)植物樣方分布Figure 1 Distribution of the sampling plots in the study area
本研究同時使用了地形數(shù)據(jù)和植被數(shù)據(jù)。地形數(shù)據(jù)采用了美國國家航空航天局(National Aeronautics and Space Administration, NASA)和日本經(jīng)濟(jì)產(chǎn)業(yè)省(Ministry of Economy, Trade andIndustry, METI)聯(lián)合發(fā)布的V3 版本30 m 分辨率ASTER GDEM 數(shù)據(jù)(https://earthdata.nasa.gov/)。研究區(qū)植被數(shù)據(jù)為《中華人民共和國植被圖(1 ? 1 000 000)》[16],下載自中國科學(xué)院植物科學(xué)數(shù)據(jù)中心(https://www.plantplus.cn/cn)。將植被類型數(shù)據(jù)轉(zhuǎn)為柵格數(shù)據(jù),并與DEM 一起重采樣至10 m,使之與遙感影像數(shù)據(jù)分辨率相統(tǒng)一。
研究區(qū)內(nèi)狼毒廣泛生長于海拔2 600~4 200 m、年均溫約0 ℃的高山及亞高山草地,分布地區(qū)植被類型主要有溫帶叢生禾草草原、溫帶叢生矮禾木、矮半灌木荒漠草原、禾草、薹草高寒草甸、蒿草和雜類草高寒草甸等。根據(jù)狼毒分布的高程范圍和植被類型對研究區(qū)Sentinel-2A 影像進(jìn)行地形及植被掩膜。
結(jié)合谷歌影像特征,對掩膜后研究區(qū)影像的耕地、道路、水體及居民地進(jìn)行采樣和監(jiān)督分類,利用分類結(jié)果進(jìn)一步剔除影像中耕地、道路、水體、居民地等非植被區(qū)域。通過逐層掩膜形成狼毒潛在發(fā)生區(qū)影像,主要包含狼毒群落與其他非狼毒草地群落,以減少其他地類對狼毒識別的干擾,提高分類精度。
根據(jù)狼毒與其他草地植被在生長期明顯的光譜差異,參考前人研究基于Sentinel-2A 影像計算相應(yīng)的狼毒敏感指數(shù),分析花期前和盛花期狼毒敏感指數(shù)的變化規(guī)律選取最優(yōu)敏感指數(shù)[9]。初步選取的狼毒敏感指數(shù)有歸一化指數(shù)NDVIb、NDVIb/r、NDVIr/b和乘積指數(shù)MIbg、MIbr,其計算公式[9]:
式中:ρred、ρgreen、ρblue、ρnir分別為紅、綠、藍(lán)、近紅外波段的地表反射率。
依據(jù)公式(1)~(5)計算狼毒敏感指數(shù),結(jié)合紅、綠、藍(lán)、近紅外波段作為狼毒分類的初選特征,并采用相關(guān)性分析和隨機(jī)森林重要性排序組合的二次降維方法對分類特征進(jìn)行優(yōu)選。首先,對所有特征兩兩之間進(jìn)行Spearman 秩相關(guān)分析,根據(jù)相關(guān)系數(shù)大小篩選出信息重復(fù)率低的特征及特征組。其次,計算每個特征在隨機(jī)森林模型中的貢獻(xiàn)值進(jìn)行隨機(jī)森林重要性排序,依據(jù)重要性大小對特征及特征組做進(jìn)一步篩選,獲取有效的分類特征。最終,通過特征逐步添加的方式構(gòu)建分類特征組合方案作為狼毒分類的信息源[17]。
本研究使用Python 語言基于PyCharm 平臺建立隨機(jī)森林分類模型,將分類特征組合方案作為信息源,46 個狼毒樣方和39 個草地樣方作為分類樣本對研究區(qū)狼毒分布進(jìn)行識別。其中80% 的樣本作為訓(xùn)練集,20%的樣本作為測試集。優(yōu)化后的隨機(jī)森林參數(shù)設(shè)置為,分類器個數(shù)為20 個,最小葉子數(shù)為4,最大深度為5。
利用地面驗證點數(shù)據(jù),在ENVI 中計算混淆矩陣對分類結(jié)果進(jìn)行精度評價,采用的參數(shù)有分類總精度(overall accuracy,OA)、Kappa 系數(shù)、制圖精度(producer accuracy,PA) 和用戶精度(user accuracy,UA)。分類總精度指正確分類樣本數(shù)與總樣本數(shù)之比,Kappa 系數(shù)是用于一致性檢驗的指標(biāo);制圖精度指某類正確分類樣本數(shù)與該類總樣本數(shù)之比;用戶精度指某類正確分類樣本數(shù)與分為該類的樣本數(shù)之比。采用以上參數(shù)對各特征組合方案的分類結(jié)果進(jìn)行精度評價[18-19]。
基于花期前及盛花期的Sentinel-2A 影像計算5 種狼毒敏感指數(shù),并對MIbg、MIbr的指數(shù)值進(jìn)行歸一化處理,在此基礎(chǔ)上計算狼毒樣方及非狼毒樣方的指數(shù)均值(表1)?;ㄆ谇芭c盛花期狼毒與非狼毒樣方的5 種敏感指數(shù)差值絕對值范圍分別為0.020~0.060、0.062~0.299。其中MIbr、MIbg指數(shù)在花期前及盛花期狼毒和非狼毒之間的差值絕對值分別為0.060、0.299 和0.039、0.248,NDVIb、NDVIb/r、NDVIr/b指數(shù)在花期前與盛花期狼毒和非狼毒之間的差值絕對值分別為0.020~0.031 及0.062~0.070。在盛花期5 種植被指數(shù)更明顯地體系了狼毒與非狼毒的差異,特別是MIbg、MIbr體現(xiàn)了兩者之間的最大差異??傮w來看,5 種指數(shù)均不同程度體現(xiàn)了花期前與盛花期狼毒與非狼毒樣方的影像光譜差異,因而初步選擇這5 種植被指數(shù)作為狼毒識別的敏感指數(shù)。
表1 花期前與盛花期敏感指數(shù)均值對比Table 1 Comparison of the sensitivity index mean in the periods before Stellera chamaejasme flowering and at full flowering
分別基于花期前及盛花期的Sentinel-2A 影像計算5 個狼毒敏感指數(shù),結(jié)合4 個多光譜波段,共提取18 項分類特征(表2)。KS (Kolmogorow-Smironov)檢驗顯示A1、A5、A6 等12 個波段不符合正態(tài)分布(P< 0.05),因此選用Spearman 秩相關(guān)分析檢驗各分類特征之間的相關(guān)性(圖3)。A4、A8 與其他各特征之間的相關(guān)性均較弱;而部分波段或植被指數(shù)間的相關(guān)性較高,如A1 與A2、A3、A9、A12、A13,A5與A6、A7,A10 與A11,A14 與A15、A16,A17 與A18,則進(jìn)一步依據(jù)隨機(jī)森林重要性排序結(jié)果從這幾組波段中各篩選出一個代表性波段。
表2 分類特征信息Table 2 Classification feature information
圖3 分類特征相關(guān)性熱力圖Figure 3 Thermodynamic chart of correlation coefficients between
隨機(jī)森林算法中決策樹每個節(jié)點的變量選擇具有隨機(jī)性,在特征變量重要性統(tǒng)計中會出現(xiàn)一定程度的不確定性[17]。為降低誤差將隨機(jī)森林算法運行10 次并計算其均值,得到各特征重要性的最終排序結(jié)果(圖4)。將各特征之間的相關(guān)性分析與其隨機(jī)森林重要性排序結(jié)果相結(jié)合,最終選擇A6、A14、A18、A8、A1、A11、A4 作為狼毒識別的最優(yōu)特征。其中,A6、A14、A18、A8 為盛花期特征,A1、A11、A4 為花期前特征。
圖4 分類特征重要性排序Figure 4 Feature importance ordering
按隨機(jī)森林重要性排序依次遞增添加特征來構(gòu)建特征組合方案,分別采用4、5、6、7 個特征構(gòu)成4 種組合方案作為分類信息源,4 個特征組合方案僅包含盛花期提取特征,5、6、7 個特征組合方案同時包含花期前與盛花期提取的特征。各方案包含的特征波段:A6、A14、A18、A8 (4 個特征),A6、A14、A18、A8、A1 (5 個特征),A6、A14、A18、A8、A1、A11 (6個特征),A6、A14、A18、A8、A1、A11、A4 (7 個特征)。
分別基于4 種特征組合方案采用隨機(jī)森林算法對掩膜的影像進(jìn)行分類,生成研究區(qū)狼毒分布圖(圖5)。4 種方案提取的狼毒群落分布趨勢基本一致,在地勢平坦的谷地附近及河流兩側(cè)狼毒群落分布較密集,總體呈現(xiàn)大面積、斑塊狀、隨機(jī)分散的分布特征。圖5 細(xì)節(jié)圖進(jìn)一步顯示,4 個特征組合方案提取的狼毒群落在整個區(qū)域的密集度及廣泛度相對較高,5、6、7 個特征組合方案提取的狼毒空間分布模式相近。
圖5 基于4 種特征組合方案的狼毒分類結(jié)果(a-d)及其細(xì)節(jié)圖(e-h(huán))Figure 5 Classification results of four feature schemes (a-d) and detailed maps (e-h(huán))
研究區(qū)總面積約4 000 km2,經(jīng)逐層掩膜形成的狼毒潛在發(fā)生區(qū)面積為2 761.51 km2。4 種組合方案分類結(jié)果的面積統(tǒng)計顯示,狼毒占適生區(qū)草地總面積的比例介于17.60%~21.86% (表3),其中,4 個特征組合提取的狼毒面積最大,占比為21.86%,5、6、7 個特征組合方案提取的狼毒面積相近,反映出4 種方案提取的狼毒面積相對較為穩(wěn)定。
表3 4 種特征組合方案提取的狼毒面積Table 3 Area information of the four feature combination schemes
4 種組合方案的分類總精度介于78.67%~84.62%(表4),Kappa 系數(shù)介于0.57~0.69,狼毒識別的生產(chǎn)者精度范圍為81.29%~87.77%,用戶精度為74.69%~81.88%。其中,6 個特征組合和7 個特征組合的分類總精度均大于80%,Kappa 均大于0.60,取得了較好的分類效果。6 個特征組合方案的狼毒分類生產(chǎn)者精度和用戶精度分別為87.77%和81.88%,取得了最優(yōu)的狼毒識別精度。通過分類特征的逐步增加,狼毒識別精度有相當(dāng)程度的提升。多時相的6 個特征組合比單時相的4 個特征組合總體分類精度提高5.25百分點,Kappa 系數(shù)增加了0.10,狼毒生產(chǎn)者精度與用戶精度分別提高了0.72 百分點和7.19 百分點,反映出多時相特征組合在狼毒群落識別中的優(yōu)勢。
表4 4 種組合方案狼毒分類精度評價Table 4 Accuracy evaluation of four feature schemes of Stellera chamaejasme classification
前人研究顯示,去云及掩膜等預(yù)處理能夠有效提高影像質(zhì)量,減少分類干擾信息[20-21]。本研究中采用時間序列去云方法去除Sentinel-2A 影像的云污染,并根據(jù)狼毒生長的地形、植被等特征逐層掩膜狼毒適生區(qū)影像,這些處理能夠有效地減少狼毒分類的干擾因素。在此基礎(chǔ)上,根據(jù)狼毒的物候特征構(gòu)建并篩選狼毒敏感指數(shù),增強狼毒與其他植物群落的影像光譜差異,進(jìn)一步為狼毒識別提供了有效的信息源[22]。Rapinel 等[23]及Noujdina 和Ustin[24]的研究顯示,與單一時相影像相比,多時相影像對草地植物群落及入侵植物旱雀麥(Bromus tectorum)的識別精度更高,體現(xiàn)出多時相遙感在物種識別研究中的優(yōu)勢。高寒草地狼毒群落的空間分布分散、覆蓋度差異大,單一時相的4 個特征組合分類總精度及狼毒用戶精度均低于80%,而包含花期前與盛花期特征的6 個特征組合各項精度參數(shù)均大于80%,分類總精度達(dá)到84.62%,表明多時相Sentinel-2A 影像在大尺度狼毒群落遙感識別中具有較優(yōu)的應(yīng)用潛力。當(dāng)狼毒呈稀疏分布(蓋度小于10%)時與其他群落光譜差異較小,利用Sentinel-2A 的10 m多光譜影像對于低覆蓋度狼毒群落的分布提取具有一定的局限性。
隨機(jī)森林重要性排序方法為蘆葦(Phragmites australis)、大蒜(Alium sativum)等植物的識別提取了有效分類特征[25],不同特征組合的構(gòu)建對比則能夠進(jìn)一步篩選出最優(yōu)的分類特征方案[26]。本研究顯示,將相關(guān)性分析與重要性排序二次降維方法相結(jié)合,可增強、提取和優(yōu)化狼毒識別特征。隨著特征數(shù)目的增加,特征組合方案的狼毒識別精度有不同程度的提高,多時相6、7 個特征組合的分類效果高于4 個特征組合。隨機(jī)森林算法在高寒草地狼毒遙感識別中具有較好的應(yīng)用精度,6 個特征組合結(jié)合隨機(jī)森林取得了最好的狼毒識別結(jié)果,本研究與前人對于入侵植物毛蕊花(Verbascum thapsus)的識別研究[27]具有一致性,進(jìn)一步驗證了該方法在生態(tài)學(xué)分類中的準(zhǔn)確性。
本研究采用花期前與盛花期狼毒入侵退化草甸的Sentinel-2A 多光譜影像,將去云及掩膜處理、敏感指數(shù)提取、特征選擇和隨機(jī)森林分類相結(jié)合在區(qū)域尺度對狼毒分布進(jìn)行遙感識別,主要結(jié)論如下:
1)影像的去云及掩膜處理能夠有效減少干擾信息,提高分類圖像質(zhì)量?;赟entinel-2A 影像計算的NDVIb、NDVIb/r、NDVIr/b、MIbg、MIbr5 種指數(shù)反映了狼毒與非狼毒群落的物候期光譜差異,可作為狼毒識別的敏感指數(shù)。
2) Spearman 秩相關(guān)性分析與隨機(jī)森林重要性排序相結(jié)合的二次降維方法有效篩選了7 個狼毒分類特征。與單時相特征組合相比,多時相特征組合有效提高了狼毒識別精度,是狼毒識別的最優(yōu)數(shù)據(jù)源。
3)結(jié)合隨機(jī)森林算法,由花期前及盛花期特征組成的6 個特征組合方案對于狼毒的識別精度最高,分類總精度達(dá)到84.62%,Kappa 系數(shù)為0.69,是區(qū)域尺度狼毒入侵調(diào)查與監(jiān)測的較優(yōu)方法。