1958-2017年太湖蒸发量年际变化趋势及主控因子*

荆思佳，肖 薇，王晶苑，郑有飞，王张 圳，胡 诚

(1：浙江省衢州市气象局，衢州 324000) (2：南京信息工程大学大气环境中心，南京 210044) (3：中国科学院地理科学与资源研究所生态系统网络观测与模拟重点实验室，北京100101) (4：无锡太湖学院苏格兰学院，无锡214063)

湖泊是地球系统各圈层相互作用的联结点，是陆地水文循环的重要组成部分. 在气候系统中，湖泊代表潮湿表面，其蒸发仅受大气条件控制，是湖泊水文循环中最直接受气候变化影响的部分，是气候变化的敏锐指标[1]. 《IPCC全球升温1.5℃特别报告》指出，相较于工业化前水平，人类活动造成全球变暖大约1.0℃[2]，全球变暖趋势明显，并且全球温室气体浓度仍在增加，2020年CO2、CH4、N2O浓度分别达到了工业化前水平的149%、262%、123%[3]. 气候变化改变了气温、风速、辐射强迫等与湖泊蒸发密切相关的气象因子，因此了解湖泊蒸发在全球变暖期间的变化趋势及其主控因子是必要的[4-6].

不同时间尺度的湖泊蒸发受不同的气象因子和物理过程控制. 在短时间尺度(小时到日尺度)上，湖泊蒸发主要由湖-气之间的水汽压差和风速协同控制[7-8]，同时天气尺度的天气系统如冷锋过境等，会在短时间内影响湖泊蒸发量[9-10]. 在季节尺度上，月蒸发量更接近净辐射的变化趋势，但是根据湖泊的深浅不同有1～5个月的相位延迟[11-12]. 而关于气候变化背景等更长时间尺度上湖泊蒸发的主控因子研究，现有的观点主要有气候变暖、辐射控制、风速变化控制3种理论. 气候变暖理论认为，气温升高会增加大气饱和水汽压，提高大气持水能力，促进开阔水面蒸发，加速全球水分循环，形成气候变暖与水分循环的正反馈[13]. 但是在过去50年间，全球气候变暖，蒸发皿观测的蒸发量却在不断减少[14-19]. 因此，Roderick等[17]提出太阳辐射的变化(即全球变“暗”或变“亮”) 导致了蒸发皿蒸发的变化，太阳辐射是水面蒸发年际变化的主控因子. 而Johnson等认为是近50年来季风减弱造成风速减小，最终影响蒸发皿蒸发[18]. 同时另有研究表明冰期缩短会引起地表反照率变化，温度升高会引起波文比变化从而改变能量分配的比例，也可能是驱动长时间尺度蒸发变化的主要因素[20].

太湖是我国的第三大淡水湖泊，是典型的富营养化亚热带浅水湖泊，以太湖为中心的太湖流域不仅是我国重要的湿地分布区，其所在的长三角地区也是我国城市化进程最快、经济最发达的区域之一. 与深水湖泊相比，浅水湖泊对气候变化的响应更为敏锐[21]. 气象观测数据表明，近50年来，太湖流域气温、太阳辐射呈先减后增的趋势. 为获得太湖蒸发量，研究者们基于间接估算和直接观测等方法量化太湖蒸发量. 毛锐、沈行毅运用扩散法、池湖蒸发差值法[22]和修正彭曼综合法[21-23]等蒸发模型，校正岸边蒸发池数据代表太湖蒸发；随后水量平衡法[24]被用于估算太湖蒸发；随着观测手段的进步，稳定同位素方法[25-26]被用于间接反算湖面蒸发. 随着涡度相关法等直接观测方法在湖泊上的成功应用，Lee等[27]、Xiao等[28]利用涡度相关法和通量梯度法对在太湖开展了原位连续观测，积累了多年数据，Xiao等[29]基于2011-2017年的涡度相关观测数据，利用能量平衡方法进行了诊断分析，发现辐射控制着太湖蒸发的年际变化，蒸发量的增加主要是由能量输入的增加造成的，减少则是由湖表发射的向上长波辐射减少造成的.

认识更长的历史时期和气候变化情景下湖泊蒸发的变化，则需要借助模型的方法. 研究者们基于CLM4.0-LISSS和E-ε等湖泊模型[30-33]研究了太湖蒸发对气候变化的响应. Hu等[33]借助CLM4.0-LISSS模型对1979-2013年太湖蒸发进行了估算，并与蒸发皿观测反算结果进行对比，两种方法均得出太湖蒸发呈上升趋势，但由逐步回归分析得出的太湖历史蒸发的主控因子不同，蒸发皿反算结果显示气温的升高是太湖历史蒸发的主控因子，模型模拟结果为向下的短波辐射. 随后，郝晓龙等[34]在Hu等[33]的基础上利用HadCEM2-ES气候模式结果驱动CLM4.0-LISSS模型模拟了2010-2100年RCP2.6、RCP4.5及RCP8.5不同温室气体排放情景下太湖蒸发量，得出各种温室气体排放情景下太湖蒸发都呈现增加的趋势，蒸发量的增加速率与辐射强迫呈正相关，风速与水汽压差的乘积可以解释64%的潜热通量变化. 但两者采用的CLM4.0-LISSS模型中使用的粗糙度参数化方案，更适用于海上强风情况，用于太湖时潜热通量会出现异常高值[35]，且模型验证时仅使用了3年的观测数据，无法表明模型可以模拟出蒸发的长期变化. 因此本文选取了2013年涡度相关观测数据与JRA-55再分析资料建立线性校正方程，通过验证2012年校正前后再分析资料与观测数据的一致性来判断校正效果，再利用校正后的JRA-55再分析资料，驱动改进了粗糙度参数化方案后的CLM4.0-LISSS模型，获得1958-2017年太湖湖面蒸发量，模拟结果与2012-2017年涡度相关观测对比验证其可用性，并分析湖面蒸发的变化趋势，寻找太湖实际蒸发年际变化的主控因子，以期为太湖地区水资源管理与气候变化研究提供科学依据.

1 材料与方法

1.1 研究区域概况

太湖(30°55′40″～31°32′58″N，119°52′32″～120°36′10″E，图1)，地处长江中下游，属于亚热带湿润季风气候区，多年平均气温15.97℃，年平均降水量1182 mm[33]，湖表面积约为2338 km2，平均水深1.94 m，是典型的亚热带大型浅水湖泊.

图1 太湖中尺度通量网观测站点位置Fig.1 Locations of the observation sites of the Lake Taihu eddy flux network

1.2 数据来源及处理

本文最终用于分析的太湖蒸发量数据均由观测或模拟所得潜热通量除以汽化潜热计算而来. 本研究利用2013年太湖中尺度通量网[27]避风港站(BFG，31°10′28″N，120°24′01″E)观测数据校正JRA-55再分析资料，并利用2012年的观测资料验证该再分析资料在太湖的适用性；再利用2012-2017年辐射及小气候观测资料驱动CLM4.0-LISSS模型，与观测所得的湖面蒸发量对比，验证该模型在太湖的适用性；最后利用1958-2017年JRA-55再分析数据驱动模型，获得太湖年累积蒸发量，分析太湖蒸发的长期变化趋势，寻找太湖实际蒸发年际变化的主控因子.

1.2.1 观测数据 本文用于校正再分析资料、验证模型适用性的相关观测数据均来自太湖中尺度通量网[27]. 前期已有研究表明，太湖月蒸发量无空间差异[11]，因此本文选择了始建于2011年12月5日、通量观测数据量较为完整且风浪区足够大的BFG站.

小气候系统架设于距湖面8.5 m处，测量空气湿度、气温(Model HMP45D/HMP155A; Vaisala, Inc, Helsinki, Finland)、风速风向(Model 03002; R. M. Young Company, Traverse City, MI, USA)及四分量辐射(Model CNR4; Kipp & Zonen B.V., Delft, the Netherlands)等气象要素. 水温梯度系统由不同深度的水温计(Model 109-L, Campbell Scientific Inc., Logan, UT, USA)构成，分别测量20、50、100、150 cm深处的水温及底泥温度.

涡度相关系统同样架设于距湖面8.5 m高度，由三维超声风速仪(CSAT3, Campbell Scientific Inc., Logan, UT, USA)、开路式红外CO2/H2O分析仪(EC150, Campbell Scientific Inc., Logan, UT, USA)组成，分别用于观测三维超声风速、大气H2O和CO2浓度，数据采样频率为10Hz；动量通量、感热通量和潜热通量数据输出步长为30 min. 通量原始数据先根据仪器工作状态和日志剔除异常值，再经过两次坐标旋转、超声虚温订正、空气密度效应订正等后处理. 对于缺失较少的通量数据(一天内少于5个数据点)进行线性内插；缺失较多的通量数据，优先使用BFG站另外一套EC数据建立线性关系插补，剩余的缺失数据，利用物质传输的相关理论插补[36].

涡度相关的观测结果存在能量平衡不闭合的问题，导致测得感热通量H、潜热通量LE偏低. 本文采用波文比(BRER)法来调整H、LE，也就是认为波文比(H/LE)的观测值是正确的，将剩余的可利用能量按波文比分配给H和LE，从而达到能量闭合. 本文的处理过程中，对3日尺度的潜热通量进行了能量闭合校正. 具体的后处理过程参见文献[11].

1.2.2 JRA-55再分析资料 用于模拟太湖历史蒸发的气象数据源于日本气象厅(Japan Meteorological Agency, JMA)的第二代大气再分析产品JRA-55再分析资料[37]. JRA-55再分析资料数据自1958年开始，在JRA-25的基础上，改进了长波辐射方案，引进了四维变分技术(4D-Var)和变分偏差校正(VarBC)，对平流层低层温度模拟更加精准，这对于本文的模拟非常重要. 本文使用的是1958-2017年JRA-55时间分辨率为3 h的气温、大气压、比湿、向下的长波辐射、向下的短波辐射数据，水平分辨率为1.25°×1.25°. 向上的短波辐射由向下的短波辐射乘以湖表反照率0.08获得.

再分析资料基于陆地站点发展而来，而湖陆下垫面在反照率、比热容、粗糙度等方面的不同导致了湖陆下垫面微气象条件存在差异. 因此本文采用线性回归的方法，利用2013年观测数据建立校正方程，对1958-2017年的JRA-55再分析资料进行简单的线性校正，通过对2012年校正前后数据的准确性验证校正效果.

1.2.3 蒸发量计算 蒸发量E可以由观测或模拟所得潜热通量LE除以汽化潜热获得：

(1)

λ=3.1475×106-2.372×103Ta

(2)

式中，λ为汽化潜热，Ta为大气温度. 年累积蒸发量由3 h潜热通量转换成蒸发量后再求和计算得到.

1.3 湖泊模型介绍

CLM4-LISSS模型(Community Land Model 4.0-Lake, Ice, Snow and Sediment Simulator)中的驱动数据包括气温、气压、湿度、风速、降水、净短波福射和向下的长波福射. 模型将湖体分为10层，利用大气强迫场，通过湖表薄层能量平衡经过4次迭代获得表面温度的唯一解Ts.在每一次迭代中，摩擦速度和空气动力学阻力将重新计算，再利用公式(3)重新计算得到湖-气间感热通量H和潜热通量LE.

湖表薄层的能量平衡方程为：

(3)

(4)

(5)

(6)

湖体各层温度由垂直温度扩散方程计算得到，公式为：

(7)

式中，cw为深度z处的比热；T为湖体不同层的水温；k为热导度；φ为到达深度z处的辐射通量.

本文选择使用CLM4.0版本模拟太湖蒸发. 目前CLM模型已更新至5.0版本，但使用CLM4.0版本并不影响年尺度太湖蒸发的模拟效果. CLM5.0模拟过程中湖泊深度为可变项，即所有的物理过程(除与深度有关的混合过程除外)均可在任意深度下进行；而在CLM4.0版本中，湖泊深度则为定值. 模式中湖泊深度主要用于计算消光系数，根据Subin等[39]的敏感性分析可知，改变湖泊深度对水热通量的模拟结果与改变消光系数获得的结果类似. 而在本研究中，消光系数已通过多组敏感性实验，以年尺度上蒸发量模拟结果最佳为目标，设置为2.5 m-1，因此使用CLM4.0版本对年尺度湖面蒸发数据影响不大.

1.4 统计方法

1.4.1 Mann-Kendall趋势检验 本文采用Mann-Kendall秩次相关检验(M-K)分析1958-2017年太湖蒸发量年际变化趋势[40-43]. 该方法对样本分布不作要求，且最终检验结果受异常值干扰较小，常应用于降水、蒸发、气温等随时间变化的要素的趋势变化分析. M-K趋势检验原假设H0是对于样本量为n的时间序列x，无时间变化趋势，趋势检验的统计量S定义为：

(8)

其中，

(9)

M-K检验利用标准化后的统计量Z值判断时间序列是否有变化趋势及是否显著，变化趋势的大小用Kendall倾斜度β表示. 计算公式为：

(10)

(11)

式中，var(S)为统计量S的方差. 当Z值为正时，时间序列呈上升趋势；反之，时间序列呈下降趋势. |Z|大于1.28、1.96和2.32分别代表变化趋势通过了90%、95%和99%的显著性检验.

1.4.2 多元逐步回归分析 各驱动变量对太湖历史蒸发的贡献率采用多元逐步回归来定量分析. 多元逐步回归方法按因变量的影响大小被依次引入回归方程中，并考察引入此因子所建立的回归方程是否显著. 当某因子的初始P值小于0.05，则此因子将被引入并重新计算P值，若重新计算的P值大于0.1，此因子将被移除. 本文利用标准化后的驱动数据和年累积蒸发量数据(变化范围0～1)建立回归方程，以此判断各气象因子对太湖蒸发的相对贡献[33,44-46]，具体如下：

Y=e+a1x1+a2x2+…aixi

(12)

式中，Y为标准化后的因变量，即年累积蒸发量；xi(i=1，2，…)为标准化后的因子，即各驱动变量；ai为变量xi在回归方程中的系数.

本文分析所用各驱动变量对年累积蒸发量的实际贡献率，由下列公式计算所得：

(13)

式中，Δxi和ΔY分别为xi和Y的变化量；μi为xi对Y的实际贡献率.

1.4.3 评价指标 本文所用评价指标为平均偏差(mean bias error,MBE)、均方根误差(root mean square error,RMSE)、决定系数(R2)、一致性指数(I)，计算公式为：

(14)

(15)

(16)

(17)

2 结果与分析

2.1 驱动数据校正结果

在驱动模型模拟蒸发量之前，需利用实测资料对再分析资料进行校正，本文使用2013年BFG站点小气候及辐射观测数据对JRA-55进行简单的线性校正，再利用2012年的实测资料验证校正结果. 表1展示了各驱动变量与实测数据的校正方程及校正前后的统计参数. 所有变量的校正方程均通过了α= 0.05的显著性检验，线性校正有效地降低了各强迫数据的系统偏差.

表1 太湖地区JRA-55再分析资料驱动数据校正方程及校正前后优度统计*

校正前JRA-55再分析资料与BFG站点实测资料对比如图2所示. 从相关性来看，再分析资料中向下短波辐射、长波辐射、气温、气压、比湿均与BFG观测值均呈现很好的相关性(P<0.01)，平均风速与观测相关性较差(r=0.29)，但仍通过了0.01的显著性检验；从一致性系数I来看，除JRA-55再分析资料的风速及气压的一致性系数较差，其余均在0.9以上；从校正前MBE及RMSE来看(如表1所示)，JRA-55再分析资料高估了太湖地区向下的短波辐射、气温，低估了向下长波辐射、气压、比湿与风速. 总体而言，JRA-55中的再分析资料中向下长波辐射、向下短波辐射、气温、比湿与实测值偏差较小，但气压存在明显的系统偏差，风速信号接近陆地观测，与湖泊风速有很大差异.

校正后，向下短波、向下长波辐射的MBE分别由23.76 和-16.74 W/m2下降至-3.10和-1.08 W/m2；气温由被高估0.030℃校正至MBE仅有0.004℃；比湿的平均偏差由-0.76 g/kg降至-0.15 g/kg；气压被低估的现象明显改善，校正后MBE仅有0.26 hPa；风速被低估的现象也得到明显改善，MBE由-2.66 m/s降低至0.14 m/s. 总体而言，线性校正方法有效地校正了JRA-55再分析资料与湖上观测数据的系统偏差.

图2 2012年校正前JRA-55再分析资料与BFG站点实测资料的对比： (a)向下短波辐射；(b)向下长波辐射；(c)气温；(d)比湿；(e)风速；(f)气压Fig.2 Comparison between JRA-55 reanalysis data without calibration and observation at the BFG site in 2012: (a) downward short-wave radiation; (b) downward long-wave radiation; (c) air temperature; (d) specific humidity; (e) wind speed; (f) air pressure

图3 2012年3日累积蒸发量时间序列 (黑色实线代表BFG站点EC观测值所得 蒸发量(经能量平衡闭合校正))Fig.3 Time series of 3-day accumulated lake evaporation in 2012 (Black line represents evaporation calculated by EC observation at BFG (forcing energy balance closure at 3-daily timescale))

图4 模拟与观测的3日累积蒸发量的1∶1图 (●代表BFG在线观测驱动模型； ○代表再分析资料驱动模型)Fig.4 Comparison between simulation and observation (●, simulation forced by in situ observation; ○, simulation forced by calibrated JRA-55)

2.2 CLM4.0-LISSS模型及再分析资料适用性验证

图3展示了2012年BFG观测蒸发量与CLM4.0-LISSS模拟蒸发量的对比结果，用于验证模型的适用性. 模型能很好地模拟出2012年蒸发量的逐日变化和季节变化趋势. 春季，模拟及观测蒸发量均迅速增加，6月因“梅雨”季太阳辐射的减少而稍有减少，7、8月达到最大值，随后蒸发量逐渐减小，并于冬季达到最小值. 图4也表明，观测与模拟值的变化范围相一致，3日平均的模拟值与观测值建立的线性回归线斜率等于0.86，MBE和RMSE分别为-0.3和2.0 mm. 这证明了CLM4.0-LISSS模型有模拟太湖蒸发的能力，可用于太湖历史蒸发量变化趋势分析研究.

校正后的JRA-55再分析资料被用于驱动模型，以验证再分析资料的适用性. 图4可以看出，相比于使用小气候系统观测数据驱动模型，利用再分析资料驱动模型时，模型表现略有下降，虽然两者拟合线接近，但模拟值更加离散(MBE=0.1 mm，RMSE=5.5 mm). 模型虽然能很好地模拟出太湖蒸发的季节变化趋势(图3)，但在秋冬季节，模型模拟值明显偏高，春夏季节蒸发量的模拟值偏低，这可能与每年只取了一个校正系数有关. 但是，季节性偏差在年尺度上相互抵消，年累积蒸发量的模拟值与观测值相近，因此本文没有再分季节/分月对再分析资料进行校正. 2012-2017年太湖年累积蒸发量的模拟值(观测值)分别为1050.3 mm(1007.0 mm )、1113.1 mm(1129.0 mm)、1005.1 mm(1011.9 mm)、977.5 mm(1036.3 mm)、983.9 mm(1079.5 mm)、1081.7 mm(1123.8 mm)，观测值和模式值的对比也表明再分析资料能模拟出太湖年累计蒸发量的趋势变化(图5). 因此本文认为校正后的JRA-55再分析资料可以驱动模型模拟年际尺度的湖泊蒸发.

2.3 历史蒸发模拟值的年际变化趋势分析

1958-2017年太湖年累积蒸发量模拟值的时间序列如图5所示. 近60年来，太湖年累积蒸发量的平均值为1018.2 mm，其中1958-1962、1965-1966、1968-1970、1972-1977、1980-1990、1999、2002、2014-2016年共31年蒸发量小于平均值，最小值出现于1977年，为870.6 mm，其余年份蒸发量高于60年平均值，最大值出现于2004年，为1120.2 mm. 需要注意的是，1977年累积蒸发量较多年平均值异常偏低约14.49%，这可能是由于1977年6月降水量较多，太阳辐射偏低(125.28 W/m2，远低于平均值155.76 W/m2)，导致1977年全年太阳辐射仅有128.65 W/m2，远低于平均水平(136.99 W/m2)约6.09%，且1977年大气比湿较高(10.7 g/kg)，比平均值(10.3 g/kg)高约6.8%，而敏感性分析结果表明潜热通量对太阳辐射及比湿极其敏感，将1977年太阳辐射增加6%，比湿减少6%后，蒸发量模拟值将增加至1022.4 mm，与多年均值相近.

线性回归和M-K趋势检验被用于分析蒸发量的变化趋势分析. 线性回归分析说明，1958-2017年，太湖蒸发量整体以1.4 mm/a的速率波动上升. M-K趋势检验也显示太湖蒸发量整体呈上升趋势，上升速率为1.4 mm/a，并通过了99%的显著性检验(Z=3.30). 但是在此期间太湖蒸发量明显被分为两个阶段：1958-1977年间，太湖蒸发量呈现了明显的下降趋势(线性回归斜率=-4.0 mm/a， Kendall倾斜度=-3.6 mm/a，Z=-1.89，通过了90%的显著性检验)；1978-2017年间，M-K趋势检验显示，太湖蒸发量以2.3 mm/a的速度增加(Z=2.69，通过了99%的显著性检验). 总体而言，1958-2017年间太湖蒸发量以1977年为界限，呈现先减后增的趋势.

图5 1958-2017年太湖蒸发量模拟值的年际变化(黑色实线代表模型模拟年累积蒸发量； 灰色虚线代表1958-2017蒸发量趋势拟合线；蓝色虚线代表1958-1977年蒸发量线性回归线； 红色虚线代表1978-2017年蒸发量线性回归线；○代表2012-2017年观测年累计蒸发量)Fig.5 Trends of simulated annual lake evaporation over Lake Taihu from 1958 to 2017 (Black solid line represents annual lake evaporation simulated by CLM-LISSS forced by JRA-55; grey dotted line represents linear fitting line for JRA-55 data from 1958 to 2017; blue dotted line represents linear fitting line for JRA-55 data from 1958 to 1977; red dotted line represents linear fitting line for JRA-55 data from 1978 to 2017；○, observed accumulated lake evaporation from 2012 to 2017)

本研究计算所得太湖蒸发量在1978-2017年间以2.3 mm/a的速度增加，这与Hu等[33]利用MERRA再分析资料驱动湖模式计算得到1979-2013年间太湖蒸发量以2.96 mm/a的速率增加，以及蒸发皿蒸发量观测值以2.54 mm/a的速率增加的结论较为一致. Rong等[45]结合Penman-Monteith方程与参考蒸散率计算了2003-2010年间距离太湖640 km的另一浅水湖泊东平湖的蒸发，其蒸发量也呈上升趋势，上升速率为4.55 mm/a. 而与太湖同纬度的青藏高原湖泊色林措蒸发量也在2007-2015年间以4.3 mm/a的速度增加[43].

2.4 驱动数据的年际变化趋势分析

1958-2017年气温(Ta)、气压(Pa)、比湿(qa)、风速(U)、向下长波辐射(L↓)、净短波辐射(S*)等强迫数据的变化趋势如图6(线性回归方法)和表2所示(M-K趋势检验法). 总体来说，1958-2017年间，Ta的增加趋势最为显著，线性回归趋势线斜率为0.015℃/a(P<0.01)，M-K趋势检验显示气温平均每年升高0.014℃，统计量Z值为3.50，通过了置信水平为99%的显著性检验. 气温最低值(15.82℃)出现于1984年，最高值(17.79℃)出现于1998年.L↓在1958-2017年也呈现显著增加的趋势，M-K方法计算得到的L↓增加量为0.053 W/(m2·a)(置信水平为95%)，线性回归计算得到的L↓增长速率为0.048 W/(m2·a)(P<0.05). M-K方法显示S↓以每年0.047 W/m2的趋势波动上升(置信水平为90%)，但线性回归结果未通过趋势检验.S↓变量范围为128.24～143.86 W/m2，其中1977年达到最小值，2007年达到最大值. 研究表明，中国的地面太阳辐射大致以1990年为界，呈现先减少后增多的趋势[47]，这与本文研究结果基本相符合. 其他变量(qa、U、Pa)在1958-2017年间无显著变化趋势. 本研究校正后的JRA-55再分析资料的风速在 1958-2017年间无显著下降趋势，这与Hu等[33]使用的MERRA再分析资料得到的结论相同，但与地面观测所得到的结论存在差异. 近50年来，地表粗糙度的增加(如城市发展、植被增加、耕地变化等)造成了内陆近地面风速持续减弱，观测中海洋和沿海地区的近地层风速大多呈现出上升趋势[48-49]. 这可能是因为JRA-55再分析资料未考虑粗糙度的变化，且同化时使用的地表观测资料与ERA-40一致，仅使用了温度、湿度、降雪资料[50-51]，未包含土地利用变化(如城市化进程等)信息，所以再分析资料中的近地层风速变化趋势弱于气象站观测数据[51]. 但是由于本文研究站点位于湖上，地表粗糙度小且变化不大，因此校正后的再分析资料因未考虑地表粗糙度变化而呈现出不显著变化趋势不影响本研究使用其作为驱动数据.

图6 1958-2017年JRA-55资料气象因子变化趋势(黑色实线代表年平均值；灰色虚线代表1958-2017年 间各气象因子线性回归拟合线；蓝色虚线代表1958-1977年间各气象因子线性回归拟合线； 红色虚线代表1978-2017年间各气象因子线性回归拟合线)Fig.6 Trends of meteorological variables of JRA-55 from 1958 to 2017 (Black lines represent annual mean meteorological variables; grey dotted lines represent linear regression lines for JRA-55 annual mean meteorological variables from 1958 to 2017; blue dotted lines represent linear regression lines for JRA-55 annual mean meteorological variables from 1958 to 1977; red dotted lines represent linear regression lines for JRA-55 annual mean meteorological variables from 1978 to 2017)

表2 1958-2017年间年均气温、气压、比湿、风速、向下长波辐射、 向下短波辐射M-K检验Z值及变化速率

驱动变量的变化趋势也被分为1958-1977年及1978-2017年两个阶段分析(表2). M-K方法显示1958-1977年间，除L↓无明显年际变化趋势外，Ta、qa、U和S↓均表现出明显的下降趋势，Pa呈显著的增加趋势；线性回归方法显示除气温以0.047℃/a的趋势下降外，其余驱动变量无显著变化特征. 1978-2017年间，各气象因子的年际变化趋势与1958-1977年存在差异.Ta、qa、S↓由显著下降趋势转变为显著上升(M-K趋势检验方法)，上升趋势分别为0.032℃/a(Z=4.46，通过了99%的显著性检验)、7.05×10-3kg/kg(Z=2.11，通过了95%的显著性检验)和0.10 W/(m2·a)(Z=-1.40，通过了90%的显著性检验).U和L↓则转变为呈明显的上升趋势，1977年后分别以9.69×10-4m/(s·a)、0.094 W/(m2·a)的速度增长，并分别通过了90%、99%的显著性检验.Pa则呈现了明显的下降趋势. 总体来说，气象因子和辐射能量在1958-1977年间与1978-2017年间呈现了不同的趋势，这对于解释蒸发量的年际变化有重要意义.

2.5 太湖蒸发年际变化主控因子分析

为确定太湖蒸发量年际变化主控因子，我们利用多元逐步回归分析了各要素对蒸发量的实际贡献率. 基于前文趋势分析结果，将研究时段分为1958-1977年和1978-2017年两个时间段分别探讨. 1958-1977年间太湖蒸发与气象因子之间的多元回归方程为E=-0.28qa+0.43L↓+0.71S↓.Ta、Pa、U被排除在方程外，仅保留了qa、L↓、S↓3个变量. 方程的R2为0.764，说明方程是可以解释76.4%太湖蒸发量的变异. 根据公式(13)计算可得，S↓对太湖蒸发变化的贡献比率最大，高达63.51%，其次来自于L↓的贡献(38.10%)，另外qa的实际贡献为-24.48%，位列第三，而Ta、Pa和U对太湖蒸发年际变化贡献很小；由1978-2017年数据建立的回归方程与1958-1977年有所不同，Pa、U、L↓被移出方程，保留了Ta、qa、S↓3个变量，回归方程为E=0.42Ta-0.31qa+0.51S↓，对太湖蒸发年际变化贡献较大的因素为Ta、qa、S↓，三者的实际贡献率分别为44.54%、-39.13%、52.51%；与此同时，1958-2017年的全部数据也被用于建立多元逐步回归方程，Ta、qa、L↓、S*的多元逐步回归系数分别为0.37、-0.39、0.17、0.55.Pa、U被排除在回归方程外. 方程的R2为0.851，说明方程可以解释的湖泊蒸发量的变异率为85.1%. 从贡献率来看，辐射是湖泊蒸发的直接能量来源，S*的增加(贡献率为54.94%)是1958-1977年太湖年累积蒸发量增加的主要原因. 第二大贡献来自qa，qa抵消了部分由辐射减少带来的正面贡献，占43.49%.Ta、L↓的贡献分别为39.59%和18.03%. 这个结果与分段建立方程的结果相吻合. 总得来说，S↓的增加是造成太湖1958-2017年蒸发量增加的主控因子，L↓、Ta、qa也会影响湖泊蒸发年际变化.

3 讨论

目前，关于控制开阔水面蒸发年际变化的主要气象因子有辐射控制、温度、风速变化3种观点. 大部分研究采用蒸发皿测算量代替(蒸发皿观测乘以系数或将其置于湖上)实际水面蒸发，或利用Penman公式及其衍生公式估算，得出蒸发量与同期气温(平均气温/最低气温/气温日较差)相关性最大[18,33,52-58]，但有研究认为蒸发皿蒸发量变化与气温关系不大，因为早有观测证明北半球冬季的总蒸发量大于夏季[59]，且部分地区气温与蒸发量变化呈相反趋势，即存在“蒸发悖论”现象，并用近50年来季风减弱造成风速减小[18-19,46,60-66]，或云层覆盖和气溶胶浓度增加造成太阳辐射(或日照时数)减少[17-19,63,67-70]来解释. 也有研究认为蒸发皿蒸发量只是陆地蒸发量增加的信号(陆地蒸发量使得蒸发皿周围湿度梯度降低，造成蒸发皿蒸发量下降)，只有在周围陆表环境水分供应充足时才可以代表开阔水面蒸发[16]. 近年来，涡度相关方法被用于观测开阔水面蒸发，并逐步发展成为国际通用的通量观测标准方法，结合各类模型模拟过去或未来气候情景下的水面蒸发用于分析蒸发量的年际变化. 对于高纬度结冰水体，大部分研究认为全球变暖导致湖泊冰期缩短，进而引发的反照率变化，是影响水面蒸发年际变化的重要因素[20,71-73]. 一部分研究利用未来气候情景的气象数据模拟蒸发时发现，当辐射输入无明显变化时蒸发量仍存在明显年际变化，因此气温加速了湖泊蒸发年际变化，且由于湖泊表面增温速度慢于空气增温速度，因此波文比发生改变，能量更多分配用于蒸发(即潜热)，造成蒸发逐年增加[20,72,74-75]. 但对于未结冰的湖泊，在不同未来气候情景下，辐射强迫存在差异，蒸发量的增加速率也会随着辐射强迫的增加而增大[34]. 大部分研究通过长年EC观测或湖模式、蒸发模式的敏感性分析、相关性分析等方法得出辐射是开阔水域蒸发的能量来源，其变化是蒸发年际变化的主控因子，辐射变化对蒸发趋势的贡献大于气温变化的贡献[1,29,33,68,76-78]，与本文的结果一致.

4 结论

本研究利用2013年BFG观测数据校正了JRA-55再分析资料并驱动CLM4.0-LISSS模型，模拟了太湖1958-2017年间湖面蒸发量，利用2012年模拟所得蒸发量与避风港站点实测数据验证模型与再分析资料的适用性，利用M-K趋势检验分析了长时间尺度上太湖蒸发量的变化趋势，再利用多元逐步回归方法寻找太湖实际蒸发的年际变化的主控因子. 主要研究成果如下：

1)利用校正后的JRA-55再分析资料驱动模型时，蒸发的模拟存在季节偏差，但季节性偏差在年尺度上相互抵消，2012年(2013年)年累积蒸发量的模拟值与观测值相近，分别为1050.3 mm(1113.1 mm)和1007.0 mm(1129.0 mm). JRA-55再分析资料可用于年际湖泊蒸发的模拟.

2)1958-2017年太湖蒸发分为两个阶段：1958-1977年间，蒸发呈现了明显的下降趋势，下降速率为-3.6 mm/a；1978-2017年间，太湖蒸发量以2.3 mm/a的速度增加.

3)通过多元逐步回归分析可知，向下的短波辐射是太湖1958-2017年间年尺度蒸发的主控因子，向下的长波辐射、气温、比湿也会影响湖泊蒸发年际变化，蒸发量对风速变化和气压变化不敏感.