基于Sentinel-2A影像干旱區棉花葉片SPAD數字制圖

唐普恩,丁建麗,3,*,葛翔宇,張振華

1 新疆大學資源與環境科學學院, 烏魯木齊 8300462 新疆大學綠洲生態教育部重點實驗室, 烏魯木齊 8300463 新疆大學智慧城市與環境建模自治區普通高校重點實驗室, 烏魯木齊 830046

新疆作為中國最大的棉花種植基地,其棉花產量能夠左右國際棉價,在我國棉花產業中占有舉足輕重的作用。而棉花作為新疆的主要經濟作物,對新疆農業經濟和社會發展也至關重要[1]。植被生長過程中,葉綠素是參與光合作用的重要色素,其既體現著植物的生理狀態,又決定著光合作用的強弱,是生態系統活力的重要體現[2]。葉綠素濃度直接影響植物代謝活動、生長周期、產量形成,是評價植物長勢的重要指標[3- 5]。葉綠素含量的高低直接或間接的反映自然環境的脅迫作用,特別是干旱區獨有的土壤高鹽背景下鹽脅迫會直接影響植物的生理生長指標和離子吸收,抑制植被光合吸收速率。植被葉片的葉綠素相對含量(SPAD)可直接用以衡量其實際葉綠素含量的高低,植物葉片的SPAD值易于與高光譜數據實現準確對應[6]。因此,快速易行的掌握葉綠素含量空間狀況對干旱區制定科學的水肥一體化管理、潛在災害脅迫 監測具有重大現實意義。

傳統實驗室測量葉綠素的方法復雜、成本高、效率低,且不具備實時、大面積監測能力,遙感技術的實施應用恰逢其時,利用遙感技術可以快速、大尺度、高效的捕獲農作物信息,很好的彌補了傳統手段的不足[7- 8]。國內外研究學者利用高光譜技術對植被葉綠素含量進行監測方面取得了很多成果。Lu等認為由紅色和近紅外反射率組合計算的可見光與植被的生物物理和生物化學特性有很好的相關性[9]。郭超凡等建立葉面尺度下不同包絡線去除衍生轉換光譜也葉綠素含量估算模型[10]。田明璐等使用多遠逐步回歸與線性回歸構建多種光譜參數建立與葉綠素相對含量的反演模型中,發現前者精度優于后者[11]。馮偉等研究了小麥葉片色素密度與冠層高光譜參數的定量關系,認為群體葉片色素密度的敏感波段主要分布在可見光區, 而紅邊區域導數光譜表現更顯著[12]。李成等[13]發現紅邊和近紅外波段構建的植被指數的蘋果樹冠層葉綠素含量估算模型更準確,更穩定。已被證明光譜與葉綠素含量間存在機理,特別是在紅邊波段,由于作物在紅光波段的強吸收以及在近紅外波段的強反射,在680—760 nm的光譜范圍內,紅光波段反射率下降,近紅外波段陡然上升,形成“吸收谷”,是綠色作物光譜最為明顯的特征。紅邊波段在實踐中光譜指數通常被視為遙感光譜數據與制備生理化學參數檢測橋梁。這些光譜指數的設計,將不良影響降低的同時具備較好簡潔性、魯棒性和準確性。然而,由于植物葉片結構和觀測角度等因素的影響,基于經驗統計的方法可能表征較差。另外,基于物理的輻射傳輸機理模型根據葉片和土壤特性,冠層結構和陽光傳感器的幾何形狀模擬冠層內的光傳播,并提供對光與植被之間相互作用的有價值的物理理解,也被廣泛用于植被重要的生理參數的反演(如葉綠素、氮含量、葉面積指數(LAI)、類胡蘿卜素等)[14- 16]。李振海等比較了僅使用SAIL+PROSPECT(PROSAIL)冠層反射率模型與農學先驗知識(APK)關聯的結果,經分析證實了PROSAIL模型反演生物物理變量的操作潛力[17]。但其容易受地域、植被類型、生長季變化等因素影響,致使其用于遙感像元尺度反演時存有較大誤差,且模型普適性較差[18- 19],但植被中的其他生理化學參數與葉綠素間具有很強的協同性和相關性,如葉片是植被光合作用的主要器官,葉面積大小和作用時間都與健康狀態有關,而葉綠素又是植被生長狀況和光合作用的重要指標。類胡蘿卜素與也葉綠素在光譜反射中存在遮蔽影響。

在干旱區復雜的下墊面和鹽分脅迫的背景下,單獨一種方法或一種屬性的協變量是無法滿足精準農業的需求。基于此,本文以囊括3個紅邊波段的Sentinel- 2A為本底,聯合最適光譜指數和基于PROSAIL物理模型的生物協變量集成估算棉花葉片SPAD值。以期獲取高精度、高分辨率區域葉綠素分布,服務以精準農業為抓手的一帶一路核心區農業信息化。

1 數據與方法

1.1 研究區概況

研究區位于新疆維吾爾自治區昌吉回族自治州(包含瑪納斯縣,沙灣縣和石河子市),地處天山北麓中段,準格爾盆地中部古爾班通古特沙漠南緣,地理坐標范圍44°20′—45°20′ N,85°20′—86°40′ E。氣候、地貌、水文條件、植被具有明顯的垂直地帶性分布[20- 21],海拔由南向北遞減,瑪納斯河流域主要形成“山地-綠洲-荒漠”系統流域,綠洲主要以片狀或帶狀分散在山前沖積平原和沖洪積扇扇緣[22]。年平均氣溫4—7℃,極端最高氣溫43℃,極端最低氣溫-42℃,多年平均降雨量150—200 mm,多年平均蒸發量達1500—2100 mm,屬于典型的溫帶大陸性干旱半干旱氣候。

1.1.1土壤特征

研究區內的土壤類型十分豐富主要有灰漠土、棕鈣土、風沙土、沼澤土、草甸土、鹽土、灌耕土和灌淤土等。由于干旱少雨、植被稀疏,大多數土壤呈現出鹽堿化、有機質少、沙性大和粘重板結等現象。但由于農業生產活動的發展,使得土壤逐漸熟化,有機質增加,土地承載力與環境容量提高,使得研究區成為了重要的農業生產基地。本研究主要以棉花為研究對象,其種植地土壤屬性如表1所示。

表1 研究區內土壤屬性

1.1.2植被類型

不同質地的土壤其理化性質差異性大,所以生長的植被也會有所不同。研究區內平原植被主要以人工種植作物棉花(Gossypiumhirssittum)、玉米(Corn)為主,天然植被主要是耐鹽耐旱植物如紅柳(Redwillow)、胡楊(Populuseuphratica)、豬毛菜(Salsola)、鹽爪爪(Kalidiumfoliatum)等。

1.2 數據采集

1.2.1棉花葉片SPAD值

為估算模型構建和驗證葉綠素精度,課題組在2019年7月進行野外調查,沿瑪納斯河流域根據地形地貌選取具有典型代表的104個樣點并采樣。利用日本美能達公司生產的SPAD快速葉綠素測定儀,在面積為10 m×10 m的樣方內,從東、西、南、北方向的四個角和中間位置上對植株冠層的上、中、下采集20片健康、大小均勻的葉片測定SPAD值,并將五個點20個測定值的算術平均值作為對應樣方內特定植被葉片的SPAD值。

1.2.2影像數據

歐洲委員會和歐空總署局聯合倡議研發的全球環境與安全監測計劃的第二顆衛星,Sentinel- 2A、B衛星分別于2015年6月23日、2017年3月7日發射成功,雙星在軌運行的重返周期為5d,最佳分辨率可達10 m[23],其遙感影像擁有13個分辨率不同的光譜波段(10、20、60 m),包括可見光、紅邊、近紅外、水汽、卷云以及短波紅外波段(表2)。根據采樣時間及云量(< 5%),擇取2019年7月26日、28日兩景影像(https：//scihub.copernicus. eu/dhus/#/home)。首先將下載的兩景影像經Sen2Cor模塊對影像實時輻射定標和大氣校正,再利用SNAP軟件使用三次內插卷積(Cubic convolution interpolation)的方法對大氣校正后的波段進行空間分辨率為10 m的重采樣,之后轉入ENVI 5.1提取本研究所用9波段(B2、B3、B4、B5、B6、B7、B8a、B11、B12),然后將其在ArcGIS 10.0使用 (Gram-Schmidt Pan Sharpening)融合方法將研究區內兩景影像拼接、融合,并按照研究區范圍使用取子區(subset)功能將影像裁剪。

表2 Sentinel- 2A影像波段光譜參數

1.3 最適光譜指數

本研究選擇Sentinel- 2A影像(B2、B3、B4、B5、B6、B7、B8a、B11、B12)9個光譜波段,其中包含3個紅邊波段,在農作物植株葉綠素含量估算中具有很大的應用潛力。因此,本研究在借鑒現有植被指數構建方法的基礎上,引入紅邊波段替代可見光波段或者近紅外波段,利用紅邊波段與近紅外波段或者紅波段的光譜反射特征表達不同葉綠素含量的差異。構建NDVI、DVI和RVI光譜指數,以采樣點的棉花葉片SPAD值來反應Sentinel- 2A影像像元的數據,計算其與棉花葉片SPAD值相關性,以期找到評估葉綠素含量的最適光譜指數,如表3所示。

表3 植被指數

1.4 PROSAIL輻射傳輸模型

PROSAIL模型是將PROSPECT與SAIL模型的集成,已成功應用于不同的生態系統(如森林和作物)常用于植被生化物理參數的反演[26-29]。將PROSPECT模型的輸出的葉片反射率和透射率,作為SAIL模型的輸入參數,再結合土壤反射率、葉面積指數、葉傾角分布、太陽天頂角,觀測天頂角等相關參數,即可得到植被冠層的反射率和透射率。遙感影像通過大氣校正,也可以獲取地表植被冠層反射率,將遙感影像與植被參數葉綠素通過物理過程聯系起來,生成一個包含所有可能的反射率的數據集,完成了模型模擬過程[30]。

本研究所用波段的光譜反射率在490—870 nm之間,在該波段范圍內LAI、葉綠素以及干物質含量均為敏感參數,將其確定為可變參數,類胡蘿卜素在模型中表現為不敏感,設為固定值,反演參數主要通過野外實測、文獻查閱、經驗值及LOPEX94數據庫獲得[31],表4為參數設置。

表4 PROSAIL模型參數設置

1.5 棉花葉片SPAD值的模型構建

在構建棉花葉片SPAD值的反演模型,是通過R- 3.6.2中sample函數k-s劃分70%(n=72)的建模集和30%(n=32)的驗證集并用set.seed函數固定選中的數據集。

隨機森林算法(Random forest,RF)作為一種較新的集成學習算法,是由Breiman[32]提出的基于決策樹分類器融合算法,算法可以綜合利用多種特征變量建模并評價,可以看成由多個弱預測器(決策樹)集成的強預測器。在隨機森林模型中[31],通過選擇原始數據集和訓練數據集(即校準數據集)的隨機樣本,使用確定性算法構建每棵樹,每一棵樹都依賴一個隨機向量,通過對數據集的列變量和行變量觀測進行隨機化,生成多個分類樹,最終將分類樹結果進行綜合。本文研究使用R- 3.6.2中隨機森林工具包對棉花葉綠素SPAD值進行預測分類,并產生重要變量性估計,研究模型的重要性分類指標由平均下降精度參數(Mean Decrease Accuracy)提供[33]。

本文將野外測得的植被葉綠素SPAD值作為模型因變量,將最適光譜指數(DVI、NDVI、RVI),Sentinel- 2A中9個波段(B2、B3、B4、B5、B6、B7、B8a、B11、B12),基于物理模型生物協變量(葉面積LAI、冠層含水量LAI-Cw、冠層葉綠素LAI-Cab、植被覆蓋度fcover、光合有效輻射fapar)共計17個特征參數分別作為自變量并構建3種方案(表5),在R- 3.6.2中利用隨機森林按照設計方案對采樣地的棉花建立植被葉片SPAD的估算模型。

表5 建模方案

1.6 模型精度評價

本研究為了評估上述3種預測模型的性能,選取了以下幾種常用的模型驗證標準。分別為決定系數(Determination coefficients,R2)、均方根誤差(Root mean squared error, RMSE)、一致性相關系數(Concordance correlation coefficient,ρc),研究思路如圖1所示。

圖1 研究流程圖Fig.1 Research flow process

圖2 樣本統計分析圖Fig.2 Sample statistical analysis chart

2 結果與分析

由圖2可知棉花葉片SPAD值全集均值為59.02,變異系數為0.092,變異系數小,說明數據離散程度小較為集中。如圖2所示,棉花葉片SPAD值全集、建模集、驗證集在統計特征中保持相似的統計形態,均值分別為：59.02、58.66、59.38,SD分別為5.43、5.25、5.43,由此可見樣本劃分合理。

2.1 棉花葉片SPAD值與高光譜植被指數的相關性

圖3為棉花葉片SPAD值與Sentinel- 2A波段(B2、B3、B4、B5、B6、B7、B8a、B11、B12)兩個隨機光譜波段的2D指數(DVI、NDVI、RVI)關系。右側的彩條表示r大小,深紅、深藍代表SPAD值與波段組合正、負相關性較高。結果表明,DVI、NDVI、RVI與棉花葉片SPAD值高度相關性的波段主要集中在紅邊(Red edge)波段(B5、B6)和近紅邊(Narrow NIR)波段B8a處。RVI與棉花葉片SPAD值相關性最高其為0.767,P**=0.195,其表達式B8a/B5。NDVI其r為0.720,DVI其r為0.692。說明紅邊波段對棉花葉綠素有較高的響應,3種指數在相應波段下對棉花葉片SPAD值均具有較高相關性。

圖3 棉花葉片SPAD值與Sentinel- 2A(B2、B3、B4、B5、B6、B7、B8a、B11、B12)中任意兩波段構建的DVI、NDVI、RVI的相關系數Fig.3 Correlation coefficient between SPAD value of cotton leaves and DVI、NDVI、RVI constructed in any two bands of sentinel- 2A (B2、B3、B4、B5、B6、B7、B8a、B11、B12) B2：藍光 Blue band；B3：綠光 Green band；B4：紅光 Red band；B5：植被紅邊Vegetation Red Edge band；B6：植被紅邊Vegetation Red Edge band；植被紅邊B7：Vegetation Red Edge band；B8a：近紅外(窄)Narrow NIR band；B11：短波紅外SWIR band；B12：短波紅外SWIR band

2.2 特征變量重要性分析

本研究通過隨機森林進行建模,并對參與不同方案建模的特征變量做出重要性分析,通過圖4可以發現在三種方案中,紅邊波段的貢獻率都大于6%,由紅邊波段參與構建的最適光譜指數在重要性占比中也較高,這也說明了紅邊波段與葉綠素響應度較強；此外三方案中RVI與LAI-Cab的貢獻始終最大其占比均超過了10%,在方案3中重要性分別10.22%和10.56%,由此可見植被指數與物理模型對于估算葉綠素含量的貢獻兩者均較為重要,皆不可忽視。

圖4 特征變量重要性Fig.4 The importance of characteristic variablesM1:方案1; M2:方案2; M3:方案3;LAI：葉面積指數 Leaf area index；LAI-Cw：冠層含水量 Canopy water content；LAI-Cab冠層葉綠素 Canopy chlorophyll；fcover：植被覆蓋度 Vegetation coverage；fapar：光合有效輻射 Photosynthetically active radiation；

2.3 棉花葉片SPAD值的估算模型評價結果

根據表6可知在所有棉花葉片SPAD估算模型中,由驗證集效果可知,方案3 >方案1 >方案2,其中方案3 (R2= 0.8258,RMSE=3.5660,ρc= 0.7205)的效果甚佳。比單一的使用最適光譜指數或物理模型生物協變量構建的模型精度提高0.05—0.07左右,這說明將Sentinel- 2A波段與植被指數、物理模型集成進行優勢互補可以提高棉花葉片SPAD值估算精度。同時一致性ρc也為三個方案中最高的。

表6 棉花葉片SPAD預測值與實測值分布

圖5 棉花葉片SPAD預測值與實測值相關性及精度驗證Fig.5 Correlation and accuracy verification of SPAD predicted value and measured value of cotton leaves

圖6 葉片SPAD遙感監測圖Fig.6 Remote sensing monitoring map of blade SPAD

3 討論

作物葉片的SPAD值反映了作物葉綠素含量的相對高低,通過了解葉綠素含量及其變化,對作物的營養診斷、長勢監測及產量估算具有重要意義[34]。本研究以Sentinel- 2A為數據源,將PROSALL模型與最適光譜指數進行結合對棉花葉片SPAD值進行估算研究。

本文將3種方案進行對比后,選擇方案3繪制圖6葉片SPAD遙感監測圖,由圖6可知研究區內葉綠素的空間分布大致為中部高四周低,其主要受水源、土壤鹽分、土壤有機質及灌溉方式影響,研究區的西北角及東北角屬于沙漠地區,土壤水分少、鹽分高、有機質含量低對棉花生長有一定抑制作用。研究區內的某些零星區域葉綠素偏高,是由于種植其他作物對其估算具有一定不準確性。同時通過表1可以發現土壤屬性存在較大差異會對所生長植被的生理參數產生較大影響。因此,今后研究要對棉花葉綠素進行空間進行分析,除了要掌握主要水源分布,還應將土壤性質(鹽分、水分、有機質等)作為重要的參考因素作為土壤特征變量參與建模,同時也要了解作物的種植分布特點。

Sentinel- 2A具有豐富的光譜波段信息,選取特征波段構建光譜指數在不同程度上抵消了由環境變化而引起的噪聲,可以用來定量估算作物的理化參數[35]。因此,本研究利用Matlab,選取最優光譜波段構建3個2D最適光譜指數(DVI、NDVI、RVI),由于選取的最優波段均為紅邊波段與葉綠素具有較強的響應,且采樣時棉花正處于花鈴期,此時棉花生長旺盛,葉片呈深綠色,對光譜的反射率較高。這與洪帥[36]等研究的結果一致,棉花的生育時期冠層葉片的葉綠素含量與紅邊指數的相關性較高。但僅構建3種光譜指數無法充分利用高光譜數據,從而在一定程度上限制估算模型的精度[10],充分利用高光譜數據構建多種光譜指數之后進行最優篩選以提高模型精度值得今后研究。

在對3種模型方案不同的特征變量進行重要性分析,易知紅邊波段與基于紅邊波段構建的光譜指數對估算模型有著關鍵作用。生物協變量中LAI-Cab貢獻最大,其在物理模型中就是表示葉綠素,同時也說明了物理模型對葉綠素估算起到一定作用。在以最適光譜指數和物理模型生物協變量構建的方案1、方案2中,RVI與LAI-Cab重要性占比為13.90%、15.24%,但在方案3中兩者占比雖依舊為最大但卻有所下降,可能是由于在特征變量增多的原因,未對變量進行篩選,往后的研究上應在優選特征變量的基礎上進行建模[37]。另外,PROSAIL模型反演過程中受到自身參數不確定性及敏感性以及植被生理參數之間的相互作用可能帶來誤差影響,本研究由于采樣條件有限,故敏感參數設置僅為LAI、Cab,沒有考慮其他參數[38],從而會影響模型精度,今后研究可多獲取參數實測值來解決。

本研究構建的三種模型方案,都是基于隨機森林回歸算法估算棉花葉片SPAD值的,出現了低值高估、高值低估的現象,同時給出不同變量對模型影響評分,在建模時進行參數優化,可以很好地學習特征變量與SPAD值關系,提高模型的預測精度。這與前人所做出的結果一致[39]。但構建的模型要想快速的監測棉花葉綠素含量,還需要獲取多年數據進行改進。其次,本文研究區以新疆人工綠洲為主,它是由綠洲農田和綠洲聚落生態系統構成的人工生態系統主要受人類控制,在人類精心維護下,綠洲人工生態系統內部生態環境趨于優化,生態穩定性較強,而生態環境的影響通常以作物物候周期為尺度表現在作物生理化學參數上。基于此,本文嘗試通過對棉花物候周期中的某一時期的SPAD高精度制圖可能,在后續研究中將棉花的而整個物候周期作為研究的重點,從而了解生態環境對于棉花生長的影響,通過整個物候周期生理參數的變化來分析生態環境的質量發展變化,以期實現對干旱區作物生長狀況的實時監測,為精準農業提供新的技術途徑。

4 結論

本文以Sentinel- 2A為高光譜數據源,將Sentinel- 2A光譜波段、植被指數和PROSAIL模型三者進行組合構建棉花葉片SPAD值的估算模型。利用Sentinel- 2A光譜波段在構建植被指數時可以發現建立的2D指數(DVI、NDVI、RVI)與棉花葉片SPAD值相關系數最好的波段在紅邊波段(B5、B6)和近紅邊波段(B8a)處,說明紅邊波段與葉綠素含量響應明顯,基于紅邊波段構建的植被指數對估算模型起到關鍵作用。在構建的3個方案模型中,利用實測數據對模型進行驗證,結果表明：將Sentinel- 2A光譜波段、植被指數和PROSAIL模型三者利用隨機森林算法構建的模型表現出了最佳SPAD預測能力,實現了棉花葉片SPAD值的遙感監測。綜上所述,單一的以經驗模型植被指數或物理模型構建的反演方法對植被葉綠素含量具有一定的預測能力,但將兩者在隨機森林中集成的模型可以提高干旱區和半干旱區棉花葉片SPAD值的預測精度,可以推廣至其他植物的葉綠素監測,為合理水肥配置,精準種植,農作物生長過程監測提供數據支撐,促進精準農業蓬勃發展。