基于GRACE陸地水儲量降尺度的塔里木河流域干旱特征及驅動因子分析

魏光輝，楊 鵬，周海鷹，夏 軍，陳 杰，桂東偉，巴音達拉

(1.新疆塔里木河流域管理局，新疆 庫爾勒 841000；2.中國地質大學(武漢)，地理與信息工程學院，武漢 430074；3.武漢大學 水資源與水電工程科學國家重點實驗室，武漢 430072；4.中國科學院地理科學與資源研究所，北京 100101;5.中國科學院新疆生態與地理研究所，烏魯木齊 830011；6.新疆農業大學水利與土木工程學院，烏魯木齊 830052)

0 引 言

塔里木河流域地處亞歐大陸腹地，遠離海洋，氣候干旱，是中國最大的內陸河流域，綠洲灌溉農業為其主要的農業發展模式，由于水資源的短缺，其間存在土壤鹽漬化、水資源流失、生態系統退化和破壞等嚴重的環境問題[1-3]。在氣候變化的背景下，水文氣象極端事件更加頻繁發生，對于靠冰雪融水補給為主的流域而言，非常不利于水資源的開發利用[3]，因此，氣候變化也加劇了塔里木河流域水安全問題[4]。雖然塔里木河流域降水和出山口徑流在近期有增加趨勢，但是從上游向中下游輸送的徑流卻明顯減少[3]。因而，分析塔里木河流域水資源量變化十分重要。2013年9月“一帶一路”經濟發展戰略的提出，以實現跨越亞非歐大陸國家和地區的區域經濟一體化為目標，推進國際合作，打造陸海內外聯動、東西雙向開放的全面開放新格局，也為我國西北干旱區經濟發展提供契機[5]。因此，分析和研究塔里木河流域水資源的變化對于該區生態環境保護和社會經濟的發展具有至關重要的作用。

陸地水儲量是全球水、能量以及生物化學循環的重要組成部分[6]，它主要包括地表水、土壤水和地下水，其中地表水包括河川徑流湖泊水、雪水當量以及地表植被冠層水，也是水資源量表示的重要方法[7]。研究陸地水儲量的時空變化對于有效管理水資源和可持續發展決策具有十分重要的意義。然而，由于地形地貌、海拔高度、地質結構的復雜，加上技術的局限，人工觀測陸地水儲量的變化存在較多困難，尤其是大尺度區域的陸地水儲量的觀測更是很難實施。然而， 2002年重力衛星GRACE的發射為大尺度陸地水儲量的觀測提供了重要工具和手段，且已經被廣泛地應用到水文與水資源領域的研究中。例如：Velicogna和Wahr[8,9]基于重力衛星GRACE數據研究發現格陵蘭島南部地區的冰蓋在加速減少；Longuevergne等[10]和Voss等[11]利用重力衛星GRACE數據分別研究了中東地表水和地下水的變化；周旭華等[12]基于重力衛星GRACE數據反演了全球陸地水儲量的變化，并結合全球陸地同化系統GLDAS模型數據進行了對比分析和研究；鐘敏等[13]利用GRACE數據發現中國華北平原、三峽庫區、青藏高原、新疆阿爾金山和青川甘交界處陸地水儲量變化幅度較大，且京津冀地區、青藏高原地區陸地水儲量呈現明顯的減少趨勢。但是，重力衛星GRACE固有的空間分辨率較低及時間序列較短的特點嚴重限制了其應用的推廣。盡管Long等[14]和Zhang等[15]應用人工智能方法結合水文氣象數據對陸地水儲量進行重建，延長了陸地水儲量的時間序列，其較低的空間分辨率依然未得到有效解決。重力衛星GRACE陸地水儲量空間分辨率的提高，需要通過科學合理的數學物理方法對其進行降尺度處理，即通過科學合理的數學物理方法在不影響重力衛星GRACE陸地水儲量信息的前提下，提高其空間分辨率，以便其在中小尺度流域的應用與推廣。除此之外，在水文領域重力衛星GRACE大多數情況都被用作分析陸地水儲量的時空變化，其水文研究應用的潛力還亟待提高和加強。

為了研究流域尺度陸地水儲量的時空變化及其影響下的流域干旱特征，本文嘗試利用重力衛星GRACE數據和全球陸地同化系統模型數據GLDAS，結合統計學方法，分別以數據格點和等效水高為權重，對塔里木河流域的陸地水儲量進行空間降尺度，分析塔里木河流域陸地水儲量的時空變化，隨后進行降尺度后的陸地水儲量分析塔里木河流域的干旱特征分析，簡要分析其驅動因素。

1 研究區概況

塔里木河流域是中國最大的內流河流域，其面積達102 萬km2，坐落于73°E～93°E，35°N～45°N之間(如圖1)，屬于典型的溫帶大陸性氣候，因此流域內降水極少，蒸散發卻極其強烈。另外流域內溫度變化幅度較大，常年冬季溫度跨度在-20和10 ℃之間，夏季溫度跨度在20和30 ℃之間[16]。就流域內的降水的空間分布而言也存在較大差異，例如：流域內山區降水較多，而平原區降水卻在較少[17-19]。另外，流域內的年溫度也表現出較大的年際差異[3]。

隨著20世紀70年代以來氣溫的上升[20]，流域內固態降水減少且冰雪融水增加。然而，對于干旱區而言山區降水及來自于山區的冰雪融水在水資源系統和水文過程起著十分重要的作用[21]。塔里木河流域的年均徑流量約為3.989 × 1010m3, 其中來自山區的冰雪融水大約占48.2%[22]。作為流域內的主要源流，阿克蘇河、和田河及葉爾羌河的徑流量分別占總徑流量的73.2%，23.2%和3.6%[23]。同時，以灌溉農業為主要農業發展方式的塔里木河流域的社會經濟的可持續發展深受流域內水資源波動和供給能力的影響[23]。

2 數據與方法

2.1 數 據

(1)GRACE數據。重力衛星 GRACE數據源自http:∥grace.jpl.nasa.gov。本文所涉及的陸地水儲量的等效水高是從第二版重力衛星球諧系數RL05轉換得到。同時，在等效水高轉換過程中還有一些其他處理細節，如用衛星激光測量數據(SLR)替代C20項[24]，重新估算1階系數[25]，冰川后回彈系數矯正[26]等。另外， 考慮到重力衛星GRACE球諧系數數據在2004年1月至2009年12月之間不缺失[27]，因此，為了得到球諧系數的距平值，本文從原始的球諧系數中去掉了該時間段的平均值。盡管遺落誤差、測量誤差等誤差存在于重力衛星GRACE的數據中，然而Chen等[25]卻證明了GRACE數據具有一定的適用性和可靠性。

(2)GLDAS數據。全球陸地同化系統模型(GLDAS)數據旨在精確反演陸地表面狀況[28]。同時，在該系統由3個陸面模型(Noah, Mosaic, Community Land Model-CLM)和一個水文模型(Variable Infiltration Capacity-VIC)組成。本文所用的數據來自于網站(http://mirador.gsfc.nasa.gov/)，文中所涉及的土壤水來自于Noah模型中的時間分辨率為月，空間分辨率為1°×1°和0.25°×0.25°兩種數據。該數據被廣泛的用到全球及各個流域的水文科學的研究領域中，具有較強的可靠性和普適性[29]。為了能夠與GRACE數據進行對比，也將2004年1月至2009年12月的平均值從原始時間序列中扣除。

(3)塔里木河流域水資源數據。本文采用的塔里木河流域水資源數據主要來自于由塔里木河流域管理局編寫的《塔里木河流域水資源公報》，其中主要涉及塔里木河流域地下水資源量數據、塔里木河流域供水量數據以及塔里木河流域用水量數據，該數據由塔里木河流域管理局出版和提供，已被生態和水文領域中的很多專家和課題組應用，具有可靠性。

(4)夏季0 ℃層高度數據及MODIS積雪產品數據。夏季(6-8月)0 ℃層高度主要基于庫車、喀什和和田3個探空站數據，來源于懷俄明大學氣象數據中心(http：// www. weather. uwyo. edu/ upperair/ sounding.html)，可免費獲取，時間分辨率為天，另外，本文是基于MODIS積雪產品提取的塔里木河流域夏季的積雪覆蓋率，該數據時間分辨率為月，空間分辨率為500 m，并且該數據可以從NASA相關網站獲取，具體為http：// modis-land. gsfc. nasa. gov/ snow. html。國內外許多專家和學者基于該數據對塔里木河流域進行過氣候變化的分析和研究[30]。

2.2 方 法

2.2.1 陸地水儲量估算方法

本文根據Warh等[9]提出的方法，將重力衛星GRACE提供的重力場球諧系數數據轉化為等效水高數據，等效水高的具體計算方法如公式(1)：

ΔSnmsin(mφ)]Pnmsinθ

(1)

式中：H為陸地水儲量的等效水高，mm；θ為緯度；φ為經度；α為赤道半徑，km；ρave為地球密度，kg/m3；kn為勒夫數(常量)；Cnm和Snm分別為重力場球諧系數的余弦和正弦項(階)；Pnmsinθ為第n階、第m次標準化后的勒讓德函數，其最高展開階數是60階。

但是，通常在計算陸地水儲量的同時還得考慮重力衛星GRACE中的誤差(遺漏誤差和測量誤差)和尺度因子，其中遺漏誤差可以通過經濾波處理前后的月均GLDAS-Noah陸地水儲量的均方根誤差(RMSD)獲得[31]，測量誤差可由研究區域與同緯度海洋區域的重力差值對應的等效水高獲得[31]，因此，經誤差和尺度因子校正后的等效水高為：

Hc=[H-(L+M)]S

(2)

式中：Hc為校正后的等效水高，mm；H為經球諧系數直接求得的等效水高，mm；L和M分別為遺漏誤差和測量誤差，mm；S為尺度因子(標量)，可以從JPL數據分發網站獲得(https:∥grace.jpl.nasa.gov/data/get-data/jpl_global_mascons/)。

值得注意的是重力衛星GARCE三級產品(即GRACE-Level 3)已經發布，它充分的考慮了GRACE數據處理成陸地水儲量過程中會出現的各種誤差和尺度因子，能很好地反映全球各區域陸地水儲量的變化[32],可以從網上免費獲得(http:∥www2.csr.utexas.edu/grace/RL05_mascons.html)。

2.2.2 陸地水儲量降尺度方法

全球陸面模式的等效水柱包含地表水和土壤水，然而這些陸面模式并未考慮地下水通量的變化，考慮到GRACE陸地水儲量的實際空間分辨率約為4°，首先對GRACE陸地水儲量數據進行升尺度處理，將空間分辨率1°×1°的GRACE陸地水儲量升尺度至4°×4°陸地水儲量，并且將0.25°的全球四個陸面模式的等效水柱升尺度至4°×4°。然后根據Wan等[20]提出的方法，以數據格點間的面積為權重因子，對GRACE陸地水儲量進行降尺度處理，具體方法如下：

sM=∑(Si×ai)/∑ai=∑(Si×ai)/∑A

(3)

式中：SM是4°×4°的陸面模式的陸地水儲量距平值的等效水高，mm；ai是0.25°×0.25°數據格點i的面積，m2；A是包含0.25°×0.25°數據格點i的4°×4°數據格點的總面積，m2。

假定GRACE陸地水儲量(mm)是可靠的，因此空間分辨率為4°×4°陸面模式水儲量距平值(SM)與4°×4°GRACE陸地水儲量距平值(SG)之間的差異可以認為是二者之間的偏差Bias(B)：

B=SM-SG

(4)

因此，基于格點等效水柱權重，進一步將4°×4°陸面模式水儲量距平值(SM)與4°×4°GRACE陸地水儲量距平值(SG)之間的總水柱差異(B×A)分解到0.25°×0.25°數據格點上：

bi={B×A[(Soi×ai)/∑(Soi×ai)]}/ai=

(B×A×Soi)/∑(Soi×ai)

(5)

式中：bi為空間分辨率0.25°×0.25°等效水高的偏差，mm；Soi為0.25°×0.25°求距平值之前的等效水高，mm。

(6)

2.2.3 基于陸地水儲量的水文干旱計算方法

首先，本文將陸地水儲量的距平值轉化為陸地水儲量變化，然后，根據陸地水儲量變化與多年同期陸地水儲量變化平均值的差值確定干旱期，具體公式表示如下：

TWSCj=TWSAj-TWSAj-1

(7)

式中：TWSCj為陸地水儲量的變化，mm；TWSAj表示整個時間序列中第j個月的陸地水儲量的距平值，mm；TWSAj-1表示整個時間序列中第j-1個月的陸地水儲量距平值，mm。

(8)

隨后，根據游程定理確定陸地水儲量干旱特征(干旱月數、干旱頻率、干旱次數及干旱強度)，TWSCA為負值確定為干旱發生，干旱月數即TWSCA為負值時的月份數量累加[公式(9)]，隨后，干旱頻率為干旱月份數量占總時間段的百分比[公式(10)]，干旱次數為所有干旱發生事件次數的累加[公式(11)]，然后，干旱強度為所有干旱事件中TWSCA累加的相反數[公式(12)]，具體表示為如下公式：

Dmonth=∑monthTWSCA<0

(9)

Dfrequency=Dmonth/Tmonth×100%

(10)

DTimes=∑eventTWSCA<0

(11)

Dintensity=-∑TWSCA, (TWSCA<0)

(12)

式中：Dmonth為干旱月數，月；Dfrequency為干旱頻率，%；DTimes為干旱次數，次；Dintensity為干旱強度，標量；monthTWSCA<0為干旱發生時的月份；Tmonth為總月數；eventTWSCA<0為干旱事件。

3 結果與討論

3.1 流域陸地水儲量變化

根據GRACE數據反演出塔里木河流域的陸地水儲量，研究發現流域北部為虧損狀態，南部為盈余狀態，虧損狀態最大達到-15 mm，盈余狀態最大達到10 mm(圖2)，這與Yang等[4]的研究結果較為一致。2002-2015年間整個時段內塔里木河流域陸地水儲量呈下降趨勢，下降速度達-1.6±1.1 mm/a，陸地水儲量最小值出現在2008年年底和2009年年初，最大值出現在2005年和2010年，在2005、2006、2010和2011年間陸地水儲量呈現盈余狀態，然而在2007、2008、2009和2014年間陸地水儲量呈現虧損狀態。另外，塔里木河流域內的GLDAS和GRACE反演的土壤水和陸地水儲量變化趨勢、相位和振幅較為一致，線性相關關系R2達到0.35。同時，研究發現由CSR、GFZ和JPL發布的球諧系數反演得到的陸地水儲量變化較為一致。

圖2 塔里木河流域陸地水儲量的時空分布

對塔里木河流域內的陸地水儲量季節變化研究發現(如圖3)，流域季節陸地水儲量分別為春季>夏季>冬季>秋季。對其空間分布研究發現，春季流域大多數地區陸地水儲量呈盈余狀態，只有在流域的東北部區域陸地水儲量為虧損狀態[圖3(a)]；夏季流域北部天山地區陸地水儲量為虧損狀態，流域南部的昆侖山地區陸地水儲量呈現盈余狀態[圖3(b)]；秋季流域內陸地水儲量幾乎全大部分區域都為出虧損狀態，而在流域東南部的局部區域陸地水儲量為盈余狀態[圖3(c)]；冬季流域東北部和西南部陸地水儲量呈現虧損狀態，流域東南部陸地水儲量呈現盈余狀態[圖3(d)]。

圖3 塔里木河流域的GRCAE陸地水儲量的季節空間分布

3.2 流域陸地水儲量變化的降尺度分析

由于重力衛星GRACE反演得到的水儲量的空間分辨存在明顯的局限性，且考慮到土壤水是陸地水儲量的重要組成部分，因此根據公式(1)～(6)及結合全球陸地水儲量GLDAS數據對重力衛星GRACE陸地水儲量進行降尺度分析(如圖4)，降尺度到0.25°×0.25°的空間分辨率。研究發現，降尺度前后的陸地水儲量在空間分布上存在一定的相似性，即在塔里木河流域北部陸地水儲量呈現虧損狀態，流域南部陸地水儲量呈現盈余狀態，這與Chen 等[33]及Deng和Chen[27]的研究結果一致，與此同時，降尺度前后的陸地水儲量在時間分布上較為一致，其振幅和相位也存在一致性，兩者的線性相關性R2達0.98，相關關系達到了顯著水平(P<0.001)。

圖4 塔里木河流域的陸地水儲量降尺度分析

為了進一步分析降尺度之后的陸地水儲量的空間分布情況，對其進行季節變化分析(如圖5)。研究發現，春季流域大多數地區陸地水儲量呈盈余狀態，只有在流域的北部和西部的部分區域陸地水儲量為虧損狀態[圖5(a)]；夏季北部天山地區陸地水儲量呈現虧損狀態，流域南部的昆侖山地區陸地水儲量呈現盈余狀態[圖5(b)]；秋季流域內陸地水儲量幾乎全大部分區域都呈現出虧損狀態，而在流域東南部的局部區域陸地水儲量呈現盈余狀態[圖5(c)]；冬季流域東北部和西南部陸地水儲量呈現虧損狀態，流域東南部陸地水儲量呈現盈余狀態[圖5(d)]。其空間分布與降尺度之前的陸地水儲量的季節空間分布較為一致。

圖5 塔里木河流域的降尺度陸地水儲量的季節空間分布

3.3 基于陸地水儲量的干旱分析

為了更科學的探究塔里木河流域的干旱情況，本文基于降尺度的陸地水儲量來分析塔里木河流域的干旱特征(如圖6)。基于公式(7)～(12)研究發現2002年4月至2015年12月間干旱總時長最長及干旱頻率最高的區域主要集中在天山南坡的東段、昆侖山北坡的東段，最長干旱時長可達112個月及相對應的干旱頻率為68.2%，相比較而言，期間干旱總時長較短及干旱頻率較低的區域主要集中在和田河流域的中下游，最短干旱時長為56個月及對應的干旱頻率為34.1%；從干旱次數方面來看，天山南坡及昆侖山北部的部分區域干旱次數發生的較少，期間干旱次數最小為2次，開都-孔雀河流域發生干旱次數較多，最多達28次；從干旱強度方面來看，天山南坡及和田河流域上游區域干旱強度較大，最大可達5 383 mm，塔里木盆地中部干旱強度最小，最小幾乎為零，這可能與高山冰雪資源消融和人類水資源開發利用有關。

圖6 基于陸地水儲量的塔里木河流域干旱情況分析

3.4 流域陸地水儲量變化的驅動因素分析

為了研究塔里木河流域陸地水儲量驅動因素的變化，本文選取了與陸地水儲量相關的水文變量進行研究(即MODIS夏季積雪覆蓋率、夏季0 ℃層高度及水資源量)(如圖7和表1)，研究發現2003-2012年間塔里木河流域年平均夏季積雪覆蓋率和0 ℃層高度分別為3.12%和4 723.7 m，二者大致呈現出一個負相關的關系，如：2003-2006年間塔里木河流域夏季積雪覆蓋率較高，但是夏季0 ℃層高度卻較低，2006-2008年間塔里木河流域夏季0 ℃層高度升高，相應的積雪覆蓋率反而降低，直到2009年夏季0 ℃層高度達到最低值，并且出現了較低的積雪覆蓋率，勢必會影響整個流域的水儲量變化，直到2010年時夏季0 ℃層高度達到次高值，積雪覆蓋率也達到了較高值，這在一定程度上可以解釋Wang等[4]發現的塔里木河流域在2009年出現極端干旱和在2010年出現極端濕潤的情況，這與Chen等[33]楊鵬等[30]的研究結果一致。

圖7 2003-2012年塔里木河流域夏季積雪覆蓋率及夏季0 ℃層高度的變化

表1 塔里木河流域水資源相關變量距平值 億m3

另外，結合GLDAS的陸地水儲量分量研究發現，塔里木河流域冠層水和雪水當量呈現不顯著的降低趨勢，與陸地水儲量的變化趨勢較為一致，而土壤水卻呈上升趨勢，因此，說明該區域的地下水呈現明顯的下降趨勢(如圖8)。同時，為了驗證該結論本研究結合《塔里木河流域水資源公報數據》進行了分析和研究，研究發現2004、2007、2008、2009、2011和2012年間流域水資源量盈余為負距平狀態，2002、2003、2005、2006、2010和2013年間流域水資源量盈余呈現正距平狀態。另外，2004、2005、2009、2010、2011和2013年間流域地下水資源量為負距平狀態，2002、2003、2006、2007、2008和2012年間流域地下水資源量呈現正距平狀態。除此之外，流域水資源盈余量的距平值和地下水資源量呈下降趨勢，上述水資源相關變量基本與流域內的陸地水儲量變化較為一致，因此，塔里木河流域陸地水儲量的變化主要是因為流域內地下水變化引起的，這一研究結果與Yang等[4]的發現一致。關于陸地水儲量時空變化研究對于流域水資源管理和水旱災害預報具有十分重要的意義。

圖8 2002-2015年塔里木河流域GLDAS陸地水儲量分量的時空分布

4 結 語

本文利用重力衛星GRACE數據、全球陸地同化系統模型GLDAS數據和塔里木河水資源數據，對塔里木河流域陸地水儲量變化、陸地水儲量變化的降尺度、塔里木河流域干旱情況及驅動因素進行了研究分析，主要結論如下。

(1)根據GRACE數據反演出塔里木河流域的陸地水儲量，研究發現流域北部為虧損狀態，南部為盈余狀態，虧損狀態最大達到-15 mm，盈余狀態最大達到10 mm。2002-2015年間整個時段內塔里木河流域陸地水儲量呈下降趨勢，下降速度達-1.6±1.1 mm/a。

(2)結合全球陸地同化模型數據和統計方法可以很好地對重力衛星GRACE反演的陸地水儲量進行降尺度分析，這使得重力衛星GRACE在反演中小尺度流域陸地水儲量變化上的應用具有十分重要的意義。

(3)盡管在氣候變化背景下塔里木河流域山區來水量有增加趨勢，但是塔里木河流域水儲量赤字情況依然很嚴峻，干旱問題突出，水資源科學合理地開發和利用十分有必要。

(4)塔里木河流域陸地水儲量受夏季0 ℃層高度及夏季積雪覆蓋率的影響較為明顯，另外，流域內水資源盈余量的距平值和地下水資源距平值的變化與塔里木河流域內的陸地水儲量變化較為一致，因此，流域地下水資源的變化是流域內陸地水儲量的變化的主要因素。