三江源2003～2020年地下水储量变化分析

张浩哲 常晓涛, 朱广彬 周 苗 刘 伟

1 山东科技大学测绘与空间信息学院，青岛市前湾港路579号，266590 2 自然资源部国土卫星遥感应用中心，北京市百胜村1号, 100048

三江源是长江、黄河、澜沧江的发源地，是长江、黄河中下游地区和东南亚国家生态安全和区域可持续发展的生态屏障[1]。三江源地下水储量影响着区域生态环境变化和长江、黄河、澜沧江中下游的地面径流量，对其进行研究有助于了解全球变暖背景下青藏高原地区水储量的变化特征。

传统的地下水监测手段主要是通过地面监测井定期监测地下水的水位变化，难以实现大空间、大尺度下的地下水变化监测[2]。自2002年GRACE卫星发射以来，利用时变重力场模型研究大范围区域地下水的时空变化成为可能。在2017年GRACE卫星停止服务后，美国国家航空航天局(NASA)与德国地球科学研究中心联合研制了GRACE-FO卫星，继续承担全球及区域水资源变化监测的任务。

目前，已有大量的研究人员利用时变重力场数据研究区域水资源的变化[3-6]，但对三江源区域地下水的分析研究较少[7]。本文利用GRACE和GRACE-FO提供的2003～2020年的时变重力场模型、GLDAS模型和湖泊水量变化数据反演三江源区域地下水的时空变化，通过奇异谱分析方法(singular spectrum analysis,SSA)提取时间序列中的趋势信号和周期信号，分别对黄河源、长江源、澜沧江源不同区域的地下水变化进行分析，并结合区域内降水数据分析三江源地下水年度补给规律。

1 研究区概况

三江源区域位于青海省南部，地理位置为31°19′～36°12′N、89°45′～102°23′E，总面积为302 500 km2，为长江、黄河和澜沧江的源头汇水区(图1)。三江源区域是青藏高原的腹地和主体，以山地地貌为主，海拔1 971～6 565 m，河流密布，湖泊、沼泽众多，是世界上海拔最高、面积最大、湿地类型最丰富的地区，也是中国淡水资源的重要补给地。

图1 三江源地区位置Fig.1 The location of the Three-River headwater region

2 数据和方法

2.1 数据与预处理

采用美国得克萨斯大学空间研究中心(CSR)提供的GRACE和GRACE-FO卫星RL06数据，模型最大阶数为96，时间范围为2003-01～2020-12。GRACE与GEACE-FO之间11个月(2017-07～2018-05)的数据缺失不作处理。此外，由于GRACE卫星在前期调试运行阶段及后期出现传感器性能下降和供能不足的问题造成20个月的模型缺失，对此，采用三次样条法进行插值[8]。

GLDAS数据选取GLDAS-NOAH水文模型，时间分辨率为1个月，空间分辨率为1°× 1°，时间范围与时变重力场数据时间范围一致。本研究主要提取水文模型中的4层土壤水(0～200 cm)、雪水当量和冠层水变化。

湖泊数据采用国家青藏高原科学数据中心的青藏高原2000～2017年高时间分辨率湖泊水位及水量变化数据集[9]，以及Cooley等[10]提供的2018～2020年全球湖泊水量变化数据集。从数据集中提取三江源地区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卫星对重力场低阶项不敏感，且模型中含有高阶噪声、奇偶阶相关误差等[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]。

依据水量平衡公式，三江源地区的地下水储量变化可表示为：

Δ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年整个三江源区域地下水储量的平均变化时间序列如图2所示。考虑到时间序列的连续性，在进行SSA分析时选择2003-01～2017-06的数据，窗口长度选为48个月(图3)。

图2 2003～2020年三江源区域地下水平均变化时间序列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)所示，可以看出，地下水变化有一个明显的年周期。对周期分量时间序列的每年相同月份相加再取平均，得到地下水储量的年度变化规律如图4所示。可以看出，在1 a中，三江源区域地下水变化基本能够维持动态平衡，随季节发生周期性波动：1～3月地下水量不断减少，3月地下水储量最少，3～9月地下水量缓慢增加，9月达到峰值，9～12月缓慢减少。

图4 三江源地下水年度变化规律Fig.4 Annual variation of groundwater in the Three-River headwater region

3.2 分时段分区域的地下水变化分析

以SSA提取的趋势信息为参考，将2003～2020年划分为2个时间段。2003-01～2012-12三江源区域地下水储量变化率如图5(a)所示，可以看出，该时间段内三江源区域地下水储量变化大部分为正增长，仅在澜沧江源头东南部分呈现缓慢负增长。地下水储量变化率空间分布整体上呈现东南向西北方向增加的趋势。长江源区域均表现为正增长，西北方向部分变化率为1.2 cm/a左右，东南方向部分变化率为0.7 cm/a左右。黄河源区域地下水增长主要在西南部分，变化率为0.6 cm/a左右，东北部分地下水变化稳定，无明显增加或减少。澜沧江源区域西北部分地下水有轻微增加，变化率为0.5 cm/a左右，东南部分表现为轻微减少，变化率为-0.4 cm/a左右。2013-01～2020-12三江源区域地下水储量变化率如图5(b)所示。该时间段内三江源区域地下水变化明显的地区有3个，分别为长江源区域西北部、黄河源区域东部、澜沧江源区域。地下水在长江源区域西北部表现为明显的正增长，变化率为1.3 cm/a左右，在黄河源区东部和澜沧江源区域表现为明显的负增长，黄河源区域东部变化率为-1 cm/a左右，澜沧江源区域变化率为-1.1 cm/a左右。

图5 三江源区域地下水储量变化率Fig.5 The change rate of the groundwater reserves of the Three-River headwater region

三江源区域地下水变化的平均时间序列如图6所示。对时间序列线性拟合后，得到区域内地下水变化的趋势，并对变化趋势进行显著性检验。结果表明，2003-01～2012-12三江源区域整体表现为显著正增长，增长率为0.7 cm/a(P<0.01)；长江源区域整体表现为显著正增长，增长率为1.0 cm/a(P<0.01)；黄河源区域整体表现为显著正增长，增长率为 0.5 cm/a(P<0.01)；澜沧江源区域整体稳定，无显著增加或减少。2013-01～2020-12三江源区域整体表现为负增长，增长率为-0.3 cm/a(P<0.01)；长江源区域整体稳定，无显著增加或减少；黄河源区域整体表现为显著负增长，增长率为-0.6 cm/a(P<0.01)；澜沧江源区域整体表现为显著负增长，增长率为-1.2 cm/a(P<0.01)。

图6 三江源区域地下水变化时间序列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年降雨量在长江源、黄河源与长江源交接处呈明显正增长，而在黄河源区域和澜沧江源区域呈现较为轻微的正增长。整体来看，这一时间段内降雨量在各个区域均呈现正增长，与这一时间段内地下水储量变化较为一致。2013～2020年降雨量在黄河源区东部和澜沧江区东南部有较为明显的正增长，在长江源区域有轻微的正增长。总体来看，区域内降水在这一时间段内依旧表现为增加，与区域内地下水变化表现为相反的趋势。对区域内降雨量和地下水变化时间序列进行相关性分析，得到相关性仅为0.35，表明三江源地下水变化的长期趋势和降雨量变化的相关性不强。

图7 三江源区域降雨量变化分布Fig.7 Distribution of rainfall changes in the Three-River headwater region

Deng等[18]研究认为，地下水变化的长期趋势受降雨量与温度的共同影响。从中国地面气温月值0.5°×0.5°格点数据集(V2.0)中提取三江源区域的温度月值，发现2003～2020年温度呈缓慢增加的趋势，18 a内平均温度上升0.76 ℃。温度的升高使区域内蒸散发量增加，和降雨量一起影响地下水变化。从图8可以看出，2013年之前三江源的降雨量多于2013年之后的降雨量，图中阴影区域为降雨量最少的年份，地下水在这些年份中都有明显的下降趋势，表明三江源地下水的短期变化主要受降雨量变化的影响。三江源地区2003～2012年地下水呈现正增长归因于降雨量多、温度低，2013～2020年地下水呈现负增长归因于降雨量少、温度高。

对降雨量进行每年相同月份取平均处理，如图9所示。可以看出，降雨量和地下水一样具有明显的年度变化特征：1月最少，1～7月逐渐增加，7月达到最大值，7～12月逐渐减少。对2个时间序列进行时间滞后相关性分析，得出地下水变化滞后降雨量变化2个月的相关性为0.97。表明地下水年度变化与降雨量年度变化规律一致，地下水变化的年度特征主要受降雨量的影响，三江源区域降雨对地下水的补给滞后期为2个月。

4 结 语

1)2003～2012年三江源区域地下水储量整体呈现增加趋势，增长率为0.7 cm/a；2013～2020年三江源地下水储量呈现减少趋势，增长率为-0.3 cm/a。

2)2003～2012年长江源地下水储量呈现增加趋势，增长率为1.0 cm/a；黄河源地下水储量呈现增加趋势，增长率为0.5 cm/a；澜沧江源地下水储量稳定，无明显增加或减少趋势。2013～2020年长江源地下水储量整体稳定；黄河源和澜沧江源地下水储量均呈现明显的减少趋势，增长率分别为-0.6 cm/a和-1.2 cm/a。

3)三江源区域地下水储量具有明显的年度变化规律：1～3月不断减少，3月达到最低值，3～9月缓慢增加，9月达到峰值，9～12月缓慢减少。地下水年度变化规律与降雨量年度变化规律一致，地下水的年度变化规律主要受降雨补给的影响，三江源区域降雨对地下水的补给滞后期为2个月。