袁煒楠，許童羽,2*，曹英麗,2，王洋，于豐華,2

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

植物通過(guò)光合作用進(jìn)行物質(zhì)和能量交換，葉片中的葉綠素作為光合作用的主要色素，是光合作用能力、葉片氮含量和作物發(fā)育階段的指示物，因此，獲取葉綠素的變化信息已成為檢測(cè)作物長(zhǎng)勢(shì)的重要手段[1]。然而，傳統(tǒng)的化學(xué)分析方法具有破壞性，且耗費(fèi)人力，時(shí)間長(zhǎng)，成本高。高光譜遙感技術(shù)是一種便捷、低成本的監(jiān)測(cè)手段，它使窄波段光譜信息的獲取成為可能。近年來(lái)發(fā)展了一系列窄波段組合運(yùn)用于葉片色素含量的估算[2]，為作物葉綠素含量的定量診斷提供了簡(jiǎn)便有效、非破壞性的數(shù)據(jù)采集和處理方法[3]。

利用高光譜技術(shù)估算葉片中的色素含量主要有2類方法模型：物理模型[4]和經(jīng)驗(yàn)?zāi)Ｐ蚚5]。近年來(lái)，國(guó)內(nèi)外學(xué)者在利用高光譜遙感數(shù)據(jù)監(jiān)測(cè)植被葉綠素含量方面做了大量的研究。HUNT等[6]和李粉玲等[7]利用高光譜特征波段組合構(gòu)建了估算作物冠層葉綠素含量的遙感模型；李媛媛等[8]利用葉片高光譜特征波段組合構(gòu)建了估算作物葉片葉綠素含量的遙感模型；丁永軍等[9]通過(guò)把原始光譜數(shù)據(jù)進(jìn)行一階微分和去除包絡(luò)線等處理后，結(jié)合波段自相關(guān)分析建立了多元線性回歸模型，并對(duì)番茄葉綠素含量進(jìn)行反演；岳學(xué)軍等[10]在主成分分析降維的基礎(chǔ)上利用支持向量機(jī)回歸算法和偏最小二乘回歸模型預(yù)測(cè)了柑橘葉片葉綠素含量；鄧小蕾等[11]應(yīng)用波段間相關(guān)系數(shù)分布對(duì)光譜數(shù)據(jù)進(jìn)行降維，建立了蘋果葉片葉綠素含量的逐步回歸模型。周冬琴等[12]選出了與水稻氮含量相關(guān)性較好的5個(gè)波段，建立回歸方程對(duì)氮素進(jìn)行反演。CHEN等[13]利用綠峰、紅邊位置、高光譜特征波段組合作為反向傳播（BP）神經(jīng)網(wǎng)絡(luò)的輸入變量對(duì)水稻葉綠素含量進(jìn)行了估算，并與線性回歸模型進(jìn)行比較分析。目前，降維方法有主成分分析、因素分析、典型相關(guān)性分析，以及偏最小二乘回歸、分段逆回歸等，這些方法的共同特點(diǎn)是依照某種最優(yōu)化原則，在原變量X1，X2，…，Xp中提取成分F1，F(xiàn)2，…，F(xiàn)m（m＜p），再利用這些成分進(jìn)行相關(guān)的分析工作。但是，由于每一個(gè)成分都是原始變量的線性組合，所以這些方法都沒(méi)有篩選變量的功能[14]。還有一些學(xué)者只選取幾個(gè)波段建立植被指數(shù)進(jìn)行降維，但這種方法只選擇了波段中的很小一部分，造成數(shù)據(jù)的大量舍棄，而且對(duì)原數(shù)據(jù)進(jìn)行了數(shù)學(xué)變形，丟失了原數(shù)據(jù)的信息。

本研究提出了一種基于主基底分析的降維方法，通過(guò)對(duì)降維后的光譜與對(duì)應(yīng)的葉綠素?cái)?shù)據(jù)建立回歸模型來(lái)進(jìn)行葉綠素含量的估算。本文選取對(duì)葉綠素敏感的400～1 000 nm[15]波段進(jìn)行Gram_Schmidt變換并找到投影空間，構(gòu)造集中波段信息的主基底來(lái)進(jìn)行建模，并與選取的植被指數(shù)降維法進(jìn)行對(duì)比分析，以期為光譜數(shù)據(jù)降維估算葉綠素含量提供一種行之有效的方法。

1 材料與方法

1.1 試驗(yàn)設(shè)計(jì)與樣品采集

試驗(yàn)于2017年6—9月在遼寧省沈陽(yáng)市沈陽(yáng)農(nóng)業(yè)大學(xué)遼中區(qū)卡力瑪村院士工作站進(jìn)行。供試水稻品種為千重浪1號(hào)，試驗(yàn)時(shí)間包含了水稻生長(zhǎng)的分蘗期、拔節(jié)孕穗期、抽穗灌漿期和成熟期。工作站近幾年種植作物均為水稻。試驗(yàn)田中共設(shè)4個(gè)施肥梯度，分別為無(wú)氮（N0）、低氮（N1）、正常氮（N2）和高氮（N3）；小區(qū)面積6 500 m2，每個(gè)水平3次重復(fù)，共12個(gè)小區(qū)，隨機(jī)分布，如圖1所示。每星期無(wú)人機(jī)飛行采集1次，共進(jìn)行13次試驗(yàn)。每次從每個(gè)小區(qū)選取3穴具有代表性的水稻作為樣本，每穴水稻隨機(jī)采集2片葉，裝入密封袋中并標(biāo)注小區(qū)的名稱及編號(hào)，采集后立即運(yùn)回實(shí)驗(yàn)室。從得到的936組樣本中隨機(jī)抽取736組作為建模樣本，剩余200組作為檢驗(yàn)樣本。

1.2 冠層光譜數(shù)據(jù)采集與葉綠素含量測(cè)定

圖1 水稻小區(qū)分布圖Fig.1 Distribution map of riceplot

圖2 無(wú)人機(jī)平臺(tái)Fig.2 Unmanned aerial vehicle(UAV)platform

本試驗(yàn)采用飛行高度50 m的中國(guó)GaiaSkymini無(wú)人機(jī)高光譜成像系統(tǒng)進(jìn)行水稻冠層光譜反射率的采集，如圖2所示。為了獲取穩(wěn)定的水稻冠層高光譜反射率數(shù)據(jù)，選擇天氣晴朗、少風(fēng)或微風(fēng)的天氣，在10：00—14：00之間飛行采集；該光譜儀的光譜范圍為400～1 000 nm，輸出間隔為1 nm。從各次試驗(yàn)的高光譜圖像中，每個(gè)小區(qū)隨機(jī)選取6個(gè)點(diǎn)作為樣本點(diǎn)，每個(gè)樣本點(diǎn)采集1條光譜曲線。

葉片的葉綠素含量采用分光光度計(jì)進(jìn)行測(cè)定。葉片用蒸餾水洗凈后吹干，去掉中脈后剪碎，混勻；稱取剪碎的新鮮葉片0.2 g，放入研缽中，加入少量石英砂和碳酸鈣粉及2～3 mL 80%丙酮，研成勻漿；再加入10 mL丙酮，繼續(xù)研磨至組織變白，靜置3～5 min。在漏斗中放置濾紙，用丙酮潤(rùn)濕，經(jīng)過(guò)濾和反復(fù)沖洗研缽和研缽棒以確保葉片色素全部進(jìn)入25 mL棕色容量瓶中；用吸管吸取丙酮，將濾紙上的葉綠素也沖洗至容量瓶中，直至濾紙上無(wú)綠色為止；最后用丙酮定容至25 mL，搖勻，于避光處?kù)o置2 h。取萃取液裝入光徑1 cm的比色杯內(nèi)，以80%丙酮為空白對(duì)照，測(cè)定在波長(zhǎng)為649、665 nm處的吸光度，利用式（1）、（2）和（3）計(jì)算葉片中葉綠素的質(zhì)量濃度。

式中：D（665 nm），D（649 nm）分別為在665 nm與649 nm波長(zhǎng)下的吸光度；Ca，Cb分別為葉綠素a與葉綠素b的質(zhì)量濃度，mg/L；C為總?cè)~綠素的質(zhì)量濃度，mg/L。

1.3 Gram_Schmidt變換原理

任意一組線性無(wú)關(guān)的變量X1，X2，…，Xs，總可以經(jīng)過(guò)Gram_Schmidt變換把它們變換成正交變量Z1，Z2，…，Zs，Gram_Schmidt變換公式見(jiàn)式（4）。

JAIN等[16]對(duì)式（4）進(jìn)行了詳細(xì)證明，且通過(guò)Gram_Schmidt變換可以得到如下推論：對(duì)于任意一組 秩為 s(s≤p)的變 量 X1，X2，…，Xp，對(duì) 它 們 做Gram_Schmidt變換后，得到 Z1，Z2，…，Zp；這其中必有Z1，Z2，…，Zs是相互直交的，且Zs+1=…=Zp=0。從推論可以看出，Gram_Schmidt變換有2個(gè)功能：1）將變量集合中的信息進(jìn)行正交分解；2）排除X1，X2，…，Xp中的冗余變量（即被變換成0的那些變量）。

1.4 高光譜降維方法

影響葉綠素的光譜波段范圍為400～1 000 nm，各波段間存在一定的相關(guān)性，如果不做降維處理，在對(duì)光譜數(shù)據(jù)進(jìn)行葉綠素反演時(shí)，需要的數(shù)據(jù)量將呈指數(shù)冪增加，給數(shù)據(jù)處理帶來(lái)困難。因此，本研究采用基于主基底分析的降維方法，其原理是通過(guò)Gram_Schmidt變換找到波段的投影空間，在投影空間下構(gòu)造出最大限度涵蓋原始變量信息的低維變量，該低維變量即為波段的主基底，利用建立的主基底變量代替原來(lái)的光譜波段數(shù)據(jù)。

投影空間的構(gòu)造原理如下：設(shè)研究的每組對(duì)象涉及p個(gè)變量，共有n組對(duì)象，則構(gòu)成一個(gè)n×p的數(shù)據(jù)矩陣，X=[x1x2…xP]。

1）選擇與健康水稻葉片最接近的一列數(shù)據(jù)xm作為h1，計(jì)算其內(nèi)積hT1h1，將其值作為初始能量，則

2）將xm從矩陣XT中去除，選擇剩余數(shù)據(jù)的第一列，與z1作Gram_Schmidt變換，得到h2，計(jì)算其內(nèi)積，則將h2作為投影空間的第二，則選擇剩余數(shù)據(jù)的下一列。

3）選擇zm列時(shí)，與z1，z2，…，zm－1進(jìn)行Gram_Schmidt變換得到hm，計(jì)算其內(nèi)積hTmhm；如果則將hm作為投影空間的第m列zm。

4）重復(fù)以上過(guò)程，直至經(jīng)Gram_Schmidt變換后得到n個(gè)相互正交的向量z1，z2，…，zn，即為投影空間。

將光譜數(shù)據(jù)與投影空間相乘，得到所求的波段主基底；具體的技術(shù)流程如圖3所示。

圖3 技術(shù)流程圖Fig.3 Techniqueflow chart

2 結(jié)果與分析

2.1 不同施氮量條件下水稻葉片高光譜特征

圖4是千重浪1號(hào)水稻在不同施氮量下的光譜特征曲線。從中可知，在不同氮素水平下，水稻光譜反射率的變化明顯，說(shuō)明可以進(jìn)行葉綠素的反演。

圖4 不同氮素水平下水稻葉片高光譜曲線Fig.4 Hyperspectral curves of rice leaves under different nitrogen levels

2.2 不同葉綠素含量的水稻葉片高光譜特征

如圖5所示：不同葉綠素含量的水稻葉片高光譜特征總體上保持一致，在400～500 nm的藍(lán)紫光波段和650～700 nm的紅光波段處形成2個(gè)吸收谷，在500～600 nm波段處形成1個(gè)反射峰，在650～750 nm的近紅外波段處光譜反射率開(kāi)始急劇上升，在750 nm波段處達(dá)到最大反射率值并趨于平穩(wěn)；同時(shí)，不同葉綠素含量會(huì)造成光譜反射率曲線的變換，即光譜曲線會(huì)隨著葉綠素含量的增大向下移動(dòng)。由此可以看出，水稻葉片光譜曲線特征與葉綠素含量有著密不可分的關(guān)系。因此，可以將光譜反射率與葉綠素含量建立相關(guān)性來(lái)估算葉綠素的含量。

圖5 不同葉綠素含量的水稻葉片高光譜曲線Fig.5 Hyperspectral curves of rice leaves with different chlorophyll(Chl)contents

2.3 基于主基底分析降維與模型的估算

選取736組光譜數(shù)據(jù)進(jìn)行主基底分析降維。在主基底分析的基礎(chǔ)上，建立波段長(zhǎng)度為5、10、15、20、25、30、35和40 nm的窗口，對(duì)波段進(jìn)行分組并找到主基底，建立主基底與葉綠素的線性函數(shù)、二次函數(shù)、指數(shù)、對(duì)數(shù)、冪函數(shù)、最小二乘回歸模型，依據(jù)決定系數(shù)（R2）、均方根誤差（root mean square error,RMSE）與主基底維數(shù)3項(xiàng)指標(biāo)，選出窗口大小正好和精度最高的回歸模型。在模型精度上，最小二乘回歸模型精度明顯高于其他回歸模型，基于最小二乘回歸模型的窗口指標(biāo)情況如表1所示。從中可以看出：隨著窗口大小的增加，R2逐漸減小；從窗口大小為30 nm開(kāi)始，R2維持在0.68左右，RMSE達(dá)到最小值，主基底在48維處趨于穩(wěn)定。因此，選擇分段長(zhǎng)度為30 nm的窗口最佳。

表1 波段窗口指標(biāo)Table1 Band window indicator

對(duì)736組400～1 000 nm波段的光譜數(shù)據(jù)以30 nm窗口分段，分成20段，應(yīng)用MATLAB R2016a編程，逐次選擇每段進(jìn)行Gram_Schmidt變換，共找到48列數(shù)據(jù)組成的投影空間。將每段光譜數(shù)據(jù)與對(duì)應(yīng)的主基底相乘，得到光譜數(shù)據(jù)的主基底。736組400～1 000 nm波段的光譜數(shù)據(jù)主基底構(gòu)成736×48的數(shù)據(jù)矩陣，光譜數(shù)據(jù)從601維降到48維。

應(yīng)用最小二乘回歸法建立葉綠素估算模型，Y=BX+A。其中：Y為葉綠素值；X為波段主基底；B為48個(gè)數(shù)據(jù)組成的數(shù)組；A為常量。因此，該模型需要48個(gè)輸入變量，主基底中的48列數(shù)據(jù)即為輸入變量。應(yīng)用剩余的200組光譜數(shù)據(jù)進(jìn)行驗(yàn)證，把波段進(jìn)行30窗口的分段，每段光譜組成的200×30數(shù)據(jù)矩陣與各組投影空間相乘，得到測(cè)試數(shù)據(jù)200×48矩陣的主基底，每行作為最小二乘回歸模型的輸入，得到200組預(yù)測(cè)葉綠素值，并與200組真實(shí)葉綠素值進(jìn)行比較，結(jié)果如圖6所示。從中可以看出，葉綠素實(shí)測(cè)值與預(yù)測(cè)值的RMSE為1.20。

2.4 植被指數(shù)與水稻葉綠素的相關(guān)性

圖6 葉綠素實(shí)測(cè)值與預(yù)測(cè)值比較Fig.6 Comparison of measured and predicted leaf chlorophyll values

依據(jù)表2中的公式計(jì)算植被指數(shù)，分析植被指數(shù)與葉綠素的相關(guān)性。從中可以看出：3種植被指數(shù)與葉綠素值均達(dá)到極顯著相關(guān)水平（P＜0.01），光化學(xué)反射指數(shù)（photochemical reflectance index,PRI）與葉綠素值呈正相關(guān)，即葉綠素值隨著植被指數(shù)的增加而增加；差值光譜指數(shù)（RD2）、葉綠素吸收比指數(shù)（modified chlorophyll absorption ratio index,MCARI）與葉綠素值呈負(fù)相關(guān)，即葉綠素值隨著植被指數(shù)的增加而降低。其中：PRI與葉綠素值相關(guān)性最高，達(dá)到0.682；MCARI與葉綠素值的相關(guān)性最低，為－0.605。根據(jù)統(tǒng)計(jì)學(xué)分析，以上3種植被指數(shù)與葉綠素值均達(dá)到中度相關(guān)，說(shuō)明這些植被指數(shù)能夠很好地反映葉綠素值的變化，能夠用來(lái)搭建估算模型。

表2 常用光譜指數(shù)Table 2 Common spectral index

2.5 基于植被指數(shù)降維的葉綠素估算模型及其檢驗(yàn)

根據(jù)植被指數(shù)構(gòu)建的線性函數(shù)、二次函數(shù)、指數(shù)、對(duì)數(shù)、冪函數(shù)、最小二乘回歸模型精度的不同，以決定系數(shù)（R2）、均方根誤差（RMSE）和相對(duì)誤差（relative error,RE）為指標(biāo)，選擇出所選植被指數(shù)擬合度較高的回歸模型作為估算模型，結(jié)果如圖7所示。擬合方程中x為PRI值，y為實(shí)測(cè)的葉綠素值。基于3種植被指數(shù)PRI、RD2和MCARI建立的回歸模型均能較好地估算葉綠素值，其中：R2分別為0.538、0.468和0.452，RMSE分別為3.295、3.761和3.872，RE為4.608%、5.407%和5.504%?；谥脖恢笖?shù)PRI所建立的指數(shù)模型精度最高，擬合方程為y=3.825e14.707x。圖8為200組檢測(cè)數(shù)據(jù)對(duì)估算模型的驗(yàn)證?？梢钥闯?，3種植被指數(shù)PRI、RD2和MCARI估算模型的RMSE分別為3.20、3.62和4.62，其中，基于植被指數(shù)PRI的估算效果最佳。

圖7 不同植被指數(shù)擬合的葉綠素含量估算模型Fig.7 Chlorophyll content estimation model fitted by different vegetation indices

3 討論

水稻冠層葉片的葉綠素含量與冠層的高光譜反射率有著密不可分的聯(lián)系，葉綠素含量越高，光譜反射率反而越低。近年來(lái)，將高光譜與葉綠素建立相關(guān)性進(jìn)而估算葉綠素的含量已經(jīng)成為研究的熱點(diǎn)。由于高光譜數(shù)據(jù)維數(shù)多，直接建立估算模型需要的數(shù)據(jù)量與維數(shù)之間呈冪指數(shù)關(guān)系，因此，對(duì)高光譜數(shù)據(jù)的降維是研究的重中之重。

圖8 不同植被指數(shù)擬合的葉綠素值驗(yàn)證模型Fig.8 Chlorophyll value verification model fitted by different vegetation indices

本研究表明，PRI、RD2和MCARI植被指數(shù)與水稻葉綠素的相關(guān)性分別為0.682、－0.645和－0.605，均達(dá)到極顯著相關(guān)水平。利用這3種植被指數(shù)建立的估算模型的決定系數(shù)（R2）、相對(duì)誤差（RE）和均方根誤差（RMSE）均達(dá)到建模精度要求，驗(yàn)證模型RMSE分別為3.20、3.62和4.62；其中基于植被指數(shù)PRI建立的估算模型反演葉綠素效果最好，擬合方程為y=3.825e14.707x，可以對(duì)水稻進(jìn)行葉綠素反演。

本文提出了一種基于主基底分析的降維方法，對(duì)降維后的光譜數(shù)據(jù)與葉綠素?cái)?shù)據(jù)建立最小二乘回歸模型，估算模型的R2為0.689，RMSE為2.20，檢驗(yàn)?zāi)Ｐ偷腞MSE為1.20，相比植被指數(shù)建立的模型精度有很大的提高，估算和驗(yàn)證模型的誤差有所降低，因而對(duì)葉綠素反演的準(zhǔn)確性更高。

本文提出的降維方法改善了傳統(tǒng)植被指數(shù)在選擇過(guò)程中丟失大量數(shù)據(jù)的弊端，在降維的同時(shí)還能夠保留光譜中起主要作用的部分，并更準(zhǔn)確地建立光譜數(shù)據(jù)與葉綠素?cái)?shù)據(jù)的估算模型，使得估算數(shù)據(jù)與實(shí)際數(shù)據(jù)偏差減小。

4 結(jié)論

本文利用高光譜數(shù)據(jù)對(duì)水稻冠層葉片的葉綠素含量進(jìn)行了估算，提出了一種基于主基底降維的方法，為水稻中葉綠素含量高精度估算提供了科學(xué)依據(jù)。該方法選取對(duì)葉綠素敏感的400～1 000 nm波段進(jìn)行Gram_Schmidt變換找到投影空間，在此投影空間下對(duì)波段進(jìn)行降維，建立估算植物葉片葉綠素值的最小二乘回歸模型，并與3種植被指數(shù)（PRI、RD2和MCARI）降維后建立的相同回歸模型進(jìn)行精度對(duì)比。結(jié)果表明：基于植被指數(shù)降維方法建立的模型建模精度和驗(yàn)?zāi)＞染陀诨谥骰追治鼋稻S法所建立的同種模型精度，說(shuō)明本研究采用的主基底分析方法對(duì)光譜數(shù)據(jù)進(jìn)行降維是可行的；同時(shí)，該方法還能更好地保留原始波段的信息，改善了傳統(tǒng)降維法的不足。因此，本研究建立的基于主基底降維的模型有較高的預(yù)測(cè)精度，對(duì)植物葉片中葉綠素含量的估算具有重要意義。