隋麗娜,房 建,郭立峰

(河北民族師范學(xué)院 數(shù)學(xué)與計(jì)算機(jī)科學(xué)學(xué)院,河北 承德 067000)

0 引言

水稻在我國廣泛種植,是我國主要的糧食作物,水稻的穩(wěn)產(chǎn)高產(chǎn)直接關(guān)系到國家糧食安全[1]。由于種植條件的特殊性,在東北等主要水稻產(chǎn)地,水稻主要為大規(guī)模的集中式種植[2]。該種植模式有助于管理,但當(dāng)某個生長時期出現(xiàn)病災(zāi)和蟲災(zāi)時,會對水稻產(chǎn)量產(chǎn)生嚴(yán)重影響[3]。因此,急需建立一套水稻全生長周期的產(chǎn)量預(yù)測模型,當(dāng)一個生長時期出現(xiàn)問題時,通過其他生長周期進(jìn)行調(diào)整,保證水稻穩(wěn)產(chǎn)。目前,水稻產(chǎn)量模型主要包括單株圖像采集法[4]、飛行器圖像分析[5]和衛(wèi)星光譜分析[6]。單株圖像分析對一株水稻不同生長時期進(jìn)行圖像處理,分析其結(jié)節(jié)數(shù)、植株高度和植株寬度等參數(shù),然后對其產(chǎn)量進(jìn)行建模分析。該方法建模精度高、針對性強(qiáng),但受限于樣本數(shù)量及樣品選擇,對于大范圍的水稻產(chǎn)量預(yù)測容易出現(xiàn)偏差[7]。飛行器圖像分析方法對水稻冠層圖像進(jìn)行分析,背景噪聲處理難度較大,且受到稻田環(huán)境因素影響,算法精度有待提高[8]。衛(wèi)星光譜分析采樣范圍過大,光譜信息龐雜,建模難度大,且獲得光譜數(shù)據(jù)成本高[9]。無人機(jī)技術(shù)的發(fā)展,使獲得低成本光譜數(shù)據(jù)成為可能。以無人機(jī)搭載光譜儀進(jìn)行40m高度均速光譜采集,通過分析水稻不同生長時期的光譜,得到水稻畝產(chǎn)預(yù)測模型。同時,進(jìn)行多組“留一法”實(shí)驗(yàn),結(jié)果表明:該方法具有良好的預(yù)測精度,性能穩(wěn)定,預(yù)測精度高。該模型數(shù)據(jù)獲得方便,預(yù)測精度較高,適于大規(guī)模一線推廣。

1 系統(tǒng)組成

系統(tǒng)包括無人機(jī)光譜系統(tǒng)、水稻產(chǎn)量模型系統(tǒng)、模型建立依據(jù)和產(chǎn)量模型建立4部分,如圖1所示。

圖1 系統(tǒng)結(jié)構(gòu)圖Fig.1 Structure of system。

無人機(jī)光譜系統(tǒng)保持無人40m高度勻速飛行,采用機(jī)載反射光譜傳感器檢測水稻冠層葉綠素光譜。模型建立依據(jù),葉綠素可以表征植物生長狀態(tài),冠層葉綠素含量可以預(yù)測水稻產(chǎn)量。通過無人機(jī)機(jī)載光譜儀檢測水稻冠層葉綠素光譜,由于葉綠素在紅邊區(qū)域、近紅外區(qū)域均有特征吸收波長,因此采用差值植被指數(shù)綜合考慮光譜包含信息。計(jì)算水稻全生命周期中的差值植被指數(shù)及其對應(yīng)的畝產(chǎn)量,對其進(jìn)行共線性分析,確定模型建立因子,并采用逐步線性擬合的方式建立模型。同時,采用“留一法”進(jìn)行畝產(chǎn)預(yù)測,對模型進(jìn)行線性擬合評價和統(tǒng)計(jì)指標(biāo)評價,測試模型可靠性。

2 建模原理

太陽與地球間距離遙遠(yuǎn),照射到地球上的太陽光是理想的白色平行光,當(dāng)白色平行光照射到水稻上時,葉綠素會吸收特定波長的光,然后將白光反射到定高飛行的無人機(jī)光譜傳觀器上。由于特定波長的光被吸收,在光譜傳感器輸出光譜曲線上形成波谷,通過分析光譜即可確定葉綠素含量。葉綠素吸收光能,將水和CO2轉(zhuǎn)化為有機(jī)物和氧氣,實(shí)現(xiàn)光合作用,可以作為衡量水稻長勢及產(chǎn)量的標(biāo)志。因此,通過無人機(jī)光譜遙感計(jì)算葉綠素含量可以預(yù)測水稻最終產(chǎn)量。

2.1 葉綠素含量與水稻產(chǎn)量

水稻不同生長周期葉綠素SPAD值與畝產(chǎn)關(guān)系如圖2所示。

圖2 水稻不同生長周期葉綠素SPAD值與畝產(chǎn)關(guān)系Fig.2 Relationship between chlorophyll SPAD and rice yield in different growth period。

葉綠素吸收光能并離子化,能量被儲存在三磷酸腺苷中,最終將水和CO2轉(zhuǎn)化為有機(jī)物和氧氣,實(shí)現(xiàn)光合作用。由于無人機(jī)反射光譜檢測重點(diǎn)針對冠層葉綠素含量,現(xiàn)采用SPAD-502型葉綠素檢測儀對不同生長時期的水稻冠層葉片進(jìn)行檢測,每葉片3點(diǎn)取樣,計(jì)算水稻冠層葉片SPAD值,即

(1)

其中,i為單株水稻冠層葉片數(shù)量,測試結(jié)果如圖2(a)所示。

葉綠素SPAD值從拔節(jié)期開始顯著增長,拔節(jié)期到孕穗期增長速率最大,到抽穗期達(dá)到最大,隨后含量開始下降。這是由于在抽穗期之前光合作用主要集中在植物生長上;從抽穗期開始,轉(zhuǎn)向生殖生長,葉片中的氮開始向稻穗集中轉(zhuǎn)移,葉綠素含量逐漸降低,此時水稻處于生殖生長階段。設(shè)冠層葉片葉綠素含量A為

(2)

其中,LAI為冠層葉面指數(shù)。不同生長時期冠層葉片葉綠素含量A與畝產(chǎn)擬合直線,相關(guān)系數(shù)如圖2(b)所示。其在孕穗期達(dá)到最高,乳熟期次之。這是由于該時期葉片發(fā)育茂盛,葉綠素含量能較好反映水稻生長狀況,抽穗期與乳熟期主要進(jìn)行稻穗生殖生長,葉片變黃,葉綠素含量降低。總體上葉綠素含量與畝產(chǎn)相關(guān)系數(shù)均達(dá)到0.72以上,擬合精度較高,冠層葉綠素含量可以作為水稻產(chǎn)量的評價對象。

2.2 冠層葉綠素?zé)o人機(jī)光譜

水稻處于不同生長時期的無人機(jī)光譜如圖3所示。其中,總體上拔節(jié)期、孕穗期、抽穗期和乳熟期光譜曲線走勢一致,葉綠素對于綠光波段吸收能力很弱,故在530～600nm形成波峰,對于藍(lán)光和紅光吸收能力較強(qiáng),在400～500nm形成較小波谷;在620～680nm形成較深波谷,同時在700～800nm紅光波段與近紅外線波段交接區(qū)域,反射率急劇增長,形成紅邊區(qū)域;在900～960nm近紅外波段形成波谷。由圖3可知:葉綠素具有3個特征吸收波段,其中紅邊區(qū)域與近紅外波段吸收效果顯著,因此采用采用紅邊區(qū)域與近紅外波段組合參數(shù)模式,可提高對于葉綠素含量檢測的穩(wěn)定性與敏感性。

圖3 不同生長時期冠層葉綠素光譜Fig.3 The rice canopy chlorophyll spectral in different growth period。

由于水稻產(chǎn)量與水稻4個生長時期的冠層葉綠素含量擬合,相關(guān)系數(shù)達(dá)到0.7以上,通過冠層葉綠素含量可以有效預(yù)測水稻產(chǎn)量;同時,冠層葉綠素含量的無人機(jī)光譜數(shù)據(jù)具有良好的特征吸收波長,通過紅邊區(qū)域與近紅外波段組合參數(shù)模式可以有效預(yù)測其含量,因此通過分析無人機(jī)光譜可以有效的預(yù)測水稻產(chǎn)量。

3 建模建立

葉綠素具有3個特征吸收波段,在紅邊區(qū)域和近紅外波段吸收效果顯著,采用單一波段不能概括光譜曲線全部信息,因此采用多波段組合進(jìn)行組合計(jì)算,得到水稻植被組合因子。在水稻生長全生命周期中,拔節(jié)期到孕穗期過程中,水稻冠層葉片處于高速生長狀態(tài);當(dāng)處于抽穗期和乳熟期時,冠層葉片由綠色向黃色發(fā)展,要充分考量冠層葉綠素含量的變化,因此選用對于背景變化敏感的差值植被指數(shù)[10]作為模型建立因子,計(jì)算公式如式(3)所示。其中,DNNIR為近紅外光譜反射率;DNR為紅邊區(qū)域光譜反射率。

DVI(NIR,R)=DNNIR-DNR

(3)

由于水稻生長分為拔節(jié)期、孕穗期、抽穗期和乳熟期,為了綜合不同時期對最終畝產(chǎn)的影響,將4個時期的差值植被指數(shù)作為模型建立的因子。

3.1 不同生長時期差值植被指數(shù)與畝產(chǎn)關(guān)系

水稻不同生長時期差值植被指數(shù)與畝產(chǎn)之間的相關(guān)系數(shù)R灰度云圖,如圖4所示。

圖4 不同時期水稻DVI(NIR,R)與產(chǎn)量相關(guān)系數(shù)云圖Fig. 4 Cloud chart of linear correlation coefficient between DVI(NIR,R)and rice yield in different growth period。

圖4中,紅邊波長為橫坐標(biāo),近紅外特征波長為縱坐標(biāo),每個坐標(biāo)相關(guān)系數(shù)R數(shù)值用該點(diǎn)的灰度表示,即坐標(biāo)點(diǎn)越亮,該點(diǎn)對于的相關(guān)系數(shù)R越大。具體數(shù)值參照右側(cè)灰度強(qiáng)度條。

圖4(a)為水稻拔節(jié)期差值植被指數(shù)與畝產(chǎn)相關(guān)系數(shù)云圖,該時期擬合相關(guān)系數(shù)R普遍在0.4以下,最大值為0.631 6。水稻處于拔穗期時,冠層葉片發(fā)育不成熟,葉片發(fā)育程度較低,土壤背景雜波處于主要因素。圖4(b)為孕穗期相關(guān)系數(shù)R云圖,總體上較拔節(jié)期有所提高,大部分處于0.5以上,最大值為0.771 3。這是由于冠層葉片在孕穗期顯著生長,葉綠素含量顯著增加,紅邊區(qū)域吸收能力增強(qiáng)。圖4(c)為抽穗期相關(guān)系數(shù)R云圖,大部分處于0.6以下,最大值為0.634 1,出現(xiàn)在1 080nm和1 089nm處。由于兩特征吸收波段相近,造成紅邊向量與近紅外向量線性相關(guān)較大,因此擬合關(guān)系系數(shù)R較低。圖4(d)為乳熟期相關(guān)系數(shù)R云圖,最大值出現(xiàn)在900nm及1 100nm處,該波段為淀粉特征吸收波段,紅光和藍(lán)光波段相關(guān)系數(shù)強(qiáng)度下降,葉綠素含量降低。在乳熟期水稻生長以生殖生長為主,稻穗中淀粉開始生成。不同生長時期差值植被指數(shù)與畝產(chǎn)相關(guān)系數(shù)最大值及出現(xiàn)波長如表1所示。

表1 不同生長時期DVI(NIR,R)特征波長及與畝產(chǎn)相關(guān)系數(shù)Table 1 Characteristic absorption wavelength and correlation coefficient between DVI(NIR,R)and rice yield in different growth period。

3.2 共線性判定

水稻在生長過程中分為4個時期,現(xiàn)對4個時期的差值植被指數(shù)與最終畝產(chǎn)量進(jìn)行線性擬合,得到水稻產(chǎn)量模型,如式(4)所示。其中,X1～X4分別為拔節(jié)期、孕穗期、抽穗期和乳熟期差值植被指數(shù)DVI(NIR,R)向量。

Y=aX1+bX2+cX3+dX4+m

(4)

在進(jìn)行多元線性擬合前必須確保X1-X4線性無關(guān),即4向量不共線。判定依據(jù)為容忍度T和主成分分析特征值。其中,容忍度T計(jì)算公式為

(4)

其中,Rij為Xj與Xi的線性相關(guān)系數(shù)。容忍度越小,該向量與其他向量的線性先關(guān)性越高。當(dāng)T<0.1時,認(rèn)為線性相關(guān)性非常嚴(yán)重;主成分分析特征值E越小,線性相關(guān)性越高?，F(xiàn)計(jì)算拔節(jié)期、孕穗期、抽穗期和乳熟期差值植被指數(shù)DVI(NIR,R)向量的容忍度和主成分特征值,如圖5所示。其中,拔節(jié)期容忍度為0.097,主成分特征值為0.14;抽穗期容忍度為0.051,主成分特征值為0.09。差值植被指數(shù)DVI(NIR,R)向量在拔節(jié)期與抽穗期高度共線,應(yīng)給予舍棄。因此,以孕穗期和乳熟期差值植被指數(shù)DVI(NIR,R)向量作為模型建立因子。

圖5 共線性判定Fig.5 The collinear test。

3.3 逐步回歸分析

由共線性分析可知:拔節(jié)期與抽穗期容忍度和主成分特征值很小,線性相關(guān)性高,在模型中剔除。因此,水稻產(chǎn)量模型簡化為

Y=aX1+bX2+m

(5)

現(xiàn)采用逐步線性擬合的方法進(jìn)行求解,主要步驟如下:

1)將孕穗期差值指數(shù)、乳熟期差值指數(shù)與畝產(chǎn)量進(jìn)行線性擬合,對回歸系數(shù)βi預(yù)測值與實(shí)際值進(jìn)行F檢驗(yàn)。

y=Xiβi+εi(i=1,2)

(6)

孕穗期與乳熟期對應(yīng)的檢驗(yàn)結(jié)果為F1、F2,取F(1)=max(F1,F2),查表的置信區(qū)間為1-a情況下,臨界值為F1,當(dāng)F(1)>F1時,引入第1個變量。

2)引入第2個量,建立畝產(chǎn)量y與X1和X2的二元線性回歸方程,對回歸系數(shù)預(yù)測值與實(shí)際值進(jìn)行F檢驗(yàn),計(jì)算結(jié)果為F(1)。查表的置信區(qū)間為1-a情況下,臨界值為F2,當(dāng)F(2)>F2引入第2個變量,完成產(chǎn)量預(yù)測模型。

現(xiàn)取8組樣本光譜數(shù)據(jù),采用“留一法”進(jìn)行交叉預(yù)測,以其中7組樣品光譜作為建模樣本,另一組作為預(yù)測樣本,進(jìn)行8次模型預(yù)測,結(jié)果表2所示。表2中,8組預(yù)測值和實(shí)際值平均相對誤差值為7.102%,表明模型具有較高的預(yù)測精度。

表2 水稻產(chǎn)量模型系數(shù)及其預(yù)測值Table 2 Coefficient and predictive value of production prediction of rice。

4 模型評價

4.1 線性擬合評價

多元線性擬合的模型評價標(biāo)準(zhǔn)為線性決定系數(shù)R2和均方根差RMSE。其中,線性決定系數(shù)R2表示模型預(yù)測值與實(shí)際值的符合程度,計(jì)算公式為

(7)

均方根差RMSE為預(yù)測值與真實(shí)值的誤差平方根的均值,計(jì)算公式為

(8)

現(xiàn)計(jì)算采用采用“留一法”進(jìn)行交叉預(yù)測,計(jì)算每一組模型R2和RMSE值,如圖6所示。

圖6 模型線性性能評價Fig. 6 Linear fitting evaluation。

其中,R2值分布在0.607～0.648之間,均值為0.630;RMSE值分布在21.612～25.145之間,均值為24.300。這說明,模型擬合線性性能良好。

4.2 統(tǒng)計(jì)指標(biāo)評價

統(tǒng)計(jì)學(xué)中常將總體相對誤差Rs和預(yù)估精度P作為模型預(yù)測精度的評價標(biāo)準(zhǔn),兩標(biāo)準(zhǔn)的計(jì)算公式為

(9)

采用式(9)計(jì)算得到8組樣品的Rs為0.0142,預(yù)估精度P為93.71。測試結(jié)果表明,該模型在統(tǒng)計(jì)學(xué)上具有良好的預(yù)測精度。

5 結(jié)論

為了實(shí)現(xiàn)大范圍的水稻產(chǎn)量預(yù)測,基于無人機(jī)水稻反射光譜分析建立了水稻產(chǎn)量模型。由于冠層葉綠素濃度可以有效表征植物生長情況,預(yù)估水稻產(chǎn)量。通過分析無人機(jī)測量的水稻冠層葉綠素光譜發(fā)現(xiàn):在620～680nm形成較深波谷;在700～800nm紅光波段與近紅外線波段交接區(qū)域,反射率急劇增長,形成紅邊區(qū)域;在900～960nm近紅外波段形成波谷,包含豐富信息。因此,采取多波段組合進(jìn)行綜合計(jì)算,以完整地提取光譜信息;選用包含紅邊區(qū)域與近紅外區(qū)域的水稻差值植被指數(shù),作為建模因子。同時,分析差值植被指數(shù)在拔節(jié)期、孕穗期、抽穗期和乳熟期的全生長周期中的特征吸收波長,將其和畝產(chǎn)進(jìn)行線性擬合,得到線性相關(guān)系數(shù)。對其在不同生長時期進(jìn)行共線性判定,發(fā)現(xiàn)差值植被指數(shù)在拔節(jié)期與抽穗期的容忍度和主成分特征值均很低,存在嚴(yán)重的線性相關(guān),因此給予舍棄,并對孕穗期與乳熟期差值植被指數(shù)與水稻畝產(chǎn)進(jìn)行逐次線性擬合。對8組樣本進(jìn)行“留一法”預(yù)測表明:模型的線性決定系數(shù)R2值分布在0.607～0.648之間,均值為0.630;RMSE值分布在21.612～25.145之間,均值為24.300。由此表明,模型精度較高。計(jì)算統(tǒng)計(jì)量總體相對誤差Rs=0.0142和預(yù)估精度P=93.71,表明該模型在統(tǒng)計(jì)學(xué)上具有較高精度。本方法數(shù)據(jù)采集方便,光譜處理簡單,適用于大規(guī)模推廣。