于豐華，馮 帥，趙依然，王定康，邢思敏，許童羽

(1 沈陽農(nóng)業(yè)大學(xué) 信息與電氣工程學(xué)院，遼寧 沈陽 110866； 2 遼寧省農(nóng)業(yè)信息化工程技術(shù)中心，遼寧 沈陽 110866)

粳稻是水稻的一種，在我國東北地區(qū)有著較大的種植面積。粳稻生育期較長，蛋白質(zhì)含量較高，食味好[1]。葉綠素含量是表征粳稻生長狀態(tài)的一個重要性狀指標(biāo)，無損、快速、大面積反演粳稻葉綠素含量一直都是國內(nèi)外精準(zhǔn)農(nóng)業(yè)學(xué)者的重要研究方向[2-3]。無人機低空遙感平臺是近幾年快速發(fā)展的一種近地遙感手段，其平臺的靈活性高、獲取數(shù)據(jù)豐富、覆蓋范圍廣，相較于其他遙感手段具有一定的優(yōu)勢，因此利用無人機高光譜遙感方式在近地尺度對粳稻冠層葉綠素含量進(jìn)行反演估算對于輔助田間精準(zhǔn)施肥、施藥具有重要的現(xiàn)實意義。

目前，國內(nèi)外在利用高光譜技術(shù)估算作物葉綠素含量已經(jīng)取得了一定的研究成果。尼加提·卡斯木等[4]研究了820～940 nm的最大反射率及反射率綜合所建立的偏最小二乘回歸法(Patial leastsquares regression，PLSR)估算春小麥葉綠素含量的效果較為理想。毛博慧等[5]研究了冬小麥冠層325～1 075 nm范圍的反射光譜，采取遺傳算法對光譜特征參數(shù)尋優(yōu)并結(jié)合相關(guān)分析結(jié)果，利用最小二乘－支持向量回歸建立冬小麥葉綠素含量預(yù)測模型，為后續(xù)施肥決策提供支持。孫紅等[6]按照葉片垂直分布位置采集馬鈴薯葉片的成像高光譜數(shù)據(jù)，選取隨機蛙跳優(yōu)化偏最小二乘(Random frog-partial least squares，RF-PLS)模型計算馬鈴薯葉片每個像素點的葉綠素含量，實現(xiàn)了不同葉位馬鈴薯葉片葉綠素含量無損檢測以及分布可視化表達(dá)。裴浩杰等[7]利用無人機低空遙感平臺，通過構(gòu)建綜合長勢指標(biāo)(Comprehensive growth index, CGI)建立小麥關(guān)鍵生育期葉片葉綠素含量反演模型，模型決定系數(shù)均大于0.7，可為小麥葉綠素含量長勢監(jiān)測提供參考。

現(xiàn)有的作物葉綠素含量反演研究中，多是通過構(gòu)建植被指數(shù)建立統(tǒng)計回歸模型來反演葉綠素含量，這對于特定地區(qū)的特定品種，會有比較理想的反演效果，但在模型的普適性等方面還存在一定的不足[8-10]。利用高光譜分析技術(shù)檢測葉綠素含量方面，前人的工作主要集中在2個方面：建立各種植被指數(shù)，利用多元線性或非線性回歸方法建立指數(shù)與葉綠素含量之間的反演模型；或者先對獲取的水稻冠層高光譜數(shù)據(jù)的全部波段進(jìn)行PCA、PLS等方法建模[11]。其中，第1種植被指數(shù)方法物理意義清晰，模型構(gòu)建簡單，但目前所構(gòu)建的光譜指數(shù)數(shù)目繁多，單一植被指數(shù)不能很好地表征整個高光譜信息；第2種方法雖然考慮了全部高光譜信息，但每次反演均需要全部波段作為模型輸入，計算量較大，模型應(yīng)用復(fù)雜[12]。

本研究針對東北粳稻冠層葉綠素含量反演開展研究。采用無人機高光譜遙感平臺，獲取粳稻關(guān)鍵生育期冠層高光譜影像數(shù)據(jù)，采用連續(xù)投影算法(Successive projection algorithm，SPA)選擇粳稻葉片葉綠素含量響應(yīng)的敏感波段，建立基于粒子群優(yōu)化的極限學(xué)習(xí)機反演模型，以期實現(xiàn)快速、精準(zhǔn)、無損診斷粳稻葉綠素含量，提高反演精度。

1 材料與方法

1.1 試驗區(qū)域與試驗設(shè)計

為覆蓋粳稻關(guān)鍵生育期如分蘗期、拔節(jié)期、抽穗期、灌漿期等，試驗于2016—2017年的6—9月，每周周三在沈陽農(nóng)業(yè)大學(xué)遼中水稻實驗站(41°47′N、122°71′E)進(jìn)行，水稻品種為‘沈稻529’。為避免陰雨多云天氣對遙感數(shù)據(jù)采集造成誤差，數(shù)據(jù)采集過程選擇天氣晴好的日期，如遇云量超過20%或不利于遙感數(shù)據(jù)采集天氣則順延采集。粳稻試驗小區(qū)設(shè)計為4個氮肥梯度處理，分別為 CK、N1、N2、N3，每種處理 3 次重復(fù)，共設(shè)立96個面積相同的長方形試驗小區(qū)。其中，CK為對照組，即不施用氮肥；N1為當(dāng)?shù)貥?biāo)準(zhǔn)氮肥施用水平，氮肥施用量為45.0 kg/hm2；N2為低氮施肥水平，施用量為22.5 kg/hm2；N3為高氮施肥水平，施用量為67.5 kg/hm2。磷肥和鉀肥的施用按照當(dāng)?shù)貥?biāo)準(zhǔn)施用量進(jìn)行，其中，磷肥標(biāo)準(zhǔn)施用量為51.75 kg/hm2，鉀肥標(biāo)準(zhǔn)施用量為18 kg/hm2。圖1為本研究試驗區(qū)域與試驗小區(qū)設(shè)計。

圖1 本研究試驗區(qū)域和試驗小區(qū)設(shè)計Fig. 1 Design of experimental district and plot in this study

1.2 無人機高光譜影像數(shù)據(jù)獲取

無人機高光譜平臺采用深圳大疆創(chuàng)新公司的M600 PRO六旋翼無人機，高光譜成像儀選用四川雙利合譜公司的GaiaSky-mini內(nèi)置推掃式機載高光譜成像系統(tǒng)，高光譜的波段范圍為400～1 000 nm，分辨率為3 nm，有效波段數(shù)為253個，高光譜成像系統(tǒng)如圖2所示。

無人機高光譜遙感平臺數(shù)據(jù)采集時間為每次試驗的上午11:00—12:00，選擇太陽光強較為穩(wěn)定的時段，無人機飛行高度為150 m。運用ENVI5.3+IDL工具軟件對獲取的高光譜遙感影像進(jìn)行小區(qū)高光譜數(shù)據(jù)提取，運用波譜角填圖方法去除干擾地物光譜的影響，對每個小區(qū)的感興趣區(qū)計算平均光譜，作為每個試驗小區(qū)的高光譜信息。

1.3 粳稻田間數(shù)據(jù)獲取與葉綠素含量測定

圖2 無人機高光譜成像系統(tǒng)Fig. 2 UAV hyperspectral imaging system

粳稻田間冠層數(shù)據(jù)獲取的關(guān)鍵生育期主要包括分蘗期、拔節(jié)期、孕穗期、抽穗期、開花期和灌漿期6個時期。粳稻田間取樣在試驗小區(qū)的中間長勢均勻區(qū)域進(jìn)行，每個小區(qū)采集4穴粳稻作為該小區(qū)的樣本點，取冠層葉片樣本50片放入密封袋，并迅速放入移動冰室進(jìn)行冷藏保鮮，以確保后期試驗數(shù)據(jù)的準(zhǔn)確性。同時利用手持差分GPS設(shè)備進(jìn)行采樣點地理坐標(biāo)位置采集，以便與無人機采集的高光譜影像數(shù)據(jù)一一對應(yīng)[7,13-14]。

待粳稻樣品送回實驗室后，采用浸提法提取樣品葉綠素含量。將丙酮、無水乙醇、蒸餾水3種溶液按9∶9∶2的體積比配置成混合溶液，選擇粳稻的完全展開葉片，剪碎后混合均勻，稱取0.4 g加入200 mL混合溶液；在實驗室遮光環(huán)境下靜置，直到樣品完全發(fā)白，利用分光光度計進(jìn)行比色，分別測定D663nm和D645nm，按公式(1)計算粳稻樣品的葉綠素含量，ρ(葉綠素)[15-16]。

1.4 數(shù)據(jù)處理

本研究共進(jìn)行了10次試驗數(shù)據(jù)采集，共采集780組有效樣本數(shù)據(jù)，選取其中的90%作為建模數(shù)據(jù)集(樣本容量為702)，另外10%為驗證數(shù)據(jù)集(樣本容量為78)。本研究采樣數(shù)據(jù)的統(tǒng)計特征如表1所示。

表1 試驗小區(qū)粳稻葉綠素含量統(tǒng)計特征Table 1 Statistical characteristics of chlorophyll content in japonica rice in experimental plot

由表1可知，本研究的建模數(shù)據(jù)集和驗證數(shù)據(jù)集2組數(shù)據(jù)除樣本量不同之外，其他統(tǒng)計參量相差不大，變異系數(shù)大于40%，說明葉綠素含量數(shù)據(jù)的離散程度較大。

如圖3所示，葉綠素吸收光的能力很強，有2個最強烈的吸收區(qū)，一個是波長為640～660 nm的紅光部分，另一個是430～450 nm的藍(lán)紫光部分。此外，在光譜的橙光、黃光和綠光部分有不明顯的吸收帶，所以本研究選擇可見光(380～780 nm)之間的光譜信息作為高光譜區(qū)間[8,17]。

1.5 高光譜特征波段提取

本研究利用Savitzky-Golay卷積平滑算法[18]對試驗小區(qū)獲取的400～800 nm的反射率光譜數(shù)據(jù)進(jìn)行平滑處理，降低噪聲對葉綠素含量反演精度的影響，在平滑處理的基礎(chǔ)上，利用雙變量相關(guān)性分析方法對葉綠素含量與高光譜波段反射率之間分別進(jìn)行相關(guān)性分析研究，同時選用皮爾森相關(guān)系數(shù)、距離相關(guān)系數(shù)、最大信息系數(shù)等3種方法進(jìn)行單波段相關(guān)性分析。

為了能夠有效的提取粳稻冠層葉綠素含量高光譜反演特征波段，本研究采用連續(xù)投影算法(SPA)提取葉綠素含量反演的特征波段。SPA是一種矢量空間共線性最小化的前向變量選擇算法，SPA的流程由3個階段構(gòu)成[19]：

第1階段，篩選出共線性最小的若干組備選波長變量子集。假設(shè)初始變量位置k(0)以及變量數(shù)N已經(jīng)給出，該階段具體步驟如下。

步驟1：第1次迭代之前(n=1)，將訓(xùn)練集光譜矩陣Xcal的第j列賦值給Xj，j=1，…，J；

步驟3：對所有的j，計算Xj在和Xk(n-1)正交的子空間上的投影，

式中：P為投影算子，

步驟 4：令xj=Pxj，j；

步驟5：令n=n+1。如果n＜N則返回步驟2結(jié)束：得到N個備選波長的位置：{k(n)；n=0，…，N－1}。

選擇過程中進(jìn)行的投影操作次數(shù)為(N－1)(J－N/2)。

第2階段，分別使用各子集中的變量建立多元線性回歸(Multiple linear regression，MLR)模型，選出均方根誤差(Root mean square error，RMSE)最小的子集。

第3階段，對第2階段選出的子集進(jìn)行逐步回歸建模，在盡量不損失預(yù)測準(zhǔn)確度的前提下得到一個變量數(shù)目較少的集合，該集合中的波長變量即是所選有效波長[20]。

1.6 葉綠素含量PSO-ELM反演建模

本研究采用基于粒子群優(yōu)化的極限學(xué)習(xí)機模型(Particle swarm optimization-extreme learning machine, PSO-ELM)反演粳稻冠層葉綠素含量。利用粒子群優(yōu)化算法優(yōu)化選擇極限學(xué)習(xí)機的輸入層權(quán)值和隱含層偏差，從而計算出輸出權(quán)值矩陣。在本研究中考慮的PSO算法的參數(shù)主要包括種群規(guī)模(p)、慣性權(quán)重(w)、學(xué)習(xí)因子C1和C2(C1、C2分別調(diào)節(jié)微粒向個體最優(yōu)或者群體最優(yōu)方向飛行的最大步長，決定微粒個體經(jīng)驗和群體經(jīng)驗對微粒自身運行軌跡的影響，C1和C2的搭配不同，將會影響到PSO算法的性能、速度位置相關(guān)系數(shù)(m)。通過與傳統(tǒng)的極限學(xué)習(xí)機反演算法比較，驗證了基于粒子群優(yōu)化的極限學(xué)習(xí)機在隱含層節(jié)點數(shù)目以及網(wǎng)絡(luò)泛化性上的優(yōu)勢[21]。

本研究使用RMSE和模型決定系數(shù)(R2)作為無人機高光譜遙感反演粳稻冠層葉綠素含量精度的評價標(biāo)準(zhǔn)。

2 結(jié)果與分析

2.1 粳稻高光譜特征波段提取結(jié)果

利用Savitzky-Golay卷積平滑算法對試驗小區(qū)獲取的400～800 nm的反射率光譜數(shù)據(jù)進(jìn)行平滑處理，處理后粳稻冠層高光譜曲線如圖4所示。

圖4 高光譜信息平滑降噪處理Fig. 4 Noise reduction of hyperspectral information smoothing

在平滑處理的基礎(chǔ)上，利用雙變量相關(guān)性分析方法對葉綠素含量與高光譜波段反射率之間分別進(jìn)行相關(guān)性分析研究，本研究選用皮爾森相關(guān)系數(shù)、距離相關(guān)系數(shù)、最大信息系數(shù)(Maximum information coefficient, MIC)等3種方法進(jìn)行單波段相關(guān)性分析(圖5)。

圖5 不同方法的單波段相關(guān)性分析性結(jié)果Fig. 5 Single band correlation analysis results of different methods

400～800 nm之間任意單一波段與粳稻冠層葉綠素含量相關(guān)性分析結(jié)果表明，基于MIC的相關(guān)性要低于基于距離相關(guān)性系數(shù)法。而利用皮爾森線性相關(guān)分析的結(jié)果在400～480 nm之間與距離相關(guān)系數(shù)的結(jié)果較為一致。而3種相關(guān)性分析方法在701～705 nm處的反射率與葉綠素含量的相關(guān)性要明顯高于其他波段與葉綠素含量之間的相關(guān)性。其中702 nm波段反射率與葉綠素含量之間的相關(guān)性最高，3種方法的相關(guān)系數(shù)分別為0.66、0.64和0.53。由圖5可知，400～800 nm波段反射率與葉綠素含量之間的相關(guān)系數(shù)介于0.50～0.65之間，在區(qū)間內(nèi)不同波段之間也存在著很強的相關(guān)關(guān)系，難以提取能夠有效建立葉綠素含量反演模型的高光譜特征波段[22]。

本研究采用SPA提取葉綠素含量反演的特征波段，結(jié)果如圖6所示，綜合分析來看，隨著特征波段提取數(shù)目的增加，交叉驗證根均方誤差 (Root mean square error of cross validation, RMSECV)現(xiàn)出逐漸變小的走勢，RMSECV值達(dá)到最低時，對應(yīng)最優(yōu)波段集合包含的波段數(shù)量，即為該時期最優(yōu)子集波段數(shù)量。根據(jù)校正集的內(nèi)部交叉驗證RMSECV值確定最佳的光譜波段數(shù)為5個。如圖7所示，最佳特征波段為 410、481、533、702和 798 nm。

圖7 基于連續(xù)投影算法(SPA)的粳稻冠層高光譜特征波段選擇Fig. 7 Selection of characteristic hyperspectral bands for japonica rice canopy based on successive projection algorithm (SPA)

2.2 葉綠素含量的反演模型

本研究通過SPA提取的粳稻冠層光譜特征波段為自變量，現(xiàn)場采樣的粳稻冠層葉片葉綠素含量為因變量，利用ELM、PSO-ELM分別建立粳稻葉綠素含量估算模型，粒子群迭代次數(shù)為100。由于PSO算法的參數(shù)設(shè)置較多，本研究設(shè)計了PSO參數(shù)正交試驗，結(jié)果如表2所示。將表2中各參數(shù)組合代入PSO-ELM算法(得到葉綠素含量預(yù)測值，并計算得出模型決定系數(shù)(R2)[22-24]。

表2 ELM和PSO-ELM建立粳稻葉綠素含量反演模型Table 2 Inversion models of chlorophyll contents in japonica rice established by ELM and PSO-ELM

表3為PSO-ELM參數(shù)正交試驗分析結(jié)果，其中Wij為參數(shù)j的i水平對應(yīng)的葉綠素含量反演模型的模型決定系數(shù)R2之和，以W11和W21為例，計算公式如下：

表3 PSO-ELM參數(shù)正交試驗分析結(jié)果1)Table 3 Analysis results of orthogonal test of PSO-ELM parameters

其余PSO參量的計算方法同理可得。

由表 3 結(jié)果可知，W11、W21、W31、W41和W51之間的差別僅由PSO種群規(guī)模的變化而印象最終模型反演精度，與PSO算法的其他參數(shù)無關(guān)。可通過對比W11、W21、W31、W41和W51的大小來選取最優(yōu)的PSO算法種群規(guī)模，其中R2最大值所對應(yīng)的種群規(guī)模即為本研究所確定的PSO-ELM算法中的種群規(guī)模。因此，本研究PSO-ELM反演粳稻葉綠素含量模型的最優(yōu)種群規(guī)模為70，同理可得最優(yōu)慣性權(quán)重w為0.9～0.3線性遞減，學(xué)習(xí)因子C1=2.80、C2=1.10，速度位置相關(guān)系數(shù)m為0.60。利用該組參數(shù)建立的PSO-ELM粳稻葉綠素含量反演模型，模型R2為0.887，RMSE為0.783，要明顯優(yōu)于ELM算法建立的葉綠素含量反演模型(R2=0.734、RMSE=1.824)，反演效果對比如圖8所示。

圖8 ELM和PSO-ELM模型的粳稻葉綠素含量反演模型檢驗結(jié)果Fig. 8 Test results of the inversion models for the chlorophyll contents in japonica rice using ELM and PSO-ELM

由反演結(jié)果可知，本研究采用的PSO優(yōu)化的極限學(xué)習(xí)機粳稻葉綠素含量反演模型的預(yù)測效果，要明顯優(yōu)于單純極限學(xué)習(xí)方法所建立的反演模型，展現(xiàn)了良好的預(yù)測能力。

3 討論與結(jié)論

本研究嘗試通過SPA方法提取無人機獲取的粳稻冠層高光譜信息的特征波段，利用PSO優(yōu)化的極限學(xué)習(xí)機算法建立粳稻葉綠素含量反演模型，以期能夠?qū)崿F(xiàn)對東北粳稻冠層葉綠素含量無人機遙感反演的快速監(jiān)測和評估。傳統(tǒng)的利用高光譜特征波段進(jìn)行一元或者多元回歸統(tǒng)計模型可能無法充分表達(dá)光譜信息與葉綠素含量之前的非線性數(shù)學(xué)關(guān)系，并在一定程度上制約了葉綠素含量反演精度[25-26]?；诖?，本研究采用PSO優(yōu)化的極限學(xué)習(xí)機方法，建立PSO-ELM葉綠素含量反演模型，選取了5個特征波段作為模型輸入，建立粳稻冠層葉綠素含量反演模型。

所構(gòu)建的PSO-ELM模型在反演精度上要高于單純ELM反演精度，這可能是因為PSO算法優(yōu)化了極限學(xué)習(xí)機的輸入權(quán)值和隱層的閾值，避免了ELM在隱層節(jié)點數(shù)較少的情況下由于隨機給定輸入權(quán)值和隱層閾值帶來的泛化能力差、標(biāo)定精度不高等缺點。本研究所選擇的特征波段與前人的結(jié)論在特征波段范圍區(qū)間上基本是相同的，但具體的特征波段并不相同[27-32]。造成這樣的原因主要是由于如下幾個因素：1)由于試驗品種不同造成的特征波段存在一定程度的差別。2)由于稻田內(nèi)存在的例如水體、雜草、土壤等干擾物，在粳稻純凈高光譜信息提取過程中有可能會參雜上述光譜信息，也會造成特征波段提取的不同。

本研究以沈陽農(nóng)業(yè)大學(xué)遼中水稻實驗站的‘沈稻529’為研究對象，利用2016—2017年的無人機高光譜粳稻冠層遙感數(shù)據(jù)，并結(jié)合響應(yīng)的地面粳稻葉綠素含量的實地測量，采用連續(xù)投影法提取能夠表征粳稻冠層高光譜信息的特征波段，在此基礎(chǔ)上通過PSO-ELM算法建立粳稻冠層葉綠素含量反演模型，得到以下結(jié)論：1)通過SPA算法提取了能夠表征400～800 nm范圍的高光譜特征波段，本研究所選取的特征波段分別為410、481、533、702和798 nm。2)通過正交試驗確定PSO算法優(yōu)化ELM模型的參數(shù)最優(yōu)值，其中最優(yōu)種群規(guī)模為80，同理可得最優(yōu)慣性權(quán)重為0.9～0.3線性遞減，學(xué)習(xí)因子1和2分別為2.80和1.10，速度位置相關(guān)系數(shù)為0.60。3)PSO-ELM粳稻葉綠素含量反演模型的R2=0.887、RMSE=0.783，明顯優(yōu)于單純極限學(xué)習(xí)機算法建立的粳稻葉綠素含量反演模型(R2=0.734、RMSE=1.824)。該模型展現(xiàn)了良好的預(yù)測能力，能夠為粳稻葉綠素含量反演研究及養(yǎng)分診斷提供一定的數(shù)據(jù)支撐和模型參考。

本研究采用無人機高光譜遙感平臺，存在一定的采集誤差，同時由于無人機平臺的限制地面樣本數(shù)量還相對有限，并且所建立的葉綠素含量反演模型僅是針對試驗粳稻品種建立的，對于其他品種葉綠素含量反演的適用性還需要進(jìn)一步完善。因此在未來的研究中將增加試驗品種，同時針對粳稻不同生育期分別建立葉綠素含量反演模型，提高模型精度和普適性。