基于增强回归树的麦-玉轮作农田蒸散影响因素分析

李 成 李兆哲 王让会 丁奠元 徐 扬

(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]。总体而言，降雨对玉米季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%。