三江源2003～2020年地下水儲量變化分析

張浩哲 常曉濤, 朱廣彬 周 苗 劉 偉

1 山東科技大學測繪與空間信息學院，青島市前灣港路579號，266590 2 自然資源部國土衛(wèi)星遙感應用中心，北京市百勝村1號, 100048

三江源是長江、黃河、瀾滄江的發(fā)源地，是長江、黃河中下游地區(qū)和東南亞國家生態(tài)安全和區(qū)域可持續(xù)發(fā)展的生態(tài)屏障[1]。三江源地下水儲量影響著區(qū)域生態(tài)環(huán)境變化和長江、黃河、瀾滄江中下游的地面徑流量，對其進行研究有助于了解全球變暖背景下青藏高原地區(qū)水儲量的變化特征。

傳統的地下水監(jiān)測手段主要是通過地面監(jiān)測井定期監(jiān)測地下水的水位變化，難以實現大空間、大尺度下的地下水變化監(jiān)測[2]。自2002年GRACE衛(wèi)星發(fā)射以來，利用時變重力場模型研究大范圍區(qū)域地下水的時空變化成為可能。在2017年GRACE衛(wèi)星停止服務后，美國國家航空航天局(NASA)與德國地球科學研究中心聯合研制了GRACE-FO衛(wèi)星，繼續(xù)承擔全球及區(qū)域水資源變化監(jiān)測的任務。

目前，已有大量的研究人員利用時變重力場數據研究區(qū)域水資源的變化[3-6]，但對三江源區(qū)域地下水的分析研究較少[7]。本文利用GRACE和GRACE-FO提供的2003～2020年的時變重力場模型、GLDAS模型和湖泊水量變化數據反演三江源區(qū)域地下水的時空變化，通過奇異譜分析方法(singular spectrum analysis,SSA)提取時間序列中的趨勢信號和周期信號，分別對黃河源、長江源、瀾滄江源不同區(qū)域的地下水變化進行分析，并結合區(qū)域內降水數據分析三江源地下水年度補給規(guī)律。

1 研究區(qū)概況

三江源區(qū)域位于青海省南部，地理位置為31°19′～36°12′N、89°45′～102°23′E，總面積為302 500 km2，為長江、黃河和瀾滄江的源頭匯水區(qū)(圖1)。三江源區(qū)域是青藏高原的腹地和主體，以山地地貌為主，海拔1 971～6 565 m，河流密布，湖泊、沼澤眾多，是世界上海拔最高、面積最大、濕地類型最豐富的地區(qū)，也是中國淡水資源的重要補給地。

圖1 三江源地區(qū)位置Fig.1 The location of the Three-River headwater region

2 數據和方法

2.1 數據與預處理

采用美國得克薩斯大學空間研究中心(CSR)提供的GRACE和GRACE-FO衛(wèi)星RL06數據，模型最大階數為96，時間范圍為2003-01～2020-12。GRACE與GEACE-FO之間11個月(2017-07～2018-05)的數據缺失不作處理。此外，由于GRACE衛(wèi)星在前期調試運行階段及后期出現傳感器性能下降和供能不足的問題造成20個月的模型缺失，對此，采用三次樣條法進行插值[8]。

GLDAS數據選取GLDAS-NOAH水文模型，時間分辨率為1個月，空間分辨率為1°× 1°，時間范圍與時變重力場數據時間范圍一致。本研究主要提取水文模型中的4層土壤水(0～200 cm)、雪水當量和冠層水變化。

湖泊數據采用國家青藏高原科學數據中心的青藏高原2000～2017年高時間分辨率湖泊水位及水量變化數據集[9]，以及Cooley等[10]提供的2018～2020年全球湖泊水量變化數據集。從數據集中提取三江源地區(qū)13個大、中型湖泊的水量數據，湖泊位置如圖1所示。

降雨和溫度數據來自國家氣象科學數據中心提供的中國地面氣溫月值0.5°×0.5°格點數據集(V2.0)和中國地面降水月值0.5°×0.5°格點數據集(V2.0)，時間范圍為2003-01～2020-12。

2.2 理論方法

2.2.1 地下水儲量變化反演

由于GRACE和GRACE-FO衛(wèi)星對重力場低階項不敏感，且模型中含有高階噪聲、奇偶階相關誤差等[11]，在計算之前需要對時變重力場數據進行預處理，步驟如下：1)用SLR(satellite laser ranging)得到的C20、C30項代替模型中的C20、C30項[12]；2)采用半徑為300 km的高斯濾波對模型系數進行處理，抑制高階噪聲[11]；3)用滑動窗進行相關濾波，消除模型系數間的相關誤差[13]。得到濾波后的陸地水反演公式為：

(ΔClmcos(mλ)+ΔSlmsin(mλ))

(1)

經式(1)處理后的等效水高值仍受到冰后回彈和誤差泄露的影響，可分別利用ICE-6G-C模型和正演模型消除冰后回彈和誤差泄露的影響[14]。

依據水量平衡公式，三江源地區(qū)的地下水儲量變化可表示為：

ΔGWS=ΔTWS-ΔSMS-

ΔSWE-ΔCWS-ΔLAKE

(2)

式中，ΔTWS為陸地水儲量變化，由時變重力場模型計算得到；ΔSMS為土壤濕度含水量變化，ΔSWE為雪水當量變化，ΔCWS為植被冠層含水量變化，這3個變量都可以從GLDAS的NOAH模型中獲得；ΔLAKE為湖泊水量變化，由湖泊水量數據集計算得到。

將13個湖泊分別看作是離散的點質量，收集到范圍內13個湖泊的經緯度位置和水量變化速率，得到三江源內13個等效點質量及其位置。根據式(3)得到球諧展開系數，并利用式(4)得到空間域的湖泊水量變化速率[15]：

(3)

(4)

式中，n為點質量總數，n=13；ΔMi為第i個點質量變化速率；θi、λi分別為第i個點的余緯和經度；N為球諧系數截斷階數。

2.2.2 奇異譜分析

奇異譜分析(SSA)是一種處理非線性時間序列的方法，通過對所研究時間序列的軌跡矩陣進行分解、重構等操作，提取出時間序列中趨勢分量、周期分量、噪聲分量等[16]。SSA的簡要處理過程為：構建軌跡矩陣、進行奇異值分解(SVD)、分組、對角平均化，具體公式可參考文獻[17]。

3 結果與分析

3.1 整體時變分析

2003～2020年整個三江源區(qū)域地下水儲量的平均變化時間序列如圖2所示?？紤]到時間序列的連續(xù)性，在進行SSA分析時選擇2003-01～2017-06的數據，窗口長度選為48個月(圖3)。

圖2 2003～2020年三江源區(qū)域地下水平均變化時間序列Fig.2 The time series of the average change of the groundwater in the Three-River headwater region from 2003 to 2020

圖3 SSA結果Fig.3 The results of SSA

對前15階重構成分(reconstruction components,RC)進行w-correlation分析，結果如圖3 (a)所示。從相關性分析結果中可以看出，第5個及之后的RC不能很好地相互分離，說明其中噪聲信息占比很大；RC1為趨勢項；RC3、RC4的相關性很好，歸為同一信號，為周期項；RC2中含有部分趨勢項和周期項，二者沒有很好地分離。

RC1的趨勢分量如圖3(b)所示，可以看出，2013-01之前三江源地下水儲量有一個緩慢增加的趨勢，2013-01之后有一個緩慢減少的趨勢，故以2013-01為界對三江源地下水變化進行空間分析。

RC3和RC4組合的周期分量如圖3(c)所示，可以看出，地下水變化有一個明顯的年周期。對周期分量時間序列的每年相同月份相加再取平均，得到地下水儲量的年度變化規(guī)律如圖4所示?？梢钥闯?，在1 a中，三江源區(qū)域地下水變化基本能夠維持動態(tài)平衡，隨季節(jié)發(fā)生周期性波動：1～3月地下水量不斷減少，3月地下水儲量最少，3～9月地下水量緩慢增加，9月達到峰值，9～12月緩慢減少。

圖4 三江源地下水年度變化規(guī)律Fig.4 Annual variation of groundwater in the Three-River headwater region

3.2 分時段分區(qū)域的地下水變化分析

以SSA提取的趨勢信息為參考，將2003～2020年劃分為2個時間段。2003-01～2012-12三江源區(qū)域地下水儲量變化率如圖5(a)所示，可以看出，該時間段內三江源區(qū)域地下水儲量變化大部分為正增長，僅在瀾滄江源頭東南部分呈現緩慢負增長。地下水儲量變化率空間分布整體上呈現東南向西北方向增加的趨勢。長江源區(qū)域均表現為正增長，西北方向部分變化率為1.2 cm/a左右，東南方向部分變化率為0.7 cm/a左右。黃河源區(qū)域地下水增長主要在西南部分，變化率為0.6 cm/a左右，東北部分地下水變化穩(wěn)定，無明顯增加或減少。瀾滄江源區(qū)域西北部分地下水有輕微增加，變化率為0.5 cm/a左右，東南部分表現為輕微減少，變化率為-0.4 cm/a左右。2013-01～2020-12三江源區(qū)域地下水儲量變化率如圖5(b)所示。該時間段內三江源區(qū)域地下水變化明顯的地區(qū)有3個，分別為長江源區(qū)域西北部、黃河源區(qū)域東部、瀾滄江源區(qū)域。地下水在長江源區(qū)域西北部表現為明顯的正增長，變化率為1.3 cm/a左右，在黃河源區(qū)東部和瀾滄江源區(qū)域表現為明顯的負增長，黃河源區(qū)域東部變化率為-1 cm/a左右，瀾滄江源區(qū)域變化率為-1.1 cm/a左右。

圖5 三江源區(qū)域地下水儲量變化率Fig.5 The change rate of the groundwater reserves of the Three-River headwater region

三江源區(qū)域地下水變化的平均時間序列如圖6所示。對時間序列線性擬合后，得到區(qū)域內地下水變化的趨勢，并對變化趨勢進行顯著性檢驗。結果表明，2003-01～2012-12三江源區(qū)域整體表現為顯著正增長，增長率為0.7 cm/a(P<0.01)；長江源區(qū)域整體表現為顯著正增長，增長率為1.0 cm/a(P<0.01)；黃河源區(qū)域整體表現為顯著正增長，增長率為 0.5 cm/a(P<0.01)；瀾滄江源區(qū)域整體穩(wěn)定，無顯著增加或減少。2013-01～2020-12三江源區(qū)域整體表現為負增長，增長率為-0.3 cm/a(P<0.01)；長江源區(qū)域整體穩(wěn)定，無顯著增加或減少；黃河源區(qū)域整體表現為顯著負增長，增長率為-0.6 cm/a(P<0.01)；瀾滄江源區(qū)域整體表現為顯著負增長，增長率為-1.2 cm/a(P<0.01)。

圖6 三江源區(qū)域地下水變化時間序列Fig.6 Time series of groundwater changes in the Three-River headwater region

3.3 降雨量分析

對中國地面降水月值0.5°×0.5°格點數據集(V2.0)進行提取，分別繪制2003～2012年和2013～2020年的降雨量變化速率圖(圖7)。由圖可見，2003～2012年降雨量在長江源、黃河源與長江源交接處呈明顯正增長，而在黃河源區(qū)域和瀾滄江源區(qū)域呈現較為輕微的正增長。整體來看，這一時間段內降雨量在各個區(qū)域均呈現正增長，與這一時間段內地下水儲量變化較為一致。2013～2020年降雨量在黃河源區(qū)東部和瀾滄江區(qū)東南部有較為明顯的正增長，在長江源區(qū)域有輕微的正增長。總體來看，區(qū)域內降水在這一時間段內依舊表現為增加，與區(qū)域內地下水變化表現為相反的趨勢。對區(qū)域內降雨量和地下水變化時間序列進行相關性分析，得到相關性僅為0.35，表明三江源地下水變化的長期趨勢和降雨量變化的相關性不強。

圖7 三江源區(qū)域降雨量變化分布Fig.7 Distribution of rainfall changes in the Three-River headwater region

Deng等[18]研究認為，地下水變化的長期趨勢受降雨量與溫度的共同影響。從中國地面氣溫月值0.5°×0.5°格點數據集(V2.0)中提取三江源區(qū)域的溫度月值，發(fā)現2003～2020年溫度呈緩慢增加的趨勢，18 a內平均溫度上升0.76 ℃。溫度的升高使區(qū)域內蒸散發(fā)量增加，和降雨量一起影響地下水變化。從圖8可以看出，2013年之前三江源的降雨量多于2013年之后的降雨量，圖中陰影區(qū)域為降雨量最少的年份，地下水在這些年份中都有明顯的下降趨勢，表明三江源地下水的短期變化主要受降雨量變化的影響。三江源地區(qū)2003～2012年地下水呈現正增長歸因于降雨量多、溫度低，2013～2020年地下水呈現負增長歸因于降雨量少、溫度高。

對降雨量進行每年相同月份取平均處理，如圖9所示?？梢钥闯觯涤炅亢偷叵滤粯泳哂忻黠@的年度變化特征：1月最少，1～7月逐漸增加，7月達到最大值，7～12月逐漸減少。對2個時間序列進行時間滯后相關性分析，得出地下水變化滯后降雨量變化2個月的相關性為0.97。表明地下水年度變化與降雨量年度變化規(guī)律一致，地下水變化的年度特征主要受降雨量的影響，三江源區(qū)域降雨對地下水的補給滯后期為2個月。

4 結 語

1)2003～2012年三江源區(qū)域地下水儲量整體呈現增加趨勢，增長率為0.7 cm/a；2013～2020年三江源地下水儲量呈現減少趨勢，增長率為-0.3 cm/a。

2)2003～2012年長江源地下水儲量呈現增加趨勢，增長率為1.0 cm/a；黃河源地下水儲量呈現增加趨勢，增長率為0.5 cm/a；瀾滄江源地下水儲量穩(wěn)定，無明顯增加或減少趨勢。2013～2020年長江源地下水儲量整體穩(wěn)定；黃河源和瀾滄江源地下水儲量均呈現明顯的減少趨勢，增長率分別為-0.6 cm/a和-1.2 cm/a。

3)三江源區(qū)域地下水儲量具有明顯的年度變化規(guī)律：1～3月不斷減少，3月達到最低值，3～9月緩慢增加，9月達到峰值，9～12月緩慢減少。地下水年度變化規(guī)律與降雨量年度變化規(guī)律一致，地下水的年度變化規(guī)律主要受降雨補給的影響，三江源區(qū)域降雨對地下水的補給滯后期為2個月。