萬(wàn)祥禹， 游為， 王海波， 范東明

西南交通大學(xué)地球科學(xué)與環(huán)境工程學(xué)院， 成都 611756

0 引言

陸地水是人類賴以生存的重要資源，研究其變化具有重要意義，但由于各種因素的制約，目前地面實(shí)測(cè)手段十分有限.GRACE(Gravity Recovery and Climate Experiment)重力衛(wèi)星為連續(xù)、大尺度監(jiān)測(cè)陸地水儲(chǔ)量變化提供了新的方法.利用GRACE衛(wèi)星重力數(shù)據(jù)反演地表質(zhì)量變化的理論最先由Wahr等提出(Wahr et al., 1998).此后Swenson和Wahr提出了利用GRACE 數(shù)據(jù)反演區(qū)域質(zhì)量變化的方法(Swenson and Wahr, 2002).國(guó)內(nèi)外已有多位學(xué)者利用GRACE數(shù)據(jù)成功探測(cè)到亞馬遜流域、印度西北地區(qū)、中東地區(qū)、加利福尼亞地區(qū)、長(zhǎng)江流域、中國(guó)華北地區(qū)和黑河流域等多個(gè)明顯的陸地水變化信號(hào)(Chen et al., 2009, 2014; Feng et al., 2013; Rodell et al., 2009; Scanlon et al., 2012; Tapley et al., 2004; Voss et al., 2013; 馮偉等, 2012; 胡小工等, 2006; 李瓊等, 2013; 尼勝楠等, 2014; 蘇曉莉等, 2012; 吳云龍等, 2015).

對(duì)維多利亞湖流域，目前也有多位學(xué)者成功用GRACE探測(cè)到其在2003—2006年?yáng)|非嚴(yán)重干旱時(shí)期的水儲(chǔ)量減少信號(hào)，并分析得出在此期間大壩排水量增大加劇了水儲(chǔ)量的減少.Swenson和Wahr(2009)利用多源衛(wèi)星數(shù)據(jù)分析了東非大裂谷地區(qū)2003—2008年間的陸地水儲(chǔ)量變化，結(jié)果表明在2003—2007年間GRACE和衛(wèi)星測(cè)高分別監(jiān)測(cè)到維多利亞湖水儲(chǔ)量以60 mm·a-1和311 mm·a-1的速率減少，并得出這一時(shí)期約有50%的水儲(chǔ)量減少由人類活動(dòng)造成.Becker等(2010)用PCA(Principal Component Analysis)方法分析了東非湖區(qū)陸地水儲(chǔ)量和降雨的變化，結(jié)果表明降雨和GRACE觀測(cè)到的陸地水儲(chǔ)量存在一個(gè)相近的變化，兩者在2005年底都處于最低值，隨后兩者在2006—2007年都迅速增長(zhǎng)，這一變化受2006年的IOD事件影響.Awange等(2013)用GRACE探測(cè)到維多利亞湖流域的水儲(chǔ)量在2003—2007年因大壩排水量增大和干旱以38.2 mm·a-1的速率減少，2007—2013年流域水儲(chǔ)量因2006—2007年和2010—2011年的大量降雨以4.5 mm·a-1的速率增加.Awange等(2014)利用ICA(Independent Component Analysis)方法分析了尼羅河流域的水儲(chǔ)量變化，結(jié)果表明維多利亞湖流域的水儲(chǔ)量變化受自然因素和人類活動(dòng)共同影響.2003-10—2006-03水儲(chǔ)量因干旱和人類活動(dòng)以84.5 mm·a-1的速率減少，2006-03—2007-04又因ENSO引起的降雨變化以145.2 mm·a-1的速率增長(zhǎng)，2007-04—2009-08間水儲(chǔ)量以25.8 mm·a-1的速率減少，2009-08—2011-01間因降雨增長(zhǎng)以49.8 mm·a-1的速率增長(zhǎng).Hassan和Jin(2014)利用GRACE等數(shù)據(jù)分析了2003-01—2012-12間東非大裂谷地區(qū)的陸地水儲(chǔ)量變化，結(jié)果表明2003-01—2005-12間維多利亞湖流域水儲(chǔ)量以41.5 mm·a-1的速率減少，2006-01—2008-12間以28.7 mm·a-1的速率增加，2003-01—2012-12間流域水儲(chǔ)量變化速率不明顯，為2.4 mm·a-1.Nanteza等(2016)的研究結(jié)果表明在2003—2010間維多利亞湖流域水儲(chǔ)量長(zhǎng)期變化趨勢(shì)不明顯，但在2003—2006年間GRACE探測(cè)到水儲(chǔ)量以4 cm·a-1的速率減少，衛(wèi)星測(cè)高與土壤水結(jié)合的結(jié)果探測(cè)到水儲(chǔ)量以6 cm·a-1的速率減少，變化速率明顯高于周邊地區(qū)(約高出80%)，這一現(xiàn)象由人類活動(dòng)引起，2006年9月至12月的大量降雨抵消了2003—2006的水儲(chǔ)量虧損 (Awange et al., 2013, 2014; Becker et al., 2010; Hassan and Jin, 2014; Nanteza et al., 2016; Swenson and Wahr, 2009).國(guó)內(nèi)外針對(duì)該流域的研究中對(duì)泄漏誤差的改正多采用尺度因子方法，Landerer和Swenson指出該方法強(qiáng)烈依賴水文模式結(jié)果(Landerer and Swenson, 2012).但水文模式在湖泊水體表面效果不佳，會(huì)降低GRACE估計(jì)結(jié)果的可靠性.由于GRACE官方機(jī)構(gòu)最新公布的RL06版本球諧位系數(shù)的處理過(guò)程相比RL05版本采用了一些新的背景模型，改善了處理方法，條帶誤差有所減小，信噪比顯著提高，這些因素都會(huì)對(duì)陸地水儲(chǔ)量變化的估算結(jié)果造成影響.此外國(guó)內(nèi)外對(duì)維多利亞湖流域水儲(chǔ)量的研究時(shí)間段多為2002—2016年，而較少關(guān)注2016—2017年的流域水儲(chǔ)量變化.針對(duì)以上這些問(wèn)題本文利用CSR發(fā)布的2003—2017年的RL06時(shí)變重力場(chǎng)模型與多種衛(wèi)星數(shù)據(jù)，采用正向建模方法改正泄漏誤差，重新估算了維多利亞湖流域水儲(chǔ)量變化速率并分析其在2016—2017年的水儲(chǔ)量變化.本文得出幾點(diǎn)與此前研究不同的結(jié)論，推測(cè)在2013-01—2016-02期間人工調(diào)控對(duì)湖泊水儲(chǔ)量造成較大影響，為內(nèi)陸湖泊流域水儲(chǔ)量變化的分析提供一種可選擇的方法和分析策略.

1 研究區(qū)域概況

維多利亞湖位于東非，是世界第二大淡水湖，面積約68600 km2，流域面積約261800 km2，流域面積足夠被GRACE衛(wèi)星(空間分辨率～300 km，等效水高精度～2 cm)以厘米級(jí)精度探測(cè)到質(zhì)量變化.流域內(nèi)約有三千萬(wàn)人口居住，湖泊為流域內(nèi)居民提供灌溉、電力等水資源，研究其水儲(chǔ)量變化對(duì)當(dāng)?shù)厮Y源合理利用具有重要意義，也可為研究我國(guó)境內(nèi)的湖泊水儲(chǔ)量變化提供參考.維多利亞湖上目前共建有Nalubaale、Kiira和Bujagali三處大壩，大壩的蓄水和放水會(huì)對(duì)湖泊水儲(chǔ)量造成較大影響.維多利亞湖的湖水僅能通過(guò)Jinja處的河流流往下游，維多利亞湖的南方還有坦噶尼喀湖 (Lake Tanganyika) 和馬拉維湖 (Lake Malawi).本文的研究范圍緯度為14°N—16°S，經(jīng)度為23°E—38°E.圖1為維多利亞湖流域概況，維多利亞湖流域范圍根據(jù)TRIP(Total Runoff Integrating Pathways)模型確定(Oki and Sud, 1998).

圖1 維多利亞湖流域概況，圖中紅線表示Total Runoff Integrating Pathways(TRIP)模型確定的維多利亞湖流域范圍(Oki and Sud, 1998)Fig.1 Overview of Lake Victoria Basin, the red line represents the boundary of the Lake Victoria Basin, which provided by the model routing Total Runoff Integrating Pathways (TRIP) (Oki and Sud, 1998)

2 數(shù)據(jù)和基本方法

2.1 GRACE

本文使用CSR最新發(fā)布的2003-01—2017-06的60階GRACE RL06時(shí)變重力場(chǎng)模型球諧位系數(shù)計(jì)算區(qū)域水儲(chǔ)量變化.同時(shí)也采用CSR RL06 Mascon產(chǎn)品(Save et al., 2016)與球諧位系數(shù)計(jì)算結(jié)果進(jìn)行對(duì)比分析，下文用SH和Mascon對(duì)兩種數(shù)據(jù)加以區(qū)分.本文用衛(wèi)星激光測(cè)距數(shù)據(jù)確定的C20項(xiàng)替代球諧位系數(shù)中的C20項(xiàng)(Cheng and Ries, 2017; Cheng and Tapley, 2004).RL05球諧位系數(shù)的一階項(xiàng)用Swenson等的計(jì)算結(jié)果替換，RL06數(shù)據(jù)則用基于Sun等的研究計(jì)算的系數(shù)替換(Sun et al., 2016; Swenson et al., 2008).GIA(glacial isostatic adjustment)用Geruo等的結(jié)果改正(Geruo et al., 2013).

球諧位系數(shù)的高階項(xiàng)在包含更多高頻信號(hào)的同時(shí)，也造成等效水高在空域分布中表現(xiàn)出很大的噪聲，因此本文用高斯濾波對(duì)其進(jìn)行平滑處理，半徑越大平滑效果越好但信號(hào)衰減或泄漏也越嚴(yán)重(Jekeli, 1981; Wahr et al., 2004).球諧位系數(shù)處理中的截?cái)嗪透咚蛊交瑫?huì)造成原始信號(hào)一定程度上的空間泄漏，本文采用Chen等提出的正向建模方法改正泄漏誤差(Chen et al., 2013, 2014, 2015, 2017).本文使用的是無(wú)空間約束的正向建模方法，處理步驟如下：

(1)將經(jīng)過(guò)截?cái)嗪透咚篂V波平滑處理后的等效水高或變化速率h0作為模型迭代初始值m0；

(2)將(1)中的模型初始值m0做球諧展開(kāi)截?cái)嘀僚c原始球諧位系數(shù)一致的階數(shù)，然后進(jìn)行相同半徑的高斯濾波，得到模型h1；

(3)將GRACE觀測(cè)值h0和模型值h1之間的差值k·(h0-h1)加入到模型初始值m0中(k為加快收斂的參數(shù))；

(4)重復(fù)步驟(2)—(3)一定次數(shù)后或者GRACE觀測(cè)值和模型之差小于一定閾值時(shí)，迭代終止.

通常GRACE數(shù)據(jù)后處理中除了高斯濾波外還需要進(jìn)行去相關(guān)濾波進(jìn)一步抑制條帶誤差.但Chen等(2014)指出去相關(guān)濾波對(duì)結(jié)果的影響是非線性的，這一部分影響正向建模方法無(wú)法恢復(fù)，故本文僅采用500 km半徑的高斯濾波來(lái)抑制噪聲.CSR的Mascon產(chǎn)品基于正則化方法僅使用GRACE信息解算而不輸入額外的模型，并且已經(jīng)利用一些約束條件改正泄漏誤差，不需要再進(jìn)行高斯濾波和去相關(guān)濾波，空間分辨率比球諧位系數(shù)高 (Save et al., 2016).

2.2 水文模式

陸地水主要由冰雪、地表水、土壤水和地下水等組成，有多種水文模式可以輸出土壤水.本文采用了兩種水文模式GLDAS(Global Land Data Assimilation System)和CPC(Climate Prediction Center)與GRACE結(jié)果對(duì)比，這兩種水文模式均不包含地表湖泊水和地下水.GLDAS水文模式包含四個(gè)模型驅(qū)動(dòng)的數(shù)據(jù)：Noah、VIC、CLM和Mosaic(D?ll et al., 2003; Rodell et al., 2004).本文采用了1°×1°的基于Noah陸地表面模型的GLDAS水文模式，取4層土壤水、植被冠層含水和冰雪之和作為陸地水.CPC水文模式來(lái)源于美國(guó)國(guó)家海洋和大氣管理局(NOAA)氣候預(yù)測(cè)中心(CPC)，每月一值，空間分辨率為0.5°×0.5°(Fan and Van Den Dool, 2004).本文采用V2版本的土壤濕度數(shù)據(jù)(0～160 cm)的月均值作為陸地水.

2.3 TRMM

降雨變化通常是造成區(qū)域陸地水變化的主要原因之一，因此本文利用TRMM(Tropical Rainfall Measuring Mission)降雨數(shù)據(jù)衡量降雨的變化.TRMM由美國(guó)國(guó)家宇航局(NASA)和日本空間發(fā)展局(JAXA)共同研發(fā)，該衛(wèi)星可提供50°S—50°N范圍內(nèi)空間分辨率0.25°×0.25°的熱帶、亞熱帶降雨數(shù)據(jù)(Kummerow et al., 1998).本文使用每月一值的3B43產(chǎn)品.

2.4 衛(wèi)星測(cè)高

湖泊水量變化除了能以質(zhì)量變化形式被重力衛(wèi)星探測(cè)到以外，還可直接利用衛(wèi)星測(cè)高監(jiān)測(cè)到水位高度變化，也可與GRACE進(jìn)行對(duì)比.衛(wèi)星測(cè)高監(jiān)測(cè)面積足夠大的水體水位變化時(shí)水面高度變化的均方根(root mean squares)能夠達(dá)到3～5 cm的水平.Hydroweb網(wǎng)站(http:∥hydroweb.theia-land.fr)提供了本文使用的湖泊水位高產(chǎn)品(Crétaux et al., 2011).

2.5 蒸散

蒸散也是影響區(qū)域陸地水儲(chǔ)量變化的因素，本文使用了蒙大拿大學(xué)NTSG(Numerical Terradynamic Simulation Group)發(fā)布的基于MODIS(Moderate Resolution Imaging Spectroradiometer)的尼羅河流域MOD16A2蒸散產(chǎn)品.該產(chǎn)品是規(guī)則的1 km2地表蒸散數(shù)據(jù)集，提供8天和1月分辨率的蒸散數(shù)據(jù).該產(chǎn)品基于Penman-Monteith方程，詳細(xì)的介紹見(jiàn)Mu等的文章(Mu et al., 2011).本文使用的數(shù)據(jù)是2003-01—2014-09的月蒸散產(chǎn)品，通過(guò)插值獲得與GRACE一致的空間分辨率.

3 計(jì)算結(jié)果與討論

3.1 RL06版本數(shù)據(jù)與RL05版本數(shù)據(jù)對(duì)比

GRACE官方機(jī)構(gòu)最新公布的RL06版本球諧位系數(shù)的處理過(guò)程相比RL05版本采用了一些新的背景模型，改善了處理方法，條帶誤差有所減小，信噪比顯著提高，有利于提取出真實(shí)變化信號(hào).現(xiàn)將60階的CSR RL05版本和RL06版本球諧位系數(shù)做相同的數(shù)據(jù)處理(500 km半徑高斯濾波)反演出2003-01—2017-06間以等效水高形式表現(xiàn)的全球質(zhì)量變化速率，兩者對(duì)比見(jiàn)圖2.

結(jié)果表明RL06版本數(shù)據(jù)在采用相同數(shù)據(jù)處理過(guò)程的情況下條帶誤差得到更好地抑制.此外還可以通過(guò)GRACE的觀測(cè)誤差來(lái)衡量反演質(zhì)量.本文采用Wahr等(2006)提出的用扣除長(zhǎng)期趨勢(shì)和周年變化后的球諧位系數(shù)殘差估算觀測(cè)誤差的方法，得到RL05版本數(shù)據(jù)和RL06版本數(shù)據(jù)進(jìn)行500 km半徑高斯平滑后的等效水高觀測(cè)誤差平均值分別為26.8 mm和22.1 mm，RL06版本數(shù)據(jù)的觀測(cè)誤差小于RL05版本.RL06版本數(shù)據(jù)在條帶誤差的抑制和觀測(cè)誤差上均優(yōu)于RL05版本，信噪比更高，因此用RL06版本數(shù)據(jù)對(duì)水儲(chǔ)量變化的估算更準(zhǔn)確.

3.2 尺度因子方法與正向建模方法對(duì)比

泄漏誤差是限制GRACE反演精度的一個(gè)重要因素.尺度因子方法廣泛應(yīng)于改正泄漏誤差，其基本原理是利用最小二乘求出尺度因子k，使得濾波前的真實(shí)信號(hào)ΔSt與濾波后的信號(hào)ΔSf的殘差平方和最小，即M=∑(ΔSt-kΔSf)2最小，可計(jì)算單一尺度因子和格網(wǎng)點(diǎn)尺度因子(Landerer and Swenson, 2012).由于Noah和CPC在湖泊水體均無(wú)數(shù)據(jù)，故本文僅計(jì)算流域單一尺度因子.利用Noah和CPC水文模式作為“真實(shí)”信號(hào)計(jì)算得到的流域單一尺度因子分別為1.34和1.27，取平均值作為流域尺度因子乘以濾波后的等效水高改正泄漏誤差.尺度因子方法的改正效果受水文模式精度影響，當(dāng)水文模式效果不佳時(shí)會(huì)降低GRACE估算結(jié)果的可靠性，正向建模方法可不依賴此類先驗(yàn)信息.為對(duì)比兩種方法在維多利亞湖流域的效果，將Mascon結(jié)果作為“真實(shí)”信號(hào)進(jìn)行與球諧位系數(shù)一致的處理，兩種方法改正泄漏誤差的效果見(jiàn)圖3.

圖3a為原始Mascon結(jié)果，圖3b為進(jìn)行與球諧位系數(shù)相同處理后的Mascon結(jié)果，圖3c為用正向建模改正泄漏誤差后的結(jié)果，圖3d為原始信號(hào)和正向建模結(jié)果的差值，圖3e為用尺度因子改正泄漏誤差后的結(jié)果.結(jié)果表明進(jìn)行與球諧位系數(shù)一致處理后的Mascon信號(hào)，信號(hào)強(qiáng)度減弱約55%，同時(shí)信號(hào)在空間上也存在泄漏，說(shuō)明高斯濾波雖然可以抑制噪聲，但是同時(shí)也會(huì)削弱真實(shí)信號(hào).經(jīng)過(guò)正向建模恢復(fù)后的信號(hào)，強(qiáng)度趨近于原始信號(hào)，與原始信號(hào)強(qiáng)度相差7%左右，而且信號(hào)空間分布上也與原始信號(hào)基本一致.相比之下尺度因子方法不論信號(hào)強(qiáng)度還是空間位置恢復(fù)結(jié)果都與原始信號(hào)有一定區(qū)別.數(shù)據(jù)模擬結(jié)果說(shuō)明正向建模方法改正維多利亞湖流域泄漏誤差的效果更好.

正向建模方法的原理是通過(guò)迭代使模型值趨近于原始信號(hào).同樣以Mascon結(jié)果作為“真實(shí)”信號(hào)，以濾波后的Mascon信號(hào)作為初始值，計(jì)算迭代100次過(guò)程中模型值流域平均值與觀測(cè)值流域平均值的差值，分析迭代次數(shù)對(duì)迭代結(jié)果的影響，結(jié)果見(jiàn)圖4.結(jié)果表明在前20次迭代過(guò)程中，模型值與觀測(cè)值之差迅速減小，之后的迭代過(guò)程中差值趨于穩(wěn)定迭代收斂.

圖2 GRACE RL05與RL06版本球諧位系數(shù)計(jì)算得到的2003-01—2017-06全球等效水高變化速率(mm·a-1)Fig.2 Global terrestrial water storage trend from January 2003 to June 2017 estimated from GRACE Release-05 spherical harmonic coefficient products and Release-06 spherical harmonic coefficient products (in mm·a-1 of equivalent water height)

圖3 以2003-01—2017-06間的維多利亞湖流域GRACE Mascon等效水高變化速率作為“真實(shí)信號(hào)”分別使用正向建模和尺度因子方法改正泄漏誤差效果對(duì)比(mm·a-1)Fig.3 Comparison between the leakage correction effect by forward modeling method and scale factor method base on using the trend of GRACE Mascon equivalent water height in Lake Victoria Basin between January 2003 and June 2017 as the simulated “real signal” (in mm·a-1 of equivalent water height)

圖4 正向建模中迭代次數(shù)的影響，差值為流域平均值Fig.4 The effect of the number of iterations on the forward modeling result, the difference is the difference of the basin average

3.3 陸地水儲(chǔ)量變化

為了與GRACE球諧位系數(shù)計(jì)算結(jié)果對(duì)比，將水文模式展開(kāi)成60階的球諧位系數(shù)，進(jìn)行500 km半徑的高斯濾波，計(jì)算成1°×1°格網(wǎng)形式的等效水高，然后用無(wú)空間約束的正向建模方法進(jìn)行泄漏誤差改正.對(duì)衛(wèi)星測(cè)高數(shù)據(jù)假設(shè)整個(gè)湖泊表面都以一致的速率變化，然后球諧展開(kāi)至60階后進(jìn)行與SH相同的濾波處理.流域陸地水儲(chǔ)量區(qū)域平均值變化時(shí)間序列如圖5，所有時(shí)間序列均扣除了2003-01—2017-06間的平均值.從圖中可得GRACE與衛(wèi)星測(cè)高探測(cè)到的趨勢(shì)較為一致， 2003-01—2017-06間流域水儲(chǔ)量均呈增加趨勢(shì)，且均探測(cè)到陸地水儲(chǔ)量在2003—2006和2016—2017有明顯的減少趨勢(shì).水文模式與GRACE和衛(wèi)星測(cè)高結(jié)果差異較大，水文模式在2003-01—2017-06間無(wú)明顯變化趨勢(shì)，雖然能探測(cè)到2003—2006的陸地水減少趨勢(shì)，但變化幅度遠(yuǎn)小于GRACE和衛(wèi)星測(cè)高結(jié)果.

為比較不同時(shí)間序列間的相似程度，分別計(jì)算了不同時(shí)間序列間的相關(guān)系數(shù)和納什效率系數(shù)，相關(guān)系數(shù)見(jiàn)表1，納什效率系數(shù)見(jiàn)表2.結(jié)合兩表中的數(shù)據(jù)得出相比于RL05 SH，RL06 SH與衛(wèi)星測(cè)高和Mascon更接近，能更準(zhǔn)確地反映水儲(chǔ)量變化.Mascon結(jié)果相比球諧位系數(shù)結(jié)果與衛(wèi)星測(cè)高更接近，認(rèn)為Mascon結(jié)果能更準(zhǔn)確地反映水儲(chǔ)量變化.Noah和CPC較為一致，但與GRACE和衛(wèi)星測(cè)高結(jié)果區(qū)別較大.

表1 2003-01—2017-06間GRACE、水文模式陸地水儲(chǔ)量和衛(wèi)星測(cè)高水位高度時(shí)間序列間的相關(guān)系數(shù)

表2 2003-01—2017-06間GRACE、水文模式陸地水儲(chǔ)量和衛(wèi)星測(cè)高水位高度時(shí)間序列間的納什效率系數(shù)

為獲取陸地水儲(chǔ)量的年際變化，對(duì)時(shí)間序列進(jìn)行13個(gè)月窗口的滑動(dòng)平均以減弱季節(jié)性的周期信號(hào)，結(jié)果見(jiàn)圖6.根據(jù)滑動(dòng)平均結(jié)果可以將本文研究的時(shí)間段分為5段，即2003-01—2006-05水儲(chǔ)量明顯減少，2006-06—2007-06水儲(chǔ)量明顯增加，2007-07—2011-03水儲(chǔ)量無(wú)明顯變化趨勢(shì)，2011-04—2016-02水儲(chǔ)量增加和2016-03—2017-06水儲(chǔ)量明顯減少.水文模式反映出的陸地水儲(chǔ)量變化幅度比GRACE和衛(wèi)星測(cè)高小，這主要因?yàn)樗哪Ｊ椒从车闹饕且驓夂蜃兓鸬耐寥浪椭脖还趯雍兓?，不包含地下水、湖泊水體和人類活動(dòng)影響.此外水文模式與GRACE和衛(wèi)星測(cè)高在2013—2016的變化趨勢(shì)存在明顯差異，GRACE和衛(wèi)星測(cè)高在此時(shí)間段探測(cè)到水儲(chǔ)量增加而水文模式探測(cè)到水儲(chǔ)量減少，推測(cè)這一現(xiàn)象由人類活動(dòng)造成.此前的研究結(jié)果表明人類活動(dòng)對(duì)維多利亞湖流域水儲(chǔ)量的影響極大，2003-01—2006-05的水儲(chǔ)量減少是由干旱和大壩排水量加大共同造成(Awange et al., 2013; Nanteza et al., 2016; Swenson and Wahr, 2009).但由于大壩排水量數(shù)據(jù)不實(shí)時(shí)公開(kāi)，因此判斷水儲(chǔ)量變化是否受到很大的人為影響有一定困難，本文將在下一節(jié)利用公開(kāi)的多源衛(wèi)星數(shù)據(jù)分析造成2013—2017年水儲(chǔ)量變化可能的原因.

為定量估算陸地水的長(zhǎng)期變化趨勢(shì)，按公式(1)擬合2003-01—2017-06間的長(zhǎng)期趨勢(shì)和周期信號(hào)的振幅和相位

ΔH=β0+β1t+β2cos(ω1t+φ1)+β3cos(ω2t+φ2)

+β4sin(ω3t),

(1)

式中，β1為長(zhǎng)期趨勢(shì)，β2、β3和β4為周年、半周年和161天周期項(xiàng)信號(hào)的振幅，φ1和φ2為周年和半周年信號(hào)的相位，ω1=2π，ω2=4π，ω3=4.5342π分別代表一年周期項(xiàng)、半年周期項(xiàng)和161天周期項(xiàng)，擬合結(jié)果見(jiàn)表3.

根據(jù)擬合結(jié)果，RL06版本數(shù)據(jù)在長(zhǎng)期趨勢(shì)上比RL05版本數(shù)據(jù)更接近衛(wèi)星測(cè)高結(jié)果(見(jiàn)表3)，RL05的結(jié)果低估了該流域的長(zhǎng)期趨勢(shì).在2003年至2017年間 RL06 SH、Mascon和測(cè)高長(zhǎng)期趨勢(shì)較接近，Mascon更接近衛(wèi)星測(cè)高的變化速率，說(shuō)明Mascon相比球諧位系數(shù)在該流域質(zhì)量變化的信號(hào)更強(qiáng)，反演結(jié)果更真實(shí).GRACE和測(cè)高得到的趨勢(shì)較一致，說(shuō)明兩種數(shù)據(jù)主要包含的質(zhì)量變化信息一致即湖泊水體.水文模式的變化相比湖泊水體明顯要小，說(shuō)明湖泊水量的變化是流域水儲(chǔ)量變化信號(hào)的主要部分.維多利亞湖流域一年存在兩個(gè)雨季，分別為3月至6月的長(zhǎng)雨季和10月至12月的短雨季，兩個(gè)雨季陸地水儲(chǔ)量變化相差較大，擬合時(shí)將一年內(nèi)兩個(gè)振幅相差較大的信號(hào)作為一個(gè)年周期信號(hào)擬合時(shí)，振幅會(huì)比長(zhǎng)雨季的小，這可能是造成流域內(nèi)半周年振幅大于周年振幅的原因.

圖6 維多利亞湖2003-01—2017-06間GRACE等效水高、衛(wèi)星測(cè)高湖泊水位高和水文模式陸地水變化時(shí)間序列13個(gè)月滑動(dòng)平均值，以mm等效水高為單位Fig.6 Smoothed time series of terrestrial water storage in Lake Victoria Basin with 13 months windows detected by GRACE, satellite altimetry and land surface models from January 2003 to June 2017 (in mm of equivalent water height)

表3 2003-01—2017-06間維多利亞湖流域GRACE、水文模式和衛(wèi)星測(cè)高水儲(chǔ)量長(zhǎng)期趨勢(shì)和周期信號(hào)擬合結(jié)果Table 3 Fitting results for trends and periodic changes of terrestrial water storage detected by GRACE, land surface models and satellite altimetry in Lake Victoria Basin from January 2003 to June 2017

2003-01—2017-06陸地水儲(chǔ)量和降雨變化速率分布見(jiàn)圖7.圖7顯示降雨在維多利亞湖流域無(wú)明顯長(zhǎng)期變化趨勢(shì).GRACE在維多利亞湖流域探測(cè)到明顯的增長(zhǎng)趨勢(shì)，在馬拉維湖則存在一個(gè)減少的趨勢(shì)，在維多利亞湖流域的西南方向也存在陸地水增長(zhǎng)信號(hào).相比之下利用水文模式均未探測(cè)到明顯的變化趨勢(shì).在水文模式效果較差的湖泊水體表面，GRACE與水文模式得到的變化趨勢(shì)差異明顯，這是由于水文模式僅包含土壤水，而GRACE是土壤水、地表水體和地下水變化之和.綜合長(zhǎng)期趨勢(shì)和周期信號(hào)的擬合結(jié)果以及空間分布，說(shuō)明在此區(qū)域水文模式與GRACE結(jié)果一致性不高.

根據(jù)滑動(dòng)平均結(jié)果將2003-01—2017-06分為5個(gè)時(shí)間段，按公式(1)擬合陸地水儲(chǔ)量在各個(gè)時(shí)間段內(nèi)的變化趨勢(shì)，結(jié)果見(jiàn)表4.

表4 2003-01—2017-06間維多利亞湖流域陸地水儲(chǔ)量5個(gè)時(shí)間段變化速率擬合結(jié)果(單位 mm·a-1)Table 4 Fitting results of five different periods for trends of terrestrial water storage in Lake Victoria Basin from January 2003 to June 2017 (unit: mm·a-1)

GRACE球諧位系數(shù)得到2003-01—2006-05維多利亞湖流域水儲(chǔ)量以79.5 mm·a-1的速率減少，2006-06—2007-06以148.9 mm·a-1的速率增加.本文結(jié)果與Awange等(2014)得出的2003-10—2006-03間84.5 mm·a-1的減少速率和2006-03—2007-06間145.2 mm·a-1的增加速率較一致，數(shù)值上的差異是由于采用了不同版本的GRACE時(shí)變重力場(chǎng)模型和不同的濾波方法.

為進(jìn)一步分析GRACE數(shù)據(jù)、土壤水和降雨數(shù)據(jù)間的關(guān)系，用PCA(Principal Component Analysis)方法分離空間模態(tài)和時(shí)間模態(tài)(Preisendorfer, 1988).圖8為2003-01—2017-06間GRACE、Noah和TRMM數(shù)據(jù)的PCA結(jié)果，圖8a為RL06 SH的時(shí)間模態(tài)和空間模態(tài)，圖8b為Mascon的時(shí)間模態(tài)和空間模態(tài)，圖8c為Noah的時(shí)間模態(tài)和空間模態(tài)，圖8d為降雨的時(shí)間模態(tài)和空間模態(tài).SH第一模態(tài)信號(hào)占比55%，第二模態(tài)信號(hào)占比19%，Mascon第一模態(tài)信號(hào)占比54%，第二模態(tài)信號(hào)占比25%，Noah第一模態(tài)信號(hào)占比57%，第二模態(tài)信號(hào)占比18%，TRMM第一模態(tài)信號(hào)占比90%，第二模態(tài)信號(hào)占比5%.四組數(shù)據(jù)第一模態(tài)呈周期信號(hào)，第二模態(tài)與年際變化信號(hào)相似.相比陸地水儲(chǔ)量，降雨的周期信號(hào)占比更高，年際變化占比較低.GRACE和Noah周期信號(hào)的波谷出現(xiàn)在每年的9月前后，波峰出現(xiàn)在3月前后， TRMM數(shù)據(jù)的波谷出現(xiàn)在6月前后，波峰出現(xiàn)在1月前后，對(duì)比陸地水和降雨波峰和波谷出現(xiàn)的時(shí)間，陸地水變化相比降雨有兩個(gè)月左右的滯后.GRACE、Noah和TRMM在2006年前后處于多年較低值，與氣象資料記錄的該地區(qū)在這一時(shí)期處于干旱狀態(tài)一致.GRACE的第二模態(tài)在維多利亞湖流域及其東南方的馬拉維湖流域呈現(xiàn)明顯的變化信號(hào).Noah和TRMM第二模態(tài)的變化信號(hào)主要出現(xiàn)在流域西側(cè)，Noah在馬拉維湖流域也呈現(xiàn)出變化信號(hào)但信號(hào)強(qiáng)度低于GRACE.GRACE與Noah和TRMM變化信號(hào)的空間分布存在一定差異，這可能是因?yàn)镚RACE探測(cè)到的水儲(chǔ)量變化受流域集水區(qū)的徑流、降雨和蒸散的共同影響，而不僅僅與流域范圍內(nèi)的降雨有關(guān).

圖8 2003-01—2017-06間GRACE等效水高、Noah陸地水和TRMM降雨主成分分析結(jié)果，圖a為GRACE SH結(jié)果，圖b為GRACE Mascon結(jié)果，圖c為GLDAS-Noah結(jié)果，圖d為TRMM降雨結(jié)果Fig.8 Principal component analysis result of GRACE TWS, GLDAS-Noah TWS and TRMM precipitation from January 2003 to June 2017, picture a is the GRACE SH result, picture b is the GRACE Mascon result, picture c is the GLDAS-Noah result, picture d is the TRMM rainfall result

圖9 2003-01—2017-06間維多利亞湖流域GRACE-SH陸地水儲(chǔ)量、降雨和蒸散變化，蒸散數(shù)據(jù)時(shí)間跨度為2003-01—2014-09，等效水高變化以mm/month為單位，累積降雨和蒸散以mm為單位Fig.9 Comparison between GRACE-SH TWS, precipitation and evapotranspiration from January 2003 to June 2017, the time span of evapotranspiration data is from January 2003 to September 2014, TWS change is measured in mm/month, the annual cumulative precipitation and evapotranspiration are measured in mm

圖10 2003-01—2017-06間Nio3.4指數(shù)與維多利亞湖流域累積降雨的對(duì)比，Nio3.4指數(shù)經(jīng)過(guò)5個(gè)月窗口滑動(dòng)平均以攝氏度為單位，累積降雨以mm為單位Fig.10 Comparison between Nio3.4 Index and annual cumulative precipitation from January 2003 to June 2017, the Nio3.4 index is smoothed with 5 months windows, and measured in degrees Celsius, the annual cumulative precipitation is measured in mm

3.4 陸地水儲(chǔ)量變化原因分析

陸地水變化是一段時(shí)間內(nèi)降雨P(guān)、蒸散E和徑流R之和，水平衡方程為ΔS=P-E-R.由于未獲取到徑流數(shù)據(jù)，且維多利亞湖流域的水儲(chǔ)量大部分通過(guò)降雨補(bǔ)充，水儲(chǔ)量大部分通過(guò)蒸散損失(Awange and Ong′ang′a, 2006)，故僅對(duì)比降雨、蒸散與等效水高.圖9為SH等效水高變化與降雨和蒸散的對(duì)比，dTWS代表的是等效水高逐月變化，時(shí)間序列均扣除了對(duì)應(yīng)時(shí)間段的平均值.圖9a是陸地水儲(chǔ)量變化ΔS和P-E的時(shí)間序列，圖9b是經(jīng)過(guò)13個(gè)月窗口滑動(dòng)平均后的ΔS和P-E的時(shí)間序列.圖9c和圖9d分別是經(jīng)過(guò)13個(gè)月窗口滑動(dòng)平均后的TWS和月降雨與月蒸散的對(duì)比.陸地水儲(chǔ)量變化受一段時(shí)間內(nèi)總的降雨和蒸散變化的影響，相比氣候變化有一定的延遲，故本文還對(duì)比了流域內(nèi)陸地水儲(chǔ)量變化和年累積降雨與蒸散，圖9e和圖9f分別是經(jīng)過(guò)13個(gè)月窗口滑動(dòng)平均后的TWS和年累積降雨與年累積蒸散的對(duì)比.

圖9a和圖9b顯示P-E和ΔS兩者間一致性較好，但兩組時(shí)間序列的變化幅度存在差異，一定程度說(shuō)明降雨和蒸散變化不是造成陸地水儲(chǔ)量變化的唯一原因.2003—2005年累積降雨低于平均水平，年累積蒸散也較多，這是造成這一時(shí)期陸地水儲(chǔ)量減少的原因之一.2006年降雨明顯增加，蒸散也減少，因此水儲(chǔ)量迅速增加.

當(dāng)然，粉絲社群內(nèi)部的這種沖突關(guān)系并不都指向破壞性，沖突關(guān)系的正功能其實(shí)也很鮮明。在爭(zhēng)端白熱化后，粉絲社群內(nèi)部之間的多種沖突關(guān)系便以一種儀式化的方式存在，粉絲們把對(duì)對(duì)方的言語(yǔ)攻擊當(dāng)成是每日實(shí)踐粉絲身份的必修課，他們認(rèn)為對(duì)立立場(chǎng)的表達(dá)是鞏固自己群體、維護(hù)自身利益重要方式。實(shí)際上，粉絲社群內(nèi)部也達(dá)成一種默契的共識(shí)，矛盾關(guān)系豐富了社群內(nèi)部溝通和交流的多樣性，保證了社群的活躍程度，維系了社群的日常。由此，粉絲內(nèi)部社群的沖突關(guān)系逐漸成為有意為之的建構(gòu)。

本文上一節(jié)得出水文模式探測(cè)到維多利亞湖流域2013—2016年水儲(chǔ)量減少，而GRACE和衛(wèi)星測(cè)高探測(cè)到這一時(shí)期流域內(nèi)水儲(chǔ)量增加.圖9顯示在這一時(shí)期流域內(nèi)降雨較少，水文模式僅包含土壤水和植被冠層含水而不包含湖泊水體，這部分水量主要受降雨和蒸散影響，這可能是造成水文模式探測(cè)到陸地水減少的原因.GRACE和衛(wèi)星測(cè)高相比水文模式能探測(cè)到湖泊水體的變化，因此推測(cè)造成這一時(shí)期水文模式和GRACE與衛(wèi)星測(cè)高得到的趨勢(shì)相反的原因是大壩蓄水造成湖泊水量增加，而干旱造成土壤水和植被冠層含水減少，人類活動(dòng)是造成這一現(xiàn)象的原因.

GRACE和衛(wèi)星測(cè)高結(jié)果在2016-03—2017-06存在一個(gè)減小趨勢(shì)，SH估算得到的變化速率為-100.3 mm·a-1，Mascon估算的速率為-129.7 mm·a-1，衛(wèi)星測(cè)高得到的速率為-105.4 mm·a-1，Noah和CPC也存在減少的趨勢(shì)但速率遠(yuǎn)小于前者(見(jiàn)表4).2016—2017的年累積降雨較少，這可能是造成2016—2017流域陸地水儲(chǔ)量減少的原因.2016-03—2017-06間 GRACE探測(cè)到的陸地水儲(chǔ)量變化速率比衛(wèi)星測(cè)高大，而2003-01—2017-06的變化速率則是衛(wèi)星測(cè)高較大，其可能的原因是GRACE在任務(wù)結(jié)束前期衛(wèi)星運(yùn)行和觀測(cè)狀態(tài)變差，誤差可能相對(duì)較大.此外GRACE和衛(wèi)星測(cè)高所觀測(cè)到的信息不同，在不同的時(shí)間段變化可能存在差異.

東非地區(qū)的降雨易受厄爾尼諾/南方濤動(dòng)ENSO(El Nio-Southern Oscillation)影響(Awange et al., 2013; Becker et al., 2010).為分析維多利亞湖流域降雨與ENSO現(xiàn)象的關(guān)系，本文采用Nio3.4指數(shù)衡量ENSO.圖10為Nio3.4指數(shù)與年累積降雨的對(duì)比，Nio3.4指數(shù)按定義進(jìn)行了5個(gè)月窗口的滑動(dòng)平均.2015年和2016年發(fā)生El Nio事件，2016年累積降雨處于多年較低水平，2008年和2011年前后發(fā)生La Nia事件，這兩年的降雨也較多.2015—2016年的El Nio事件在近幾年強(qiáng)度較大持續(xù)時(shí)間也較長(zhǎng)，累積降雨量減少較多，因此推測(cè)2016年后區(qū)域陸地水儲(chǔ)量變化與El Nio事件引起的降雨變化有關(guān).

4 結(jié)論

本文利用2003-01—2017-06間的GRACE RL06版本時(shí)變重力場(chǎng)模型并結(jié)合水文模式和衛(wèi)星測(cè)高等數(shù)據(jù)分析維多利亞湖流域的陸地水儲(chǔ)量變化.受限于GRACE重力衛(wèi)星以及其他數(shù)據(jù)的精度，目前利用衛(wèi)星數(shù)據(jù)高精度地監(jiān)測(cè)區(qū)域水儲(chǔ)量變化有一定困難.

數(shù)據(jù)處理過(guò)程中產(chǎn)生的泄漏誤差是影響GRACE結(jié)果可靠性的一個(gè)重要因素，因此在GRACE現(xiàn)有空間分辨率下較高精度地監(jiān)測(cè)水儲(chǔ)量變化，需要采用合適的泄漏誤差改正方法.為對(duì)比不同方法的泄漏誤差改正效果，本文將 Mascon數(shù)據(jù)作為“真實(shí)”信號(hào)分別用正向建模方法和單一尺度因子方法改正泄漏誤差，兩種方法的結(jié)果對(duì)比表明在缺乏足夠先驗(yàn)信息的情況下正向建模方法能更好地改正維多利亞湖流域的泄漏誤差.但是不論尺度因子方法還是正向建模方法對(duì)泄漏誤差的改正均有局限性.尺度因子方法嚴(yán)重依賴先驗(yàn)的水文模式精度，當(dāng)先驗(yàn)?zāi)Ｐ途炔蛔銜r(shí)改正后的GRACE結(jié)果的可靠性也會(huì)降低，而且多數(shù)水文模式都不含地下水和極地的信息，研究地下水變化和極地質(zhì)量變化時(shí)尺度因子方法效果較差.正向建模方法雖然不需要先驗(yàn)的水文信息，但是此方法在質(zhì)量變化信號(hào)分布明確的地區(qū)恢復(fù)效果才具有明顯優(yōu)勢(shì).國(guó)內(nèi)外針對(duì)泄漏誤差的改正還有多種不同方法，因此為了獲得較好的改正效果對(duì)不同的研究區(qū)域要選用合適的泄漏誤差改正方法.

本文還在缺乏直接觀測(cè)數(shù)據(jù)的情況下，利用公開(kāi)的衛(wèi)星觀測(cè)數(shù)據(jù)分析人工調(diào)控是否對(duì)近期維多利亞湖流域的水儲(chǔ)量變化造成影響，為采用相似方法研究我國(guó)湖泊水儲(chǔ)量變化提供參考.對(duì)比GRACE、衛(wèi)星測(cè)高和水文模式的時(shí)間序列，發(fā)現(xiàn)在2013-01—2016-02降雨減少的情況下，水文模式探測(cè)到流域內(nèi)水儲(chǔ)量減少而GRACE和測(cè)高探測(cè)到湖泊水量增長(zhǎng)，推測(cè)這一現(xiàn)象由大壩蓄水造成.根據(jù)本文的計(jì)算結(jié)果在2017年該流域的水儲(chǔ)量已有明顯減少趨勢(shì)，而且歷史上此湖泊的水量虧損也曾因?yàn)榧哟笈潘慷鴩?yán)重加劇，若仍不加控制地對(duì)湖泊水量進(jìn)行人工調(diào)控，該地區(qū)的水平衡極易被打破.此外由于我國(guó)湖泊流域面積大多小于維多利亞湖流域，因此只有在GRACE的空間分辨率和數(shù)據(jù)處理方法進(jìn)一步改善后，利用GRACE重力衛(wèi)星高精度地監(jiān)測(cè)我國(guó)湖泊水量變化才更具有可行性.GRACE衛(wèi)星已終止運(yùn)行， GRACE Follow-On于2018年5月發(fā)射升空，該任務(wù)可繼續(xù)為大尺度地監(jiān)測(cè)陸地水儲(chǔ)量變化提供衛(wèi)星重力數(shù)據(jù).

致謝感謝CSR提供的GRACE衛(wèi)星數(shù)據(jù)，NOAA提供的CPC水文模式，LEGOS/CNES提供的衛(wèi)星測(cè)高產(chǎn)品，NASA提供的GLDAS和TRMM數(shù)據(jù)，蒙大拿大學(xué)NTSG提供的蒸散產(chǎn)品.感謝審稿專家對(duì)本文提出的寶貴修改意見(jiàn).