地球物理学报  2021, Vol. 64 Issue (4): 1167-1180   PDF    
应用GPS数据和Slepian基函数反演川云渝地区陆地水储量变化
成帅, 袁林果, 姜中山, 刘中冠, 张迪, 徐小凤     
西南交通大学地球科学与环境工程学院, 成都 611756
摘要:局部Slepian函数是将局部区域内的地球物理信号转化为空间谱的一种方法,其可以保证在球面上局部范围内获得最优谱平滑解,非常适用于局部范围地球物理信号的研究.本文利用中国陆态网西南地区72个测站的连续GPS观测资料分析川云渝地区陆地水负荷形变特征,并基于Slepian函数方法解算60阶的空间谱基函数,结合弹性质量负荷理论研究了川云渝地区2011年至2015年陆地水储量变化的时空分布模式.针对Slepian函数的边界效应问题,本文使用GLDAS格网数据计算得到站点处垂直负荷位移时间序列,然后利用该位移数据来进行水储量变化恢复实验,结果表明当边界扩充为3°时能较好地恢复GLDAS模型输出的陆地水储量变化.通过对比区域内GPS、GRACE、GLDAS得到的等效水高以及降雨数据,发现季节性降水是陆地水变化的一个重要驱动因子,GPS反演结果与GRACE和GLDAS数据具有较强的空间一致性.云南地区周年变化要强于川渝地区,其中云南西部的山区陆地水变化最大,约为30 cm,最小为川北以及重庆地区仅为7 cm.相较于GPS反演结果,GRACE与GLDAS明显低估了陆地水储量的季节性变化,分别达到24%和47%.比较分析地区内平均等效水高时间序列的相位发现,GPS得到的陆地水变化与降雨数据一致性较好,而GRACE与GLDAS存在一到两个月左右的时延.同时GPS能较好的探测出2015年1月左右南方地区大范围的强降水,而GRACE与GLDAS并没有体现出该现象,说明GPS能更为灵敏地探测到局部地区陆地水的变化.在站点等效水高时间序列上,GPS与GRACE的相关性总体上要优于GPS与GLDAS,陆地水周年变化较大的云南和四川西部地区站点三种数据间相关性较好,而其他季节性信号不明显的地区则相关性较差.本文的研究表明运用GPS-Slepian方法能够独立地监测高时空分辨率的陆地水储量变化,是作为当前补充GRACE观测资料空缺期的有益尝试.
关键词: Slepian基函数      GPS垂向位移      陆地水储量变化      负荷形变     
Investigating terrestrial water storage change in Sichuan, Yunnan and Chongqing using Slepian basis functions
CHENG Shuai, YUAN LinGuo, JIANG ZhongShan, LIU ZhongGuan, ZHANG Di, XU XiaoFeng     
Faculty of Geosciences and Environmental Engineering, Southwest Jiaotong University, Chengdu 611756, China
Abstract: The Slepian basis function is the approach that is capable of transforming the geophysical signals into spatial spectrum information and obtaining the optimal spectral smoothing solution on the spherical surface, which is suitable to study the geophysical signal at a regional spatial scale. In this study, the characteristics of the hydrological loading deformation in Sichuan, Yunnan, and Chongqing were analyzed by using long-time GPS coordinate time series at 72 Crustal Movement Observing Network of China (CMONOC) stations. Hereafter, based on the elastic mass loading theory and the Slepian basis functions up to 60 order, we further investigated the spatiotemporal patterns of terrestrial water storage change in Southwest China throughout the period of 2011-2015. For testing the effect of different expanded boundary on the quality of terrestrial water storage, the gridded GLDAS data were employed to calculate the vertical loading displacement time series at each site, then they were used to recover the water storage change. Results show that the terrestrial water storage change can be well reproduced by GLDAS data when the boundary is expanded up to 3 degrees. Comparison of the equivalent water heights (EWH) determined by GPS, GRACE, and GLDAS data indicates that the GPS-derived EWH is consistent with those from GRACE and GLDAS in terms of spatial pattern and that seasonal precipitation is one of the major driving factors for terrestrial water storage change. The annual variation of terrestrial water storage is large in Yunnan in comparison with that in Sichuan and Chongqing, which peaks at the amplitude of about 30 cm in Western Yunnan. Whereas it presents small fluctuations of 7 cm in northern Sichuan and Chongqing. Compared with the GPS-derived values, the GRACE and GLDAS significantly underestimate the seasonal amplitude of terrestrial water storage change by 24% and 47%, respectively. The average EWH time series inferred by GPS data show a good correlation with precipitation data, while temporal variations between precipitation and GRACE- and GLDAS-derived EWH present a time delay of one or two months. In contrast to the GRACE and GLDAS, GPS-derived terrestrial water storage product has a good capability of tracing the regional-scale extreme rainfall as illustrated in the strong rainfall event occurred in South China on January 2015. In general, GNSS-inverted time-varying water height has a better correlation with GRACE results than that with GLDAS. We could obtain good consistency between these three datasets where there are large annual water variations, in contrast to poor correlation in regions without conspicuous seasonal water changes. Our results imply that the GPS inversion for water changes based on Slepian basis function succeeds in tracking the high spatial-temporal variations of terrestrial water storage, and it is a useful tool to fill the data gap when GRACE measurements are not available.
Keywords: Slepian basis function    GPS vertical displacement    Terrestrial water storage variation    Loading deformation    
0 引言

陆地水包括了地表水(湖泊、河流和水库蓄水)、土壤水、地下水、冰雪以及树冠水等(Famiglietti et al., 2011),它们在人类的日常生活、农业灌溉、工业生产等领域都扮演着至关重要的角色.中国是一个极度缺水的国家,水资源短缺已经成为我们国家社会经济发展面临的一个重大资源问题(孙冬营等,2020).川云渝地区作为南水北调西线工程的重要源头(张金良等,2020),监测和量化当地陆地水的变化对地方政府进行干旱与洪涝的预防、气候变化分析甚至国家水资源管理政策的制定等都着有重大意义.目前多种现代空间大地测量技术和水文同化建模已广泛应用于监测陆地水的变化.

全球尺度下的水文模型主要包含陆面模型和水文水资源模型,当前应用最为广泛的是美国国家航空航天局NASA (Nation Aeronautics and Space Administration)发布的全球陆面数据同化系统GLDAS(Global Land Data Assimilation System),该模型的优点是能提供全球高时空分辨率的多种陆地表状态和通量场,例如降水、蒸散、地表径流和地下径流等(Rodell et al., 2004).Wang和Zeng(2012)利用青藏高原地区的63个气象站点观测资料验证了GLDAS能提供较为准确降水数据,同时气温的描述也较好.李霞等(2014)利用水文站观测数据验证了模型在黄河源区的适用性,发现GLDAS数据能较好地体现黄河源区的流量特征.但GLDAS数据也存在着一些缺陷,例如Noah模型缺乏部分地表水和深层地下水数据(Jin and Zhang, 2016),不能完整地监测区域陆地水变化.同时部分地区存在径流与蒸散结果不够精确.例如,王文等(2017)研究表明GLDAS低估了中国外流区的径流以及高估了蒸散作用的影响,且该模型简化了物理过程,只考虑了气候因素作用,没有考虑流域内人为因素对水储量造成的影响(Tatsumi and Yamashiki, 2015).因此,GLDAS数据很难以精确反映实际的陆地水变化.

重力恢复与气候试验GRACE(Gravity Recovery and Climate Experiment) 是美国国家航空航天局和德国航空航天中心于2002年3月17日发射的重力卫星.其主要目标是获取每月的时变地球重力场来探究水圈、冰层和海洋中的质量变化(Tapley et al., 2004).当前GRACE提供的时变重力场数据已广泛应用于研究全球范围内季节和年际尺度上的陆地水储量变化、地下水储量变化、冰层质量再分布等(Seo et al., 2006; Tiwari et al., 2009; Velicogna, 2009).除此之外,GRACE对局部地区较为明显的质量变化也有良好的响应,例如有学者研究利用GRACE探测三峡水坝以及巴贡大坝的蓄水量变化(Wang et al., 2011; Tangdamrongsub et al., 2019). 但GRACE一方面由于条带误差和高频误差的原因解算时需要进行滤波处理(Swenson and Wahr, 2006),这在一定程度上抑制真实信号特征,另一方面,GRACE在全球范围内只能保持一个月的时间分辨率和约300 km的空间分辨率,导致小范围内的高频质量变化的监测受到限制.此外,在GRACE与后期的GRACE Follow-on任务衔接期间有着接近一年的观测空缺,因此需要一种代替手段来弥补这部分空缺.

地表季节性的陆地水变化会导致弹性地球由于负荷作用产生周期性的形变.Van Dam and Wahr(1998)发现单个GPS站点的位移时间序列包含了地表负荷对地球造成的显著位移,Bevis等(2005)认为GPS观测技术能对周围100 km内的负荷变化有着较好的灵敏度,后续也有许多学者将GPS位移时间序列与GRACE计算的位移数据进行对比,发现GPS时间序列中的季节性成分主要与陆地水变化有关(Van Dam et al., 2001; Davis et al., 2004; 王林松等,2014; Gu et al., 2017; 丁一航等,2018). 于是有学者探索利用GPS垂直位移时间序列来研究陆地水储量变化,例如Argus等(2014)Fu等(2015)基于GPS站垂直位移时间序列,应用圆盘负荷格林函数的方法分别反演计算了美国加利福利亚地区以及美国西部地区的陆地水储量变化,Borsa等(2014)基于GPS监测的地表隆起计算了美国西部2013年干旱时期流失的水质量.国内也有学者根据GPS监测的地表形变来研究了云南境内的水储量变化以及干旱特征(何思源等,2018; Jiang et al., 2017).Hsu等(2020)基于GPS垂直形变反演了台湾地区2005至2016年间的陆地水储量变化,并分析了降雨模式、入渗率、土壤饱和度和径流等多种因素对陆地蓄水的影响.他们的结论表明除了GRACE与水文模型以外,GPS可以作为一种独立的数据源来监测陆地水变化.且相较于GRACE而言,GPS时间序列时间分辨率高(达到天或天以下),监测时间长(部分站点观测时长达二十余年),使得GPS在小尺度上的反演计算上有着独特的优势.

然而格林函数方法对区域内GPS站点的分布和密度有着较高的要求,不适合用于大尺度陆地水储量的反演工作(Han and Razeghi, 2017).为了能够使GNSS用于大尺度水储量的研究工作,将其转为与GRACE类似的空间谱是一种有效的分析方法.本文探究应用局部Slepian函数将负荷位移信号转换为Slepian空间谱,并以此反演地表质量变化.Slepian问题最早的提出,主要用于解决有限信号截断问题和提高信号功率谱密度的估算等应用(Slepian, 1983),后来进一步完善后逐渐开始应用于地球科学等各个领域.Slepian函数与GRACE所使用的球谐系数法类似,是用一定阶数的Slepian基函数表示局部区域的地球物理信号场.Simons等(1997)基于该方法研究了金星的重力地形构造.Han等用该方法研究了2004年苏门答腊地震的震后局部重力场变化以及局部月球重力场的优化(Han et al., 2008; Han, 2008).Harig和Simons(2012)将该方法应用于监测冰川消融现象,利用重力卫星跟踪数据计算了格陵兰岛地区2002到2011年间的质量变化趋势.在国内该方法也已普遍应用于地月重力场方面的相关研究(朱广彬等,2012孙雪梅等,2015陈石等,2017).在反演陆地水方面,Han和Razeghi(2017)使用Slepian函数反演了澳大利亚大陆地区质量变化,结果显示该方法可以适当降低区域内站点空间覆盖率的要求,实现更大尺度的陆地水储量变化的反演计算.

川云渝三省位于我国西南地区,地处青藏高原东南部,地势西高东低,是青藏高原向东挤压的构造单元中活动性最强的一个,其内部结构具有明显的空间变异性(Jin et al., 2019).气候由西部青藏高原的高原季风气候向东逐渐过渡为亚热带季风气候,区域内水资源丰富,有着较强的季节性现象且空间差异明显.本文目标是开发基于Slepian基函数的GNSS陆地水反演方法并验证该方法在川云渝地区的适用性.文中将GPS垂直位移信号应用局部Slepian基函数结合地表质量负荷理论反演得到2011—2015年川云渝三省陆地水储量变化,并在空间和时间上与GRACE、水文模型和降雨数据进行对比来分析该方法的可行性.

1 数据处理 1.1 GPS位移时间序列及预处理

本文采用国家地震科学数据共享中心提供的单天解坐标时间序列(http://www.eqdsc.com/data/pgv-sjxl.htm),选取了川云渝地区分布较均匀的72个GPS站,站点分布如图 1所示.本文选取的数据时间跨度为2011—2015年,且在研究时段内具有较好的完整性.数据解算使用麻省理工学院提供的GAMIT软件,采用IGS提供的轨道和钟差信息,先验对流层延迟采用GMF模型改正,固体潮和极潮改正使用IERS03模型,海潮使用FES2004模型进行改正.为了能精确地提取陆地水变化所造成的地表形变信号,本实验在后续的时间序列处理中,使用国际质量负荷服务(http://massloading.net/)所提供的大气负荷(Merra2模型)和非潮汐海洋负荷(OMTC05模型)的垂向位移改正产品.图 2为两者产生的垂向位移在站点处的周年振幅大小.由图可知,川云渝地区受到大气负荷的影响较为明显,垂向周年振幅可达6.5 mm,受到非潮汐海洋负荷的影响较小,最大约为0.4 mm.在站点原始坐标时间序列中将这两种负荷信号扣除得到主要由陆地水变化所产生的位移信号.最后为了准确获取季节性位移,本文拟合并去除了时间序列中的线性趋势以及跳变,以减小构造运动,地震以及仪器更换等其他因素的影响.拟合公式如下:

(1)

图 1 西南地区GPS站分布图 Fig. 1 Distribution of GPS stations in Southwest China
图 2 站点处大气负荷(Non-tidal atmospheric loading: NTAL)与非潮汐海洋负荷(Non-tidal ocean loading: NTOL) 垂直位移周年振幅 Fig. 2 Annual vertical displacement amplitude of atmospheric loading and non-tidal oceanic loading at GPS station

式中,ti为观测时间,y为观测序列,a0为线性项,H(ti- Tgj) 为地震与仪器更换所造成的阶跃项,CkSk表示待估的正余弦振幅,e(ti)为残差.当k为1,2时,时间序列分别拟合周年变化与半周年变化.

为了保证反演计算时数据的完整性,本文采用了Kriging-Kalman-Filter模型补齐缺失的GPS数据(Liu et al., 2018),该方法插值时考虑了GPS站之间的空间相关性,可以更好地恢复GPS时间序列中的地球物理信号.

1.2 GRACE和GLDAS数据处理

本文使用的GRACE时变重力场数据来自美国德克萨斯大学空间研究中心(Center of Space Research,CSR) 提供的RL06 v02 Mascon数据,本文选取数据的时间跨度为2011年至2015年,共47个月.这些Mascon解纯粹由GRACE与GRACE-FO信息进行正则化约束计算得到,为了便于分析重采样为0.25°×0.25°的格网形式.数据解算时采用卫星激光测距的估计值替代C20项和C30项(仅GRACE-FO解),并根据Swenson等(2008)计算的估计值校正了地心改正项以及基于ICE6G-D模型修正冰川均衡调整.Mascon数据较好地解决了条带误差,不需要再附加滤波以及缩放因子等处理过程,相较于GRACE球谐解能更好地保留真实的质量变化特征.

水文模型使用美国航天航空局提供的全球陆地数据同化模型GLDAS-Noah(http://daac.gsfc.nasa.gocv),数据为格网形式.文中选择空间分辨率为1°×1°,时间分辨率为一个月的L4表面陆地模型,GLDAS中主要体现了浅层土壤水以及冰雪等成分的变化,将模型中的积雪量、树冠总蓄水和四层(0~2 m)土壤含水量累加就可得到陆地总水储量. 再去均值后即可得到研究时段内的陆地水储量变化.

2 原理和方法 2.1 Slepian函数原理

球谐函数法是在保证全局正交性的条件下将球面上全球分布的地球物理信号用一系列不同阶次的球谐系数表示.但研究时常有目标信号只集中于某一特定区域内的情况,Slepian函数原理是寻找带边界限制的局部球谐函数,函数的信号能最佳地集中在单位球表面上的某个封闭部分,进而可以表示特定区域内的局部物理场(Simons and Dahlen, 2006),Slepian函数只与研究区域与选取的最大阶数有关.基于最大阶数L解算第a项Slepian基函数ga的公式为

(2)

θ为余纬,λ为经度,Ylmlm次的表面球谐函数,gα, lm为第α项Slepian函数的lm次球谐膨胀系数,它可以由球面Ω上某个封闭区域R内的能量集中度γ来确定.

(3)

Slepian函数的目的在于得到γ接近于1的一系列gα, lm值来求得信号集中于R内的基函数gα.方法是用各阶次两个球谐函数积的积分构造矩阵Dlm, lm

(4)

解上式中的特征值与特征向量即可得到能量集中度γgα, lm.本文选择封闭积分区域为西南地区川云渝三省边界,最大阶数由积分区域内GPS站点密度决定,试验后本文选择60阶,计算得到特征值和特征向量分别为γgα, lm.图 3为60阶Slepian函数下求得的3721个能量集中度γ.

图 3 各阶次下Slepian函数能量集中度 Fig. 3 Concentration ratio of each Slepian function

图 3可以得出越靠前的项Slepian函数能量集中度越接近于1,即信号在指定区域所占的比重越大.图 4为前六项基函数信号在积分区域的空间分布.其中红色表示信号为正,蓝色表示信号为负,颜色的深浅表示信号的大小,内轮廓为实际研究区域边界,外轮廓是为了能让信号能充分覆盖区域周围,将研究区域向外扩充一定缓冲区后的积分边界.结合图 3图 4也可以得到γ越接近于1时,信号越集中于积分区域内.随着项数α的增大,基函数信号逐渐往积分区域外扩散,当能量集中度减少至0时,表明信号主要集中于实验区域以外.因此一般只需要选择γ大于一定阈值的前J项基函数,使其能充分覆盖积分区域即可有效转换内部所有地球物理信号.本文中的实验选择γ大于0.1以上的Slepian基函数.

图 4 前6项Slepian基函数信号分布 α为基函数项数,γ表示能量集中度.内轮廓为川云渝三省边界,外轮廓为扩充后积分边界. Fig. 4 Signal distribution of the first 6 Slepian basis functions α is the number of base function items, γ represents the energy concentration.The inner contour is the boundary of Sichuan, Yunnan and Chongqing, and the outer contour is the extended boundary for integral calculation.
2.2 Slepian函数计算等效水高

根据上节中Slepian函数的特性,我们可以将地球表面局部位移场信号用一组Slepian基函数来表示.各个站点处的垂向位移信号u可以近似为

(5)

sa为前J项Slepian函数下的位移基函数系数.又根据地表弹性负荷理论(Farrell, 1972),可将区域内位移基函数系数转化为质量负荷基函数系数sβewh

(6)

ρe为地球密度,ρw为水密度,hl为负荷勒夫数,最后得到区域内(θ, λ)处等效水高形式的质量负荷为

(7)

3 反演结果与分析 3.1 边界敏感性测试

由于Slepian基函数是从中心向外平滑扩散的,边界附近远离积分区域中心的Slepian基函数灵敏度会有所降低(Harig and Simons, 2012).因此为了保证有足够的灵敏度来响应边界附近的站点信号,也为了验证该方法反演水储量变化的可靠性.本文进行了恢复川云渝三省质量变化趋势的实验来确定缓冲区的大小.实验首先将2012年的GLDAS水储量变化转化为类似于GRACE的60阶球谐系数,并以此计算了区域内站点处相应的垂直位移时间序列,用来进行不同缓冲区下的信号恢复实验.实验选择最大阶数为60,缓冲区边界在0°至3°之间,分别计算了各缓冲区下水储量变化恢复情况如图 5所示.图中第一幅(Simulation)为GLDAS数据2012年水储量周年振幅,后面为不同缓冲区下恢复的区域内水储量变化图像.可以发现当缓冲区小于2°时,由于边界附近Slepian基函数信号灵敏度较低,边界会有明显的信号缺失现象,同时信号大量集中于区域中央部分.随着缓冲区的增大,边界信号浓度逐渐提升,当缓冲区为2°至3°时趋势逐渐稳定且较好地恢复了模拟数据的水储量变化,当缓冲区为3°时可恢复90%以上的模拟水储量变化.分析模拟实验中的误差来源一是Slepian函数的截断误差的影响,二是由于站点分布不均匀,导致无法完全恢复部分地区的质量变化信号.针对这个问题,陈石等(2017)的相关实验表明在站点覆盖的空白区容易产生较大的剩余误差,而均匀分布的站点对Slepian函数恢复信号效果有着较好的提升.

图 5 不同缓冲区下模拟信号恢复实验 (a) GLDAS模拟的陆地水变化; (b)—(h) 不同扩充缓冲区下恢复的信号变化. Fig. 5 Experiment of signal recovery under different buffers The first sub-picture is the land water change simulated by GLDAS, and the follow pictures are the signal changes recovered under different expansion buffers.
3.2 反演计算结果

基于积分区域大小内站点分布密度,本文最终选择最大阶数60,缓冲区为3°,γ>0.1的前31项Slepian函数将处理后GPS垂向位移时间序列通过式(5)转换为区域内的Slepian位移场,转换后的GPS空间位移场的周年振幅如图 6所示.图 6a为川云渝区域内空间插值后的垂直位移周年振幅,图 6b为用Slepian函数转化的空间位移场周年振幅.图中显示川云渝地区陆地水造成的垂向季节性变化有着明显的从西南至东北减少的趋势,具体表现为云南以及四川西南部地区位移周年变化较为明显,往四川盆.地方向开始变化逐渐减小,最大位于云南西部山地可达13 mm,最小位于重庆地区约为4 mm.图 6a图 6b比较发现Slepian函数能有效地将位移信号展开到了整个区域,图 6b相比图 6a信号过渡更加平滑.分析认为这是由于两种方法恢复的信号成分不同,空间插值侧重对空间局部点值变化的拟合,而Slepian函数是从空间谱频率上合成地球物理信号,该方法对局部的高频干扰具有明显的压制作用.

图 6 (a) 空间插值后的周年振幅与(b)Slepian函数合成的周年振幅 Fig. 6 (a) Annual amplitude obtained by spatial interpolation and (b) Slepian function synthesis

将时变Slepian位移谱场通过式(6)、(7)转换为时变质量负荷分布,在计算过程中用高斯滤波抑制不稳定的高频噪声信号,得到区域等效水高形式的地表质量变化.为了更好分析陆地水的季节性变化,本文利用中国气象数据网(http://data.cma.cn/) 提供的中国地面降水数据计算了区域平均降水量.采用该降水量与GPS-Slepian、GRACE和GLDAS计算的川云渝地区陆地水储量周年振幅进行比较,结果如图 7所示.从图中可以发现在空间分布上GPS-Slepian,GRACE和GLDAS结果有较好的一致性,三种数据在区域内分布都表现为南北分化的趋势,质量变化最大均为云南西部地区,GPS结果周年振幅达300 mm左右,大于GRACE的230 mm和GLDAS的160 mm.往川渝地区开始逐渐减少,最小在四川北部和重庆地区GPS结果为70 mm左右,GRACE和GLDAS分别为40 mm和25 mm,结合降雨数据可以发现川云渝地区降雨最大位于云南西南部山区,最大高达160 mm,这点与前面数据结果保持一致,表明在云南地区季节性降水是陆地水变化的一个重要驱动因子.而在川西高原处降水丰富却没有造成明显的陆地水周年变化,推测原因是区域内地形及地质条件不利于降水的存储.

图 7 川云渝地区2011—2015年陆地水储量和降雨数据周年振幅 Fig. 7 Annual amplitude of terrestrial water storage and precipitation data in Southwest China from 2011 to 2015

为了体现川云渝地区整体陆地水随时间的变化,本文分别计算了GPS-Slepian、GRACE和GLDAS结果在地区内的平均等效水高,并结合区域内平均降水量进行了综合对比(图 8).图中可以看出,在随时间变化特征上GPS与GRACE和GLDAS都有着很强的一致性,表现为在每年夏季的8至10月份出现高峰,在冬季的1月至3月达到低谷.对比降雨数据可知川云渝地区降水集中于夏季,冬季降雨全年最少,这与前面三种数据结论一致.但在峰值出现时间上,GPS结果要比GRACE与GLDAS要早一定时间,这点在冬季特别明显,分析原因为西南地区冬季降雨较少,陆地水变化主要受区域内径流和蒸散的影响,而GPS对周围轻微的水文变化更为灵敏,因此造成GPS结果与GRACE和GLDAS之间存在延迟.姚朝龙等(2019)在西南地区干旱的研究中也有类似现象.在时间序列大小上,GPS得到的陆地水变化与降雨数据相关性较好,同时发现在2015年1—2月GPS-Slepian结果表现出了小幅度的上升,降雨数据在2015年1月也出现了异常增大,经查询资料发现这是由于该时期南方地区大范围强降水过程的影响(李妍和李春,2019).但这点在GRACE与GLDAS上没有明显体现,表明GPS-Slepian结果能反映更为真实的局部陆地水变化,而GRACE与GLDAS都在一定程度上抑制了时空信号特征.

图 8 区域内平均陆地水储量变化及降雨时间序列 Fig. 8 Variation of average terrestrial water storage and precipitation time series
3.3 站点时间序列分析

为了观察GPS反演结果在区域内各个部分的时间变化趋势,本文选取了位于四川炉霍、巴中、越西以及云南东川、腾冲、通海的六个站进行分析.各站反演的等效水高时间序列如图 9所示.比较各个站等效水高时间序列,位于云南西部山区的腾冲站(YNTC)变化幅度最大,云南东部云贵高原处的东川站(YNDC)与通海站(YNTH)次之,位于四川盆地的巴中站(SCBZ)振幅最小,总体上站点处GPS反演得到的等效水高变化与GRACE和GLDAS之间较为一致,GPS振幅比其余两种数据要大. SCYX、YNDC和YNTH站点处2011年5月至2012年1月之间GLDAS结果与GRACE和GPS结果出现较大差异,存在明显的水储量减少现象,对比图 8并查阅资料认为可能是由于该时期四川南部云贵高原附近出现的局部干旱,降水相较于往年有所减少(周莉等,2018).而降雨量大小对GLDAS计算陆地水储量的影响更为明显,且王文等学者分析过GLDAS高估了中国区域内流域蒸散作用,综合导致了GLDAS数据在发生极端气候时期精度较低.

图 9 站点处等效水高时间序列 Fig. 9 Time series of equivalent water height at GPS stations

为了体现各站点处GPS结果与其他数据时间序列上的一致性,本文计算了站点处互相关系数来判断三种数据间的相关性.由于地球物理信号类似于红噪声,具有显著的自相关现象(Grinsted et al., 2004),因此文中采用了高春春等(2019)分析南极冰盖地区所使用的Monte Carlo模拟法进行了互相关函数的显著性检验.

本文将2011—2015年内GPS反演的等效水高进行月平均,GRACE缺失的数据采用三次多项式插值方法补齐后,分别计算各个站点处GPS与GRACE,GPS与GLDAS,GRACE与GLDAS三种等效水高时间序列间的互相关系数,同时拟合每个站点三种数据等效水高的一阶自相关系数,采用一阶自回归模型(Autoregressive model:AR)来模拟地球物理信号中的红噪声信号,模拟的时间序列与站点处等效水高的长度一致,Monte Carlo模拟次数为5000.如图 10,本文计算了在不同AR1系数下90%、95%和99%置信水平的相关显著临界值进行Monte Carlo模拟来计算相关显著的临界值,比较互相关系数与临界值来判断各个站点处三种数据的相关程度.

图 10 Monte Carlo模拟的不同AR1模型在90%、95%和99%置信度检验时相关显著的临界值 Fig. 10 The cross-correlation 90%, 95% and 99% confidence threshold from different AR1 model derived by Monte Carlo experiments

文中选取90%置信水平的临界值作为标准,为了对比明显,将未达到临界值的站点互相关系数置为0.由此得到各站点处三种数据间的互相关系数如图 11所示.比较三幅子图发现,GPS、GRACE和GLDAS之间互相关系数在云南地区和四川西部地区大于四川东部以及重庆地区,分析是因为川西和云南地区站点密集且有较明显的季节性信号.比较GPS与其他两种数据间的相关程度,发现GPS与GRACE在整体上相关性更强,站点互相关系数平均为0.77.而GPS与GLDAS之间的相关系数普遍要小,平均互相关系数为0.72.其中GPS和GRACE主要在四川东部以及重庆地区相关性较差,区域内的有10个站没有通过90%置信度检验,推测是由于该区域陆地水季节性较弱,GPS结果反映的是站点附近局部范围的陆地水变化,而GRACE反映的是分辨率300 km左右的整体变化特征,因此在季节性变化较小的地区两种数据反映的站点处陆地水变化会有不同.在该区域GPS与空间分辨率较高的GLDAS相关性明显更好也可以验证这一点.子图三中GRACE与GLDAS之间在云南以及四川西部地区相关性最好,根据图 9分析,这是因为两者的陆地水变化时间序列上都有着相似的相位以及振幅,使得两者的等效水高变化特征更为相近.比较分析三幅子图,发现GLDAS与GPS和GRACE之间在四川以北的陕甘青地区相关性都较差,分析原因可能是这些地区水资源匮乏,径流与降雨较少导致GLDAS反映的水文变化不够精确,且GLDAS不能探测到完整的陆地水变化,导致在部分地区与另外两者有差异.

图 11 GPS,GRACE与GLDAS数据之间的互相关系数,图中值为0表示站点处互相关系数小于90%置信水平的相关显著临界值 Fig. 11 The cross-correlation coefficient between GPS, GRACE and GLDAS, and value is 0 means that the cross-correlation coefficient is less than 90% confidence level
4 结论

本文应用Slepian函数反演了西南地区川云渝三省内的陆地水储量变化,基于GPS坐标时间序列数据反演得到区域内2011年至2015年的等效水高时间序列.实验中选择川云渝三省边界解算60阶的Slepian函数,为了保证对边界站点有足够的灵敏性,本文实验模拟了将边界扩充0°至3°的条件下对Slepian函数恢复质量变化的影响.最后通过对比区域内GPS、GRACE、GLDAS得到的等效水高和降雨数据的空间分布和时间序列变化来分析反演结果.本文的主要研究结论如下:

(1) 局部Slepian函数可以将空间上离散的GPS位移信号转化为位移空间谱.实验表明Slepian函数合成的空间位移场能很好的反映实际的位移分布,与插值得到的结果相比之间的空间过渡更加平滑.川云渝地区位移周年振幅有着从西南至东北减少的趋势,表现为四川南部青藏高原、云南地区位移周年变化较为明显,往四川平原方向开始变化逐渐减小,最大位于云南西部山区可达13 mm,最小位于重庆地区约为4 mm.

(2) 川云渝地区陆地水储量周年变化在空间分布上体现为南部的云南地区大于北部的川渝地区.其中云南西部山地陆地水周年变化最大,四川北部以及重庆地区总体变化最小.GPS得到的陆地水周年振幅约为7~30 cm,而GRACE和GLDAS仅为GPS的76%和53%左右,这是因为GPS结果分辨率高且在监测局部区域环境质量变化上更为灵敏,而Mascon数据实际分辨率仅为300 km左右,解算时的正则化约束也与200 km高斯滤波效果相当(Save et al., 2016),GLDAS则缺乏了部分地表水以及深层地下水信息,因此两者振幅相较于GPS偏小.在平均水储量变化和站点处变化趋势中可以发现,GPS结果与降雨数据在时间上更为接近,GPS能得到近实时的监测陆地水变化,而GRACE与GLDAS在冬季存在一定的延迟.在极端气候的监测上,GPS结果与降雨数据均能体现出2015年1月的强降雨事件,而GRACE与GLDAS则没有表现出明显特征.在时间分辨率上,GPS能监测以天为单位的陆地水储量变化,而GRACE与GLDAS只能反映月变化,所以GPS在监测高频的水文特征上有着一定优势.

(3) GPS反演的站点处等效水高时间序列与GRACE和GLDAS结果具有较好的互相关性,均有80%以上的站点通过90%的置信度检验.在空间上相关程度为南部云南地区强于北部川渝地区,四川西部地区强于四川东部以及重庆地区,原因是云南和四川西部地区站点更为密集且季节性现象较明显.整体来看GPS结果在时空上与GRACE更为相似,GPS-GRACE之间的互相关系数平均为0.77,强于GPS-GLDAS之间的0.72.GPS与GRACE之间主要在四川东部和重庆地区相关性较差,分析是因为地区内陆地水的季节性变化较小且GPS结果主要反映站点周围局部的陆地水变化,而GRACE反映的是300 km左右的整体变化特征.GRACE-GLDAS整体相关性最优,但在四川以北的陕甘青地区效果较差,推测是由于地区内降雨与径流较少,且GLDAS缺乏地下水等因素综合导致GLDAS反映的水文变化不精确.

实验证明GPS站网也可以作为独立观测量来有效的监测陆地水储量变化,是GRACE与水文模型之外的另一种可行方案.在反演计算中,Slepian函数在GPS站分布相对均匀的条件下可以适当降低站点密集程度且能应用于面积更大的地区.反演结果可以得到区域内随时间变化的等效水高时间序列,结合学者利用GRACE探测干旱的方法(Thomas et al., 2014),本文认为此结果可同样用于监测西南区域内的干旱事件.在目前GRACE和GRACE Follow-on之间任务衔接导致出现大量数据空缺的情况下,该方法可以对应用GPS手段填补数据空缺以及联合多种大地测量手段监测全球各地区陆地水储量变化提供一定的研究和参考.

致谢  感谢中国陆态网提供的GPS坐标时间序列数据,感谢美国德克萨斯大学空间研究中心提供的时变重力场模型,感谢美国航天航空局提供的全球陆地数据同化模型数据,感谢中国气象数据网提供的中国地面降水数据,同时感谢Simons及其团队提供的Slepian谱分析函数相关源码.以及本文实验中画图所用的GMT软件及其相关开发者们在此一并感谢.
References
Argus D F, Fu Y N, Landerer F W. 2014. Seasonal variation in total water storage in California inferred from GPS observations of vertical land motion. Geophysical Research Letters, 41(6): 1971-1980. DOI:10.1002/2014GL059570
Bevis M, Alsdorf D, Kendrick E, et al. 2005. Seasonal fluctuations in the mass of the Amazon River system and Earth's elastic response. Geophysical Research Letters, 32(16): L16308. DOI:10.1029/2005GL023491
Borsa A A, Agnew D C, Cayan D R. 2014. Ongoing drought-induced uplift in the western United States. Science, 345(6204): 1587-1590. DOI:10.1126/science.1260279
Chen S, Xu W M, Wang Q S. 2017. The spherical harmonic model of gravity field in mainland China by Slepian local spectrum method. Acta Geodaetica et Cartographica Sinica (in Chinese), 46(8): 952-960. DOI:10.11947/j.AGCS.2017.20150542
Davis J L, Elósegui P, Mitrovica J X, et al. 2004. Climate-driven deformation of the solid Earth from GRACE and GPS. Geophysical Research Letters, 31(24): L24605. DOI:10.1029/2004GL021435
Ding Y H, Huang D F, Shi Y L, et al. 2018. Determination of vertical surface displacements in Sichuan using GPS and GRACE measurements. Chinese Journal of Geophysics (in Chinese), 61(12): 4777-4788. DOI:10.6038/cjg2018L0654
Famiglietti J S, Lo M, Ho S L, et al. 2011. Satellites measure recent rates of groundwater depletion in California's Central Valley. Geophysical Research Letters, 38(3): L03403. DOI:10.1029/2010GL046442
Farrell W E. 1972. Deformation of the earth by surface loads. Reviews of Geophysics, 10(3): 761-797. DOI:10.1029/rg010i003p00761
Fu Y N, Argus D F, Landerer F W. 2015. GPS as an independent measurement to estimate terrestrial water storage variations in Washington and Oregon. Journal of Geophysical Research Solid Earth, 120(1): 552-566. DOI:10.1002/2014JB011415
Gao C C, Lu Y, Shi H L, et al. 2019. Detection and analysis of ice sheet mass changes over 27 Antarctic drainage systems from GRACE RL06 data. Chinese Journal of Geophysics (in Chinese), 62(3): 864-882. DOI:10.6038/cjg2019M0586
Grinsted A, Moore J C, Jevrejeva S. 2004. Application of the cross wavelet transform and wavelet coherence to geophysical time series. Nonlinear Processes in Geophysics, 11(5-6): 561-556. DOI:10.5194/npg-11-561-2004
Gu Y C, Yuan L G, Fan D M, et al. 2017. Seasonal crustal vertical deformation induced by environmental mass loading in mainland China derived from GPS, GRACE and surface loading models. Advances in Space Research, 59(1): 88-102. DOI:10.1016/j.asr.2016.09.008
Han S C, Sauber J, Luthcke S B, et al. 2008. Implications of postseismic gravity change following the great 2004 Sumatra-Andaman earthquake from the regional harmonic analysis of GRACE intersatellite tracking data. Journal of Geophysical Research: Solid Earth, 113(B11): B11413. DOI:10.1029/2008JB005705
Han S C. 2008. Improved regional gravity fields on the Moon from Lunar Prospector tracking data by means of localized spherical harmonic functions. Journal of Geophysical Research: Planets, 113(E11): E11012. DOI:10.1029/2008JE003166
Han S C, Razeghi S M. 2017. GPS recovery of daily hydrologic and atmospheric mass variation: A methodology and results from the Australian continent. Journal of Geophysical Research: Solid Earth, 122(11): 9328-9343. DOI:10.1002/2017JB014603
Harig C, Simons F J. 2012. Mapping Greenland's mass loss in space and time. Proceedings of the National Academy of Sciences of the United States of America, 109(49): 19934-19937. DOI:10.1073/pnas.1206785109
He S Y, Gu Y C, Fan D M, et al. 2018. Seasonal variation of terrestrial water storage in Yunnan province inferred from GPS vertical observations. Acta Geodaetica et Cartographica Sinica (in Chinese), 47(3): 332-340.
Hsu Y J, Fu Y N, Bürgmann R, et al. 2020. Assessing seasonal and interannual water storage variations in Taiwan using geodetic and hydrological data. Earth and Planetary Science Letters, 550: 116532. DOI:10.1016/j.epsl.2020.116532
Jiang W P, Yuan P, Chen H, et al. 2017. Annual variations of monsoon and drought detected by GPS: a case study in Yunnan, China. Scientific Reports, 7: 5874. DOI:10.1038/s41598-017-06095-1
Jin H L, Gao Y, Su X N, et al. 2019. Contemporary crustal tectonic movement in the southern Sichuan-Yunnan block based on dense GPS observation data. Earth and Planetary Physics, 3(1): 53-61. DOI:10.26464/epp2019006
Jin S G, Zhang T Y. 2016. Terrestrial water storage anomalies associated with drought in Southwestern USA from GPS observations. Surveys in Geophysics, 37(6): 1139-1156. DOI:10.1007/s10712-016-9385-z
Li X, Gao Y H, Wang W Z, et al. 2014. Climate change and applicability of GLDAS in the headwater of the Yellow River Basin. Advances in Earth Science (in Chinese), 29(4): 531-540. DOI:10.11867/j.issn.1001-8166.2014.04.0531
Li Y, Li C. 2019. An analysis on a heavy rain weather process during 8-13 January, 2015 in Southern China. Periodical of Ocean University of China (in Chinese), 49(7): 11-19. DOI:10.16441/j.cnki.hdxb.20180157
Liu N, Dai W J, Santerre R, et al. 2018. A MATLAB-based Kriged Kalman Filter software for interpolating missing data in GNSS coordinate time series. GPS Solutions, 22: 25. DOI:10.1007/s10291-017-0689-3
Rodell M, Houser P R, Jambor U, et al. 2004. The global land data assimilation system. Bulletin of the American Meteorological Society, 85(3): 381-394. DOI:10.1175/BAMS-85-3-381
Save H, Bettadpur S, Tapley B D. 2016. High-resolution CSR GRACE RL05 Mascons. Journal of Geophysical Research: Solid Earth, 121(10): 7547-7569. DOI:10.1002/2016JB013007
Seo K W, Wilson C R, Famiglietti J S, et al. 2006. Terrestrial water mass load changes from Gravity Recovery and Climate Experiment (GRACE). Water Resources Research, 42(5): W05417. DOI:10.1029/2005WR004255
Simons F J, Dahlen F A. 2006. Spherical Slepian functions and the polar gap in geodesy. Geophysical Journal of the Royal Astronomical Society, 167(1): 1039-1061. DOI:10.1111/j.1365-246X.2006.03065.x
Simons M, Solomon S C, Hager B H. 1997. Localization of gravity and topography: constraints on the tectonics and mantle dynamics of Venus. Geophysical Journal of the Royal Astronomical Society, 131(1): 24-44. DOI:10.1111/j.1365-246X.1997.tb00593.x
Slepian D. 1983. Some comments on Fourier analysis, uncertainty and modeling. SIAM Review, 25(3): 379-393. DOI:10.1137/1025078
Sun D Y, Liu X B, Long X L, et al. 2020. Spatial and temporal heterogeneity of water poverty in China: based on WPI model. Journal of Nanjing Tech University (Social Science Edition) (in Chinese), 19(5): 104-114.
Sun X M, Li F, Yan J G, et al. 2015. An analysis of the applicability of Slepian function in analyzing lunar local gravity field. Acta Geodaetica et Cartographica Sinica, 44(3): 264-273. DOI:10.11947/j.AGCS.2015.20130728
Swenson S, Wahr J. 2006. Post-processing removal of correlated errors in GRACE data. Geophysical Research Letters, 33(8): L08402. DOI:10.1029/2005GL025285
Swenson S, Chambers D, Wahr J. 2008. Estimating geocenter variations from a combination of GRACE and ocean model output. Journal of Geophysical Research: Solid Earth, 113(B8): B08410. DOI:10.1029/2007JB005338
Tangdamrongsub N, Han S C, Jasinski M F, et al. 2019. Quantifying water storage change and land subsidence induced by reservoir impoundment using GRACE, Landsat, and GPS data. Remote Sensing of Environment, 233: 111385. DOI:10.1016/j.rse.2019.111385
Tapley B D, Bettadpur S, Watkins M, et al. 2004. The gravity recovery and climate experiment: mission overview and early results. Geophysical Research Letters, 31(9): L09607. DOI:10.1029/2004GL019920
Tatsumi K, Yamashiki Y. 2015. Effect of irrigation water withdrawals on water and energy balance in the Mekong river basin using an improved VIC land surface model with fewer calibration parameters. Agricultural Water Management, 159: 92-106. DOI:10.1016/j.agwat.2015.05.011
Thomas A C, Reager J T, Famiglietti J S, et al. 2014. A GRACE-based water storage deficit approach for hydrological drought characterization. Geophysical Research Letters, 41(5): 1537-1545. DOI:10.1002/2014GL059323
Tiwari V M, Wahr J, Swenson S C. 2009. Dwindling groundwater resources in northern India, from satellite gravity observations. Geophysical Research Letters, 36(18): L18401. DOI:10.1029/2009gl039401
Van Dam T, Wahr J. 1998. Modeling environment loading effects: a review. Physics and Chemistry of the Earth, 23(9-10): 1077-1087. DOI:10.1016/S0079-1946(98)00147-5
Van Dam T, Wahr J, Milly P C D, et al. 2001. Crustal displacements due to continental water loading. Geophysical Research Letters, 28(4): 651-654. DOI:10.1029/2000GL012120
Velicogna I. 2009. Increasing rates of ice mass loss from the Greenland and antarctic ice sheets revealed by GRACE. Geophysical Research Letters, 36(19): L19503. DOI:10.1029/2009GL040222
Wang A H, Zeng X B. 2012. Evaluation of multireanalysis products with in situ observations over the Tibetan Plateau. Journal of Geophysical Research: Atmospheres, 117(D5): D05102. DOI:10.1029/2011JD016553
Wang L S, Chen C, Zou R, et al. 2014. Using GPS and GRACE to detect seasonal horizontal deformation caused by loading of terrestrial water: A case study in the Himalayas. Chinese Journal of Geophysics (in Chinese), 57(6): 1792-1804. DOI:10.6038/cjg20140611
Wang W, Cui W, Wang P. 2017. Comparison of GLDAS Noah model hydrological outputs with ground observations and satellite observations in China. Water Resources and Power, 35(5): 1-6.
Wang X W, De Linage C, Famiglietti J S, et al. 2011. Gravity recovery and climate experiment (GRACE) detection of water storage changes in the three gorges reservoir of China and comparison with in situ measurements. Water Resources Research, 47(12): W12502. DOI:10.1029/2011WR010534
Yao C L, Luo Z C, Hu Y M, et al. 2019. Detecting droughts in Southwest China from GPS vertical position displacements. Acta Geodaetica et Cartographica Sinica (in Chinese), 48(5): 547-554. DOI:10.11947/j.AGCS.2019.20180308
Zhang J L, Ma X Z, Jing L H, et al. 2020. Plan optimization of the West Route of South-to-North Water Diversion Project. South-to-North Water Transfers and Water Science & Technology (in Chinese), 18(5): 109-114. DOI:10.13476/j.cnki.nsbdqk.2020.0098
Zhou L, Cai R H, Lan M C, et al. 2018. Temporal and spatial distribution characteristics of drought in southwest China during 2011-2012 and its diagnostic analysis of dynamic causes. Hubei Agricultural Sciences (in Chinese), 57(13): 28-33, 45.
Zhu G B, Li J C, Wen H J, et al. 2012. Slepian localized spectral analysis of the determination of the earth's gravity field using satellite gravity gradiometry data. Acta Geodaetica et Cartographica Sinica (in Chinese), 41(1): 1-7.
陈石, 徐伟民, 王谦身. 2017. 应用Slepian局部谱方法解算中国大陆重力场球谐模型. 测绘学报, 46(8): 952-960. DOI:10.11947/j.AGCS.2017.20150542
丁一航, 黄丁发, 师悦龄, 等. 2018. 利用GPS和GRACE分析四川地表垂向位移变化. 地球物理学报, 61(12): 4777-4788. DOI:10.6038/cjg2018L0654
高春春, 陆洋, 史红岭, 等. 2019. 基于GRACE RL06数据监测和分析南极冰盖27个流域质量变化. 地球物理学报, 62(3): 864-882. DOI:10.6038/cjg2019M0586
何思源, 谷延超, 范东明, 等. 2018. 利用GPS垂直位移反演云南省陆地水储量变化. 测绘学报, 47(3): 332-340.
李霞, 高艳红, 王婉昭, 等. 2014. 黄河源区气候变化与GLDAS数据适用性评估. 地球科学进展, 29(4): 531-540. DOI:10.11867/j.issn.1001-8166.2014.04.0531
李妍, 李春. 2019. 2015年1月中国南方大范围持续暴雨的过程分析. 中国海洋大学学报(自然科学版), 49(7): 11-19. DOI:10.16441/j.cnki.hdxb.20180157
孙冬营, 刘新波, 龙兴乐, 等. 2020. 基于WPI模型的中国水贫困时空异质性研究. 南京工业大学学报(社会科学版), 19(5): 104-114.
孙雪梅, 李斐, 鄢建国, 等. 2015. Slepian函数在月球局部重力场分析中的适用性分析. 测绘学报, 44(3): 264-273. DOI:10.11947/j.AGCS.2015.20130728
王林松, 陈超, 邹蓉, 等. 2014. 利用GPS与GRACE监测陆地水负荷导致的季节性水平形变: 以喜马拉雅山地区为例. 地球物理学报, 57(6): 1792-1804. DOI:10.6038/cjg20140611
王文, 崔巍, 王鹏. 2017. GLDAS Noah模型水文产品与中国地面观测及卫星观测数据的对比. 水电能源科学, 35(5): 1-6.
姚朝龙, 罗志才, 胡月明, 等. 2019. 利用GPS垂向位移监测西南地区干旱事件. 测绘学报, 48(5): 547-554. DOI:10.11947/j.AGCS.2019.20180308
张金良, 马新忠, 景来红, 等. 2020. 南水北调西线工程方案优化. 南水北调与水利科技(中英文), 18(5): 109-114. DOI:10.13476/j.cnki.nsbdqk.2020.0098
周莉, 蔡荣辉, 兰明才, 等. 2018. 2011-2012年冬季西南地区干旱的时空分布特征及其动力成因诊断. 湖北农业科学, 57(13): 28-33, 45.
朱广彬, 李建成, 文汉江, 等. 2012. 卫星重力梯度数据确定地球重力场的Slepian局部谱分析方法. 测绘学报, 41(1): 1-7.