玉米农田生态系统蒸散发模型参数优化

邱中齐，周琳琳，刘红娟，田强龙，赵子敬，张晓梅，魏国孝*

（1.兰州大学 资源环境学院，兰州 730000;2.会宁县太平店镇人民政府 农业农村综合服务中心，甘肃 白银 730799）

0 引 言

【研究意义】黑河流域是我国西北地区第二大内陆流域，是以水循环过程为纽带的冰川/冻土、绿洲和荒漠为主要特征的多元生态系统。黑河中下游的绿洲是甘肃重要的农业生产基地，需要大量的灌溉用水[1]。蒸散发模型能够反映蒸散发过程进而揭示水能循环过程，对于灌溉制度及高效利用水资源有重要意义。然而由于技术设备落后等原因，并不是所有区域都能根据实测值得到具体的模型参数，需根据条件相似获得该区域的经验参数[2]。因此，如何减小输出参数的不确定性，提高模型模拟的精度是一个重要的研究课题。【研究进展】大量研究者对Shuttleworth-Wallace（SW）蒸散发模型进行了验证，Stannard等[3]利用Penman-Monteith（PM）模型、Shuttleworth-Wallace模型以及Priestley–Taylor（PT）模型模拟了科罗拉南部半干旱牧场的天然植被的蒸散发，发现SW原模型模拟精度最好。Wei等[4]将PM模型、SW原模型以及PT模型应用于黑河流域的玉米农田生态系统，发现SW原模型效果最好。目前，对于能量闭合问题的研究，前人利用感热法、潜热法、可利用能量闭合法、波文比法和平均闭合法，使观测的能量强制闭合，研究结果表明5种闭合方法得到的能量闭合结果都不理想[5]。对于没有实测资料的地域，目前最好的方法就是根据地域的相似性获得该区域的经验参数[2]。【切入点】本文利用黑河流域大满超级站玉米农田生态系统的涡动相关系统数据和气象数据，基于蒸散模型中潜热通量和感热通量，引入能量闭合因子，结合贝叶斯参数优化方法，得到更为合理的蒸散发模型参数，提高蒸散发模型的精度。【拟解决的关键问题】本文主要通过研究2014年大满超级站水热通量，气象观测数据及叶面积指数数据，分别采用SW原模型以及经过引入能量闭合因子的多目标优化模型对大满超级站数据进行模拟及验证，绘制出大满超级站校准期和验证期的SW原模型和优化后模型的感热通量和潜热通量模拟曲线图，分别计算出2种方案对于大满站数据模拟的评价指标，对比2种方案的模拟效果。

1 材料与方法

1.1 研究区概况

黑河流域是我国西北重要的内陆河流域之一，处于河西走廊的中部，是甘蒙西部最大的内陆河流域。黑河全长为821 km，流域面积为14.29万km2，黑河流域内包括3种类型的自然环境单元，有不同的地质地貌及生态环境。其北部与蒙古接壤，东以大黄山与武威盆地相连，西部以黑山与疏勒河流域毗邻。同时，黑河流域处于重要的地理位置，它对多个地区产生影响，并且涉及3个省级行政区，上游位于青海省祁连县，中游位于甘肃的部分市县，下游位于甘肃金塔和内蒙古自治区额济纳旗。流域上游为主要产流区，中下游分别为流域的耗水区和径流消失区。本研究主要选取黑河流域的大满超级站，其地理位置如图1所示。

大满超级站（100°24′37″E，38°51′26″N）地处黑河流域中部（以下简称大满站），位于张掖市，以绿洲和灌溉农业为特征，是一个生活和农业用水量较大的地区。该地区年平均降水量为125 mm，年平均气温为7.20 ℃。年平均潜在蒸发量约为2 290 mm。年平均日照时间为310 h，无霜期148 d。土壤类型以粉质黏壤土为主。研究区为典型的农业灌溉区，当地主要种植玉米和小麦。玉米的生长季为5—9月，行距40 cm，株距30 cm。

1.2 数据来源及处理

本文数据全部来源于“黑河流域生态-水文过程综合遥感试验（HiWATER）水文气象数据观测集”，由数字黑河（http://heihe.tpdc.ac.cn/zh-hans/）提供[6]。包括2014年大满站的水热通量、气象观测数据及叶面积指数数据[7-9]。

1.2.1 涡动（EC）数据

大满站涡动相关仪的架设高度为4.5 m，超声风速仪（CSAT3）与CO2/H2O分析仪（Li7500A）之间的距离是20 cm左右，超声朝向为正北方向。

涡动数据包含潜热通量和感热通量，均为0.5 h时间尺度，数据处理的主要步骤包括：野点值剔除，延迟时间校正，坐标旋转（二次坐标旋转），频率响应修正，超声虚温修正和密度（WPL）修正等[10-11]。对各通量值进行质量评价，主要是大气平稳性和湍流相似性特征（ITC）的检验，对Eddypro软件输出的30 min通量值进行质量控制：①剔除仪器出错时的数据；②剔除降水前后1 h的数据；③剔除原始数据中每30 min内缺失率大于10%的数据；④剔除夜间弱湍流的观测数据（u*小于0.1 m/s），对于缺失程度为2 h以内的数据进行线性插补，对于缺失少于14 d的用平均日变化放大进行插补[12]。

1.2.2 气象观测数据

大满站仪器所测量的数据主要有降水量、气温、风速、大气压强、土壤含水率、土壤热通量、净辐射等数据。数据均为0.5 h尺度，由于天气问题及仪器故障会产生一定的缺失数据，需要进行数据插补，对于缺失程度为2 h以内的数据进行线性插补，对于缺失少于14 d的用平均日变化放大进行插补。

1.2.3 其他数据

叶面积指数（LAI）数据来自2015年黑河生态水文遥感试验：黑河流域1 km分辨率5 d合成叶面积指数（LAI）数据集(http://www.tpdc.ac.cn/zh-hans/data)，模型中需要的日叶面积指数（LAI）采用三次样条函数插值所得[13-14]。

1.3 模型简介

1.3.1 Shuttleworth-Wallace原模型

Shuttleworth-Wallace原模型是将植物冠层和裸露土壤分为2层，将潜热通量分为来自植物的潜热通量和土壤的潜热通量。因此，太阳净辐射被分为植物冠层上方接受的太阳净辐射（Rn）和地面接受的太阳净辐射（Rns）。叶片中储存的水分通过气孔蒸腾到植物冠层时，受到植物群体冠层阻力和冠层边界阻力的影响；而土壤水通过汽化从土壤表面蒸发到冠层高度的过程中，在裸露的土壤表层受到土壤表面阻力影响[15-18]。

1.3.2 模型参数优化方案及方法

本文采取差分进化自适应算法（Differential Evolution Adaptive Metropolis），简称DREAM算法，进行模型参数优化，其核心思想为贝叶斯理论。DREAM算法的思想是以贝叶斯理论为基础，通过构造多条马尔科夫链来估计参数的后验信息[19]。贝叶斯理论参数估计的思想是将参数视为符合某种先验概率分布的随机变量。根据贝叶斯理论，后验概率密度的求解可转化为似然函数的求解。由于对数似然函数计算简单。因此，以对数函数作为似然函数。基于贝叶斯不确定性方法，结合潜热和感热数据，对模型进行参数优化，能量闭合因子α的先验区间基于实测计算得到的能量不闭合率得到。

1.3.3 模型评价指标

为了更好地判断和理解模型的模拟效果，引入传统模型评价指标：决定系数（R2）、线性回归率、均方根误差（RMSE）、一致性指数（IA）、纳什系数（NSE）。

2 结果与分析

2.1 能量闭合分析

能量不闭合是指涡动相关系统测量的潜热通量（LE）和感热通量（H）之和低于可用能量（A），目前常用能量闭合率（Energy Balance Ratio，EBR）或能量不闭合率（I）来衡量能量的闭合情况。

图2为大满站30 min尺度的能量平衡闭合状况，纵坐标为潜热通量和感热通量之和，横坐标为净辐射减去土壤通量，即可用能量。线性拟合回归直线的斜率为能量闭合率EBR。大满站玉米农田生态系统能量闭合率为0.83，即能量不闭合率为0.17，R2为0.89，表明拟合结果的可信度高，说明其涡动系统测量的相关数据可靠，可以用于模型模拟数据对比分析。

图2 大满站能量闭合散点图Fig.2 Energy closed scatter diagram of Daman station

2.2 Shuttleworth-Wallace原模型参数分析

根据相关文献，可以选择适当的模型参考值，具体数值如表1所示。其中gmax根据Kelliher等[20]对不同植被类型的最大气孔导度研究得出，Q50、D50、Kq和Ka根据Leuning等[21]对15个不同植被站点的PM模型计算求得，b1和b2根据Sellers等[16]对土壤阻力的研究得出，进而分析黑河流域玉米生态系统的SW原模型模拟效果。

表1 SW原模型各参数经验值Table 1 Empirical values of parameters of SW model

2.3 模型优化后参数分析

根据上文可知，计算出大满站能量闭合率为0.83，即能量不闭合率为0.17，所以能量闭合因子的先验区间为[0,0.2]。结合前人研究成果，大满站的SW原模型参数先验区间通过文献调查可得，见表2。

对比表2及表3可知，经过DREAM算法优化后，参数不确定性降低，能量闭合因子α也接近于最初计算的0.17。

表2 模型各参数的先验区间Table 2 A priori interval of model parameters

表3 模型参数优化后的最大似然函数值、平均值及95%置信区间Table 3 Maximum likelihood function value, mean value and 95% confidence interval after model parameter optimization

2.4 模型评价

为检验参数优化后的模型对潜热通量和感热通量的模拟性能是否提高，基于大满站的数据，将玉米生长季分为模型校准期（6月3日—8月22日）和模型验证期（8月23日—9月27日），对不同时期模型的潜热通量及感热通量进行综合评价（图3—图5）。由于能量闭合中潜热通量是用LE表示，SW原模型中的潜热通量是用λET表示，二者表示的都是潜热通量，为了方便表达，下文中统一使用λET来表示潜热通量。

图3 校准期模型潜热通量和感热通量模拟曲线Fig.3 Simulation curve of latent heat flux and sensible heat flux of model in calibration period

图4 验证期模型潜热通量和感热通量模拟曲线Fig.4 Simulation curves of latent heat flux and sensible heat flux of the model in the validation period

图5 校准期模型潜热通量和感热通量线性回归拟合Fig.5 Linear regression fitting diagram of latent heat flux and sensible heat flux of model in calibration period

SW原模型和优化方案在校准期和验证期对潜热通量和感热通量的模拟效果不一致。校准期SW原模型模拟潜热通量较实测值偏低；优化方案模拟值与实测值基本一致，表明优化方案较原模型对于潜热通量模拟性能有明显提高。

验证期优化方案对潜热通量模拟效果较好，模拟值与实测值基本一致，而SW原模型值低于实测值，说明优化方案相对于SW原模型性能有较大的提高。但是在感热通量模拟方面，优化方案模拟数据较SW原模型有所提高，但是达不到模型校准期的效果。

由表4—表5和图5—图6可知，SW原模型和优化方案在校准期和验证期对潜热通量和感热通量的模拟评价指标有差异。

表4 校准期模型潜热通量和感热通量模拟评价指标Table 4 Simulation evaluation index of latent heat flux and sensible heat flux of model in calibration period

表5 模型验证期潜热通量和感热通量模拟评价指标Table 5 Simulation and evaluation indexes of latent heat flux and sensible heat flux during model validation

图6 验证期模型模拟潜热通量和感热通量的线性回归拟合Fig.6 Linear regression fitting diagram of latent heat flux and sensible heat flux simulated by the model in the validation period

在模型校准期，在潜热通量方面，优化方案的评价指标都好于SW原模型，SW原模型和优化方案的NSE分别为0.19和0.82，RMSE分别为156.20和74.26，IA分别为0.81和0.95。优化方案相对于SW原模型，RMSE降低54.46%，IA提高17.3%；在感热通量模拟方面，SW原模型和优化方案后模型的NSE分别为-27.07、-2.67，RMSE分别为87.62、68.52；可见，在模型校准期，经过优化方案参数优化的模型性能较原模型有较大的提高（表4，图5）。

在模型验证期，关于潜热通量模拟方面，SW原模型和优化方案的NSE分别为0.23、0.80，RMSE分别为125.76和62.24，IA分别为0.83和0.95。优化方案相对于SW原模型，RMSE降低50.51%，IA提高14.46%。在感热通量模拟方面，SW原模型和优化方案的NSE分别为-5.61、-4.28，RMSE分别为120.47、107.71；可以得出，优化后对于潜热通量和感热通量的模拟都优于原模型，但是在感热通量模拟方面，模型性能提升不是特别的明显（表5，图6）。

3 讨 论

由于全球气温升高，农田灌溉用水将受到影响。如何制定更为合理的农田灌溉制度，提高农田生态系统的水分利用效率是今后应对全球变暖的一个重要研究方向[22]，可靠水文资料在研究中不可或缺[23-25]。水文模型作为一种研究流域内复杂水文现象的有效工具，一直在无资料区域预测的研究中发挥着重要作用。由于在无资料地区没有或缺乏历史观测资料，不能进行模型参数的率定与检验，因此，参数的确定问题是水文模型在无资料地区应用的关键所在。目前，常用的模型参数确定方法包括参数移植、参数估计以及参数的区域回归方法等[2]。但是以这些方式获取的模型参数的参考值存在不确定性，本文以SW原模型为基础，采用差分自适应算法，以贝叶斯理论为核心思想，对SW原模型中参数的文献参考值进行优化，以潜热通量和感热通量的模拟值与实测值进行数据同化，同时引入能量闭合因子，从而得到模型的最优参数，并将参数代入模型中，对模型及其参数进行诊断评估，将经过优化参数得到的模拟结果与SW原模型进行对比分析。结果表明优化方案降低了参数的不确定性，同时较原模型有更好的模拟效果。

SW原模型作为复杂的蒸散发模型，描述蒸散发及其能量过程，包括潜热通量和感热通量的变化过程，2个过程同时发生并且互相影响，如果使用单一目标的实测数据对模型进行同化，会造成参数和模型性能趋于单一目标的数据匹配，导致其他过程出现过度不匹配[22]。本文在SW原模型的数据同化过程及其评估过程中，以多目标优化为目标，结合大满站的实测数据，考虑潜热通量和感热通量对模型的约束，提取实测数据中的可用信息，分析参数的不确定，评估模型模拟的性能。但是SW原模型在模拟潜热通量方面有着详细的描述，对于感热通量只有简单的能量平衡方程，因此基于潜热通量和感热通量的多目标优化SW原模型，也并不能全面提高SW原模型的模拟效果。

优化方案是基于感热通量和潜热通量对SW原模型的参数进行优化，但是优化后的效果却不甚相同。在优化模型校准期，优化方案对感热通量和潜热通量有着很好的模拟效果，与实测数值非常接近；在模型的验证期，优化方案对潜热通量依旧有着很好的模拟效果，但对于感热通量的模拟效果没有明显的提高。原因可能是：①由于校准期和验证期，植物的生长状态及环境条件不一致，而参数会随着植物的生长状态及环境条件而变化，从而导致模拟效果改进不一致；②感热通量数据仅通过简单的能量平衡公式计算而得到，缺乏完善的机理研究；③玉米农田生态系统冠层和地面间的能量交换作用较为复杂，模型中未考虑能量交换作用，会对模型模拟造成一定的影响[26]。

4 结 论

1）在无资料地区，根据区域条件相似通过参数移植或参数估计等方法获取的参数存在较大的不确定性，而经过DREAM算法优化后，参数的不确定性减小。

2）利用多种有效信息，约束平衡模型的参数，能够大大降低模型参数的不确定性；针对模型的能量不闭合问题，定量计算能量不闭合率，将能量闭合因子α引入模型中，不仅能够降低关键参数的不确定性，而且能使参数受到更好的约束。

3）参数优化后的模型与SW原模型同时对2014年大满站玉米农田生态系统的潜热通量和感热通量进行模拟。对比2个方案的模拟效果，可以证实经过DREAM算法优化参数后的模型对潜热通量和感热通量的模拟效果有了显著的提升，验证期模拟潜热通量R2从0.72提升至0.84，RMSE从125.76减少至62.24。

4）优化后模型较SW原模型潜热通量模拟效果好于感热通量的模拟，模型校准期和模型验证期的效果也存在差异。