蟲害脅迫下毛竹葉綠素含量高光譜估算方法

李 凱，陳蕓芝*，許章華, ，黃旭影，胡新宇，汪小欽

1. 福州大學(xué)衛(wèi)星空間信息技術(shù)綜合應(yīng)用國家地方聯(lián)合工程研究中心，福建 福州 350116 2. 空間數(shù)據(jù)挖掘與信息共享教育部重點實驗室，福建 福州 350116 3. 福州大學(xué)環(huán)境與資源學(xué)院，福建 福州 350116 4. 南京大學(xué)國際地球系統(tǒng)科學(xué)研究所，江蘇 南京 210093

引 言

毛竹作為我國竹林產(chǎn)業(yè)的重要組成部分，在我國林業(yè)的經(jīng)濟發(fā)展中起到關(guān)鍵作用。剛竹毒蛾是毛竹的主要蟲害，嚴重威脅毛竹的生長，制約了竹林產(chǎn)業(yè)的可持續(xù)發(fā)展[1]。葉綠素是綠色植被吸收光能的物質(zhì)，作為參與植被光合作用最重要的色素，直接影響到綠色植被對光能的利用[2]。當(dāng)植被受蟲害脅迫而“失綠”時，常伴隨著內(nèi)部葉綠素含量的變化，并在光譜特性上發(fā)生變化[3]。因此，植被葉綠素含量的狀態(tài)已成為監(jiān)測植被蟲害的一項重要指標。

高光譜數(shù)據(jù)能夠提供豐富的光譜信息，可以更加精確地檢測植被色素光譜特征，估計植被的生化參數(shù)進而診斷植被的健康狀態(tài)。實際應(yīng)用中，高光譜數(shù)據(jù)在提供豐富信息的同時，也存在信息冗余量大的問題，導(dǎo)致統(tǒng)計參數(shù)的估計偏差增大[3]。特征波長的提取可有效地對高光譜數(shù)據(jù)進行降維，從而達到優(yōu)化模型的效果。通過對高光譜數(shù)據(jù)進行光譜變換能夠增強特征波長信息，進而提高模型的估算精度[4]。國內(nèi)外學(xué)者通過利用不同光譜的特征光譜位置、 多波段的光譜指數(shù)、 特征光譜參數(shù)等光譜信息開展了大量的植被生化參數(shù)研究。如王雪等[5]利用競爭性自適應(yīng)重加權(quán)采樣算法進行波長篩選，并結(jié)合SPXY樣本劃分法建立玉米含水率的偏最小二乘估算模型。宮兆寧等[6]系統(tǒng)分析了典型濕地植被的多項光譜特征參數(shù)，確定了“紅邊”位置、 565和735 nm歸一化光譜指數(shù)為估算葉綠素含量的最佳光譜特征參數(shù)。劉文雅等[7]利用相關(guān)系數(shù)法篩選出與馬尾松葉綠素含量高度相關(guān)的植被指數(shù)，并對馬尾松冠層葉綠素含量進行估算。毛博慧等[8]采取遺傳算法進行波長篩選，利用優(yōu)選波長組合植被指數(shù)，建立基于敏感植被指數(shù)的冬小麥葉綠素含量估算模型。

為提高竹葉葉綠素含量估算精度，通過對不同蟲害程度的毛竹葉片光譜和葉綠素含量進行測定，采用連續(xù)投影算法(successive projections algorithm，SPA)對光譜數(shù)據(jù)進行壓縮，提取特征波長，對比不同樣本集劃分方法對葉綠素含量建模精度的影響，建立蟲害脅迫下竹葉葉綠素含量的多元逐步回歸(multiple stepwise regression，MSR)模型，以期為今后毛竹蟲害監(jiān)測提供理論依據(jù)。

1 實驗部分

試驗于2019年4月在福建省順昌縣內(nèi)毛竹林進行；該縣位于福建省西北部，季風(fēng)氣候明顯，全年雨量充沛，森林資源豐富，是福建省的毛竹生產(chǎn)基地縣。根據(jù)剛竹毒蛾蟲害發(fā)生的特點，隨機選取6個不同樣方，參照《植物病理學(xué)》的蟲害劃分標準，將毛竹樣本分為健康、 輕度蟲害、 中度蟲害及重度蟲害4個等級[9]，對采集的不同蟲害等級毛竹葉片進行光譜及葉綠素含量的現(xiàn)場測定。

1.1 數(shù)據(jù)采集

使用美國ASD FieldSpec 3便攜式光譜儀進行光譜數(shù)據(jù)測量，光譜儀波長范圍為350～2 500 nm，其中350～1 000 nm波長的采樣間隔為1.4 nm，1 000～2 500 nm波長的采用間隔為2 nm，共2 151個波段。選擇不同蟲害程度的葉片進行測定，測定時間為10:00—15:00，探頭距葉片15 cm處垂直向下測量，共采集102個葉片樣本，每個樣本測得10條光譜，取平均值作為該樣本的光譜，每組數(shù)據(jù)采樣前均進行標準白板校正。為保證數(shù)據(jù)的精度，每個樣本光譜測量后采用日本SPAD-502葉綠素計測定葉片相對葉綠素含量(soil and plant analyzer development，SPAD)。

1.2 光譜數(shù)據(jù)處理

利用ASD FieldSpec 3光譜儀配套軟件ViewSpecPro對采集的光譜數(shù)據(jù)進行處理，計算葉片光譜反射率，保留400～900 nm范圍內(nèi)的光譜數(shù)據(jù)，并通過Savitzky-Golay濾波算法對光譜數(shù)據(jù)進行平滑去噪處理。對原始光譜(original spectrum，OS)進行一階導(dǎo)數(shù)(first derivative，F(xiàn)D)處理，可消除與波長無關(guān)的漂移[10]。包絡(luò)線去除(continuum removal，CR)是將原始光譜曲線與包絡(luò)線曲線進行比值運算后得到，該方法能夠增強各生化參量敏感波段的光譜吸收特征[11]，并計算包絡(luò)線去除后光譜的一階導(dǎo)數(shù)(continuum removal-first derivative，CR-FD)。

1.3 連續(xù)投影算法

連續(xù)投影算法是一種使矢量空間共線性最小化的前向變量選擇方法，通過對光譜向量進行投影操作，能夠找尋到最具代表性的變量組，獲取具有最小共線性的波長子集。該算法的步驟如下[12]：

(1)在光譜矩陣中任選一條光譜列向量作為起始向量；

(2)分別計算剩余列向量在起始向量的正交平面上的投影向量；

(3)挑選出最大投影作為下一次投影的起始向量，直到挑選變量個數(shù)達到所需個數(shù)；

(4)將提取的所有波長組合進行多元線性回歸，最小的RMSE所對應(yīng)的組合即為最優(yōu)波長組合。

1.4 SPXY算法

樣本數(shù)據(jù)集的劃分會在一定程度上影響模型估算的精度，SPXY(sample set partitioning based on joint x-y distance)法[13]是一種同時考慮樣本光譜以及生化參量歐氏距離的樣本集選擇方法，使選擇的樣本集更具有代表性，從而改善所建模型的估算能力。計算公式如式(1)和式(2)

(1)

(2)

式中，J為光譜中的波長數(shù)，xp(j)和xq(j)分別為樣本p和q的光譜值，dx(p，q)和dy(p，q)分別表示樣本p和q的x變量和y變量的歐式距離。該算法的步驟如下：

(1)計算樣本間的兩兩距離，選擇距離最大的兩個樣本作為起始樣本；

(2)分別計算剩余樣本與已選樣本間的距離，對于每個剩余樣本而言，選擇與起始樣本間的最短距離；

(3)在所選擇的最短距離中選擇最長距離所對應(yīng)的樣本，作為第三個樣本，以此逐步選擇，直至所選樣本數(shù)達到目標個數(shù)。

1.5 模型精度檢驗

采用決定系數(shù)(R2)，均方根誤差(RMSE)和相對分析誤差(RPD)3個指標對葉綠素含量估算模型進行精度評價。其中，R2越大模型估算效果越好；RMSE越小模型估算效果越好；當(dāng)RPD>3時，表明模型的估算效果極好，當(dāng)RPD為3～2.25時，表示模型的估算效果較好，當(dāng)RPD為2.25～1.75時，模型的估算效果一般，當(dāng)RPD<1.75時，模型的估算效果較差[14]。各指標計算公式分別如式(3)—式(5)

(3)

(4)

(5)

2 結(jié)果與討論

2.1 蟲害毛竹葉片光譜特征

圖1為不同蟲害程度竹葉光譜曲線，可以看出，健康竹葉在520～570 nm波長范圍因葉綠素吸收較少而形成了一個綠色波段的小反射峰，稱之為“綠峰”。在640～700 nm波長范圍因葉綠素的強吸收形成了一個紅色波段吸收谷，稱之為“紅谷”。由于葉片的多孔薄壁細胞結(jié)構(gòu)特性對近紅外光強烈的反射，形成光譜上的最高峰。當(dāng)毛竹遭受蟲害侵蝕后光譜發(fā)生明顯變化，主要表現(xiàn)為可見光波段范圍內(nèi)的“綠峰”和“紅谷”的逐漸消失，“紅邊”斜率減小，近紅外波長反射率降低，隨著蟲害程度加重，這一現(xiàn)象更加明顯。導(dǎo)致這一現(xiàn)象發(fā)生的主要原因是毛竹遭受蟲害侵蝕后，竹葉葉綠素含量相對減少，細胞結(jié)構(gòu)遭到破壞。因此，根據(jù)光譜間的差異，進行竹葉葉綠素含量的估算是可行的。

圖1 不同蟲害程度的竹葉光譜曲線Fig.1 Spectra of bamboo leaves under pests stress

2.2 葉綠素含量與光譜數(shù)據(jù)之間的相關(guān)性分析

利用Pearson相關(guān)分析法計算不同光譜與竹葉葉綠素含量的相關(guān)性，如圖2所示。從圖2可以看出，原始光譜(original spectrum，OS)在可見光波段與葉綠素含量呈負相關(guān)，在波長644 nm處相關(guān)性最高，相關(guān)系數(shù)為-0.759，在近紅外波長基本不相關(guān)；包絡(luò)線去除(continuum removal，CR)與葉綠素含量的相關(guān)性曲線與原始光譜相似，相關(guān)系數(shù)有明顯提升，其中在波長702 nm處相關(guān)性最高，相關(guān)系數(shù)為-0.852；一階導(dǎo)數(shù)(first derivative，F(xiàn)D)與葉綠素含量的相關(guān)系數(shù)呈現(xiàn)正負起伏變化，在波長629 nm處相關(guān)性最高，相關(guān)系數(shù)為-0.776；包絡(luò)線去除一階導(dǎo)數(shù)(continuum removal-first derivative，CR-FD)相關(guān)性曲線與一階導(dǎo)數(shù)相似，在波長724 nm處的相關(guān)系數(shù)高達0.862?？梢姡袢~葉綠素含量與不同光譜相關(guān)性較高的波長均處于600～750 nm，通過光譜變換可以有效提升與葉綠素含量的相關(guān)性。

圖2 竹葉葉綠素含量與光譜的相關(guān)系數(shù)Fig.2 Correlation coefficient between chlorophyll content and spectra of bamboo leaves

2.3 特征波長提取

使用連續(xù)投影算法對不同光譜數(shù)據(jù)集進行特征波長選擇，結(jié)果如表1所示。OS所選波長個數(shù)為7個，集中分布在“綠峰”、 “紅谷”、 “紅邊”位置；CR有5個波長被選擇，所選波長主要分布在“紅邊”位置；FD有4個波長被選擇，所選波長位置較為分散；CR-FD有5個波長被選擇，所選波長大多分布在“紅邊”位置。從各數(shù)據(jù)集波長選擇結(jié)果可以看出，不同光譜數(shù)據(jù)集選擇的波長位置及個數(shù)有所差異，但所選波長位置基本上集中在綠光、 紅光、 “紅邊”3個區(qū)域，各光譜數(shù)據(jù)集都有多個與葉綠素含量相關(guān)性較高的波長(600～750 nm)被選擇。經(jīng)包絡(luò)線去除后的近紅外波長光譜值較為接近，且與可見光波長的光譜值差異變大，從而導(dǎo)致在近紅外范圍有1個波長被選擇。由于算法本身的原因?qū)е乱自诓ㄩL范圍邊緣處選取波長[15]。

表1 SPA波長選擇結(jié)果Table 1 SPA wavelengths selection result

2.4 模型的建立和預(yù)測

共采集不同蟲害程度竹葉樣本102個，利用SPA算法對各光譜數(shù)據(jù)集特征波長進行提取后，分別利用SPXY法和隨機法進行樣本集的劃分，建模集與驗證集的劃分比例為7∶3，然后建立基于樣本劃分及特征波長的多元逐步回歸模型，并計算每種模型在驗證集中的葉綠素含量估算精度，結(jié)果如表2所示。從表2可以看出，利用SPXY法對樣本進行劃分后的模型對葉綠素含量的估算精度明顯優(yōu)于基于隨機法樣本劃分的模型估算精度?；赟PXY樣本劃分法的OS，F(xiàn)D，CR以及CR-FD模型的R2均高于0.8，其中OS和CR-FD模型R2高達0.840和0.835；OS和CR-FD模型的RMSE分別為2.925和2.604；CR和CR-FD模型的RPD值分別達到2.294和2.364，表明這2種模型具有較好的預(yù)測能力。圖3為基于SPXY樣本劃分法的竹葉葉綠素含量預(yù)測結(jié)果，CR和CR-FD模型的葉綠素含量預(yù)測效果較好。綜合考慮，選CR-FD模型為竹葉葉綠素含量的最優(yōu)估算模型。表3為各模型的回歸方程式。

表2 竹葉葉綠素含量估算模型參數(shù)Table 2 Estimation model parameters of chlorophyll content in bamboo leaves

表3 竹葉葉綠素含量多元逐步回歸模型Table 3 Multiple Stepwise Regression model of chlorophyll content in bamboo leaves

圖3 葉綠素含量預(yù)測值與實測值比較Fig.3 Comparison of predicted and measured chlorophyll contents

3 結(jié) 論

(1)蟲害脅迫下竹葉光譜特征發(fā)生明顯變化，其可見光波段范圍內(nèi)的“綠峰”和“紅谷”逐漸消失，“紅邊”斜率減小，近紅外波長反射率降低。

(2)對剛竹毒蛾蟲害脅迫下竹葉原始光譜進行光譜變換可以有效提升與葉綠素含量的相關(guān)性，可見光區(qū)域是估算竹葉葉綠素含量的敏感區(qū)域，包絡(luò)線去除一階導(dǎo)數(shù)光譜在724 nm處與竹葉葉綠素含量的相關(guān)系數(shù)最高，達到0.826。

(3)采用連續(xù)投影算法提取的特征波長主要分布于綠光、 紅光、 “紅邊”位置，且每種光譜數(shù)據(jù)集都有多個與葉綠素含量相關(guān)性較高的波長(600～750 nm)被選擇。

(4)基于SPXY樣本劃分法對建模樣本進行選擇，建立的葉綠素含量估算模型的R2，RMSE，RPD相比于隨機樣本劃分法都有明顯的提升，其中R2和RPD平均提高0.1和0.5，RMSE平均降低0.7。

(5)多元逐步回歸可對連續(xù)投影算法提取的特征波長進一步壓縮，結(jié)合連續(xù)投影算法與SPXY樣本劃分法建立的包絡(luò)線去除和包絡(luò)線去除一階導(dǎo)數(shù)多元逐步回歸模型可有效地對竹葉葉綠素含量進行估算，其中包絡(luò)線去除一階導(dǎo)數(shù)模型的估算精度最高，為剛竹毒蛾蟲害脅迫下竹葉葉綠素含量的最優(yōu)估算模型。