荊思佳，肖 薇，王晶苑，鄭有飛，王張 圳，胡 誠

(1：浙江省衢州市氣象局，衢州 324000) (2：南京信息工程大學(xué)大氣環(huán)境中心，南京 210044) (3：中國科學(xué)院地理科學(xué)與資源研究所生態(tài)系統(tǒng)網(wǎng)絡(luò)觀測與模擬重點(diǎn)實(shí)驗(yàn)室，北京100101) (4：無錫太湖學(xué)院蘇格蘭學(xué)院，無錫214063)

湖泊是地球系統(tǒng)各圈層相互作用的聯(lián)結(jié)點(diǎn)，是陸地水文循環(huán)的重要組成部分. 在氣候系統(tǒng)中，湖泊代表潮濕表面，其蒸發(fā)僅受大氣條件控制，是湖泊水文循環(huán)中最直接受氣候變化影響的部分，是氣候變化的敏銳指標(biāo)[1]. 《IPCC全球升溫1.5℃特別報(bào)告》指出，相較于工業(yè)化前水平，人類活動(dòng)造成全球變暖大約1.0℃[2]，全球變暖趨勢明顯，并且全球溫室氣體濃度仍在增加，2020年CO2、CH4、N2O濃度分別達(dá)到了工業(yè)化前水平的149%、262%、123%[3]. 氣候變化改變了氣溫、風(fēng)速、輻射強(qiáng)迫等與湖泊蒸發(fā)密切相關(guān)的氣象因子，因此了解湖泊蒸發(fā)在全球變暖期間的變化趨勢及其主控因子是必要的[4-6].

不同時(shí)間尺度的湖泊蒸發(fā)受不同的氣象因子和物理過程控制. 在短時(shí)間尺度(小時(shí)到日尺度)上，湖泊蒸發(fā)主要由湖-氣之間的水汽壓差和風(fēng)速協(xié)同控制[7-8]，同時(shí)天氣尺度的天氣系統(tǒng)如冷鋒過境等，會(huì)在短時(shí)間內(nèi)影響湖泊蒸發(fā)量[9-10]. 在季節(jié)尺度上，月蒸發(fā)量更接近凈輻射的變化趨勢，但是根據(jù)湖泊的深淺不同有1～5個(gè)月的相位延遲[11-12]. 而關(guān)于氣候變化背景等更長時(shí)間尺度上湖泊蒸發(fā)的主控因子研究，現(xiàn)有的觀點(diǎn)主要有氣候變暖、輻射控制、風(fēng)速變化控制3種理論. 氣候變暖理論認(rèn)為，氣溫升高會(huì)增加大氣飽和水汽壓，提高大氣持水能力，促進(jìn)開闊水面蒸發(fā)，加速全球水分循環(huán)，形成氣候變暖與水分循環(huán)的正反饋[13]. 但是在過去50年間，全球氣候變暖，蒸發(fā)皿觀測的蒸發(fā)量卻在不斷減少[14-19]. 因此，Roderick等[17]提出太陽輻射的變化(即全球變“暗”或變“亮”) 導(dǎo)致了蒸發(fā)皿蒸發(fā)的變化，太陽輻射是水面蒸發(fā)年際變化的主控因子. 而Johnson等認(rèn)為是近50年來季風(fēng)減弱造成風(fēng)速減小，最終影響蒸發(fā)皿蒸發(fā)[18]. 同時(shí)另有研究表明冰期縮短會(huì)引起地表反照率變化，溫度升高會(huì)引起波文比變化從而改變能量分配的比例，也可能是驅(qū)動(dòng)長時(shí)間尺度蒸發(fā)變化的主要因素[20].

太湖是我國的第三大淡水湖泊，是典型的富營養(yǎng)化亞熱帶淺水湖泊，以太湖為中心的太湖流域不僅是我國重要的濕地分布區(qū)，其所在的長三角地區(qū)也是我國城市化進(jìn)程最快、經(jīng)濟(jì)最發(fā)達(dá)的區(qū)域之一. 與深水湖泊相比，淺水湖泊對氣候變化的響應(yīng)更為敏銳[21]. 氣象觀測數(shù)據(jù)表明，近50年來，太湖流域氣溫、太陽輻射呈先減后增的趨勢. 為獲得太湖蒸發(fā)量，研究者們基于間接估算和直接觀測等方法量化太湖蒸發(fā)量. 毛銳、沈行毅運(yùn)用擴(kuò)散法、池湖蒸發(fā)差值法[22]和修正彭曼綜合法[21-23]等蒸發(fā)模型，校正岸邊蒸發(fā)池?cái)?shù)據(jù)代表太湖蒸發(fā)；隨后水量平衡法[24]被用于估算太湖蒸發(fā)；隨著觀測手段的進(jìn)步，穩(wěn)定同位素方法[25-26]被用于間接反算湖面蒸發(fā). 隨著渦度相關(guān)法等直接觀測方法在湖泊上的成功應(yīng)用，Lee等[27]、Xiao等[28]利用渦度相關(guān)法和通量梯度法對在太湖開展了原位連續(xù)觀測，積累了多年數(shù)據(jù)，Xiao等[29]基于2011-2017年的渦度相關(guān)觀測數(shù)據(jù)，利用能量平衡方法進(jìn)行了診斷分析，發(fā)現(xiàn)輻射控制著太湖蒸發(fā)的年際變化，蒸發(fā)量的增加主要是由能量輸入的增加造成的，減少則是由湖表發(fā)射的向上長波輻射減少造成的.

認(rèn)識(shí)更長的歷史時(shí)期和氣候變化情景下湖泊蒸發(fā)的變化，則需要借助模型的方法. 研究者們基于CLM4.0-LISSS和E-ε等湖泊模型[30-33]研究了太湖蒸發(fā)對氣候變化的響應(yīng). Hu等[33]借助CLM4.0-LISSS模型對1979-2013年太湖蒸發(fā)進(jìn)行了估算，并與蒸發(fā)皿觀測反算結(jié)果進(jìn)行對比，兩種方法均得出太湖蒸發(fā)呈上升趨勢，但由逐步回歸分析得出的太湖歷史蒸發(fā)的主控因子不同，蒸發(fā)皿反算結(jié)果顯示氣溫的升高是太湖歷史蒸發(fā)的主控因子，模型模擬結(jié)果為向下的短波輻射. 隨后，郝曉龍等[34]在Hu等[33]的基礎(chǔ)上利用HadCEM2-ES氣候模式結(jié)果驅(qū)動(dòng)CLM4.0-LISSS模型模擬了2010-2100年RCP2.6、RCP4.5及RCP8.5不同溫室氣體排放情景下太湖蒸發(fā)量，得出各種溫室氣體排放情景下太湖蒸發(fā)都呈現(xiàn)增加的趨勢，蒸發(fā)量的增加速率與輻射強(qiáng)迫呈正相關(guān)，風(fēng)速與水汽壓差的乘積可以解釋64%的潛熱通量變化. 但兩者采用的CLM4.0-LISSS模型中使用的粗糙度參數(shù)化方案，更適用于海上強(qiáng)風(fēng)情況，用于太湖時(shí)潛熱通量會(huì)出現(xiàn)異常高值[35]，且模型驗(yàn)證時(shí)僅使用了3年的觀測數(shù)據(jù)，無法表明模型可以模擬出蒸發(fā)的長期變化. 因此本文選取了2013年渦度相關(guān)觀測數(shù)據(jù)與JRA-55再分析資料建立線性校正方程，通過驗(yàn)證2012年校正前后再分析資料與觀測數(shù)據(jù)的一致性來判斷校正效果，再利用校正后的JRA-55再分析資料，驅(qū)動(dòng)改進(jìn)了粗糙度參數(shù)化方案后的CLM4.0-LISSS模型，獲得1958-2017年太湖湖面蒸發(fā)量，模擬結(jié)果與2012-2017年渦度相關(guān)觀測對比驗(yàn)證其可用性，并分析湖面蒸發(fā)的變化趨勢，尋找太湖實(shí)際蒸發(fā)年際變化的主控因子，以期為太湖地區(qū)水資源管理與氣候變化研究提供科學(xué)依據(jù).

1 材料與方法

1.1 研究區(qū)域概況

太湖(30°55′40″～31°32′58″N，119°52′32″～120°36′10″E，圖1)，地處長江中下游，屬于亞熱帶濕潤季風(fēng)氣候區(qū)，多年平均氣溫15.97℃，年平均降水量1182 mm[33]，湖表面積約為2338 km2，平均水深1.94 m，是典型的亞熱帶大型淺水湖泊.

圖1 太湖中尺度通量網(wǎng)觀測站點(diǎn)位置Fig.1 Locations of the observation sites of the Lake Taihu eddy flux network

1.2 數(shù)據(jù)來源及處理

本文最終用于分析的太湖蒸發(fā)量數(shù)據(jù)均由觀測或模擬所得潛熱通量除以汽化潛熱計(jì)算而來. 本研究利用2013年太湖中尺度通量網(wǎng)[27]避風(fēng)港站(BFG，31°10′28″N，120°24′01″E)觀測數(shù)據(jù)校正JRA-55再分析資料，并利用2012年的觀測資料驗(yàn)證該再分析資料在太湖的適用性；再利用2012-2017年輻射及小氣候觀測資料驅(qū)動(dòng)CLM4.0-LISSS模型，與觀測所得的湖面蒸發(fā)量對比，驗(yàn)證該模型在太湖的適用性；最后利用1958-2017年JRA-55再分析數(shù)據(jù)驅(qū)動(dòng)模型，獲得太湖年累積蒸發(fā)量，分析太湖蒸發(fā)的長期變化趨勢，尋找太湖實(shí)際蒸發(fā)年際變化的主控因子.

1.2.1 觀測數(shù)據(jù) 本文用于校正再分析資料、驗(yàn)證模型適用性的相關(guān)觀測數(shù)據(jù)均來自太湖中尺度通量網(wǎng)[27]. 前期已有研究表明，太湖月蒸發(fā)量無空間差異[11]，因此本文選擇了始建于2011年12月5日、通量觀測數(shù)據(jù)量較為完整且風(fēng)浪區(qū)足夠大的BFG站.

小氣候系統(tǒng)架設(shè)于距湖面8.5 m處，測量空氣濕度、氣溫(Model HMP45D/HMP155A; Vaisala, Inc, Helsinki, Finland)、風(fēng)速風(fēng)向(Model 03002; R. M. Young Company, Traverse City, MI, USA)及四分量輻射(Model CNR4; Kipp & Zonen B.V., Delft, the Netherlands)等氣象要素. 水溫梯度系統(tǒng)由不同深度的水溫計(jì)(Model 109-L, Campbell Scientific Inc., Logan, UT, USA)構(gòu)成，分別測量20、50、100、150 cm深處的水溫及底泥溫度.

渦度相關(guān)系統(tǒng)同樣架設(shè)于距湖面8.5 m高度，由三維超聲風(fēng)速儀(CSAT3, Campbell Scientific Inc., Logan, UT, USA)、開路式紅外CO2/H2O分析儀(EC150, Campbell Scientific Inc., Logan, UT, USA)組成，分別用于觀測三維超聲風(fēng)速、大氣H2O和CO2濃度，數(shù)據(jù)采樣頻率為10Hz；動(dòng)量通量、感熱通量和潛熱通量數(shù)據(jù)輸出步長為30 min. 通量原始數(shù)據(jù)先根據(jù)儀器工作狀態(tài)和日志剔除異常值，再經(jīng)過兩次坐標(biāo)旋轉(zhuǎn)、超聲虛溫訂正、空氣密度效應(yīng)訂正等后處理. 對于缺失較少的通量數(shù)據(jù)(一天內(nèi)少于5個(gè)數(shù)據(jù)點(diǎn))進(jìn)行線性內(nèi)插；缺失較多的通量數(shù)據(jù)，優(yōu)先使用BFG站另外一套EC數(shù)據(jù)建立線性關(guān)系插補(bǔ)，剩余的缺失數(shù)據(jù)，利用物質(zhì)傳輸?shù)南嚓P(guān)理論插補(bǔ)[36].

渦度相關(guān)的觀測結(jié)果存在能量平衡不閉合的問題，導(dǎo)致測得感熱通量H、潛熱通量LE偏低. 本文采用波文比(BRER)法來調(diào)整H、LE，也就是認(rèn)為波文比(H/LE)的觀測值是正確的，將剩余的可利用能量按波文比分配給H和LE，從而達(dá)到能量閉合. 本文的處理過程中，對3日尺度的潛熱通量進(jìn)行了能量閉合校正. 具體的后處理過程參見文獻(xiàn)[11].

1.2.2 JRA-55再分析資料 用于模擬太湖歷史蒸發(fā)的氣象數(shù)據(jù)源于日本氣象廳(Japan Meteorological Agency, JMA)的第二代大氣再分析產(chǎn)品JRA-55再分析資料[37]. JRA-55再分析資料數(shù)據(jù)自1958年開始，在JRA-25的基礎(chǔ)上，改進(jìn)了長波輻射方案，引進(jìn)了四維變分技術(shù)(4D-Var)和變分偏差校正(VarBC)，對平流層低層溫度模擬更加精準(zhǔn)，這對于本文的模擬非常重要. 本文使用的是1958-2017年JRA-55時(shí)間分辨率為3 h的氣溫、大氣壓、比濕、向下的長波輻射、向下的短波輻射數(shù)據(jù)，水平分辨率為1.25°×1.25°. 向上的短波輻射由向下的短波輻射乘以湖表反照率0.08獲得.

再分析資料基于陸地站點(diǎn)發(fā)展而來，而湖陸下墊面在反照率、比熱容、粗糙度等方面的不同導(dǎo)致了湖陸下墊面微氣象條件存在差異. 因此本文采用線性回歸的方法，利用2013年觀測數(shù)據(jù)建立校正方程，對1958-2017年的JRA-55再分析資料進(jìn)行簡單的線性校正，通過對2012年校正前后數(shù)據(jù)的準(zhǔn)確性驗(yàn)證校正效果.

1.2.3 蒸發(fā)量計(jì)算 蒸發(fā)量E可以由觀測或模擬所得潛熱通量LE除以汽化潛熱獲得：

(1)

λ=3.1475×106-2.372×103Ta

(2)

式中，λ為汽化潛熱，Ta為大氣溫度. 年累積蒸發(fā)量由3 h潛熱通量轉(zhuǎn)換成蒸發(fā)量后再求和計(jì)算得到.

1.3 湖泊模型介紹

CLM4-LISSS模型(Community Land Model 4.0-Lake, Ice, Snow and Sediment Simulator)中的驅(qū)動(dòng)數(shù)據(jù)包括氣溫、氣壓、濕度、風(fēng)速、降水、凈短波福射和向下的長波福射. 模型將湖體分為10層，利用大氣強(qiáng)迫場，通過湖表薄層能量平衡經(jīng)過4次迭代獲得表面溫度的唯一解Ts.在每一次迭代中，摩擦速度和空氣動(dòng)力學(xué)阻力將重新計(jì)算，再利用公式(3)重新計(jì)算得到湖-氣間感熱通量H和潛熱通量LE.

湖表薄層的能量平衡方程為：

(3)

(4)

(5)

(6)

湖體各層溫度由垂直溫度擴(kuò)散方程計(jì)算得到，公式為：

(7)

式中，cw為深度z處的比熱；T為湖體不同層的水溫；k為熱導(dǎo)度；φ為到達(dá)深度z處的輻射通量.

本文選擇使用CLM4.0版本模擬太湖蒸發(fā). 目前CLM模型已更新至5.0版本，但使用CLM4.0版本并不影響年尺度太湖蒸發(fā)的模擬效果. CLM5.0模擬過程中湖泊深度為可變項(xiàng)，即所有的物理過程(除與深度有關(guān)的混合過程除外)均可在任意深度下進(jìn)行；而在CLM4.0版本中，湖泊深度則為定值. 模式中湖泊深度主要用于計(jì)算消光系數(shù)，根據(jù)Subin等[39]的敏感性分析可知，改變湖泊深度對水熱通量的模擬結(jié)果與改變消光系數(shù)獲得的結(jié)果類似. 而在本研究中，消光系數(shù)已通過多組敏感性實(shí)驗(yàn)，以年尺度上蒸發(fā)量模擬結(jié)果最佳為目標(biāo)，設(shè)置為2.5 m-1，因此使用CLM4.0版本對年尺度湖面蒸發(fā)數(shù)據(jù)影響不大.

1.4 統(tǒng)計(jì)方法

1.4.1 Mann-Kendall趨勢檢驗(yàn) 本文采用Mann-Kendall秩次相關(guān)檢驗(yàn)(M-K)分析1958-2017年太湖蒸發(fā)量年際變化趨勢[40-43]. 該方法對樣本分布不作要求，且最終檢驗(yàn)結(jié)果受異常值干擾較小，常應(yīng)用于降水、蒸發(fā)、氣溫等隨時(shí)間變化的要素的趨勢變化分析. M-K趨勢檢驗(yàn)原假設(shè)H0是對于樣本量為n的時(shí)間序列x，無時(shí)間變化趨勢，趨勢檢驗(yàn)的統(tǒng)計(jì)量S定義為：

(8)

其中，

(9)

M-K檢驗(yàn)利用標(biāo)準(zhǔn)化后的統(tǒng)計(jì)量Z值判斷時(shí)間序列是否有變化趨勢及是否顯著，變化趨勢的大小用Kendall傾斜度β表示. 計(jì)算公式為：

(10)

(11)

式中，var(S)為統(tǒng)計(jì)量S的方差. 當(dāng)Z值為正時(shí)，時(shí)間序列呈上升趨勢；反之，時(shí)間序列呈下降趨勢. |Z|大于1.28、1.96和2.32分別代表變化趨勢通過了90%、95%和99%的顯著性檢驗(yàn).

1.4.2 多元逐步回歸分析 各驅(qū)動(dòng)變量對太湖歷史蒸發(fā)的貢獻(xiàn)率采用多元逐步回歸來定量分析. 多元逐步回歸方法按因變量的影響大小被依次引入回歸方程中，并考察引入此因子所建立的回歸方程是否顯著. 當(dāng)某因子的初始P值小于0.05，則此因子將被引入并重新計(jì)算P值，若重新計(jì)算的P值大于0.1，此因子將被移除. 本文利用標(biāo)準(zhǔn)化后的驅(qū)動(dòng)數(shù)據(jù)和年累積蒸發(fā)量數(shù)據(jù)(變化范圍0～1)建立回歸方程，以此判斷各氣象因子對太湖蒸發(fā)的相對貢獻(xiàn)[33,44-46]，具體如下：

Y=e+a1x1+a2x2+…aixi

(12)

式中，Y為標(biāo)準(zhǔn)化后的因變量，即年累積蒸發(fā)量；xi(i=1，2，…)為標(biāo)準(zhǔn)化后的因子，即各驅(qū)動(dòng)變量；ai為變量xi在回歸方程中的系數(shù).

本文分析所用各驅(qū)動(dòng)變量對年累積蒸發(fā)量的實(shí)際貢獻(xiàn)率，由下列公式計(jì)算所得：

(13)

式中，Δxi和ΔY分別為xi和Y的變化量；μi為xi對Y的實(shí)際貢獻(xiàn)率.

1.4.3 評(píng)價(jià)指標(biāo) 本文所用評(píng)價(jià)指標(biāo)為平均偏差(mean bias error,MBE)、均方根誤差(root mean square error,RMSE)、決定系數(shù)(R2)、一致性指數(shù)(I)，計(jì)算公式為：

(14)

(15)

(16)

(17)

2 結(jié)果與分析

2.1 驅(qū)動(dòng)數(shù)據(jù)校正結(jié)果

在驅(qū)動(dòng)模型模擬蒸發(fā)量之前，需利用實(shí)測資料對再分析資料進(jìn)行校正，本文使用2013年BFG站點(diǎn)小氣候及輻射觀測數(shù)據(jù)對JRA-55進(jìn)行簡單的線性校正，再利用2012年的實(shí)測資料驗(yàn)證校正結(jié)果. 表1展示了各驅(qū)動(dòng)變量與實(shí)測數(shù)據(jù)的校正方程及校正前后的統(tǒng)計(jì)參數(shù). 所有變量的校正方程均通過了α= 0.05的顯著性檢驗(yàn)，線性校正有效地降低了各強(qiáng)迫數(shù)據(jù)的系統(tǒng)偏差.

表1 太湖地區(qū)JRA-55再分析資料驅(qū)動(dòng)數(shù)據(jù)校正方程及校正前后優(yōu)度統(tǒng)計(jì)*

校正前JRA-55再分析資料與BFG站點(diǎn)實(shí)測資料對比如圖2所示. 從相關(guān)性來看，再分析資料中向下短波輻射、長波輻射、氣溫、氣壓、比濕均與BFG觀測值均呈現(xiàn)很好的相關(guān)性(P<0.01)，平均風(fēng)速與觀測相關(guān)性較差(r=0.29)，但仍通過了0.01的顯著性檢驗(yàn)；從一致性系數(shù)I來看，除JRA-55再分析資料的風(fēng)速及氣壓的一致性系數(shù)較差，其余均在0.9以上；從校正前MBE及RMSE來看(如表1所示)，JRA-55再分析資料高估了太湖地區(qū)向下的短波輻射、氣溫，低估了向下長波輻射、氣壓、比濕與風(fēng)速. 總體而言，JRA-55中的再分析資料中向下長波輻射、向下短波輻射、氣溫、比濕與實(shí)測值偏差較小，但氣壓存在明顯的系統(tǒng)偏差，風(fēng)速信號(hào)接近陸地觀測，與湖泊風(fēng)速有很大差異.

校正后，向下短波、向下長波輻射的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；氣壓被低估的現(xiàn)象明顯改善，校正后MBE僅有0.26 hPa；風(fēng)速被低估的現(xiàn)象也得到明顯改善，MBE由-2.66 m/s降低至0.14 m/s. 總體而言，線性校正方法有效地校正了JRA-55再分析資料與湖上觀測數(shù)據(jù)的系統(tǒng)偏差.

圖2 2012年校正前JRA-55再分析資料與BFG站點(diǎn)實(shí)測資料的對比： (a)向下短波輻射；(b)向下長波輻射；(c)氣溫；(d)比濕；(e)風(fēng)速；(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日累積蒸發(fā)量時(shí)間序列 (黑色實(shí)線代表BFG站點(diǎn)EC觀測值所得 蒸發(fā)量(經(jīng)能量平衡閉合校正))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日累積蒸發(fā)量的1∶1圖 (●代表BFG在線觀測驅(qū)動(dòng)模型； ○代表再分析資料驅(qū)動(dòng)模型)Fig.4 Comparison between simulation and observation (●, simulation forced by in situ observation; ○, simulation forced by calibrated JRA-55)

2.2 CLM4.0-LISSS模型及再分析資料適用性驗(yàn)證

圖3展示了2012年BFG觀測蒸發(fā)量與CLM4.0-LISSS模擬蒸發(fā)量的對比結(jié)果，用于驗(yàn)證模型的適用性. 模型能很好地模擬出2012年蒸發(fā)量的逐日變化和季節(jié)變化趨勢. 春季，模擬及觀測蒸發(fā)量均迅速增加，6月因“梅雨”季太陽輻射的減少而稍有減少，7、8月達(dá)到最大值，隨后蒸發(fā)量逐漸減小，并于冬季達(dá)到最小值. 圖4也表明，觀測與模擬值的變化范圍相一致，3日平均的模擬值與觀測值建立的線性回歸線斜率等于0.86，MBE和RMSE分別為-0.3和2.0 mm. 這證明了CLM4.0-LISSS模型有模擬太湖蒸發(fā)的能力，可用于太湖歷史蒸發(fā)量變化趨勢分析研究.

校正后的JRA-55再分析資料被用于驅(qū)動(dòng)模型，以驗(yàn)證再分析資料的適用性. 圖4可以看出，相比于使用小氣候系統(tǒng)觀測數(shù)據(jù)驅(qū)動(dòng)模型，利用再分析資料驅(qū)動(dòng)模型時(shí)，模型表現(xiàn)略有下降，雖然兩者擬合線接近，但模擬值更加離散(MBE=0.1 mm，RMSE=5.5 mm). 模型雖然能很好地模擬出太湖蒸發(fā)的季節(jié)變化趨勢(圖3)，但在秋冬季節(jié)，模型模擬值明顯偏高，春夏季節(jié)蒸發(fā)量的模擬值偏低，這可能與每年只取了一個(gè)校正系數(shù)有關(guān). 但是，季節(jié)性偏差在年尺度上相互抵消，年累積蒸發(fā)量的模擬值與觀測值相近，因此本文沒有再分季節(jié)/分月對再分析資料進(jìn)行校正. 2012-2017年太湖年累積蒸發(fā)量的模擬值(觀測值)分別為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)，觀測值和模式值的對比也表明再分析資料能模擬出太湖年累計(jì)蒸發(fā)量的趨勢變化(圖5). 因此本文認(rèn)為校正后的JRA-55再分析資料可以驅(qū)動(dòng)模型模擬年際尺度的湖泊蒸發(fā).

2.3 歷史蒸發(fā)模擬值的年際變化趨勢分析

1958-2017年太湖年累積蒸發(fā)量模擬值的時(shí)間序列如圖5所示. 近60年來，太湖年累積蒸發(fā)量的平均值為1018.2 mm，其中1958-1962、1965-1966、1968-1970、1972-1977、1980-1990、1999、2002、2014-2016年共31年蒸發(fā)量小于平均值，最小值出現(xiàn)于1977年，為870.6 mm，其余年份蒸發(fā)量高于60年平均值，最大值出現(xiàn)于2004年，為1120.2 mm. 需要注意的是，1977年累積蒸發(fā)量較多年平均值異常偏低約14.49%，這可能是由于1977年6月降水量較多，太陽輻射偏低(125.28 W/m2，遠(yuǎn)低于平均值155.76 W/m2)，導(dǎo)致1977年全年太陽輻射僅有128.65 W/m2，遠(yuǎn)低于平均水平(136.99 W/m2)約6.09%，且1977年大氣比濕較高(10.7 g/kg)，比平均值(10.3 g/kg)高約6.8%，而敏感性分析結(jié)果表明潛熱通量對太陽輻射及比濕極其敏感，將1977年太陽輻射增加6%，比濕減少6%后，蒸發(fā)量模擬值將增加至1022.4 mm，與多年均值相近.

線性回歸和M-K趨勢檢驗(yàn)被用于分析蒸發(fā)量的變化趨勢分析. 線性回歸分析說明，1958-2017年，太湖蒸發(fā)量整體以1.4 mm/a的速率波動(dòng)上升. M-K趨勢檢驗(yàn)也顯示太湖蒸發(fā)量整體呈上升趨勢，上升速率為1.4 mm/a，并通過了99%的顯著性檢驗(yàn)(Z=3.30). 但是在此期間太湖蒸發(fā)量明顯被分為兩個(gè)階段：1958-1977年間，太湖蒸發(fā)量呈現(xiàn)了明顯的下降趨勢(線性回歸斜率=-4.0 mm/a， Kendall傾斜度=-3.6 mm/a，Z=-1.89，通過了90%的顯著性檢驗(yàn))；1978-2017年間，M-K趨勢檢驗(yàn)顯示，太湖蒸發(fā)量以2.3 mm/a的速度增加(Z=2.69，通過了99%的顯著性檢驗(yàn)). 總體而言，1958-2017年間太湖蒸發(fā)量以1977年為界限，呈現(xiàn)先減后增的趨勢.

圖5 1958-2017年太湖蒸發(fā)量模擬值的年際變化(黑色實(shí)線代表模型模擬年累積蒸發(fā)量； 灰色虛線代表1958-2017蒸發(fā)量趨勢擬合線；藍(lán)色虛線代表1958-1977年蒸發(fā)量線性回歸線； 紅色虛線代表1978-2017年蒸發(fā)量線性回歸線；○代表2012-2017年觀測年累計(jì)蒸發(fā)量)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)

本研究計(jì)算所得太湖蒸發(fā)量在1978-2017年間以2.3 mm/a的速度增加，這與Hu等[33]利用MERRA再分析資料驅(qū)動(dòng)湖模式計(jì)算得到1979-2013年間太湖蒸發(fā)量以2.96 mm/a的速率增加，以及蒸發(fā)皿蒸發(fā)量觀測值以2.54 mm/a的速率增加的結(jié)論較為一致. Rong等[45]結(jié)合Penman-Monteith方程與參考蒸散率計(jì)算了2003-2010年間距離太湖640 km的另一淺水湖泊東平湖的蒸發(fā)，其蒸發(fā)量也呈上升趨勢，上升速率為4.55 mm/a. 而與太湖同緯度的青藏高原湖泊色林措蒸發(fā)量也在2007-2015年間以4.3 mm/a的速度增加[43].

2.4 驅(qū)動(dòng)數(shù)據(jù)的年際變化趨勢分析

1958-2017年氣溫(Ta)、氣壓(Pa)、比濕(qa)、風(fēng)速(U)、向下長波輻射(L↓)、凈短波輻射(S*)等強(qiáng)迫數(shù)據(jù)的變化趨勢如圖6(線性回歸方法)和表2所示(M-K趨勢檢驗(yàn)法). 總體來說，1958-2017年間，Ta的增加趨勢最為顯著，線性回歸趨勢線斜率為0.015℃/a(P<0.01)，M-K趨勢檢驗(yàn)顯示氣溫平均每年升高0.014℃，統(tǒng)計(jì)量Z值為3.50，通過了置信水平為99%的顯著性檢驗(yàn). 氣溫最低值(15.82℃)出現(xiàn)于1984年，最高值(17.79℃)出現(xiàn)于1998年.L↓在1958-2017年也呈現(xiàn)顯著增加的趨勢，M-K方法計(jì)算得到的L↓增加量為0.053 W/(m2·a)(置信水平為95%)，線性回歸計(jì)算得到的L↓增長速率為0.048 W/(m2·a)(P<0.05). M-K方法顯示S↓以每年0.047 W/m2的趨勢波動(dòng)上升(置信水平為90%)，但線性回歸結(jié)果未通過趨勢檢驗(yàn).S↓變量范圍為128.24～143.86 W/m2，其中1977年達(dá)到最小值，2007年達(dá)到最大值. 研究表明，中國的地面太陽輻射大致以1990年為界，呈現(xiàn)先減少后增多的趨勢[47]，這與本文研究結(jié)果基本相符合. 其他變量(qa、U、Pa)在1958-2017年間無顯著變化趨勢. 本研究校正后的JRA-55再分析資料的風(fēng)速在 1958-2017年間無顯著下降趨勢，這與Hu等[33]使用的MERRA再分析資料得到的結(jié)論相同，但與地面觀測所得到的結(jié)論存在差異. 近50年來，地表粗糙度的增加(如城市發(fā)展、植被增加、耕地變化等)造成了內(nèi)陸近地面風(fēng)速持續(xù)減弱，觀測中海洋和沿海地區(qū)的近地層風(fēng)速大多呈現(xiàn)出上升趨勢[48-49]. 這可能是因?yàn)镴RA-55再分析資料未考慮粗糙度的變化，且同化時(shí)使用的地表觀測資料與ERA-40一致，僅使用了溫度、濕度、降雪資料[50-51]，未包含土地利用變化(如城市化進(jìn)程等)信息，所以再分析資料中的近地層風(fēng)速變化趨勢弱于氣象站觀測數(shù)據(jù)[51]. 但是由于本文研究站點(diǎn)位于湖上，地表粗糙度小且變化不大，因此校正后的再分析資料因未考慮地表粗糙度變化而呈現(xiàn)出不顯著變化趨勢不影響本研究使用其作為驅(qū)動(dòng)數(shù)據(jù).

圖6 1958-2017年JRA-55資料氣象因子變化趨勢(黑色實(shí)線代表年平均值；灰色虛線代表1958-2017年 間各氣象因子線性回歸擬合線；藍(lán)色虛線代表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年間年均氣溫、氣壓、比濕、風(fēng)速、向下長波輻射、 向下短波輻射M-K檢驗(yàn)Z值及變化速率

驅(qū)動(dòng)變量的變化趨勢也被分為1958-1977年及1978-2017年兩個(gè)階段分析(表2). M-K方法顯示1958-1977年間，除L↓無明顯年際變化趨勢外，Ta、qa、U和S↓均表現(xiàn)出明顯的下降趨勢，Pa呈顯著的增加趨勢；線性回歸方法顯示除氣溫以0.047℃/a的趨勢下降外，其余驅(qū)動(dòng)變量無顯著變化特征. 1978-2017年間，各氣象因子的年際變化趨勢與1958-1977年存在差異.Ta、qa、S↓由顯著下降趨勢轉(zhuǎn)變?yōu)轱@著上升(M-K趨勢檢驗(yàn)方法)，上升趨勢分別為0.032℃/a(Z=4.46，通過了99%的顯著性檢驗(yàn))、7.05×10-3kg/kg(Z=2.11，通過了95%的顯著性檢驗(yàn))和0.10 W/(m2·a)(Z=-1.40，通過了90%的顯著性檢驗(yàn)).U和L↓則轉(zhuǎn)變?yōu)槌拭黠@的上升趨勢，1977年后分別以9.69×10-4m/(s·a)、0.094 W/(m2·a)的速度增長，并分別通過了90%、99%的顯著性檢驗(yàn).Pa則呈現(xiàn)了明顯的下降趨勢. 總體來說，氣象因子和輻射能量在1958-1977年間與1978-2017年間呈現(xiàn)了不同的趨勢，這對于解釋蒸發(fā)量的年際變化有重要意義.

2.5 太湖蒸發(fā)年際變化主控因子分析

為確定太湖蒸發(fā)量年際變化主控因子，我們利用多元逐步回歸分析了各要素對蒸發(fā)量的實(shí)際貢獻(xiàn)率. 基于前文趨勢分析結(jié)果，將研究時(shí)段分為1958-1977年和1978-2017年兩個(gè)時(shí)間段分別探討. 1958-1977年間太湖蒸發(fā)與氣象因子之間的多元回歸方程為E=-0.28qa+0.43L↓+0.71S↓.Ta、Pa、U被排除在方程外，僅保留了qa、L↓、S↓3個(gè)變量. 方程的R2為0.764，說明方程是可以解釋76.4%太湖蒸發(fā)量的變異. 根據(jù)公式(13)計(jì)算可得，S↓對太湖蒸發(fā)變化的貢獻(xiàn)比率最大，高達(dá)63.51%，其次來自于L↓的貢獻(xiàn)(38.10%)，另外qa的實(shí)際貢獻(xiàn)為-24.48%，位列第三，而Ta、Pa和U對太湖蒸發(fā)年際變化貢獻(xiàn)很??；由1978-2017年數(shù)據(jù)建立的回歸方程與1958-1977年有所不同，Pa、U、L↓被移出方程，保留了Ta、qa、S↓3個(gè)變量，回歸方程為E=0.42Ta-0.31qa+0.51S↓，對太湖蒸發(fā)年際變化貢獻(xiàn)較大的因素為Ta、qa、S↓，三者的實(shí)際貢獻(xiàn)率分別為44.54%、-39.13%、52.51%；與此同時(shí)，1958-2017年的全部數(shù)據(jù)也被用于建立多元逐步回歸方程，Ta、qa、L↓、S*的多元逐步回歸系數(shù)分別為0.37、-0.39、0.17、0.55.Pa、U被排除在回歸方程外. 方程的R2為0.851，說明方程可以解釋的湖泊蒸發(fā)量的變異率為85.1%. 從貢獻(xiàn)率來看，輻射是湖泊蒸發(fā)的直接能量來源，S*的增加(貢獻(xiàn)率為54.94%)是1958-1977年太湖年累積蒸發(fā)量增加的主要原因. 第二大貢獻(xiàn)來自qa，qa抵消了部分由輻射減少帶來的正面貢獻(xiàn)，占43.49%.Ta、L↓的貢獻(xiàn)分別為39.59%和18.03%. 這個(gè)結(jié)果與分段建立方程的結(jié)果相吻合. 總得來說，S↓的增加是造成太湖1958-2017年蒸發(fā)量增加的主控因子，L↓、Ta、qa也會(huì)影響湖泊蒸發(fā)年際變化.

3 討論

目前，關(guān)于控制開闊水面蒸發(fā)年際變化的主要?dú)庀笠蜃佑休椛淇刂?、溫度、風(fēng)速變化3種觀點(diǎn). 大部分研究采用蒸發(fā)皿測算量代替(蒸發(fā)皿觀測乘以系數(shù)或?qū)⑵渲糜诤?實(shí)際水面蒸發(fā)，或利用Penman公式及其衍生公式估算，得出蒸發(fā)量與同期氣溫(平均氣溫/最低氣溫/氣溫日較差)相關(guān)性最大[18,33,52-58]，但有研究認(rèn)為蒸發(fā)皿蒸發(fā)量變化與氣溫關(guān)系不大，因?yàn)樵缬杏^測證明北半球冬季的總蒸發(fā)量大于夏季[59]，且部分地區(qū)氣溫與蒸發(fā)量變化呈相反趨勢，即存在“蒸發(fā)悖論”現(xiàn)象，并用近50年來季風(fēng)減弱造成風(fēng)速減小[18-19,46,60-66]，或云層覆蓋和氣溶膠濃度增加造成太陽輻射(或日照時(shí)數(shù))減少[17-19,63,67-70]來解釋. 也有研究認(rèn)為蒸發(fā)皿蒸發(fā)量只是陸地蒸發(fā)量增加的信號(hào)(陸地蒸發(fā)量使得蒸發(fā)皿周圍濕度梯度降低，造成蒸發(fā)皿蒸發(fā)量下降)，只有在周圍陸表環(huán)境水分供應(yīng)充足時(shí)才可以代表開闊水面蒸發(fā)[16]. 近年來，渦度相關(guān)方法被用于觀測開闊水面蒸發(fā)，并逐步發(fā)展成為國際通用的通量觀測標(biāo)準(zhǔn)方法，結(jié)合各類模型模擬過去或未來氣候情景下的水面蒸發(fā)用于分析蒸發(fā)量的年際變化. 對于高緯度結(jié)冰水體，大部分研究認(rèn)為全球變暖導(dǎo)致湖泊冰期縮短，進(jìn)而引發(fā)的反照率變化，是影響水面蒸發(fā)年際變化的重要因素[20,71-73]. 一部分研究利用未來氣候情景的氣象數(shù)據(jù)模擬蒸發(fā)時(shí)發(fā)現(xiàn)，當(dāng)輻射輸入無明顯變化時(shí)蒸發(fā)量仍存在明顯年際變化，因此氣溫加速了湖泊蒸發(fā)年際變化，且由于湖泊表面增溫速度慢于空氣增溫速度，因此波文比發(fā)生改變，能量更多分配用于蒸發(fā)(即潛熱)，造成蒸發(fā)逐年增加[20,72,74-75]. 但對于未結(jié)冰的湖泊，在不同未來氣候情景下，輻射強(qiáng)迫存在差異，蒸發(fā)量的增加速率也會(huì)隨著輻射強(qiáng)迫的增加而增大[34]. 大部分研究通過長年EC觀測或湖模式、蒸發(fā)模式的敏感性分析、相關(guān)性分析等方法得出輻射是開闊水域蒸發(fā)的能量來源，其變化是蒸發(fā)年際變化的主控因子，輻射變化對蒸發(fā)趨勢的貢獻(xiàn)大于氣溫變化的貢獻(xiàn)[1,29,33,68,76-78]，與本文的結(jié)果一致.

4 結(jié)論

本研究利用2013年BFG觀測數(shù)據(jù)校正了JRA-55再分析資料并驅(qū)動(dòng)CLM4.0-LISSS模型，模擬了太湖1958-2017年間湖面蒸發(fā)量，利用2012年模擬所得蒸發(fā)量與避風(fēng)港站點(diǎn)實(shí)測數(shù)據(jù)驗(yàn)證模型與再分析資料的適用性，利用M-K趨勢檢驗(yàn)分析了長時(shí)間尺度上太湖蒸發(fā)量的變化趨勢，再利用多元逐步回歸方法尋找太湖實(shí)際蒸發(fā)的年際變化的主控因子. 主要研究成果如下：

1)利用校正后的JRA-55再分析資料驅(qū)動(dòng)模型時(shí)，蒸發(fā)的模擬存在季節(jié)偏差，但季節(jié)性偏差在年尺度上相互抵消，2012年(2013年)年累積蒸發(fā)量的模擬值與觀測值相近，分別為1050.3 mm(1113.1 mm)和1007.0 mm(1129.0 mm). JRA-55再分析資料可用于年際湖泊蒸發(fā)的模擬.

2)1958-2017年太湖蒸發(fā)分為兩個(gè)階段：1958-1977年間，蒸發(fā)呈現(xiàn)了明顯的下降趨勢，下降速率為-3.6 mm/a；1978-2017年間，太湖蒸發(fā)量以2.3 mm/a的速度增加.

3)通過多元逐步回歸分析可知，向下的短波輻射是太湖1958-2017年間年尺度蒸發(fā)的主控因子，向下的長波輻射、氣溫、比濕也會(huì)影響湖泊蒸發(fā)年際變化，蒸發(fā)量對風(fēng)速變化和氣壓變化不敏感.