洪苗，柳平增，張艷，馬學(xué)文，鄭勇，柳建增

(1. 山東農(nóng)業(yè)大學(xué)信息科學(xué)與工程學(xué)院，山東泰安，271000； 2. 山東農(nóng)業(yè)肥業(yè)科技有限公司，山東肥城，271600；3. 山東省現(xiàn)代農(nóng)業(yè)農(nóng)村發(fā)展研究中心，濟南市，250100； 4. 山東省萬興食品有限公司，濟南市，250000)

0 引言

黃瓜是我國主栽蔬菜之一，規(guī)模和產(chǎn)量都位居世界第一[1-2]。環(huán)境條件是實現(xiàn)作物高產(chǎn)、高效、優(yōu)質(zhì)生產(chǎn)的關(guān)鍵因素[3-4]，作物生長與長期所處的環(huán)境有著緊密的聯(lián)系。研究基于主要環(huán)境因子的設(shè)施黃瓜生長模型對提高黃瓜生產(chǎn)水平和管理能力具有十分重要的意義。

目前國內(nèi)外學(xué)者基于溫室環(huán)境的作物生長模型主要從兩個研究角度展開多項研究工作。針對環(huán)境因素的選擇，其研究內(nèi)容以光溫為主。Heuelink[5]用有效積溫法建立番茄的生長模擬模型。陳永快等[6]以有效積溫為變量，采用Logistic曲線方程，建立NFT下小白菜生長模型。劉志剛等[7]以光溫因子—輻熱積(TEP)為驅(qū)動變量，建立溫室基質(zhì)栽培生菜營養(yǎng)生長模擬模型。倪紀(jì)恒等[8]根據(jù)溫室黃瓜果實對溫度和輻射的響應(yīng)，建立以輻熱積為尺度的溫室黃瓜果實模型，能較準(zhǔn)確地預(yù)測溫室黃瓜各節(jié)位的果實生長。還有一些學(xué)者增加了變量的內(nèi)容，在方法的選擇上具有多樣化。張傳帥等[9]采用偏最小二乘回歸研究日光溫室內(nèi)空氣溫度、光照強度、CO2濃度等多個環(huán)境因素的番茄生長模型。員玉良等[10]采用主成分回歸研究了向日葵基于空氣溫度、相對濕度及光合輻射等環(huán)境變量的莖直徑生長模型。Kurtar等[11]采用多元回歸分析的方法研究基于溫度、光照以及SPAD值的黃瓜生長模型。趙亞威等[12]利用BP神經(jīng)網(wǎng)絡(luò)構(gòu)建了葉類蔬菜如小白菜、油菜的環(huán)境數(shù)據(jù)與長勢特征葉面積的生長模型，為溫室植物的種植提供參考依據(jù)。

綜上所述，學(xué)者們多以溫度、光照等環(huán)境因素為尺度研究黃瓜生長模型，不僅環(huán)境變量的選擇具有局限性，并且研究方法上忽略了變量之間多重共線性問題。本文擬利用主成分分析提取能夠代表黃瓜生長狀況的指標(biāo)，利用嶺回歸解決環(huán)境因素與生長指標(biāo)之間的多重共線性問題，并為試驗增加環(huán)境變量個數(shù)，建立基于溫室內(nèi)主要環(huán)境因子的設(shè)施黃瓜不同時期生長模型，為設(shè)施黃瓜種植應(yīng)對環(huán)境因素的變化提供決策依據(jù)。

1 材料與方法

1.1 試驗材料

試驗于2020年7月20日—11月13日在山東省德州市陵城區(qū)麋鎮(zhèn)創(chuàng)新園區(qū)進行。試驗黃瓜為嫁接品種“德瑞特2號”，砧木為南瓜種，種植方式采用無土基質(zhì)盆栽法，2020年7月20日種植第一批，2020年8月30日種植第五批，每批次種植160顆。定植時黃瓜苗已兩葉一心，定植行間距為70 cm，灌溉采用雙行滴管帶，滴箭插置于幼苗根部附近。各批次黃瓜定植后統(tǒng)一進行水肥管理，用量按照常規(guī)園區(qū)常規(guī)管理方式，本試驗不做特殊處理。

1.2 數(shù)據(jù)采集類型與方法

試驗園區(qū)溫室內(nèi)環(huán)境數(shù)據(jù)采用了實驗室自主研發(fā)的“神農(nóng)物聯(lián)”系列物聯(lián)網(wǎng)設(shè)備連續(xù)自動采集完成，采集時間間隔30 min，采集數(shù)據(jù)種類為空氣溫度、空氣濕度、二氧化碳濃度以及光照強度等指標(biāo)，并從數(shù)據(jù)挖掘的角度提煉出晝夜溫差字段，其測量所需設(shè)備信息如表1所示。黃瓜地上部生長數(shù)據(jù)每周采集一次，每次取10株長勢相同樣本采集株高、莖粗、葉片數(shù)、葉面積以及單片葉面積等指標(biāo)，其中葉面積指標(biāo)使用國產(chǎn)LA-S植物圖像分析儀所得，莖粗作物莖粗測量儀測量得到，株高采用直尺或卷尺所得。

表1 采集設(shè)備Tab. 1 Acquisition tools

1.3 數(shù)據(jù)處理方法

1.3.1 主成分分析原理及方法

主成分分析(principal component analysis，PCA)，又名主分量分析，是通過降維技術(shù)，把具有一定相關(guān)性的變量，按照貢獻率的大小，重新組合成一組新的、且相互獨立、少數(shù)能充分反映原始數(shù)據(jù)信息的綜合指標(biāo)以替代原來的指標(biāo)，進而在保留重要信息的前提下，避開原變量之間的線性相關(guān)關(guān)系，便于進行下一步的數(shù)據(jù)分析[13]。

設(shè)施黃瓜地上部生長指標(biāo)有5個樣本，每個樣本有18個數(shù)據(jù)，則構(gòu)成5×18維矩陣

(1)

當(dāng)數(shù)據(jù)類型不一致時，通常具有不同的量綱和量綱單位，此時需要對數(shù)據(jù)進行標(biāo)準(zhǔn)化處理以消除量綱的影響，從而解決數(shù)據(jù)之間的可比性。首先將黃瓜生長指標(biāo)進行標(biāo)準(zhǔn)化處理，將矩陣X標(biāo)準(zhǔn)化為Z，其中

(2)

式中：μ——所有生長指標(biāo)數(shù)據(jù)的均值；

σ——所有生長指標(biāo)數(shù)據(jù)的標(biāo)準(zhǔn)差。

將式(1)經(jīng)過標(biāo)準(zhǔn)化處理后得到設(shè)施黃瓜生長指標(biāo)的相關(guān)系數(shù)矩陣

(3)

計算生長指標(biāo)的相關(guān)矩陣R的特征根λi及其特征向量αi，并將特征根由大到小排列。λi越大則所對應(yīng)的主成分反映的信息就越大。第1個主成分F1可表示為

F1=α1，1Z1+α2，1Z2+…+α18，1Z18

(4)

求解貢獻率并確定主成分個數(shù)，第i個主成分的貢獻率

(5)

m個主成分的累計貢獻率

(6)

當(dāng)累積貢獻率大于80%時，則可以選取m個因子作為主成分，這m個因子能較大限度的保留原數(shù)據(jù)的信息。

1.3.2 多重共線性分析

在多元回歸中，當(dāng)兩個或更多的自變量相關(guān)時，就有可能出現(xiàn)多重共線性(multicollinearity)的情況[14]。當(dāng)多重共線性嚴(yán)重時，模型或數(shù)據(jù)的微小變化可能造成系數(shù)估計的較大變化，這將使得結(jié)果模型不穩(wěn)定，也不容易解釋。一般來說，只要不存在嚴(yán)重共線性，對預(yù)測不會有較大影響，但高度的多重共線性會造成計算困難。

方差膨脹因子(VIF)是多重共線性的度量之一，其通過檢查指定的解釋變量能夠被回歸方程中其他全部解釋變量所解釋的程度來檢測多重共線性。其定義為

(7)

式中：Rj2——第j個變量在其他所有變量回歸時的確定系數(shù)。

VIF取值接近于1，多重共線性越輕，反之越重。通常以10作為共線性判斷邊界。當(dāng)VIF<10，可認(rèn)為不存在多重共線性；當(dāng)10≤VIF<100，則認(rèn)為存在較強的多重共線性；當(dāng)VIF≥100，則認(rèn)為多重共線性非常嚴(yán)重。設(shè)施黃瓜不同生長期生長指標(biāo)和環(huán)境因子的VIF檢查，其數(shù)值多大于10，說明變量之間存在較強共線性，模型建立時應(yīng)選擇適用于多重共線性的方法。

1.3.3 嶺回歸分析原理及方法

針對多重共線性問題，1962年，霍爾(A.E.Hoerl)首先提出改進最小二乘估計的方法，稱為嶺估計。1970年，霍爾與肯納德(Kennard)合作，對嶺估計做了進一步的改進[15]。嶺估計通過嶺參數(shù)k的不同取值獲得回歸參數(shù)的估計值，當(dāng)k=0時就是普通最小二乘估計。

嶺估計通過損失部分信息，降低精度，放棄最小二乘法的無偏性的代價得到更加穩(wěn)定的估計結(jié)果，相較于最小二乘法，嶺估計對病態(tài)數(shù)據(jù)的擬合度更好，且均方誤差更小，因此專用于共線性數(shù)據(jù)分析[16]。當(dāng)變量存在多重共線性，即數(shù)學(xué)描述為|X′X|≈0，那么在原本的相關(guān)矩陣X′X基礎(chǔ)上加上一個正常數(shù)矩陣kI(k>0)，使得X′X+kI接近奇異的程度將會大大降低，即X′X+kI=0的可能性比X′X=0的可能性小很多，從而使嶺估計具有更小的均方誤差，設(shè)嶺估計

(8)

2 結(jié)果與分析

2.1 設(shè)施黃瓜地上部生長指標(biāo)主成分分析

提取出表現(xiàn)設(shè)施黃瓜苗期和初花期長勢好壞的主要特征，其中X1、X2、X3、X4、X5分別表示株高、莖粗、葉面積、葉片數(shù)以及平均單葉面積。首先利用R語言中scale函數(shù)對生長數(shù)據(jù)進行標(biāo)準(zhǔn)化處理，將其轉(zhuǎn)化為無量綱的純數(shù)值，便于不同單位或者量級的指標(biāo)能夠進行和加權(quán)，之后選用主成分分析法進行長勢特征的提取。苗期結(jié)果如表2所示。

表2 苗期生長指標(biāo)相關(guān)性分析Tab. 2 Correlation analysis of seedling growth indexes

顯然，苗期各生長指標(biāo)之間顯著相關(guān)，故此選擇主成分分析作為數(shù)據(jù)分析的主要手段，對生長數(shù)據(jù)進行主成分提取。由表3苗期總方差解釋和圖1碎石圖可見，第1個主成分因子的貢獻率大于0.8，即可得到m=1，即選定一個主成分F1。

表3 苗期總方差解釋Tab. 3 Explanation of total variance at seedling stage

圖1 苗期碎石圖Fig. 1 Seedling lithotripsy

表4 苗期初始因子成分矩陣Tab. 4 Initial factor component matrix of seedling stage

從表4可以看出，葉面積、葉片數(shù)、莖粗、平均單葉面積、株高在第一主成分上均有較高的載荷，說明第一主成分F1完全包含了這5個指標(biāo)的信息。即選取這一個主成分可以完全反映設(shè)施黃瓜苗期長勢的指標(biāo)信息。所以根據(jù)初始因子載荷矩陣及其特征值計算出對應(yīng)的單位特征向量，并結(jié)合式(4)，可以得到

F1=0.48×X3+0.45×X4+0.44×X2+0.44×X5+0.42×X1

(9)

由以上分析可知，葉面積在該主成分上占據(jù)最高載荷，則葉面積可作為評估黃瓜苗期長勢的主要指標(biāo)。

同理，對初花期生長數(shù)據(jù)進行主成分分析挑選出最能代表初花期長勢的指標(biāo)。

表5 初花期總方差解釋Tab. 5 Interpretation of total variance at initial flowering

圖2 初花期碎石圖Fig. 2 Early flowering lithotripsy

表6 初花期初始因子成分矩陣Tab. 6 Initial factor component matrix of early flowering stage

由表5初花期總方差解釋和圖2初花期碎石圖可見，第一個和第二個主成分因子的貢獻率之和大于0.8，即可得到m=2，即選定兩個主成分F1、F2。從表6可以看出，葉面積、葉片數(shù)、莖粗在第一主成分上具有較高的載荷，株高、單葉面積在第二主成分上具有較高的載荷，說明主成分F1、F2基本上包含了這5個指標(biāo)的信息。即選取這2個主成分完全可以反映設(shè)施黃瓜初花期長勢的指標(biāo)信息。因此得到初花期主成分函數(shù)

F1=0.59×X3+0.58×X1+0.50×X4+0.09×X2+0.22×X5

(10)

F2=0.12×X3-0.15×X1-0.37×X4+0.66×X2+0.63×X5

(11)

由以上分析可知，葉面積在第一主成分上具有最高載荷，則葉面積可作為評估黃瓜初花期長勢的主要指標(biāo)。

綜上所述，最能表現(xiàn)設(shè)施黃瓜苗期和初花期生長的指標(biāo)皆是葉面積。因此在黃瓜生長時，葉面積越大，黃瓜的長勢越好，結(jié)瓜越多，瓜農(nóng)收益越高。與現(xiàn)場調(diào)研以及與專家咨詢所得結(jié)果一致，葉面積越大，光合作用越強，植株產(chǎn)生營養(yǎng)物質(zhì)的能力越強，結(jié)瓜越多。所以，若想要黃瓜長勢良好，則應(yīng)適度增大葉面積。

2.2 基于主要環(huán)境因子的設(shè)施黃瓜苗期生長模型

將環(huán)境數(shù)據(jù)以及長勢數(shù)據(jù)進行標(biāo)準(zhǔn)化處理后，使用R程序包car包的函數(shù)vif()將所采集苗期葉面積數(shù)據(jù)及環(huán)境數(shù)據(jù)進行共線性檢查，結(jié)果顯示VIF溫度=182.926 091，VIF光照=426.808 457，VIFCO2=3.538 106，VIF晝夜溫差=63.243 900，VIF濕度=155.156 345。因其數(shù)值除CO2外皆大于10，則說明環(huán)境數(shù)據(jù)之間存在多重共線性問題，故考慮嶺回歸算法建立設(shè)施黃瓜苗期生長模型。

在嶺回歸中，最困難的是確定嶺參數(shù)k*。使用程序包ridge中l(wèi)inearRidge()函數(shù)分析苗期葉面積和環(huán)境變量之間的關(guān)系，由linearRidge()函數(shù)得嶺跡k*=0.133 520 1，滿足嶺跡選擇k的條件。所得結(jié)果如表7所示。

表7 苗期嶺回歸結(jié)果Tab. 7 Ridge regression results at seedling stage

經(jīng)過回歸后擬合得嶺回歸方程

Y葉面積=1.419×10-7+0.212 4×Z溫度+0.260 6×Z光照-0.222 0×ZCO2+0.184 5×Z晝夜溫差+0.284 0×Z濕度

(12)

分析嶺回歸方程可得：當(dāng)Z光照、ZCO2、Z晝夜溫差、Z濕度不變時，Z溫度每增加1個單位，Y葉面積則平均增加0.212 4個單位；當(dāng)ZCO2、Z晝夜溫差、Z溫度、Z濕度不變時，Z光照每增加1個單位，Y葉面積則平均增加0.260 6個單位；當(dāng)Z溫度、Z光照、Z晝夜溫差、Z濕度不變時，ZCO2每增加1個單位，Y葉面積則平均減少0.222 0個單位；當(dāng)Z光照、ZCO2、Z溫度、Z濕度不變時，Z晝夜溫差每增加1個單位，Y葉面積則平均增加0.184 5個單位；當(dāng)Z光照、ZCO2、Z溫度、Z晝夜溫差不變時，Z濕度每增加1個單位，Y葉面積則平均增加0.284 0個單位。

2.3 基于主要環(huán)境因子的設(shè)施黃瓜初花期生長模型

表8 初花期嶺回歸結(jié)果Tab. 8 Early flowering ridge regression results

經(jīng)過回歸后擬合得嶺回歸方程

Y葉面積=2.584×10-6+0.291 8×Z溫度+0.337 1×Z光照+0.271 7×ZCO2

(13)

分析擬合方程可得：當(dāng)Z光照、ZCO2不變時，Z溫度每增加1個單位，Y葉面積則平均增加0.291 8個單位；當(dāng)ZCO2、Z溫度不變時，Z光照每增加1個單位，Y葉面積則平均增加0.337 1個單位；當(dāng)Z光照、Z溫度不變時，ZCO2每增加1個單位，Y葉面積則平均增加0.271 7個單位。

2.4 模型驗證分析

對所得的嶺回歸生長模型進行顯著性檢驗與方差分析。從表9中可以看出苗期該回歸模型的決定系數(shù)大于90%，方差為1.965，殘差為2.839，檢驗結(jié)果為極顯著，可用于黃瓜苗期長勢與環(huán)境因素的相關(guān)性分析。

表9 模型檢驗Tab. 9 Model test

對黃瓜苗期長勢影響排名前三的環(huán)境因素依次是光照、溫度和濕度，且影響情況為正相關(guān)，其他變量影響較小。初花期該回歸模型的決定系數(shù)大于90%，方差為1.43，殘差為1.972，檢驗結(jié)果極顯著，可用于黃瓜初花期長勢與環(huán)境因素的相關(guān)性分析。以此證明，對黃瓜初花期長勢影響最大的是溫度和光照，其次是CO2濃度。

3 結(jié)論

本研究以山東省德州市陵城區(qū)麋鎮(zhèn)創(chuàng)新園區(qū)秋冬茬設(shè)施黃瓜為試驗對象，對基于溫室主要環(huán)境因子設(shè)施黃瓜生長模型進行研究，主要結(jié)論如下。

1) 采用PCA提取黃瓜地上部分多個長勢指標(biāo)序列特征發(fā)現(xiàn)，苗期和初花期的葉面積載荷值在第一主成分中分別為0.994、0.980，說明葉面積能夠較明顯的表現(xiàn)黃瓜的苗期和初花期長勢。在實際生產(chǎn)中可以幫助用戶更便捷的判斷黃瓜長勢好壞，葉面積越大，農(nóng)戶后期收益越高，因此在生長發(fā)育期適度增大葉面積對設(shè)施黃瓜種植管理至關(guān)重要。

2) 由嶺回歸建立的設(shè)施黃瓜生長模型發(fā)現(xiàn)，苗期主要環(huán)境影響因素為光照、溫度和濕度，初花期主要影響因素為溫度和光照，且都為正向顯著影響，模型的決定系數(shù)皆大于90%，檢驗結(jié)果顯著性較高。

3) 在實際生產(chǎn)種植中環(huán)境參數(shù)變化較大，不同參數(shù)之間關(guān)系復(fù)雜，將變量進行關(guān)聯(lián)性分析更有利于指導(dǎo)實際應(yīng)用，可為后期溫室環(huán)境調(diào)控提供參考。

作物生長是一個復(fù)雜的過程，各長勢指標(biāo)之間互有關(guān)聯(lián)，本研究僅探尋了環(huán)境因素對黃瓜長勢主要指標(biāo)的影響，未對其他長勢指標(biāo)的作用進行研究，并且種植管理等因素對黃瓜長勢的影響也有待進一步的探討。