基于EFAST的不同生产水平下WOFOST模型参数敏感性分析

兴 安 卓志清 赵云泽 李 勇 黄元仿,2

(1.中国农业大学土地科学与技术学院, 北京 100193; 2.自然资源部农用地质量与监控重点实验室, 北京 100035)

0 引言

目前，作物模型在生产潜力评价、气候变化影响评估、作物生长监测、作物产量预报、农业风险评估、农业管理决策等诸多领域的应用越来越广泛[1-6]，其中模型参数校正和优化仍然是模型本地化应用的关键过程，而敏感性分析是进行参数校正和优化的重要前提[7-9]。敏感性分析能够将不确定性结果分配至不同的模型参数，量化模型参数对输出结果的影响，识别并筛选出模型中的关键参数，从而降低参数优化过程中的工作量和输出结果的不确定性[10]。敏感性分析通常分为局部敏感性分析方法和全局敏感性分析方法[11]，其中，全局敏感性分析方法有EFAST方法[12-13]、Morris方法[14]、Sobol’方法[15]等，均在作物模型敏感性分析中得到了广泛的应用[16-20]。

作物模型中参数的重要性不仅与模型结构有关，还受到其他参数和输入数据的影响。特别是在区域尺度上，气候、土壤以及管理条件的复杂性导致区域尺度上参数敏感性分析结果存在一定的不确定性[21-22]。因此，许多研究针对不同区域年份的气候、土壤、管理条件的差异和作物参数取值范围进行了参数敏感性的影响分析[22-24]。如文献[25]选取不同气候类型下的4个冬小麦代表性站点，对APSIM-Wheat模型品种参数、土壤参数在不同气候类型及不同产量水平下的敏感性差异进行分析，结果表明，不同气候区间产量和生育期蒸散对参数的敏感性分析结果存在差异，不同产量水平间差异不大；文献[26]采用Morris和Sobol’方法，分别计算了欧洲5个国家的主要水稻种植区及代表不同大陆性气候条件的3个年份下的WARM模型参数敏感性，结果表明，模型对作物参数的敏感性在不同气候条件下并不总是一致的；文献[27]采用Morris和EFAST方法对欧洲不同气候条件(包括未来气候变化)下WOFOST模型不同作物参数的敏感性进行了分析，结果表明，敏感性分析受气候条件，尤其是极端气候条件的影响，从而表现出不同的结果，但整体上最敏感参数趋于一致。以上研究均强调了在区域尺度和长时间序列上进行参数敏感性分析的必要性，从结果上看，在较大时空尺度上敏感性分析结果因气候类型的多样性及复杂性而存在较大差异，而在较小尺度相同气候类型下，不同气候条件对参数敏感性的影响还有待进一步探讨。此外，WOFOST模型作为应用最广泛的机理作物模型之一，其参数敏感性分析也是目前的研究热点。然而，大部分研究在进行参数敏感性分析[24,27-29]和参数校验[3-4,30-32]时未考虑蒸散速率修正因子(CFET)和土壤水分消耗作物群数量(DEPNR)等作物水分利用相关品种参数，而这些参数能够表征作物对水分胁迫的响应，这会影响模型模拟精度，尤其是水分限制生产水平下的作物生长模拟。同时WOFOST模型对所有作物生长的描述相同，在同一种植制度、不同作物间的参数敏感性差异分析有助于进一步理解模型运作过程。

本文选取温带季风气候类型的黄淮海平原旱作区不同积温区内的3个站点，基于农业气象试验站点数据以及实地采样数据，采用EFAST方法分析WOFOST模型冬小麦和夏玉米产量在不同生产水平、不同气候条件下的参数敏感性及其排序一致性，并探讨敏感性差异形成的原因，为后续区域尺度模型参数标定和验证提供基础。

1 材料与方法

1.1 研究区概况

黄淮海平原旱作区[33-34]是我国重要旱作作物种植区及商品粮基地，其大部分属温带季风气候，年均温度14～16℃，年均降雨量500～1 200 mm，土壤类型主要为潮土、褐土和砂浆黑土，主要种植粮食作物为冬小麦和夏玉米，冬小麦播种时间一般为10月，次年5月下旬或6月上、中旬成熟，夏玉米常年播种在6月，9月上、中旬成熟。

本研究基于黄淮海旱作区2000—2015年日平均温度大于等于10℃的年均积温数据，以200℃为间隔划分不同积温区[35]，在积温为4 600～4 800℃、4 800～5 000℃、5 200℃以上的积温区中选取黄骅、商丘和驻马店3个农业气象试验站点作为研究地区(图1)，其基本信息和气候特征见表1。

图1 黄淮海平原旱作区分布图

表1 3个代表站点的气候基本特征

1.2 WOFOST模型

WOFOST模型由荷兰瓦赫宁根大学和世界粮食研究中心(CWFS)共同开发研制，以逐日气象数据为驱动，通过土壤、管理和作物参数限制与调整作物生长过程的动态解释性模型，可以模拟潜在、水分限制和养分限制3种水平下的作物生长。WOFOST模型中物候期用无量纲的状态变量——生长阶段(Development stage, DVS)来描述，DVS为0表示出苗期、1表示开花期、2表示成熟期[36]。

模型运行所需的气象数据(逐日最高气温、最低气温、降雨量)、管理数据(播种日期、出苗期、成熟期)从中国气象数据网(http:∥data.cma.cn/site/index.html)获取，逐日辐射、早晨水汽压、地上2 m处平均风速采用FAO提出的Angstrom公式将从中国气象数据网获取的日照时数、平均相对湿度以及平均风速进行转换获取。土壤数据(土壤含水率、田间持水量、萎蔫系数、水分特征曲线、饱和导水率)由旱作区取样数据和农业气象站数据获取，土样取于2017年5—6月，采用重铬酸钾外加热法、环刀法测定土壤有机质含量和容重，采用激光粒度仪测定土壤颗粒组成，并采用Soil Water Characteristics软件[37]计算转换为所需土壤参数。

1.3 研究方法

1.3.1参数选取及取值范围

根据冬小麦和夏玉米不同生长发育特点以及前人关于WOFOST模型作物参数敏感性分析结果[24,27,29]，选取初始、绿叶、同化、同化物转化、呼吸作用、干物质分配、死亡、根系以及水分利用相关参数进行敏感性分析，2种作物参数各50个，其中生育期积温参数未考虑，运行模型时各生育期参数根据农业气象站点数据中的实际生育期计算获取。模型作物参数缺省值基础上浮动±10%并考虑模型参数默认范围作为参数的取值范围，敏感性分析中将所有参数的分布假定为均匀分布，详见表 2。表2中TDWI为初始总干物质质量；LAIEM为出苗时叶面积指数；RGRLAI为叶面积指数最大增长率；SLATB为不同生育期的比叶面积；SPAN为35℃时叶片生命周期；TBASE为叶龄的低温阈值；KDIFTB为不同生育期散射光消光系数；EFFTB为日平均气温为0℃或40℃时单叶同化CO2的光能利用效率；AMAXTB为单叶最大CO2同化速率；TMPFTB为最大同化速率在不同均温条件下的校正因子；TMNFTB为总同化速率低温校正因子；CVL为叶同化物转化效率；CVO为储存器官同化物转化效率；CVR为根同化物转化效率；CVS为茎同化物转化效率；Q10为温度变化10℃时呼吸作用变化率；RML为叶维持呼吸速率；RMO为储存器官的维持呼吸速率；RMR为根的维持呼吸速率；RMS为茎的维持呼吸速率；RFSETB为不同生育期衰老校正因子；FRTB为不同生育期茎干物质的分配系数；FLTB为不同生育期叶干物质的分配系数；FOTB为不同生育期储存器官干物质的分配系数；RDRRTB为不同生育期根相对死亡速率；RDRSTB为不同生育期茎相对死亡速率；CFET为蒸散速率修正因子；DEPNR为土壤水分消耗作物群数量；RDI为初始根深；RRI为根深最大日增加量；RDMCR为最大根深。

表2 WOFOST模型冬小麦-夏玉米参数取值范围

1.3.2敏感性分析方法

EFAST法是SALTELLI等[13]在傅里叶幅度敏感性检验法(Fourier amplitude sensitivity test, FAST)的基础上，结合Sobol’法的优点所提出的一种基于方差分解的全局敏感性分析方法，具有稳定、计算快速等特点。EFAST法原理是选取合适的搜索曲线在多维参数空间中运行，将一组非线性相关的整数频率分配给模型的所有输入参数，并对模型中选取的参数引入一个具有共同独立参数的函数，使模型作为独立参数的周期函数，将多维积分降为一维积分。通过把目标函数转换成傅里叶级数，得到各频率的傅里叶频谱曲线，由频谱曲线可计算由参数xi引起的模型输出方差和总方差，其比值即为该参数的敏感度，具体算法详见文献[13,38-39]，其简单算法公式为

(1)

(2)

(3)

STi=Si+Sij+…+Sij…k

(4)

式中V——模型总方差

Vi——xi输入变化单独引起的模型方差

Vij——xi通过xj作用贡献的耦合方差

Vij…k——xi通过与xj、…、xk相互作用贡献的耦合方差

Si——xi的一阶敏感性指数

Sij——xi的二阶敏感性指数

Sij…k——xi的多阶敏感性指数

STi——xi的总敏感性指数

1.3.3敏感性分析方案

基于农业气象试验站点多年降雨量数据，选取2000年(丰水)、2007年(平水)、2013年(枯水)为研究年份，采用EFAST方法，分别对黄骅、商丘、驻马店等站点冬小麦和夏玉米不同生产水平下的产量进行参数敏感性分析，旨在分析不同气候条件下和不同生产水平下参数敏感性的差异。模型产量的输出结果为总储存器官干物质量(Total dry weight of living storage organs，TWSO)，不同气候条件由不同站点(不同积温区)和不同年份(不同降水年型)下的气候条件来表征，即本研究中共考虑9个不同气候条件下的参数敏感性，不同生产水平考虑潜在和水分限制生产水平。全局敏感性分析试验借助软件Simlab2.2的EFAST模块实现。EFAST方法认为采样个数大于等于65倍的参数个数为有效，本研究采样数取73，总采样数为3 650，能够满足EFAST方法基本要求。模型3年3站点2种作物共需运行3 650×3×3×2=65 700次，本研究全部的模型运算通过编写Python程序实现。

通过对不同气候条件下各参数敏感性指数进行平均，计算出各参数总敏感指数，采用参数敏感性指数大于或等于0.10作为界定主要敏感参数的标准[40]。

1.3.4一致性检验方法

不同生产水平、气候条件的敏感性排序的一致性由TDCC(Top-down concordance coefficient)系数来衡量[41]。

假设有排序矩阵Rij，其中有m个因子，n个观测值。rnm为变量m在所有n个观测值的排序序号。在本研究中，m为不同生产水平或不同气候条件，n为待分析的参数个数。

(5)

计算排序矩阵Rij的savage得分

(6)

并代入TDCC系数计算公式

(7)

式中CTDCC——TDCC系数

CTDCC的显著性P值由统计值T计算，计算式为

T=m(n-1)CTDCC

(8)

T为自由度n-1的χ2分布。一般认为，CTDCC越接近1，且P值小于0.05，表明各变量排序具有显著的一致性。

2 结果分析

2.1 冬小麦-夏玉米参数敏感性

2.1.1冬小麦参数敏感性

WOFOST模型冬小麦产量在2种不同生产水平下的参数敏感性分析结果如图2所示。潜在生产水平下产量主要敏感参数有TBASE、CVO、TMNFTB3、AMAXTB1.50、FLTB0.25、FOTB1.00等，其中TBASE和CVO的全局敏感性指数最高，为0.21和0.17。水分限制生产水平下对产量敏感的冬小麦参数有CFET、CVO、TBASE、KDIFTB0.00、TMNFTB3、SLATB0.00、DEPNR等，其中CFET为最敏感的参数，其全局敏感性指数为0.56。由此可发现，2种生产水平下的主要敏感参数有所区别，潜在生产水平下，冬小麦的生长发育主要受温度和太阳辐射影响，因此温度相关的参数如叶龄低温阈值与总同化速率低温校正因子有较高的敏感性；而水分限制生产水平下，冬小麦所需水分仅由降水提供，当降水不足导致水分亏缺会影响冬小麦的正常生长发育，因此产量对蒸散速率修正因子等作物水分利用相关参数高度敏感。

图2 冬小麦参数敏感性分析结果

2.1.2夏玉米参数敏感性

图3 夏玉米参数敏感性分析结果

夏玉米产量在不同生产水平下的参数敏感性分析结果如图3所示。潜在生产水平下，夏玉米产量主要敏感参数有TMNFTB10、EFFTB40、FLTB0.33、SPAN、CVL、SLATB0.00、AMAXTB1.75以及FLTB0.00等，其中最敏感的参数为TMNFTB10和EFFTB40，其全局敏感性指数分别为0.88和0.54。水分限制生产水平下，夏玉米作物参数敏感性分析结果与潜在生产水平差异不大，除了上述潜在生产水平主要敏感参数以外，RDI和SLATB0.78也表现较高的敏感性，其中TMNFTB10依然是最敏感的参数，其全局敏感性指数为0.69。该结果表明，在潜在和水分限制生产水平下，受温度影响的同化速率相关参数是影响夏玉米产量的主要敏感参数，而水分则没有对夏玉米生长发育起到太大的限制作用。

2.2 冬小麦-夏玉米参数敏感性排序一致性

2.2.1不同生产水平

采用TDCC系数量化冬小麦和夏玉米产量在不同生产水平作物参数敏感性的排序一致性，结果表明冬小麦和夏玉米不同生产水平TDCC系数分别为0.82和0.98，P均小于0.01。可见冬小麦和夏玉米在不同生产水平下参数敏感性排序均有较高的显著一致性，但相较于夏玉米，冬小麦不同生产水平TDCC系数较低，其不同生产水平间敏感性差异较夏玉米大，这与上述冬小麦和夏玉米参数敏感性分析结果相一致。

2.2.2不同气候条件

不同气候条件冬小麦作物参数敏感排序一致性检验结果如表3所示，潜在生产水平下由不同站点和年份组成的不同气候条件的总排序一致性较好，TDCC系数为0.92，不同站点和不同年份上也均表现出较高的一致性，TDCC系数最低达0.91，最高达0.97，P均小于0.01，说明不同气候条件下潜在生产水平参数敏感性排序具有显著的一致性，不同气候条件对参数敏感性的影响较低。而水分限制生产水平下不同气候条件排序一致性较差，TDCC系数仅为0.61，同时在不同站点和不同年份参数敏感排序一致性检验结果中可发现，2000、2007、2013年各年份不同站点间TDCC系数分别为0.84、0.94和0.77，黄骅、商丘、驻马店各站点不同年份间TDCC系数为0.60、0.70和0.78，表明各站点在不同年份间的一致性较差，而同一年份在不同站点间的一致性较好，可见不同气候条件对水分限制生产水平下的参数敏感性有较大的影响，且不同降水年型气候条件差异造成的敏感性差异较大。

表3 不同气候条件下冬小麦参数敏感性的排序一致性

不同气候条件夏玉米参数敏感性的排序一致性检验结果如表 4所示，潜在生产水平下夏玉米参数敏感性TDCC系数为0.98，在不同年份和不同站点间TDCC系数均在0.97以上，P均小于0.01，表明潜在生产水平下夏玉米参数敏感性排序在不同气候条件下具有显著的一致性。而在水分限制生产水平下，其TDCC系数为0.86，与潜在生产水平相比，不同气候条件下夏玉米水分限制生产水平参数敏感性排序一致性相对较差。不同年份和不同站点的一致性检验表明商丘和驻马店站点在3个不同年份和2013年3个不同站点TDCC系数较低，参数敏感性排序一致性较差，可见不同气候条件对水分限制生产水平下的参数敏感性有影响，其中不同降水年型间气候条件差异是造成总体一致性较差的主要原因。

表4 不同气候条件下夏玉米参数敏感性的排序一致性

2.3 不同气候条件下水分限制生产水平参数敏感性差异

对于冬小麦产量而言，2000、2013年所有站点产量对CFET、DEPNR等水分利用相关参数高度敏感，而2007年所有站点产量对CFET、DEPER等参数不敏感，相反TBASE等受温度影响的同化速率相关参数是其主要敏感参数。夏玉米相较于冬小麦，在不同气候条件下的参数敏感性差异不大，在2013年商丘和驻马店站点夏玉米主要敏感参数为SPAN、CFET、DEPNR、RDI等叶片和水分利用相关参数，而其余站点、年份主要敏感的参数均为TMNFTB10、EFFTB40等(图4，图中站点年份编号中，HH表示黄骅，SQ表示商丘，ZMD表示驻马店，后面数字表示年份)。总体而言，水分限制生产水平下不同气候条件对参数敏感性的影响有以下2种表现形式：降水不足、水分胁迫程度较大时(枯水年)，作物产量对CFET、DEPER、RDI等水分利用和根系相关参数敏感；降水较充足、水分胁迫程度较低时(丰水年和平水年)，作物产量表现出对TBASE、TMNFTB10及EFFTB40等受温度影响的同化速率和光能利用效率相关参数敏感，这与潜在和水分限制生产水平下的参数敏感性分析结果相似。此外，2000年虽为丰水年，但因1999年(枯水年)降雨量仅为517.3 mm，冬小麦生育期前期水分胁迫程度较大，因此也会对CFET和DEPER等参数有较大的敏感性。

3 讨论

3.1 参数全局敏感性

上述研究表明冬小麦潜在和水分限制生产水平下产量主要敏感参数有所区别。潜在生产水平下主要对TBASE、TMNFTB3、CVO等参数敏感，其中前2个参数表征冬小麦低温条件下进入越冬期同化速率减慢或停止的过程，CVO表示的是同化物转化为储存器官干物质的效率，这些过程均对作物生长有着重要的影响，文献[27,30,42]等研究中也表明了这些参数的敏感性和重要性。冬小麦在水分限制生产水平下主要敏感参数为CFET、CVO、TBASE、TMNFTB3等，其中敏感性指数最高的CFET参数表征作物面对水分胁迫时的响应[43]。目前大部分研究在进行敏感性分析时对参数CFET的考虑较少，但该参数在本研究中表现出对水分限制生产水平下产量有较大的影响，该参数准确校验可能对黄淮海旱作区作物产量模拟有重要意义。

夏玉米潜在和水分限制生产水平下的产量主要敏感参数相一致，主要为TMNFTB10、EFFTB40、SPAN，当水分充足或水分胁迫较小时，作物产量形成主要由光温影响，因此在低温或高温下产量对光能利用效率和最大同化速率有较高的敏感性。但文献[24,27-28]关于夏玉米参数敏感性的研究表明，除EFFTB40、SPAN 等参数外，夏玉米产量还对CVO、TBASE、FOTB1.10、FLTB0.88、SLATB0.00和SLATB0.5等参数有着较高的敏感性，其中部分参数在本研究中并没有表现出很强的敏感性，可能的原因是参数的取值范围、区域尺度、模拟的气象条件、环境或田间管理措施等因素不同[44]，这也说明了在特定操作环境下运用模型之前进行敏感性分析的重要性[27]。

3.2 参数敏感性一致性检验

冬小麦和夏玉米不同生产水平参数敏感性分析结果表明冬小麦不同生产水平主要敏感参数有所差别，而夏玉米则差异不大，冬小麦不同生产水平一致性相对于夏玉米较差。这主要是由2种作物生育期内降雨量的不同所导致，冬小麦生长期主要为10月至次年6月，在黄淮海旱作地区正是降水较少时期，水分缺乏情况下会提高作物水分利用相关参数的敏感性，而夏玉米生育期为6月至10月，这期间雨水较充足，基本能够满足作物生长发育需求，夏玉米生长发育主要受光温条件限制，因此其不同生产水平的参数敏感排序一致性较高，主要敏感参数间差异不大。

冬小麦和夏玉米参数敏感在不同气候条件下的排序一致性检验结果表明：在不同气候条件下WOFOST模型潜在生产水平下的产量对作物参数的敏感性有较高的一致性和稳定性，而在水分限制生产水平下一致性较差。这可能是由于在黄淮海平原旱作区处于同一气候类型区，光温条件差异并不大，模型对外界光温条件的响应较为一致，因此会导致不同气候条件下潜在生产水平敏感性参数有较高的一致性。相反，研究区内无论年际间或是不同地区降水差异均较大，从而增加了不同气候条件下水分限制生产水平参数敏感性差异。本研究的TDCC系数相对低于文献[17]的研究结果，其主要原因可能在于本研究计算一致性时考虑了全部的参数，而文献[17]选取前9个敏感参数进行一致性检验；而不同气候条件参数敏感性TDCC系数高于文献[25]的研究结果，这主要是由于文献[25]的研究主要是针对大区域尺度，选取的4个站点分属不同的气候类型，气候条件相差较大，因此各气候区间敏感参数TDCC系数较低。同时，本研究中不同生产水平间参数敏感性TDCC系数高于前人的研究结果，其原因可能是由于本研究只考虑了潜在和水分限制生产水平，如果考虑养分限制甚至是实际产量，其影响因素更加复杂，相应的敏感性分析结果可能差异更大。由此可进一步明确作物模型是通过作物参数来表达作物对不同外界条件的不同响应，不同外界条件会导致不同的参数敏感性，考虑的外界条件越复杂，参数敏感性的差异会越大。因此，在特定环境下运用模型时进行敏感性分析是不可缺少的过程，同时，多种生产水平及更为复杂的影响因素下的模型敏感性分析可能是进一步了解模型并提高区域上应用的有效措施之一。

4 结论

(1)WOFOST模型冬小麦2种生产水平的TDCC系数为0.82，主要敏感参数间有所差别，潜在生产水平下的主要敏感参数有叶龄的低温阈值(TBASE)、储存器官同化物转化效率(CVO)、总同化速率在低温3℃时的校正因子(TMNFTB3)等，而冬小麦水分限制生产水平下蒸散速率修正因子(CFET)和土壤水分消耗作物群数量(DEPNR)等水分利用相关参数以及储存器官同化物转化效率(CVO)等参数表现出较大的敏感性。

(2)WOFOST模型夏玉米2种生产水平TDCC系数为0.98，主要敏感参数间差异不大，主要为总同化速率低温10℃时的校正因子(TMNFTB10)、每日温度为40℃时单叶片同化CO2的初始光能利用效率(EFFTB40)、35℃时叶片生命周期(SPAN)等。

(3)相同气候类型下不同气候条件对不同生产水平下的参数敏感性的影响不同。潜在生产水平下冬小麦和夏玉米不同气候条件TDCC系数分别为0.92和0.98，不同气候条件对潜在生产水平下影响较小，而水分限制生产水平下2种作物不同气候条件TDCC系数分别为0.61和0.86，不同气候条件对水分限制生产水平影响较大，其中冬小麦水分限制生产水平下由不同气候条件引起的参数敏感性差异相对较大，这均与气候条件在不同时空上的差异有关。

(4)作物模型是通过作物参数来表达作物对不同外界条件的响应，不同外界条件会导致不同的参数敏感性，外界条件越复杂，参数敏感性差异越大。因此，在特定环境下运用模型时进行敏感性分析是不可缺少的过程。