基于增強回歸樹的麥-玉輪作農田蒸散影響因素分析

李 成 李兆哲 王讓會 丁奠元 徐 揚

(1.揚州大學植物保護學院， 揚州 225009； 2.南京信息工程大學江蘇省農業氣象重點實驗室， 南京 210044；3.揚州大學水利科學與工程學院, 揚州 225009； 4.揚州大學農學院， 揚州 225009；5.揚州大學植物功能基因組學教育部重點實驗室， 揚州 225009)

0 引言

農田水分循環過程是土壤-植物-大氣連續體研究的重要內容[1]，蒸散(Evapotranspiration，ET)作為其中的關鍵變量之一，在農田與大氣間的水分和能量交換中扮演了重要角色[2]。農田ET包括植物蒸騰與土壤蒸發，它與農田灌溉、作物生長等過程關系密切[3-4]。因此，在淡水資源日益趨緊的背景下，理解農田ET的動態變化對提高農業水資源利用效率、發展氣候智慧型農業，具有重要的現實意義。

農田ET的實測方法有很多，主要包括渦度相關法、波文比法、穩定同位素法、蒸滲儀法、大孔徑閃爍儀法等[5]。其中，渦度相關法在測定時不會損壞植物且測量精度較高，是近年來應用較為廣泛的方法之一[6]。基于渦度相關法的農田ET研究表明：ET在作物生長季內有明顯的日、季節變化特征，但由于地理位置、水熱條件、栽培管理方式的差異，使不同農田ET的變化幅度、峰值及其出現時間也有所不同[7-8]。另一方面，目前基于該方法的研究大多采用短期觀測資料進行分析(小于3年)，在一定程度上易受天氣條件的影響[9]，不利于全面認識農田ET的變化規律。已有研究表明：不同環境因子，如凈輻射(Net radiation，Rn)、氣溫(Air temperature，T)、飽和水汽壓差(Vapor pressure deficit，VPD)、風速(Wind speed，WS)、土壤含水率(Soil water content，SWC)等，對作物生長和農田ET變化產生重要影響[10-11]。但這些環境因子之間往往存在共線性問題[12]，使一些傳統分析方法(如多元回歸分析、結構方程模型等)在量化農田ET的影響因子和影響程度方面，可能存在一定偏差。為厘清這些環境因子對農田ET變化的影響，近年來興起的機器學習方法，如隨機森林方法、增強回歸樹(Boosted regression tree，BRT)方法等，在探究多因子相互關系方面有重要潛力[13-14]。以BRT方法為例，它能在不預先設定數學表達式的前提下，模擬輸入變量和輸出變量之間的非線性關系，并通過提取各輸入變量的邊際效應，厘清它們對輸出變量的影響[15]；同時BRT方法對各變量的數據類型、概率分布和共線性等具有較大的容忍度[16]，因而將該方法用于農田ET的多因子影響分析，可為闡明不同因子對農田ET變化的影響貢獻與影響方式等提供方法參考。然而，目前基于該方法對輪作田不同作物生長季及農閑期階段ET變化及其影響因子的對比研究，還相對缺乏。

為此，本研究以華北典型麥-玉(冬小麥-夏玉米)輪作田為例，基于逐日觀測數據(包括7個小麥生長季和8個玉米生長季)，分析不同作物生長季和農閑期階段ET的變化特征；并利用BRT方法，闡明農田ET對環境因子的響應特征，以期為變化環境下農業水資源高效利用以及作物模型優化等提供科學依據。

1 材料與方法

1.1 研究站點概況

研究站點位于山東省禹城市西南部(36.83°N，116.56°E)，是中國科學院生態系統研究網絡(Chinese Ecosystem Research Network，CERN)農田觀測站之一。該站年平均氣溫為13.2℃，年均降雨量為582 mm，太陽總輻射達5 225 MJ/m2，無霜期約 200 d[17]。 地貌類型為黃河沖積平原，土壤母質為黃河沖積物，耕層土壤有機質質量比在10～12 mg/kg之間，pH 值介于7.8～8.0之間[18]。種植模式為小麥與玉米輪作，二者的播種與收獲日期如表1所示。

1.2 數據獲取與處理

本研究所用到的逐日觀測數據，如潛熱通量、Rn、T、VPD、WS及表層(0～5 cm)SWC等，由CSAT3型三維超聲風速儀(Campbell Scientific公司，美國)、LI-7500型H2O/CO2紅外氣體分析儀(LICOR Biosciences公司，美國)和配套的自動微氣象觀測系統實測獲得。將儀器設立在面積約0.13 km2的小麥-玉米輪作田中。其中，三維超聲風速儀和H2O/CO2紅外氣體分析儀安裝在植被冠層上方2.5 m處，盛行風向的風浪區長度為5 km，滿足通量觀測試驗的基本要求。并且相關實測數據已通過質量控制，具有較高的可信度[19]。觀測時段為2003年6月至2010年10月，逐日ET通過潛熱通量與汽化潛熱的比值計算得到[9]。根據表1中作物播種/收獲日期，將逐日ET及環境要素數據劃分為小麥季、玉米季和農閑期階段，分別進行統計分析。

表1 禹城站小麥與玉米的播種/收獲日期Tab.1 Harvest date and sowing date of wheat and maize at Yucheng site

1.3 BRT方法

BRT方法是一種在傳統的分類回歸樹基礎上，結合提升算法(Boosting)開發而來的自學習方法。它使用二元遞歸分割算法消除不同因子之間的相互作用，并通過對多個回歸樹進行加權平均，從而提高模型的準確性[14]。具體計算步驟如下[15-16]：①隨機抽取一定比例的數據樣本，產生大量具有固定節點數的回歸樹模型，并利用剩余數據對擬合結果進行檢驗。②采用逐步迭代的方法，按順序學習每一個回歸樹。在每個回歸樹前乘以一個較小學習速率后加入總預測模型。③利用交叉驗證的方法，確定最優回歸樹的數量，以提高模型的穩定性和準確性，且避免出現過擬合現象。④根據輸入變量在分類回歸樹中對輸出變量離差平方和的減少量，計算各輸入變量的重要性程度，該重要性程度值是模型中多個回歸樹離差平方和減少量的平均值。由于BRT方法在處理非線性觀測數據方面具有較高的預測精度，并且計算結果穩定、易于解釋等優點，目前已被廣泛用于歸因分析和邊際效應研究等領域[20]。

本研究基于R Studio軟件和gbm包(版本：2.1.8)，分析不同環境因子(T、Rn、VPD、WS和SWC)與ET之間的關系，進而識別影響ET的關鍵因子及其貢獻率。BRT涉及的相關參數如下：決策樹復雜度為8，學習速率為0.001，分割比例為0.75。具體的計算流程見文獻[15]。

2 結果與分析

2.1 麥-玉輪作田小氣候特征

麥-玉輪作田T、Rn、VPD、WS、SWC和ET的逐日變化曲線，如圖1所示。T表現出冬季低、夏季高的單峰型變化特征，研究時段內T的日最高值和最低值分別出現在2005年6月下旬(32.7℃)和2003年1月上旬(-11.3℃)，年平均T為13.0～13.8℃，多年平均值為13.4℃(表2)。對于不同作物生長季而言，小麥季平均T低于年平均水平，變化范圍在7.1～9.2℃之間；玉米季平均T變化范圍為23.6～25.0℃，多年平均值為24.2℃。農閑期平均T介于小麥季和玉米季之間，達19.2℃。

圖1 研究時段內麥-玉輪作田T、Rn、VPD、WS和SWC的季節變化曲線Fig.1 Seasonal changes of T, Rn, VPD, WS and SWC in wheat-maize rotation field during study period

與T類似，Rn也呈單峰型變化特征。研究時段內逐日Rn在-47.2～217.8 W/m2波動，年平均Rn介于56.3～67.0 W/m2，多年平均值為63.2 W/m2(表2)。除2003年和2008年外，其他年份Rn的年平均值大多為65 W/m2左右。對于不同階段而言，小麥季Rn為(47.4±2.2) W/m2，明顯低于玉米季(100.6±9.1) W/m2，而農閑期Rn則介于二者之間，為(60.3±13.2) W/m2。

VPD在年內呈現出雙峰型特征，兩個峰值分別出現在6月和10月，通常每年第1個峰值的VPD明顯高于第2個峰值。研究時段內逐日VPD變化范圍在0～3.0 kPa之間，多年平均值為0.56 kPa(表2)。就不同作物生長季而言，小麥季VPD的多年平均值較玉米季略低，二者分別為0.48 kPa和0.71 kPa。

WS和SWC在年內變化較為復雜，具有明顯的波動特征。逐日WS變化范圍在0～8.2 m/s之間，多年平均值為1.8 m/s。其中WS超過4 m/s以上的天數主要集中在小麥季，約有135 d；而玉米季和農閑期WS超過4 m/s以上的天數不足小麥季的10%。研究時段內SWC基本在0.2 m3/m3以上，其中小麥季SWC為(0.33±0.04) m3/m3，較玉米季SWC整體偏低(表2)。

2.2 麥-玉輪作田ET變化特征

麥-玉輪作田ET的逐日變化曲線如圖1所示。ET在年內表現出雙峰型變化特征，兩個峰值分別出現在小麥季和玉米季，通常每年第1個高峰的ET明顯高于第2個高峰。研究時段內逐日ET變化范圍在0～9.6 mm/d之間，多年平均值為682.6 mm(表2)。

就不同階段而言，小麥季ET的日最高值為9.6 mm/d，出現在2009年5月8日。研究時段內小麥季單日ET超過5 mm/d共有109 d，分布在4—6月。其中，出現在5月的天數最多，約有70 d；而出現在4月和6月分別有26 d和13 d。小麥季ET變化范圍在315.1～503.9 mm之間，多年平均值為408.5 mm。

研究時段內玉米季ET的日最高值為7.0 mm/d，出現在2010年7月31日。玉米季單日ET超過 5 mm/d 共有18 d，分布在7—9月，其中以出現在 8月 的天數最多(11 d)，7月次之(5 d)。玉米季ET變化范圍在221.7～315.0 mm之間，多年平均值為261.2 mm。

農閑期ET的日最高值為3.5 mm/d，出現在2006年6月16日。與小麥季和玉米季相比，農閑期ET相對較低，變化范圍在14.1～43.9 mm之間，多年平均值為25.8 mm。

2.3 麥-玉輪作田ET變化影響因子

考慮到環境因子之間可能存在較高的相關關系，本研究利用BRT方法分析不同作物生長季和農閑期階段ET變化的多環境因子影響。對BRT模型運行結果進行交叉驗證：小麥季、玉米季和農閑期階段觀測值與模擬值的R2分別為0.92、0.81和0.89；均方根誤差分別為0.5、0.51、0.24 mm。上述結果反映了模型的擬合效果整體較好，證明BRT模型結果可信。各環境因子對麥-玉輪作田ET的影響貢獻率，如表3所示。

表3 環境因子對麥-玉輪作田ET的影響貢獻率Tab.3 Contribution of environmental factors to ET in wheat-maize rotation field %

在影響小麥季ET變化的環境因子中(表3)，Rn貢獻率最大(81.4%)，VPD次之(5.9%)，其余因子按貢獻率排序分別為SWC(5.6%)、T(4.2%)、WS(2.9%)。圖2為貢獻率排在前3位的環境因子與小麥季ET之間的邊際效應，進而揭示主要環境因子(Rn、VPD和SWC)對小麥季ET的影響。Rn對小麥季ET的影響隨Rn增大而呈正相關關系(圖2a)，特別是Rn在30～180 W/m2時，影響程度陡然上升，并達到最大值；但當Rn超過180 W/m2后，Rn對ET的影響程度基本保持不變。VPD與小麥季ET之間的影響曲線呈波動特征(圖2b)。當VPD在 0～1.2 kPa 時，隨著VPD的增加，曲線呈上升趨勢；并在1.4 kPa處曲線達到最高點，產生最高正影響力；此后隨VPD的增加，曲線呈下降趨勢，并逐步達到穩定狀態。當SWC小于0.4 m3/m3時，對小麥季ET具有一定程度的影響；當SWC高于0.4 m3/m3時，則影響程度保持相對平穩(圖2c)。

圖2 主要環境因子對小麥季ET影響的邊際效應Fig.2 Marginal effects of main environmental factors on ET in wheat seasons

玉米季ET受各環境因子的影響也存在明顯差異。Rn是影響玉米季ET變化的關鍵因子，貢獻率達52.7%(表3)。如圖3a所示，當Rn在10～160 W/m2之間時，影響曲線隨Rn的增大而快速上升，并在Rn達到180 W/m2時，正影響力達到最大；當Rn超過180 W/m2后，Rn對ET的影響程度基本保持不變。VPD對玉米季ET的影響僅次于Rn，其貢獻率為17.6%(表3)。從VPD對ET的影響曲線可知(圖3b)，整體呈倒“U”形特征；當VPD約為0.7 kPa時，影響曲線出現最高值，然后開始出現局部下降趨勢，并在VPD達到1.7 kPa時，曲線逐步達到穩定狀態。SWC也是影響玉米季ET變化的重要因子之一，其貢獻率為12.7%(表3)。當SWC處于0.35～0.60 m3/m3時，影響曲線呈現出一定的波動增加趨勢(圖3c)。

圖3 主要環境因子對玉米季ET影響的邊際效應Fig.3 Marginal effects of main environmental factors on ET in maize seasons

農閑期ET受環境因子的影響如表3和圖4所示。與小麥季和玉米季不同，影響農閑期ET變化的關鍵環境因子依次為Rn、SWC和WS，貢獻率分別為36.8%、22.5%和17.0%(表3)。圖4為貢獻率排在前3位的環境因子與ET之間的邊際效應。影響曲線分別在一定范圍內隨Rn、SWC和WS的增大而增大，但當Rn超過135 W/m2時，影響曲線不再隨Rn的增加而繼續增大(圖4a)。SWC的影響曲線具有波動變化特征(圖4b)。當SWC在0.36～0.43 m3/m3之間時，曲線隨SWC增加而呈上升趨勢；當SWC增加至高于0.50 m3/m3時，影響曲線不再隨SWC的增大而變化。與SWC相似，當WS小于1.8 m/s時，隨著WS的增加，曲線呈上升趨勢；而當WS高于2.7 m/s時，其影響程度保持相對平穩(圖4c)。

圖4 主要環境因子對農閑期ET影響的邊際效應Fig.4 Marginal effects of main environmental factors on ET in non-growing seasons

3 討論

為進一步理解禹城站農田ET特征，表4列出了華北平原不同站點農田ET的監測結果。如表4所示，基于渦度相關法的華北麥-玉輪作田ET大多介于640～791 mm之間，并表現出一定的空間異質性[21-26]，如封丘站和欒城站ET較高，分別為791、693 mm，而館陶站ET較低，為608.7 mm(表4)。由于各站點均位于華北平原，氣候條件較為相似，受作物生長和農業管理(灌溉、施肥、作物品種等)等因素的影響，因而造成不同站點間ET數值上的差異[27]。

表4 研究時段內華北平原不同站點農田ET結果對比Tab.4 Comparisons of ET in different sites during study period across North China Plain

就不同階段而言，禹城站小麥季ET的多年平均值明顯高于玉米季，這除了與二者生長季環境條件不同外，也與它們分屬于不同作物類型有關。一般而言，C4植物(玉米)比C3植物(小麥)擁有更高的碳同化速率且氣孔導度較小，使蒸騰耗水較低[26,28]。而農閑期階段地表近似“裸地”的狀況，使該階段ET占全年ET總量的比例不足4%?？紤]到機器學習方法在探究多因子相互關系方面的潛力，本研究利用BRT方法分析了不同作物生長季和農閑期階段ET變化的多環境因子影響。結果表明：Rn是影響麥-玉輪作田不同階段ET變化的重要因子，這是因為無論是土壤蒸發還是植物蒸騰，輻射提供的能量可以改變系統中水熱分配狀況，使植物體內水分和土壤水分發生相變汽化過程，進而影響ET變化[29]。但同時Rn對不同階段ET變化的影響也存在明顯差異，其貢獻率由高到低依次為小麥季(81.4%)、玉米季(52.7%)、農閑期(36.8%)。除Rn外，VPD對作物生長季ET也存在一定的影響。一方面，VPD作為衡量大氣干燥程度的重要指標，在土壤-植物-大氣連續體中，能直接影響氣孔開閉，控制植物蒸騰等生理過程[30]，對農田ET產生重要影響；另一方面，ET不隨VPD增大而無限增大，存在一個閾值范圍。FLETCHER等[31]發現VPD閾值受到植物耐旱程度的影響，一般不超過2.0 kPa，而本研究發現：小麥季和玉米季VPD的這一閾值雖然在FLETCHER等[31]的研究結果范圍之內，但二者并不相同，分別為1.4 kPa和0.7 kPa。SWC作為衡量土壤干旱的重要指標，一般通過影響木質部的導水率，從而控制氣孔導度[32]。雖然SWC和VPD對農田碳-水交換的調節效應可能是獨立的，但部分研究證實VPD對ET的影響可能受到SWC的調節[33]。在地表近似“裸地”的農閑期時段，ET主要以土壤蒸發為主，此時SWC和WS對該時段ET變化的影響不可忽略。因此從模型建模的角度考慮，今后仍需要在作物模型中，針對不同作物類型和不同階段，建立相應的參數化方案，從而有效降低ET模擬結果的不確定性[34]。

由于研究區玉米季屬于多降雨期，因而降雨對ET及各環境因子也可能產生一定的影響。為此，本研究選擇7月下旬至8月上旬兩次典型降雨時段，對比分析降雨時段前后ET與各環境因子的變化差異(表5)。由表5可知，由于降雨期間太陽輻射減少、濕度增大，使Rn、T、VPD較降雨前出現不同程度地下降，而降雨后土壤水分得到補充，使SWC出現增大；對ET而言，其在降雨期間出現下降，但降雨結束后的晴天階段，ET會大幅增加，這與前人的研究結果相似[35-36]。但郭映等[29]基于黃土高原丘陵溝壑區玉米田ET資料，發現8月中上旬白天氣溫較高，若降雨發生在傍晚時段，則對玉米田ET不會產生較大影響；而8月下旬后，由于玉米處于生長晚期，此時氣溫驟降和陰雨天氣會使ET減少[29]?？傮w而言，降雨對玉米季ET影響具有一定的復雜性。一方面，全球變暖背景下降雨的時空異質性日趨增加，且有向極端化發展的傾向[37]。因此，降雨的強度、發生時間、持續期和極值變化可能會對不同生育階段玉米ET的影響存在明顯差異。另一方面，北方地區ET受能量或水分的限制隨降水類型不同而變化明顯[38]，因而厘清降水對玉米季ET的影響還需要考慮種植區域的氣候背景條件。為此，在后續的研究中，將通過加入更多不同地域、不同時間尺度(如逐小時、逐日等)的觀測數據定量分析這一影響。

除環境因子外，作物生長和農業管理等可能也會影響農田ET的變化。例如在作物生長初期，葉面積指數(LAI)較小，Rn主要轉化為顯熱通量；隨著作物生長，LAI逐步增加，光合速率增大，植物蒸騰作用加強，Rn更多地轉化為潛熱通量，從而影響ET變化[39]。對于灌溉而言，它使SWC升高，從而導致土壤蒸發量增加；同時SWC升高也使得植物根系吸收水量增加，進而使得植物蒸騰增加，促使ET升高[40]。但受限于相關數據資料缺乏，目前對這方面的定量分析仍有待于深入探討。以LAI為例，目前長時序LAI變化通常是利用遙感數據進行表征，但是受時空分辨率的影響，遙感LAI數據并不一定能反映地面渦度通量塔附近農田LAI的真實情況[11]，因此為了降低本研究的不確定性，暫未考慮LAI的定量影響。后續將圍繞LAI、灌溉和施肥等因子，盡可能多地獲得相關特征變量的實測數據，并融入已有的分析結果中，以便進一步提高對麥-玉輪作田ET變化的解釋能力。

4 結論

(1)麥-玉輪作田ET具有明顯的雙峰型季節變化特征，通常每年小麥季ET的峰值明顯高于玉米季。就不同階段而言，小麥季和玉米季多年ET平均值分別為408.5 mm和261.2 mm，而農閑期ET多年平均值為25.8 mm。

(2)Rn是不同作物生長季ET季節變化的重要影響因素，特別是對小麥季ET影響的貢獻率達81.4%，遠高于玉米季。此外，VPD對小麥季和玉米季ET也存在一定的影響，且這種影響具有閾值差異。對于農閑期而言，影響該階段ET變化的關鍵環境因子依次為Rn、SWC和WS，貢獻率分別為36.7%、22.5%和17.0%。