通徑分析結合BP 神經網絡方法估算夏玉米作物系數及蒸散量

王怡寧，張曉萌，路 璐，顧 南，王振龍，劉 猛，王國慶

（1. 南京水利科學研究院，南京 210029；2. 河北工程大學水利水電學院，邯鄲 056000；3. 安徽省（水利部淮委）水利科學研究院，蚌埠 233000）

0 引 言

蒸散是地表能量平衡和生態水循環的重要組成部分[1-3]。準確估算蒸散量對指導農田灌溉，提高田間水利用率具有重要意義。目前，蒸散量主要有2 種獲取方式：1）通過儀器監測，2）利用以理論模型為基礎的計算方法進行估算，通過調整模型參數計算蒸散量[4]。通過儀器監測蒸散發的技術已相對成熟[5]，且對儀器設備、安裝環境及操作方法均有較高要求，難以廣泛應用[6]。因此，國內外學者提出了諸多理論估算方法，包括水量平衡法、Penman-Monteith 模型（P-M 模型）、雙層模型、Priestley-Taylor 模型（P-T 模型）等[7]。在總結了現有蒸散量模型優缺點的基礎上，FAO 提出了計算蒸散的標準方法，即作物系數法[8]，該方法目前被廣泛應用。

作物系數作為作物系數法的重要參數，受作物類型、土壤濕度和氣象等要素影響。作物系數的計算方法包括雙作物系數法和單作物系數法2 種[9]，雙作物系數法主要用于區分土壤蒸發和植物蒸騰，所需參數較多，計算復雜；僅研究作物需水量而不需區分蒸發和蒸騰情況下，單作物系數法則簡單實用。但不同地區的作物系數也存有差異，FAO 僅推薦了作物特定生育階段的作物系數值，簡單地采用推薦值計算所得蒸散量與實際值存在差異，且無法獲取連續值[10]。對此，有關學者利用作物生長天數、累積溫度和葉面積指數等指標分別與作物系數建立關系，以此動態地估算作物系數[11-15]。然而，相關研究多通過分析單變量和作物系數間數值關系構建函數，以此估算作物系數。本研究分析多氣象要素對作物系數的影響，并結合葉面積指數構建作物系數估算模型計算蒸散量，以此了解作物系數的動態特征，準確估算作物需水量。

1 材料與方法

1.1 研究區概況

試驗地位于淮北平原南部具有代表性的大型綜合實驗站五道溝水文水資源實驗站內。該站地處蚌埠市固鎮縣新馬橋鎮，屬南北氣候過渡帶。研究區年內降雨分配不均，冬春降雨較少，夏秋降雨較多，其中6—9 月屬于汛期，全年60%以上降雨發生在該時期。砂姜黑土為研究區主要土壤類型，垂直裂隙發育。0～0.1、 ＞0.1～0.2、＞0.2～0.3 和 ＞0.3～0.4 m 各土層土壤平均容重分別為1.25、1.45、1.41 和1.46 g/cm3。站內設有10 套大型稱重式蒸滲儀，可控制水位埋深，自動觀測作物蒸散、降雨入滲、潛水蒸發和土壤水分變化，作物布局以冬小麥和夏玉米輪作為主。

1.2 試驗設計

以淮北平原分布廣泛的砂姜黑土為研究對象，在大型稱重式蒸滲儀內種植夏玉米。蒸滲儀型號為FR101A，分辨率0.025 mm，口徑面積4.0 m2，土柱高4.0 m。考慮該區為地下水淺埋區(地下水埋深為1～3 m)，故蒸滲儀內埋深分別控制為1 和3 m。夏玉米于2018 年6 月22 日播種，10 月9 日收割。玉米生育期內未灌溉，降雨補給水分。

1.3 測定和計算指標

獲取夏玉米全生育期內非雨期葉面積、蒸滲儀及高精度氣象站數據。其中，3 m 埋深情景下，蒸滲儀在6 月30 日、7 月5 日、7 月6 日、8 月8 日、8 月14 日、8 月18 日受人為活動干擾，數據表現異常，予以剔除。

蒸滲儀實測蒸散量：通過大型稱重式蒸滲儀實測夏玉米生育期內蒸散量變化。試驗期間蒸滲儀自動采集數據，10 min/次。

作物生長圖像觀測：蒸滲儀上方設有攝影機采集筒內作物生長圖像，1 h/次。

土壤含水率：蒸滲儀內10、30、50、100、130、180 cm埋深處設有時域反射儀，自動監測包氣帶土壤水分體積含水率，1 h/次。

氣象數據獲取：蒸滲儀南側高2 m 處設有高精度氣象站，可獲取氣溫（T）、空氣濕度（Relative Humidity，RH）、太陽輻射（Q）、凈輻射（Rn）、土壤熱通量（G），大氣壓力（P0）、水面蒸發（E）、風速（U2）、風向（W）和雨量（Pr），10 min/次。

玉米葉面積指數計算：選取每個蒸滲儀測筒內具有代表性的3 株玉米，采用Yaxin-1242 的葉面積儀（北京雅欣理儀科技有限公司）測定玉米葉片面積，測定頻率為7 d/次，并且計算單位土地面積上玉米葉片總面積與土地面積的比值，得出葉面積指數（leaf area index，LAI）。

1.4 夏玉米生長階段劃分

參考FAO 建議的作物生育期階段的劃分，初期為作物生長早期，土壤基本無作物覆蓋；發育期從初期結束到作物基本覆蓋土壤表面；中期從充分覆蓋到開始成熟；后期從中期結束到收獲[16-17]，結合每日采集的夏玉米生長圖像，將夏玉米生長階段劃分為初期、發育期、中期和后期共4 個階段，時間分別為6 月22 日—7 月15 日、7 月16 日—8 月19 日，8 月20 日—9 月12 日和9 月13日—10 月8 日。

1.5 作物系數計算及模擬方法

1.5.1 基于作物系數法計算作物系數（簡稱實測值）

實際蒸散量可利用作物系數法[17]計算，即

式中ET 表示實際蒸散量，mm；ET0表示參考作物蒸散量，mm；Ks表示土壤水分脅迫系數；Kc表示作物系數。本文所用日蒸散量為每小時累加值。

根據式（1）得到作物系數的計算公式，即

實際蒸散量由大型稱重式蒸滲儀測得，前后采集蒸滲儀總質量之差為該時段ET。其基本原理如下：

式中Pr表示降水量，mm；I 表示灌溉量，mm；Eg表示潛水蒸發量，mm；Pa表示深層土層滲漏量，mm；R 表示徑流量，mm；ΔS 表示土壤水蓄變量，mm。資料選取為非雨期，故R 為0。Eg和Pa根據馬里奧特瓶原理測得，灌溉量I 為0。故式（2）可簡化為

式（1）中參考日作物蒸散量ET0采用FAO P-M 公式計算，表達式為[18]

式中ET0表示參考作物蒸散量，mm；Rn表示作物表面凈輻射，MJ/(m2·d)；G 表示土壤熱通量，MJ/(m2·d)；γ 表示濕度計常數，kPa/℃；T 表示平均氣溫，℃；es表示飽和水氣壓，kPa；ea表示實際水汽壓，kPa；Δ 表示飽和水氣壓-溫度曲線斜率，kPa/℃；U2表示距地面2 m 高處風速，m/s。

土壤水分脅迫系數Ks通過式（5）計算[19]，其反映了土壤水分對蒸散發的影響。

式中Ksi表示第i 天土壤脅迫系數；θi表示第i 天0～40 cm土層的平均含水率，%；θ10i和θ30i分別表示第i 天地面以下10 和30 cm 處平均體積含水率，%；θw表示凋萎含水率，%；θf表示田間持水率，%；ρb表示0～40 cm 的土壤容重，g/cm3。根據研究區域砂姜黑土的土壤特性和實測數據，θw為10%，θf為28%，ρb為1.4 g/cm3。

1.5.2 基于構建模型估算作物系數（簡稱模擬值）

由于作物系數受環境因子和作物自身生長狀況影響，可設作物系數為關于葉面積指數和氣象要素的函數

式中LAI 表示葉面積指數；w1...wn表示系列氣象因子；z1...zn表示系列水文因子；其余變量意義同上。

氣象因子影響作物的生長環境，而葉面積指數則影響作物自身的生理生態指標，這些都與作物系數有關?，F有研究表明[4]，葉面積指數是影響作物系數的重要指標，尤其在生長前期，葉面積指數變化較快，響應較大。考慮到氣象要素較多，作物系數對不同因子響應程度不一，故本文利用通徑分析法篩選出對作物系數影響關鍵的因子，再結合BP 神經網絡構建作物系數模型。具體方法如下：

1）在研究多個自變量與因變量間關系、變量間關系復雜以及自變量如何通過對其他自變量影響間接作用于因變量時，多采用通徑分析法[20]。自變量對因變量的作用大小可用通徑系數來表征，通徑系數可表示為

式中rjk表示自變量間相關系數，Rj表示自變量與因變量間相關系數。當j=k 時，稱為直接通徑系數，表示自變量y 對因變量的xj直接作用；當j≠k 時，稱為間接通徑系數，表示自變量xj通過xk對因變量y 的間接作用。通經系數正值表示正向影響，負值表示負向影響。

2）以通徑分析法篩選所得關鍵因子作為模型輸入變量，作物系數為輸出變量。試驗期間共129 組數據，其中1 m 埋深共67 組數據，3 m 埋深共62 組數據，隨機選取70%作為訓練集，驗證集和檢驗集各為15%。模型構建時先確定隱含層神經元數個數，通過MATLAB 神經網絡工具箱實現。隱含層神經元個數計算式為

式中a 為輸入神經元數目；b 為輸出神經元數目；c 為常數，

1.6 評價指標

對實測值和模擬值進行評價，評價指標有相關系數（Correlation Coefficient，r）、均方根誤差（Root Mean Square Error，RMSE）、平均絕對誤差（Mean Absolute Error，MAE）和一致性指數（Concordance Index，dIA）。通常MAE 和RMSE 越小、預報準確率越大，相關系數r和一致性指數dIA越接近于1，表明實際值與模擬值越接近，模型的預測能力越強。

2 結果與分析

2.1 不同埋深水平下葉面積指數隨時間的變化規律

計算得出LAI，并線性插補得到逐日LAI。圖1 展現了夏玉米生育期內不同埋深下LAI 的變化過程：2 種埋深下LAI 變化趨勢基本一致，但數值大小存在差異。生長初期LAI 變化較小，發育期LAI 迅速增大，中期達到峰值，后期隨著葉片的枯黃，LAI 快速減小。初期2 種埋深下LAI 相差不大，其余階段內總體上3 m 埋深的LAI 略小于1 m 埋深。根據五道溝地中蒸滲儀群研究成果得出：砂姜黑土區適宜玉米生長的埋深為0.8～1.5 m，埋深較淺時，玉米生長能有效利用地下水補給量，LAI 相應較大。由上表明，不同的種植條件下，LAI 趨勢雖相似，但大小不同。

2.2 作物系數模型構建及驗證

2.2.1 作物系數主要影響因子篩選

利用通徑分析法篩選氣象要素對作物系數的關鍵影響因子。首先，對夏玉米生育期內2 種埋深下Kc、LAI、T、RH、Q、Rn、G、W、E、U2、P0和各層土壤含水率（10、30、50、100、130、180 cm）進行正態檢驗，P0未通過檢驗。將剩余因子逐步回歸，2 種埋深下Kc均受T、U2和LAI 的影響顯著，3 m 埋深下還受130 cm 處含水率θ130影響。計算各因子對Kc的通徑系數時，將不顯著因子剔除，結果如表1 所示。1 m 埋深下，各因子對Kc的直接作用由大到小為T、U2、LAI，除LAI 的直接通徑系數小于0 外，其他均大于0；間接作用由大到小為T、U2、LAI，U2和T 通過LAI 對Kc的負向作用較為明顯；從決定系數可以看出，各因子對Kc的綜合決定能力LAI＞T＞U2。θ130、T、U2和LAI 對3 m 埋深Kc影響較大，直接作用由大到小為T、θ130、U2、LAI；間接作用由大到小為T、θ130、U2、LAI，除LAI 的間接通徑系數小于0 外，其他均大于0，說明LAI 對Kc的間接影響為負向作用；2 種埋深LAI 通徑系數均為負值，對Kc的影響均為負向，葉面積指數影響地表覆蓋度、冠層蒸騰表面積和通風狀況，葉面積指數在前期日尺度效應中對作物系數的影響較大，后期隨著葉面積指數的減少，在大尺度效應中綜合效應為負效應。從決定系數可以看出T 對Kc的綜合決定能力最大，其他因子對Kc決定能力的大小依次是θ130、LAI 和U2。2 種埋深的Kc均受T、U2和LAI 影響，3 m 埋深Kc同時受θ130影響顯著。這從另一方面驗證了葉面積是影響作物系數的關鍵因子，并且氣溫和風速對作物系數至關重要，這與文獻[5]的研究結果相似。同時，埋深較深時，作物系數還受地表以下130 cm 處土壤水分含水率的影響。130 cm 處于包氣帶的中間位置，可以綜合反映整個包氣帶內水量變化情況。

圖1 夏玉米生育期葉面積指數變化過程 Fig.1 Change processes of leaf area index (LAI) during growth period of summer maize

2.2.2 作物系數模型驗證

依據表1 中Kc的影響因子分析結果，選取影響Kc的關鍵因子（1 m：T、U2和LAI；3 m：T、U2、LAI 和θ130）作為輸入變量，Kc為輸出變量，利用BP 神經網絡分別建立不同埋深下作物系數模型（式（8）），并將模擬值與實測值比較。圖2 展現了夏玉米生育期作物系數模擬值與實際值變化過程。不同埋深水平下，Kc的模擬值與實際值趨勢相似，呈單峰變化。玉米生長初期，田間基本為裸土狀態，無植被覆蓋，Kc較??；發育期間，玉米快速生長，Kc迅速增大；中期及后期隨著玉米快速成熟，Kc也迅速減小。1 和3 m 2 種埋深水平下生育期模擬值與實際值相關系數r 分別為0.94 和0.92，MAE 分別為0.04 和0.08 mm/d，其中1 m 埋深Kc初期、發育期、中期和后期MAE 分別為0.06、0.09、0.05 和0.03 mm/d，3 m 埋深Kc各生育期MAE 分別為0.11、0.10、0.07 和0.03 mm/d，各生長階段模擬精度均較高。將各生育階段每日的Kc取平均值，得到4 個生長階段的Kc。模擬的地下水位1 m 埋深的夏玉米生長初期、發育期、中期和后期4個生長階段Kc平均值分別為0.21、0.63、0.60 和0.29，3 m埋深的夏玉米Kc平均值則分別為0.42、0.67、0.64 和0.19。

表1 不同埋深作物系數與葉面積指數和氣象文因子的通徑分析 Table 1 Path analysis of crop coefficient and leaf area index and meteorological factors at different burial depths

圖2 不同埋深夏玉米生育期作物系數模擬值與實測值 Fig.2 Simulated and actual values of crop coefficient during growth period of summer maize under different burial depths

2.3 不同埋深蒸散量模擬變化

2.3.1 夏玉米全生育期水分脅迫系數變化

Ks是計算作物需水量的一個重要指標，圖3 列出了夏玉米生育期Ks變化過程，圖3a 為1 m 埋深的夏玉米全生育期Ks變化過程，圖3b 為3 m 埋深的夏玉米全生育期Ks變化過程。由圖3 可知，前期作物不受土壤水分脅迫影響，中期、發育期和后期夏玉米均受水分脅迫影響，3 m 埋深的作物受水分脅迫影響較1 m 埋深的大，其中，中期表現最為明顯。

圖3 不同埋深夏玉米生育期土壤水分脅迫系數變化 Fig.3 Changes of soil water stress coefficients during growth period of summer maize at different burial depths

2.3.2 夏玉米全生育期參考蒸散量變化

由式（5）計算1 和3 m 埋深夏玉米生育期ET0，如圖4 所示。由于ET0計算采用實驗站氣象資料，3 m 埋深夏玉米ET0與1 m 埋深夏玉米ET0在數值上相等，2 種不同埋深水平下夏玉米各生長階段ET0均表現出相同的變化特征，初期、生育期、中期和后期ET0呈逐漸波動上升趨勢。由于Rn在生長階段變化較大并出現負值，增強了ET0的波動。2 種埋深的夏玉米生長初期、發育期、中期和后期ET0各生長階段平均值分別為8.9、8.9、13.2 和14.4 mm。

2.3.3 不同埋深水平下蒸散量模擬

夏玉米生育期模擬與實際蒸散變化過程如圖5 所示，ET 模擬值與實際值的變化趨勢相似，初期玉米需水量較小，ET 減少；發育期和中期需水增大，ET 較大；后期需水減少，ET 逐漸變小。1 和3 m 埋深日均ET 分別為4.75 和4.54 mm。隨埋深增大，毛管水難以上升到地表，土壤蒸發減少。同時，作物根系所在土層含水率減少，作物蒸騰減少，故ET減少。全生育期1 m 埋深ET 模擬與實際比值為0.98，相關系數r 為0.93，不同生長階段ET 實際值（模擬值）分別為：2.35(2.19)、5.01(5.45)、7.20(7.16)、3.61(3.67)mm，MAE 分別為0.56、0.59、0.66 和0.45 mm/d；3 m 埋深ET 模擬與實際比值為0.87，r 為0.89，不同生長階段ET 實際值(模擬值)分別為3.15(3.80)、5.38(5.58)、6.42(6.55)、2.24(2.28) mm，MAE 分別為0.82、0.98、0.68、0.29 mm/d。總體而言，2 種埋深ET 模擬精度均較高，且不同生長階段日預報能力較好。

圖4 夏玉米生育期參考蒸散變化 Fig.4 Change of reference evapotranspiration during growth period of summer maize

圖5 不同埋深夏玉米生育期模擬與實際蒸散變化過程 Fig.5 Change of simulated and actual value of evapotranspiration in growth period of summer maize at different burial depths

2.4 不同時間尺度下模型精度評價

基于模型式（8）獲得的作物系數代入式（1）獲得不同時間尺度下2 種埋深下蒸散模擬值，與實際值進行比較，結果見表2。對于日尺度蒸散，1 m 埋深全生育期MAE 為0.72 mm/d，RMSE 為1.00 mm/d，預報準確率（MAE＜1 mm/d）為77.6%。生長初期、發育期、中期和后期MAE 分別為0.56、0.59、0.66 和0.45 mm/d，RMSE分別為0.48、0.47、0.60 和0.30，各生育階段模擬精度均較高。3 m 埋深全生育期MAE 為0.73 mm/d，RMSE 為0.95，預報準確為72.3%，誤差值較小，預報精度高。各生育階段MAE 分別為0.82、0.98、0.68 和0.29 mm/d，RMSE 分別為0.99、0.87、0.76 和0.43，。2 種埋深下不同時間尺度各生育階段dIA均接近1，絕對誤差小于1.0 mm/d，預報準確率（＜1 mm/d）達70%，表明基于葉面積指數和氣象要素的蒸散量模擬模型擬合優度較高，有很好的預報能力。隨著時間尺度由1 d 增至5 d，2 種埋深的dIA及預報準確率逐漸增加，MAE 和RMSE 基本減小，模型擬合優度升高。這可能是由于隨著時間尺度的增加，作物需水量中部分誤差正負值相互抵消，使得整體來看時間尺度越大、誤差越小，并且實際值和模擬值的相關性也隨之增大[21-22]。

表2 夏玉米不同時間尺度下蒸散模擬精度 Table 2 Evapotranspiration simulation accuracy at different time scales of summer maize

3 討 論

對Kc模擬的研究很多，不同的地域，Kc的影響因素也不同。雷志棟等[15]研究結果表明Kc受太陽輻射、相對濕度和風速影響較大，并指出Kc在生育期內的變化與LAI 具有良好的相關關系，建立了兩者的經驗公式。王宇等[22]研究結果表明，夏玉米Kc主要受葉面積指數、氣溫、凈輻射和土壤含水率影響。本研究綜合考慮了葉面積指數和氣象要素對Kc的影響，篩選出影響Kc的關鍵因子模擬Kc，為進一步提高蒸散估算的精度提供了理論依據。

本文預測的Kc和ET 與王振龍等[23]基于三基點溫度估算1 m 埋深的Kc和ET 相比，預報精度較好。1 m 埋深夏玉米Kc的MAE 為0.06 mm/d，而王振龍等[23]得出的MAE 為0.12 mm/d。由于FAO 推薦的Kc是在特定條件下，根據五道溝實驗站的氣象資料進行修訂[24]，修訂后的Kc與模擬的1 m 埋深夏玉米各生長階段Kc相比基本一致。模擬的3 m 埋深夏玉米各生長階段Kc與修訂值相比，初期Kc偏高，發育期和中期Kc基本一致，后期Kc偏低?？傮w來說，Kc的模擬精度較高。利用模擬的Kc估算ET，1 m 埋深的夏玉米ET 全生育期MAE 為0.72 mm/d，與王振龍等[23]估算得出的MAE 0.83 mm/d 相比，前者模擬精度較高。

1 m 埋深的Kc和ET 模擬精度較3m 埋深效果好，這主要是因為埋深較深時Kc和ET 受氣象影響減弱，地下水位埋深對其影響增強。由于模型主要模擬非雨期作物系數和蒸散量變化，對于雨期模擬仍需進一步研究。本研究利用通徑分析比較了氣象多要素對作物系數的影響，并結合BP 神經網絡利用篩選所得關鍵因子模擬作物系數。但神經網絡為黑箱模型，如何依據通徑分析結果提出一種包含關鍵因子的概念性模型，以及與大田尺度實際蒸散對比有待進一步研究。

4 結 論

考慮到氣象因子和葉面積指數對作物系數的影響，基于通徑分析法篩選出影響作物系數的關鍵因子，將葉面積指數和氣象要素對作物系數進行模擬，再利用作物系數法估算作物蒸散量，獲得的結論主要有：

1）利用通徑分析法篩選出對作物系數影響關鍵的因子。結果表明：2 種埋深的作物系數均受氣溫、風速和葉面積指數影響，3 m 埋深作物系數同時受130 cm 處土壤含水率影響顯著。這從另一方面驗證了葉面積是影響作物系數的關鍵因子，并且溫度和風速對作物系數至關重要。同時，埋深較深時，作物系數還受地表以下130 cm 處土壤含水率的影響。

2）以篩選所得關鍵因子為輸入，作物系數為輸出，再結合BP 神經網絡構建作物系數模型。結果表明2 種埋深下蒸散模型的擬合優度較高，蒸散量各生育階段模擬值與實際值一致性指數均接近1.00，驗證了利用葉面積指數和氣象要素模擬作物系數及計算蒸散量的可行性。

3）模擬夏玉米蒸散量時，2 種埋深下不同時間尺度（1、3、5 d）的預報精度均較高，且隨預報時間尺度的增大，絕對誤差減小，一致性指數增大，即預報精度和預報能力提高。本研究利用氣象要素模擬作物系數及計算蒸散量的2 種模型均具有較高的預報能力，能夠滿足灌溉計劃制定和農業生產用水管理的預報要求。