农田生态系统在全球生态环境、能量平衡和区域气候研究中具有重要地位, 对农田地表蒸散量的研究一直是国际上关注的热点问题^[1]，特别是提高农田蒸散量估算精度不仅对于研究全球气候变化和水资源评价具有重要作用，而且对于农作物需水生产管理、旱情监测预测、水资源有效开发利用等方面亦有十分重要的应用价值^[2-6]。Priestley-Taylor (简称P-T) 模型在诸多蒸散模型中具有较强的机理性，其忽略了湍流输送条件对水分蒸散的影响，所需气象变量较少，在区域水分蒸散模拟研究中得到了广泛应用。Nurit等^[7]在双源地表能量平衡模型中基于P-T模型估算了农作物冠层蒸散量，刘绍民^[8]在农田蒸散试验资料的基础上用P-T模型计算了新疆棉田实际蒸散量，吕厚荃等^[9]基于P-T模型估算了华北地区玉米田实际蒸散量，米娜等^[10]比较分析了Hargreaves，FAO Penman-Monteith模型和P-T模型对玉米农田蒸散量的模拟精度，研究表明，作物系数采用分段单值法时P-T模型模拟结果最优。基于P-T模型修正式，Sumner等^[11]对美国佛罗里达州非灌溉牧草地、Pereira^[12]对巴西湿热气候和美国半干旱气候条件下的野生植被以及Stannard^[13]对美国圣路易斯山谷半干旱牧草地野生植被的蒸散量均获得了较好的模拟效果。然而P-T模型参数α的取值对蒸散量的估算精度具有决定性作用。在下垫面湿润条件下参数α取常数1.26能够取得较好的模拟效果^[14]，但一些学者分析实际观测资料发现，参数α不是常数，具有日变化和季节变化特征^[15-16]，同时许多研究也表明：受研究区地理条件、海陆间平流输送和年降水量变化等影响，参数α变异性显著^{[3, 12-13, 17-19]}，与土壤湿度^{[5, 8-9, 17]}、饱和水汽压差^{[7, 20]}、太阳辐射^[11]、叶面积指数^{[5, 8-9, 11]}、耦合系数^[12]等因子密切相关。因此，应用P-T模型估算蒸散量时，要依据研究区下垫面和气候特征等因素确定参数α合适的数值。

三江平原是我国重要的商品粮生产基地，大豆田为典型雨养旱作农田。近年来，在全球气候异常背景下，该区域生长季内和年际间降水量差异比较显著，降水量直接控制着土壤含水量和大豆田蒸散量，进而影响大豆生长发育和产量。研究表明：土壤水分含量充足时，蒸散量主要受大气蒸发需求的影响，P-T模型参数α随大气蒸发需求的增加而增加^[14]，一旦土壤含水量下降，参数α会因水分蒸散表层阻力的明显增加而减小^[17]。就三江平原气候特征而言，P-T模型对大豆田蒸散量的模拟效果如何、参数α是否需要修正以及修正后模拟效果如何，这些问题目前尚不明确，有待深入探究与评价。本研究以涡度相关系统测量值为标准，通过分析P-T模型对三江平原大豆田蒸散量的模拟效果，明确估算P-T模式对于大豆田蒸散量的适用性，为评价该地区大豆田水分状况、提高水资源利用率以及预测大豆产量提供依据。

中国科学院三江平原沼泽湿地生态试验站 (47°35′ N，133°31′E) 位于东北三江平原腹地黑龙江省同江市辖区，属温带湿润大陆性季风气候，夏季温暖湿润，冬季严寒漫长，冻结期长达7~8月。年平均气温为1.9℃，1月平均气温为-21℃，7月平均气温为22℃，年平均降水量约为600 mm，60%以上降水集中于6—8月。本研究试验点位于旱地试验场内的大豆田，地表水分主要来源于大气降水，在大豆生长期间没有任何灌溉措施，属典型雨养农田。大豆田土壤类型为潜育白浆土，地下水埋深6~7 m。由于白浆土具有难以透水的紧实黏土层，因此，降水量较大时地表有短期积水。观测期间大豆种植品种为绥化14-3，行距为0.65 m，南北陇向。5月初开始观测，5月中、下旬播种，播种前一次性施入肥料 (N, P₂O₅和K₂O施肥量分别为40 kg·hm^-1, 90 kg·hm^-1和15 kg·hm^-1)，5月下旬6月上旬出苗，9月中、下旬至10月初成熟收割，10月下旬结束观测。观测期间大豆发育期状况如表 1所示。

表 1

表 1 大豆发育期状况 Table 1 Developmental stages of soybean 发育阶段 开始日期 2005年 2006年 2007年 播种前 05-08 05-01 05-04 播种期 05-25 05-17 05-20 出苗期 06-07 05-26 06-10 开花期 07-03 07-05 07-12 结荚期 07-28 08-10 07-30 鼓粒期 08-17 08-25 08-10 成熟期 09-12 09-15 09-25 收割期 09-28 10-04 10-05 收割后 10-23 10-23 10-23   注：播种前日期为涡度相关系统观测开始日期，收割后日期为涡度相关系统观测结束日期。

表 1 大豆发育期状况 Table 1 Developmental stages of soybean

2005—2007年5—10月利用布置在本研究试验点内的涡度相关系统对大豆田蒸散量进行连续观测。涡度相关系统主要由LI-7500开路式红外气体分析仪和CSAT3三维超声风速仪组成，分别用于测定CO₂/H₂O密度和三维风速，采样频率为10 Hz, 实时数据和在线计算的30 min平均通量均存储于PC卡内。同时还利用小气候观测系统对主要环境因子进行观测，观测仪器主要包括CNR1净辐射传感器、HMP45C空气温湿度传感器和52203降水量传感器。此外，每隔10 d采用CI-203型叶面积仪测定大豆叶面积指数，用以反映植被发育的季节变化。

涡度相关技术作为测定地气间水热通量的标准方法在全球得到广泛应用，目前存在的主要问题是在涡度相关系统通量观测中能量不闭合现象较为普遍，国际通量网 (FLUXNET) 平均能量不闭合度约为20%^[21]，能量平衡闭合程度已被作为评价涡度相关数据可靠性的方法，FLUXNET许多站点将能量平衡闭合状况分析作为一种标准程序用于通量数据的质量评价。该研究中由于缺乏必要的辅助观测数据，无法进行能量平衡闭合状况分析，但现有研究表明：作物生长季能量闭合状况好于非生长季^[21-23]，且植株低矮的农作物能量闭合程度明显较好^[23-24]。该研究主要观测了大豆发育期间地表蒸散量，且大豆发育期间最大高度约为75 cm，赵晓松^[25]对该研究中的涡度相关数据进行了大气湍流统计特征分析、谱分析、通量贡献区评价、倾斜校正和超声湿度校正，结果表明：总体而言，该观测所获得的通量数据可信度较高，能够代表观测区域通量大小。该研究中，大豆田观测点下垫面平坦均质，有足够的通量贡献区，因此，主要考虑夜间通量低估带来的系统误差。目前多数研究者利用夜间湍流通量与表征湍流强度大小的摩擦风速 (u^*) 的关系剔除夜间弱湍流交换下的通量值^[26]，该研究中剔除标准为u^* < 0.1 m·s^-1，被剔除的数据按缺测数据处理。此外，在涡度相关系统实际观测中，由于传感器故障、降水、感应器维护和断电等原因，不可避免地会出现缺测。2005—2007年观测期间数据缺测率 (包括被剔除的异常数据) 分别为16.7%，31.9%和21.4%，明显低于FLUXNET 35%的平均值^[27]。缺测通量数据依据下述方法进行插值：小于2 h的数据空缺，采用线性外推法插补；不小于2 h的数据空缺，采用查表法插补，查表法是生态系统能量通量标准插补方案的首选方案^[28]。

Priestley和Taylor以平衡蒸发为基础，通过引进参数α, 导出最小平流条件下蒸发力的估算公式^[4]，即

(1)

式 (1) 中，E₀为参考作物蒸散量 (单位：mm·d^-1)，Δ为饱和水汽压曲线斜率 (单位：kPa·℃^-1)，γ为干湿球湿度方程常数 (取值0.066 kPa·℃^-1)，R_n为日总净辐射，G为土壤热通量，R_n和G单位需换算为mm·d^-1, α为量纲为1的经验系数，当取值为1.26时，该公式估算值为湿润条件下的参考作物蒸散量。依据参考作物蒸散量估算实际蒸散量需考虑植被发育和地表水分供给条件的影响。为了简化计算过程，本研究利用P-T修正式估算实际蒸散量E_a，其表达形式为

(2)

式 (2) 主要考虑了非理想条件下的蒸散量，参数α是一个关于环境因子x₁，x₂，…, x_i的经验函数，要依据研究区域的实际环境条件确定^[11]。日尺度上计算蒸散量时土壤热通量G与净辐射R_n相比很小，通常可以忽略^[29]。

本研究以涡度相关系统测量值为标准，通过分析修正参数α前后P-T模型对大豆田蒸散量的模拟效果，明确该模型在三江平原大豆田的适用性。首先评估参数α取常规值1.26时P-T模型对2006年大豆田的模拟效果，然后依据2005年和2007年实测值反推出参数，利用SPSS19.0软件分析其与环境因子的关系，以2006年实测值对P-T修正式进行验证。本研究中评价模型模拟效果的方法除了对实测值与模拟值进行线性回归分析外，还引入平均偏差、均方根误差和模拟效率，多角度统计分析模拟值的准确性^[30]。

图 1为P-T模型参数α取常规值1.26时，依据小气候观测资料和式 (1) 估算大豆田蒸散量，并与涡度相关系统实测值进行比较得出的不同时期模拟结果。由图 1a可见，大豆出苗前模拟值明显大于实测值，且散点分布相当分散，平均偏差为1.65 mm·d^-1，均方根误差为2.10 mm·d^-1。由图 1b可见，大豆生长期内 (从出苗至成熟期) 散点绝大多数位于1:1线上侧，偏差明显减小，平均偏差和均方根误差分别为1.22 mm·d^-1和1.49 mm·d^-1。由图 1c可见，大豆收割后模拟值明显小于实测值，散点分布也较分散，平均偏差和均方根误差分别为-0.74 mm·d^-1和0.86 mm·d^-1。图 1给出了蒸散量实测值与模拟值之间的线性回归方程，显然在大豆生长期和收割后回归方程可靠性较高，其决定系数R²较大，且达到0.01显著性水平 (p < 0.01)，说明模拟值变化特征与实测值一致，但不同时期的模拟效率均为负值 (大豆出苗前、生长期和收割后模拟效率分别为-3.00，-0.90，-1.43)，表明模拟值与实测值差距较大，P-T模型不能用于估算大豆田蒸散量。

图 1

图 1. 基于P-T模型参数α常规值的大豆田不同阶段蒸散量模拟值与实测值比较 (a) 出苗前，(b) 生长期，(c) 收割后 Fig 1. Comparison of observed evapotranspirations and those simulated by P-T model based on conventional values of α during different period (a) before emergence, (b) growing season, (c) after harvest

由图 1可见，P-T模型估算值有大有小，在大豆出苗前明显高于实测值，在大豆收割后又明显低于实测值，在大豆生长期间系统偏大。偏差的这种变化可能与参数α取常数有关。依据2006年实测资料反推出参数α的季节变化趋势如图 2所示。大豆出苗前和生长期参数α平均值分别为0.76(±0.42) 和0.86(±0.27)，而在大豆收割后参数α明显增大，平均值为2.20(±0.74)。Priestley等^[14]推荐的参数α(平均值为1.26) 适用于湿润下垫面，而三江平原大豆田水分主要来源于大气降水，降水量的季节和年际变化使地表水分条件变化很大，从而使参数α变异性较大，具有很强的经验性，因此，P-T模型参数α采用常规值1.26必然会产生较大误差，有必要依据三江平原大豆田实际条件修正参数α。

图 2

图 2. 不同时期大豆田P-T模型参数α变化特征 (a) 出苗前，(b) 生长期，(c) 收割后 Fig 2. Variation characteristics parameter α of P-T model over soybean field during different period (a) before emergence, (b) growing season, (c) after harvest

利用2005年和2007年5—10月涡度相关系统蒸散量实测值和同期观测的小气候资料，依据式 (2) 反推计算得到各日参数α，分析其与生物和气候因子的相关性。作物生长期内蒸腾作用是水分蒸散的重要组成部分，在植被生长旺盛期尤为明显^[31]，因此大豆叶面积指数是农田蒸散量的一个重要影响因子。如图 3所示，参数α随叶面积指数的增加呈增加趋势，二者呈对数显著正相关关系，表明叶面积指数是影响参数α的一个重要因子。这与Burba等^[5]对美国俄克拉荷马州冬小麦、刘绍明^[8]对中国新疆棉田P-T模型α参数的研究结论一致。饱和水汽压差是水汽从蒸发面输送到周围大气的直接驱动力，但饱和水汽压差对参数α的影响方式与饱和水汽压差变化范围有关。如图 3b所示，参数α随饱和水汽压差的增加表现出先增加、后降低的变化趋势。利用二者建立的一元二次回归方程得出转折点为饱和水汽压差D=5.05 hPa，当D < 5.05 hPa时，参数α与D呈正相关关系；当D> 5.05 hPa时，参数α与D呈负相关关系。

图 3

图 3. 大豆生长期内叶面积指数 (a) 和饱和水汽压差 (b) 对P-T模型参数α的影响 Fig 3. Effects of LAI (a) and D(b) on parameter α of P-T model during growing season of soybean

进一步分析饱和水汽压差D不同变化范围内二者的关系，结果如图 4所示。当D < 5.05 hPa时，参数α与D呈幂函数正相关关系 (幂函数的指数为正值)；当D>5.05 hPa时，参数α与D呈幂函数负相关关系 (幂函数指数为负值)。通常情况下，饱和水汽压差D增加能够促进水分输送，使蒸散量增加，因此，在其他条件相同的情况下，特别是地表水分供给充足的条件下，参数α与D呈正相关关系^[7]。但地表水分条件是控制蒸散的首要因子，三江平原大豆田属于雨养农田，农田土壤水分仅来源于自然降水，久旱少雨使土壤含水量下降，地表蒸发和作物蒸腾作用逐渐减弱，导致近地层水汽含量下降而气温升高，干旱条件下蒸散量与近地表温度呈负相关关系^[32]，进而表现为参数α与D呈负相关关系^[20]。由此可见，久旱少雨时饱和水汽压差D的大小在一定程度上反映了地表的干旱程度。

图 4

图 4. 不同饱和水汽压差 (D) 对P-T模型参数α的影响 (a)D＜5.05 hPa，(b)D＞5.05 hPa Fig 4. Effects of D on parameter α of P-T model under different conditions (a)D < 5.05 hPa, (b)D > 5.05 hPa

为确定P-T模型修正式的表达式，以参数α为因变量，以叶面积指数L与饱和水汽压差D为自变量，依据图 3和图 4所揭示的函数关系，确定函数α的非线性回归方程表达式为α=alnL+bD^c。a, b和c为方程参数，大小反映了叶面积指数L与饱和水汽压差D对参数α影响程度，其值越大，表明影响程度越大，且参数c的正负决定了参数α与D是正相关还是负相关。利用SPSS19.0统计软件进行迭代计算，以确定各参数值大小。结果如表 2所示，非线性回归方程中参数α与叶面积指数L与饱和水汽压差D的关系，与图 3和图 4所揭示的规律完全一致，但回归方程参数a和b的标准差均较大，说明此估算值的显著性水平不显著。相比之下，参数c的标准差均较小，说明此估算值的显著性水平较显著。非线性回归方程的决定系数R²较高，叶面积指数L与饱和水汽压差D能够解释参数α变化的78.6%和67.9%，表明方程的拟合度较高。

表 2

表 2 大豆生长期P-T模型参数α的非线性统计模型 Table 2 Non-linear statistical models of parameter α of P-T during soybean growing season 条件 参数a(±标准差) 参数b(±标准差) 参数c(±标准差) R2 样本量 D < 5.05 hPa 0.083(±0.087) 0.195 (±0.105) 1.294 (±0.376) 0.786 80 D > 5.05 hPa 0.108 (±0.128) 5.246 (±6.102) -0.953 (±0.596) 0.679 150

表 2 大豆生长期P-T模型参数α的非线性统计模型 Table 2 Non-linear statistical models of parameter α of P-T during soybean growing season

大豆出苗前和收割后，气象因子是控制水分蒸散的主要因子，为进一步定量描述参数α与环境因子的关系，在剔除模型中已包含的变量后，以参数α和环境因子分别为因变量和自变量，利用SPSS19.0软件采用多元逐步回归方法建立多元一次统计模型，入选模型的影响因子、回归系数和显著性检验结果如表 3所示。大豆出苗前，参数α主要受太阳辐射R_a与饱和水汽压差D控制，参数α与R_a呈正相关关系，而与D呈负相关关系。饱和水汽压差D标准回归系数的绝对值明显大于太阳辐射R_a的标准回归系数，说明D对参数α的影响作用要强于R_a的影响。参数α与R_a和D建立的回归方程能够解释参数α变化的39.1%。大豆收割后参数α主要受风速W_s控制，W_s能够解释参数α变化的39.2%。D和W_s是维持地表水分蒸散的必要条件，二者增加均能够促进水分蒸散，因此与P-T模式参数α呈正相关关系。三江平原地区属春旱多发区，大豆出苗前D因子回归系数为负值，反映了干旱对水分蒸散的抑制作用。太阳辐射是地表能量的直接来源，因此, 参数α与R_a呈正相关关系，Sumner等^[11]针对牧草地蒸散量建立的参数α经验函数也证明了这个结论。

表 3

表 3 大豆出苗前和收割后P-T模型参数α的统计模型 Table 3 Statistical models of parameter α of P-T model before soybean emergence and after soybean harvest 发育阶段 常数项及回归因子 回归系数 (±标准差) 标准回归系数 R2 出苗前 C0 0.965*(±0.137) 0.391 Ra 0.028*(±0.008) 0.450 D -0.100*(±0.016) -0.786 收割后 C0 0.278 (±0.494) 0.392 Ws 0.550*(±0.134) 0.626   注：*表示达到0.001显著性水平。

表 3 大豆出苗前和收割后P-T模型参数α的统计模型 Table 3 Statistical models of parameter α of P-T model before soybean emergence and after soybean harvest

依据修正后的α函数代入式 (2) 对P-T模型修正式进行验证，结果如图 5所示。由图 5a可见，大豆出苗前P-T修正式模拟精度明显提高，散点较均匀地分布于1:1线两侧，平均偏差和均方根误差分别为0.06 mm· d^-1和0.60 mm· d^-1，较参数修正前分别下降了96.4%和71.4%，明显接近于理想值 (0 mm·d^-1)；模拟效率由负值转变为正值 (0.57)，接近于理想值1，表明模型模拟精度较高。由图 5b可见，大豆生长期内P-T修正式模拟精度也明显提高，散点较均匀地分布于1:1线两侧，尽管散点离散程度相对偏大，但平均偏差和均方根误差较参数修正前分别下降87.7%和38.3%，达到0.15 mm·d^-1和0.92 mm·d^-1，模拟效率由负值转变为正值 (0.28)，表明模型能够用于估算大豆田生长期蒸散量。由图 5c可见，大豆收割后P-T修正式模拟精度也明显提高，尽管大部分散点分布于1:1线下侧，但散点离散程度明显减小，平均偏差和均方根误差分别为-0.21 mm·d^-1和0.41 mm·d^-1，较参数修正前分别下降了71.6%和52.3%；模拟效率由负值转变为正值 (0.42)，表明模型能够用于估算大豆收割后的地表蒸散量。蒸散量实测值与模拟值之间线性回归方程的截距、斜率和决定系数R²如图 5中公式所示，其数值大小较参数α修正前明显接近于其理想值 (分别为0，1和1)，这也表明，P-T修正式估算精度明显提高。P-T模型参数α是关于环境因子的经验函数，因此，依据研究区域环境条件确定参数α亦能提高模型估算精度^[33]，其他学者修正P-T模型参数α后也得到了相同的结论。Sumner^[11]拟合建立了关于太阳辐射和叶面积指数的参数α函数后，很好地估算了牧草地半小时蒸散量；Pereira^[12]研究了湿热气候和半干旱气候条件下P-T模型的模拟精度，结果表明:P-T模型参数α与耦合系数呈倒数关系，据此修正参数α后P-T模型的模拟值与蒸渗仪测量值具有很好的一致性。Stannard^[13]利用叶面积指数和土壤蒸发量修正参数α后，P-T模型很好地模拟了半干旱牧草地野生植被的蒸散量。刘绍民^[8]和吕厚荃等^[9]分别基于P-T模型，结合作物叶面积指数和土壤湿度, 建立了适用于干旱地区棉田和华北地区玉米田的实际蒸散量计算模式。

图 5

图 5. 基于P-T模型参数α修正值的大豆田不同阶段蒸散量模拟值与实测值比较 (a) 出苗前，(b) 生长期，(c) 收割后 Fig 5. Comparison of observed evapotranspirations and those simulated by P-T model based on modified α value during different periods (a) before emergence, (b) growing season, (c) after harvest

整体而言，参数α修正后，P-T模型对模拟精度均明显提高，能够用于估算雨养大豆田不同阶段的蒸散量。但在大豆收割后模拟值明显小于实测值，这可能与大豆收割后农田翻耕、气温显著下降有关。另外，本研究在分析P-T修正式对大豆田蒸散量模拟效果时，仅用1年数据进行验证，其结果缺乏普遍性，今后将利用三江平原其他区域、其他年份的大豆田蒸散量实测资料进一步评估其模拟效果。

本文以2005—2007年实测蒸散量和气象因子资料，分析了P-T模型对雨养大豆田蒸散量的模拟效果，主要结论如下：

1) P-T模型参数α采用常规值1.26模拟效果较差，在大豆出苗前和生长期模拟值明显大于实测值，而在大豆收割后模拟值明显小于实测值，不能用于估算大豆田蒸散量。

2) P-T模型参数α与环境因子密切相关。大豆生长期内参数α与叶面积指数呈对数正相关关系；当饱和水汽压差D较大时参数α与D呈幂函数负相关关系，而当饱和水汽压差D较小时参数α与D呈幂函数正相关关系；大豆出苗前参数α与太阳辐射呈正相关关系，与D呈负相关关系；大豆收割后α与风速呈显著正相关关系。

3) 依据回归分析结果修正参数α后，P-T模型模拟精度明显提高，对不同时期大豆田蒸散量均取得了较好的模拟效果，表明参数α修正后P-T模型才可用于估算三江平原雨养大豆田蒸散量。