土壤墑情預(yù)測相關(guān)性數(shù)據(jù)耦合性試驗分析*

安小宇，趙復(fù)興，柳海濤

(鄭州輕工業(yè)大學(xué)電氣信息工程學(xué)院，鄭州市，450002)

0 引言

我國是一個水資源非常缺乏的國家，同時農(nóng)業(yè)用水又占據(jù)了很大的比例[1]。在有限的資源下，改善農(nóng)作物的生長狀態(tài)，就需要不斷提高灌溉技術(shù)[2]。提高墑情預(yù)測的精確度，對農(nóng)作物的生長起到至關(guān)重要的作用，并對我國的水資源做出了巨大的貢獻(xiàn)[3]。因此，對墑情的研究是近些年相關(guān)科研機構(gòu)的重要課題之一[4]。

段浩等[5]提出在遙感Penman-Monteith模型中土壤含水量與土壤蒸發(fā)量的關(guān)系，通過運用遙感的方法對蒸發(fā)量進(jìn)行試驗，驗證了土壤含水量的變化與土壤蒸發(fā)量之間的關(guān)系；Scott等[6]通過飽和的含水量和蒸騰系數(shù)建立了土壤墑情變化的模型，驗證了蒸騰因素對土壤含水量的影響；劉小剛等[7]將時間序列分析法與ArcGIS普通克里金差值法相結(jié)合的方法對墑情變化情況進(jìn)行試驗，得出了墑情變化與大氣環(huán)境之間的聯(lián)系；Schmer等[8]通過監(jiān)測到的環(huán)境數(shù)據(jù)進(jìn)行研究，得出了墑情與地表溫度間的耦合關(guān)系；裴源生等[9]以滿足水資源實時調(diào)度模型系統(tǒng)的需求為目標(biāo)，結(jié)合地表因素、大氣環(huán)境因素和人工灌溉因素的空間分布情況，得出對墑情變化造成的實時影響；Wilson等[10]通過對溫度與墑情之間相互關(guān)系的研究，得出墑情的季節(jié)性變化規(guī)律和每日溫度變化的關(guān)系；Lin等[11]采用土壤水動力的研究方法，得出農(nóng)作物種植條件的變化對墑情的影響。

傳統(tǒng)的方法對墑情進(jìn)行預(yù)測時，都會考慮到影響墑情變化的各種變量，這樣就會增加墑情預(yù)測的難度、增加預(yù)測的時間同時降低預(yù)測的精確度。本文首先采用相關(guān)性分析法，來得到墑情的變化與其它變量的耦合大小，然后通過ROC曲線分析得出在有無降水狀態(tài)下的閾值，進(jìn)而通過卡方分析，分別得出在有無降水狀態(tài)下墑情與其它變量之間的耦合關(guān)系，最后通過線性回歸分析和BP神經(jīng)網(wǎng)絡(luò)進(jìn)行試驗對比，為提高墑情預(yù)測的準(zhǔn)確率提供技術(shù)參考。

1 材料與方法

1.1 試驗材料

本文是針對河南省平頂山市的農(nóng)田環(huán)境進(jìn)行分析試驗，該城市處于中國華東地區(qū)，屬于溫帶大陸季風(fēng)性氣候，平均海拔高度為400 m。從該地區(qū)采集到的農(nóng)田環(huán)境變量有：蒸發(fā)量、地溫、降水、氣壓、土壤相對濕度(墑情)、日照時數(shù)、氣溫、風(fēng)速。實測數(shù)據(jù)為2019年3—9月的數(shù)據(jù)，其中3—9月有降水為57 d(26.64%)，無降水為157 d(73.36%)，總數(shù)據(jù)(N)為214組。

1.2 試驗工具

本試驗所使用的數(shù)據(jù)庫借助Excel錄入各種數(shù)據(jù)，然后使用SPSS17.0統(tǒng)計軟件對所錄入的數(shù)據(jù)進(jìn)行試驗分析，得出最終的試驗結(jié)論。在錄入墑情和大氣環(huán)境的數(shù)據(jù)后，針對有無降水的狀態(tài)進(jìn)行分析檢驗，檢驗水平為P<0.05(P為顯著性程度)。

1.3 試驗方法

本文首先采用Pearson相關(guān)系數(shù)法，得出墑情與各變量之間的相關(guān)系數(shù)，然后使用ROC曲線來對有無降水兩種情況進(jìn)行分析，分別得出在這兩種情況下各個變量的閾值，而后通過卡方分析得出在有無降水狀態(tài)下的卡方值，最后利用線性回歸分析和BP神經(jīng)網(wǎng)絡(luò)，來對有無剔除與墑情變化相關(guān)性較小的變量進(jìn)行試驗對比，得出在剔除相關(guān)性較小的變量后，對墑情的預(yù)測更加顯著。

1.3.1 Pearson相關(guān)分析

(1)

式中：EX——X的平均值；

EY——Y的平均值；

N——數(shù)據(jù)的總量；

RX,Y——相關(guān)系數(shù)。

1.3.2 ROC曲線分析

受試者工作特性曲線(ROC曲線)分析主要對診斷的閾值進(jìn)行修正，獲得多對真(假)陽性率值。在ROC曲線下面的面積稱做AUC值，用來判斷診斷試驗的可靠性[13]。通常認(rèn)為AUC的值在0.5～0.7之間時，診斷的可靠性較低；在0.7～0.9之間時，診斷的可靠性中等；大于0.9時診斷的可靠性較高[14]。

本文提出ROC曲線分析是來區(qū)分在有無降水兩種情況下，各變量對墑情變化的影響程度。即在有雨狀態(tài)和無雨狀態(tài)兩種情況下，來分析影響墑情變化的主要因素。根據(jù)試驗的有關(guān)數(shù)據(jù)，對有雨組和無雨組進(jìn)行判定分析，得出在不同狀態(tài)下的臨界值，并且對不同狀態(tài)下各變量的坐標(biāo)分布進(jìn)行計算，同時以此狀態(tài)的敏感性做為縱坐標(biāo)代表真陽性率，特異性作為橫坐標(biāo)代表假陽性率，來對ROC曲線進(jìn)行繪制。

1.3.3 卡方檢驗

卡方檢驗是一種假設(shè)性檢驗的方法之一，用來檢驗兩組或兩組以上的樣本率之間的差別、變量與變量之間有無相關(guān)性等方面的問題[15]。卡方檢驗是在ROC曲線分析后的第一次檢驗，根據(jù)ROC得出的閾值來劃分土壤的高低墑情，然后對所有試驗數(shù)據(jù)進(jìn)行分類，從而得出期望值E(r)；測得的實際值為O(r)，從而根據(jù)式(2)得出各變量與墑情之間的卡方值。

(2)

當(dāng)r為1項集時，E(r)=O(r)，則式(2)轉(zhuǎn)為式(3)。

E(r)=N×E(r1)/N×…×E(rk)/N

(3)

卡方檢驗就是判斷土壤數(shù)據(jù)實際的測試值與理論值的偏差程度，卡方值的大小就表示測試值與理論值的相關(guān)性程度[16]。

1.3.4 線性回歸分析

線性回歸分析在統(tǒng)計學(xué)中，使用線性回歸方程的最小平方函數(shù)對一個或多個自變量的變化情況，來預(yù)測與之相關(guān)的某變量的未來值，而建立的一種分析方法[17-18]。本文以墑情作為因變量Y，影響墑情變化的各個變量作為自變量X，建立線性回歸模型如式(4)所示。通過線性回歸分析，來對比在有無剔除與墑情變化相關(guān)性較小變量時的預(yù)測誤差，從而得出在墑情預(yù)測時的更優(yōu)預(yù)測模型。

Y=β0+β1X1+β2X2+…+βn·Xn+ε

(4)

式中：Y——墑情值；

β0——常量；

ε——殘差；

Xi——影響墑情變化的各個變量，i=1，2,3…n；

βi——影響墑情的變量Xi的系數(shù)。

1.3.5 BP人工神經(jīng)網(wǎng)絡(luò)

截止到目前為止，使用最為廣泛的神經(jīng)網(wǎng)絡(luò)是BP算法的多層前饋網(wǎng)絡(luò)，在許多非線性的模型中，BP算法可以較為準(zhǔn)確的反應(yīng)出變量與變量之間的關(guān)系[19]。常用到的BP結(jié)構(gòu)如圖1所示，由輸入層、隱含層、輸出層組成。

圖1 BP神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)Fig. 1 BP neural network structure

本文將影響墑情變化的7個變量作為輸入層，墑情作為輸出層，隱含層的節(jié)點數(shù)目由式(5)進(jìn)行計算，由式(5)得出本試驗隱含層節(jié)點數(shù)范圍為4～13。

(5)

式中：M——輸入層節(jié)點數(shù)；

L——輸出層節(jié)點數(shù)；

a——整數(shù)，1～10。

本文采用試算法來對隱含層節(jié)點數(shù)進(jìn)行確定。表1列出了隱含層不同節(jié)點個數(shù)下的訓(xùn)練情況，可以看出，相對誤差的變化不大，而決定系數(shù)R2的變化較大，而隱含層節(jié)點數(shù)為9時，R2最小，從而確定隱含層節(jié)點數(shù)為9。最終確定BP神經(jīng)網(wǎng)絡(luò)的拓?fù)浣Y(jié)構(gòu)為7∶9∶1。

表1 不同隱含層節(jié)點數(shù)模型訓(xùn)練的表現(xiàn)Tab. 1 Performance of different number of neurons in the hidden layer under training

本文將影響墑情變化的各變量作為輸入層，隨后在隱含層進(jìn)行算法處理，將經(jīng)過處理后的數(shù)據(jù)與真實數(shù)據(jù)進(jìn)行對比，若其差值不滿足設(shè)定的精度要求，則進(jìn)入反向傳播，此時各層神經(jīng)元的權(quán)值、閾值將進(jìn)行改變，通過以上過程循環(huán)不斷對權(quán)值、閾值進(jìn)行改變，直到滿足設(shè)置的最小差值或者預(yù)設(shè)的訓(xùn)練次數(shù)[20-21]，其過程如下。

1) 參數(shù)初始化。首先設(shè)置網(wǎng)絡(luò)輸入層、隱含層、輸出層節(jié)點數(shù)為M、q、L，其次初始化各個神經(jīng)元層之間的權(quán)值、閾值。

2) 正向傳播算法。隱含層輸出如式(6)所示。

(6)

式中：Hi——隱含層輸出；

j——輸入層；

i——隱含層；

ωij——輸入層與隱含層連接權(quán)值；

a——隱含層閾值；

f——隱含層激勵函數(shù)。

把得出的Hi作為輸出層的輸入值，進(jìn)而對輸出層進(jìn)行計算。輸出層的輸出如式(7)所示。

(7)

式中：Ok——輸出層輸出；

bk——輸出層閾值；

ωki——隱含層與輸出層連接權(quán)值。

3) 輸出O與期望y′k的誤差如式(8)所示。

ek=y′k-Ok,k=1,2,…,L

(8)

式中：ek——第k個節(jié)點輸出與期望的誤差。

4) 權(quán)值、閾值的更新如式(9)～式(12)所示。

(9)

更新輸入層與隱含層之間的連接權(quán)值，得出下一次神經(jīng)元的連接權(quán)值。

ωki=ωki+ηHiek,i=1,2,…q;k=1,2,…,L

(10)

更新隱含層與輸出層之間的連接權(quán)值，得出下一次神經(jīng)元的連接權(quán)值。

(11)

更新隱含層的閾值，得出下一次神經(jīng)元的閾值。

bk=bk+ηek,k=1,2,…L

(12)

更新輸出層的閾值，得出下一次神經(jīng)元的閾值。

式中：η——學(xué)習(xí)效率。

5) 返回步驟(2)，更新過權(quán)值與閾值后重新計算各神經(jīng)元輸出值，當(dāng)輸出誤差小于設(shè)定誤差時，保留當(dāng)前權(quán)值與閾值，否則重復(fù)以上過程直到輸出誤差滿足為止。

通過BP神經(jīng)網(wǎng)絡(luò)，來對比在有無剔除與土壤墑情變化相關(guān)性較小變量時的預(yù)測誤差，從而得出在墑情預(yù)測時的更優(yōu)預(yù)測模型，同時與線性回歸分析的預(yù)測結(jié)果進(jìn)行對比，來驗證本試驗的準(zhǔn)確性。

2 試驗結(jié)果與分析

2.1 Pearson結(jié)果分析

本文首先對傳感器采集的所有數(shù)據(jù)進(jìn)行相關(guān)分析，得出墑情與各變量之間的相關(guān)性數(shù)值，如表2所示。從表2中可以看出，相對濕度與蒸發(fā)量的相關(guān)系數(shù)最大為-0.534，呈負(fù)相關(guān)，并且顯著性程度P<0.05，相對濕度與降水、地溫、氣壓、日照時數(shù)、風(fēng)速呈相關(guān)性，同時顯著性程度P<0.05。氣溫與相對濕度的相關(guān)性較弱，為0.018，把該變量作為墑情預(yù)測時擬剔除的變量。

表2 Pearson相關(guān)性分析結(jié)果Tab. 2 Pearson correlation analysis results

2.2 ROC曲線結(jié)果分析

在有無降水兩種狀態(tài)下，對土壤墑情做進(jìn)一步的ROC分析。圖2為在有雨狀態(tài)下的ROC曲線，從圖2中可以看出在有雨狀態(tài)下，降水的ROC曲線下的面積最大。

在表3中可以看出，在有雨狀態(tài)下，降水量的AUC值為1，相對濕度的AUC值為0.851，而其余幾個變量的AUC值都小于0.5，故在有雨狀態(tài)下，相對濕度的變化只與降水量有較大的相關(guān)性。在表3中還列出了各變量在有雨狀態(tài)下ROC曲線的標(biāo)準(zhǔn)誤差和顯著性的具體數(shù)值。在無雨狀態(tài)下的ROC曲線如圖3所示。

圖2 有雨狀態(tài)下的ROC曲線Fig. 2 ROC curve under rain condition

表3 有雨狀態(tài)曲線下的面積Tab. 3 Area under rain state curve

圖3 無雨狀態(tài)下的ROC曲線Fig. 3 ROC curve under non rain condition

從圖3中可以看出在無雨狀態(tài)下，日照時數(shù)、蒸發(fā)量、地溫、氣壓的ROC曲線下的面積相對較大。在表4中可以看出，在無雨狀態(tài)下，蒸發(fā)量的AUC為0.656，P<0.001；地溫的AUC為0.637，P<0.01；日照時數(shù)的AUC為0.814，P<0.001；氣溫的AUC為0.589，P<0.05，故在無雨狀態(tài)下，相對濕度的變化與蒸發(fā)量、地溫、日照時數(shù)有較大的相關(guān)性。

表4 無雨狀態(tài)曲線下的面積Tab. 4 Area under rain-water condition curve

最終根據(jù)ROC曲線可以得出在有無降水兩種狀態(tài)下每個變量的閾值，如表5所示，進(jìn)而來為卡方分析提供一個標(biāo)準(zhǔn)閾值，進(jìn)一步得出與墑情變化的相關(guān)性較大的變量。

表5 分類后的閾值Tab. 5 Thresholds for classification

2.3 卡方檢驗結(jié)果分析

根據(jù)ROC曲線得出的每個變量的閾值進(jìn)行卡方分析，其中，降水量1表示有降水，0表示無降水。分別對每個變量在高墑值和低墑值兩種情況下進(jìn)行分類分析，即在高墑情況下，大于或者小于各個變量的閾值劃分情況如表6所示，從表6中可以看出蒸發(fā)量的卡方值最大，為67.608，則顯著性最強，說明蒸發(fā)量越大對墑情的影響就越大；其次是日照時數(shù)，卡方值為39.731，顯著性為0.001，說明日照時間越大，對墑情的影響越大；氣溫的卡方值最小，為1.385，顯著性為0.275，表示氣溫與土壤墑情不顯著相關(guān)。經(jīng)過表6分析可得，與墑情變化相關(guān)性較大的變量有蒸發(fā)量、地溫、氣壓、日照時數(shù)和風(fēng)速。

表6 卡方分析結(jié)果Tab. 6 Chi-square analysis result

2.4 線性回歸結(jié)果分析

通過線性回歸來對土壤墑情進(jìn)行預(yù)測，表7～表9為在未剔除與墑情變化相關(guān)性較小的變量時的分析結(jié)果，以此來與本文第2.5節(jié)進(jìn)行對比。從表7中可以看出，回歸分析的Sig值為1.557 2e-32，從表8中可以看出，線性回歸模型預(yù)測的標(biāo)準(zhǔn)偏差為13.229 5，殘差為11.895 9，表9為各變量的回歸系數(shù)。圖4為墑情的回歸標(biāo)準(zhǔn)化的正態(tài)分布圖，從圖4中可以看出試驗所用的數(shù)據(jù)呈正態(tài)分布，說明試驗分析有效。

表7 線性回歸方差分析Tab. 7 Variance analysis of linear regression

表8 線性回歸殘差統(tǒng)計Tab. 8 Residual statistics of linear regression

表9 線性回歸回歸系數(shù)Tab. 9 Regression coefficient of linear regression

圖4 相對濕度的正態(tài)分布圖Fig. 4 Normal distribution diagram of relative humidity

2.5 降維后的線性回歸對比分析

經(jīng)過剔除相關(guān)性較小變量的方差分析和線性回歸系數(shù)的分析如表10所示，從表中可以看出，線性回歸的顯著性值為9.606e-33，而未剔除相關(guān)性較小的變量的顯著性值為1.557 2e-32，說明在剔除相關(guān)性較小變量后的顯著性更好。表11為剔除相關(guān)性較小變量的殘差統(tǒng)計表，從表中可以看出經(jīng)過剔除相關(guān)性較小變量之后的預(yù)測值的標(biāo)準(zhǔn)偏差為12.888 3，殘差為12.264 7，而未經(jīng)過剔除相關(guān)性較小變量的預(yù)測的標(biāo)準(zhǔn)偏差為13.229 5，殘差為11.895 9。由此可見，經(jīng)過剔除相關(guān)性較小的變量之后，對墑情預(yù)測的準(zhǔn)確性更為有利。表12為剔除相關(guān)性較小變量后的回歸系數(shù)。

表10 剔除相關(guān)性較小變量后的方差分析Tab. 10 Variance analysis after the elimination of

表11 剔除相關(guān)性較小變量后的殘差統(tǒng)計Tab. 11 Residual statistics after the elimination of

表12 剔除相關(guān)性較小變量后的回歸系數(shù)Tab. 12 Regression coefficient after the elimination of

2.6 BP神經(jīng)網(wǎng)絡(luò)預(yù)測結(jié)果分析

為了綜合體現(xiàn)BP預(yù)測模型的性能，通過絕對誤差的平均值(PING)、方差(VAR)、相對誤差的平均值(Relativeerror)三個指標(biāo)來評價預(yù)測結(jié)果，在無剔除相關(guān)性較小的變量時，BP神經(jīng)網(wǎng)絡(luò)預(yù)測后的相對誤差為0.004 9，平均誤差為0.024 4，方差為3.248 9×10-4，如表13所示。在刪除相關(guān)性較小的變量時，BP神經(jīng)網(wǎng)絡(luò)預(yù)測后的相對誤差為0.004 7，平均誤差為0.023 2，方差為2.708 0×10-4，如表13所示。經(jīng)過對比可以說明，在剔除相關(guān)性較小的變量后，墑情的預(yù)測精度要優(yōu)于未剔除相關(guān)性較小的變量的預(yù)測精度。

表13 BP神經(jīng)網(wǎng)絡(luò)結(jié)果對比Tab. 13 Comparison of BP neural network results

3 結(jié)論

本文針對影響土壤墑情變化的因素，通過運用相關(guān)性分析、ROC曲線分析、卡方分析相結(jié)合的手段得出與墑情變化的相關(guān)性較大的變量，然后通過線性回歸分析與BP神經(jīng)網(wǎng)絡(luò)，來對有無剔除相關(guān)性較小變量時的墑情進(jìn)行預(yù)測對比分析。

1) 經(jīng)過相關(guān)性分析可以得出，影響墑情變化的相關(guān)性較大的變量為蒸發(fā)量、降水、氣壓、地溫、日照時數(shù)和風(fēng)速。

2) 經(jīng)過ROC曲線分析，得出在有降水和無降水狀態(tài)下的閾值，且對相關(guān)性分析得出的結(jié)果進(jìn)行了進(jìn)一步的驗證。

3) 根據(jù)ROC曲線得出的閾值，然后通過卡方分析可以得出，在墑情較高或者較低的情況下，影響墑情變化最主要的因素是蒸發(fā)量、地溫、氣壓、日照時數(shù)和風(fēng)速。

4) 本文對比了在有無剔除與墑情變化相關(guān)性較小變量的預(yù)測分析，通過分析結(jié)果可以得出，經(jīng)過剔除與墑情變化相關(guān)性較小的變量后，線性回歸分析預(yù)測的標(biāo)準(zhǔn)偏差為12.888 3，BP神經(jīng)網(wǎng)絡(luò)在剔除相關(guān)性較小變量時預(yù)測值的相對誤差為0.004 7，二者預(yù)測的結(jié)果均優(yōu)于未剔除相關(guān)性較小變量時的預(yù)測結(jié)果。

5) 在剔除與墑情變化相關(guān)性較小的變量后，不僅可以降低預(yù)測模型的維數(shù)，其對墑情的預(yù)測值也會更加準(zhǔn)確，為今后農(nóng)場更加精準(zhǔn)的獲悉土壤墑情提供可能。