基于GRACE/GRACE-FO和Swarm卫星研究2002—2020年南极和格陵兰岛冰盖质量时空变化

陈威,钟敏,冯伟*，王长青，李文浩，梁磊

1 中国科学院精密测量科学与技术创新研究院，大地测量与地球动力学国家重点实验室,武汉 430077 2 中国科学院大学地球与行星科学学院,北京 100049 3 中山大学测绘科学与技术学院，广东珠海 519082 4 南京工业大学测绘科学与技术学院，南京 211800

0 引言

在全球气候变暖的背景下，冰川消融不仅会导致全球海平面上升，还会使得极端天气事件频发(Church et al.,2013;Meredith et al.,2019).南极和格陵兰岛作为全球最大的两个大陆冰盖，约占地球总冰川面积97%、总冰量99%(Grinsted,2013).若两极冰盖完全消融，全球海平面大约上升68 m(Alley et al.,2005;Church et al.,2013).准确估计南极和格陵兰岛冰盖的质量变化不仅有助于解释全球海平面变化的成因以及对全球水循环的贡献，还有助于理解海洋-大气之间的相互作用和全球气候变化.

自2002年以来，GRACE(Gravity Recovery and Climate Experiment)重力卫星在监测南极和格陵兰岛冰盖质量变化中发挥了不可替代的作用.Ramillien等(2006)利用GRACE较早地给出了2002—2005年期间南极和格陵兰岛冰盖质量变化，指出由于研究时间段较短很难从年际项中分离出趋势信号.Barletta等(2008)顾及到冰川均衡调整(Glacial Isostatic Adjustment,GIA)的影响后，估计的南极和格陵兰岛质量变化速率分别为-171 Gt·a-1、-101 Gt·a-1.Velicogna(2009)发现2002—2009年期间南极和格陵兰岛的冰盖消融速率并非均匀变化，它们分别以30 Gt·a-2、26 Gt·a-2的加速度项在加速消融.Nerem和Wahr(2011)通过比较去除GIA的J2项与南极、格陵兰岛冰盖的质量变化，发现2002—2010年期间地球的扁率变化主要与南极和格陵兰岛冰盖消融有关.Velicogna等(2014)详细分析南极和格陵兰岛冰盖子流域的质量变化后，发现格陵兰岛70%的冰川融化来自于其西南和东北区域，南极81%的冰川消融来自于阿蒙森海域地区和南极半岛.Forsberg等(2017)通过分析2002—2015年期间南极和格陵兰岛冰盖的质量变化，推断出它们分别等效于全球平均海平面每年上升约0.72 mm和0.26 mm.Groh等(2019)采取多种GRACE数据后处理方法分别估计了南极和格陵兰岛冰盖的质量变化，认为2003—2013年期间南极冰盖的消融速率范围为99～108 Gt·a-1，格陵兰岛的消融速率则为252～274 Gt·a-1.

1 数据和方法

1.1 Swarm卫星数据

(1)

(2)

(3)

1.2 GRACE和GRACE-FO数据

本文使用的GRACE/GRACE-FO level-2数据产品来自于CSR、GFZ和JPL三家公布的GSM 60阶球谐系数(http:∥icgem.gfz-potsdam.de/series).其中GRACE数据的时间范围选用2002年4月至2017年6月，GRACE-FO数据选用2018年6月到2020年3月.GRACE和GRACE-FO均为60阶球谐系数，其后处理方法除采用300 km半径的高斯滤波和nmax=60外，其余后处理方法与Swarm数据的处理步骤一致.本文将CSR、GFZ和JPL三家机构的算术平均值作为最终估计结果.此外，我们还选用了RL06 Version 2版本的CSR Mascon和 JPL Mascon数据来验证本文的最终结果.

1.3 ARIMA-MC预测方法

为弥补GRACE/GRACE-FO卫星之间近一年的间断期，本文除使用Swarm时变重力场数据外，还引入了ARIMA-MC(Autoregressive Integrated Moving Average Model-Monte Carlo)方法来预测间断期两极冰盖的质量变化.ARIMA即自回归移动平均模型，它是由Box等(2016)提出的一种基于统计学的时间序列预测分析方法.ARIMA包括自回归模型(Autoregressive,AR)，差分模型(Integrated,I)和移动平均模型(Moving Average,MA)三个部分.首先，为保证样本时间序列具有平稳性需要进行d阶差分，差分后的时间序列采用自相关函数来判断其是否具有平稳性.进行d阶差分后，我们采用最小化信息量准则AIC(Akaike Information Criterion)来确定自回归移动平均模型ARMA的阶数p和q，最小化AIC公式为：

AIC=-2logL+2(p+q+1),(4)

这里的L为d阶差分后时间序列的似然函数值.确定d、p和q之后，利用Matlab中的ARIMA函数获取最优的拟合模型.ARIMA模型公式如下：

(5)

其中L为滞后算子，Øi为自回归模型系数，θi为移动平均模型系数，εt为白噪声序列.确定最优的ARIMA模型后，依据ARIMA模型的概率分布随机数进行蒙特卡罗模拟(Monte Carlo,MC)，并以此来预测2017年7月—2018年5月期间的南极和格陵兰岛冰盖质量变化.本文的MC预测方法选用Matlab中的Simulate函数实现.

2 结果

2.1 两极冰盖质量变化的时空分布

图1表示的是2002年4月至2020年3月期间，南极和格陵兰岛冰盖质量变化趋势的空间分布.南极半岛(Antarctic Peninsula)和西南极阿蒙森海域(Amundsen Sea)的冰川消融信号较大，威尔克斯地区(Wilkes Land,WL)存在较小的冰川消融信号，而莫德女王地区(Queen Maud Land)和玛丽伯德地区(Marie Byrd Land)则呈现出冰川累积.格陵兰岛冰盖的主要消融区域位于其西北、西南、南部和东南的临海区域，其北部和东北沿海地区也存在一定的冰川消融，但是格陵兰岛中央区域存在微弱的冰川累积.

图1 2002年4月至2020年3月期间(a)南极和(b)格陵兰岛冰盖质量变化趋势的空间分布Fig.1 The spatial patterns of ice mass rate in (a)Antarctica and (b)Greenland from April 2002 to March 2020

2002年4月至2020年3月，南极和格陵兰岛冰盖质量变化(等效海平面高变化)的时间序列如图2所示.GRACE/GRACE-FO卫星探测到2002—2020年期间南极和格陵兰岛冰盖质量损失速率分别为-119±23 Gt·a-1和-259±20 Gt·a-1，等效于全球平均海平面每年约上升0.33 mm和0.72 mm.这一结果与CSR RL06 Mascon v2和JPL RL06 Mascon v2结果十分接近.整个南极冰盖质量变化的时间序列显示2007年以前南极冰盖质量处于净平衡状态(即该时间段内整个南极区域的冰盖消融与累积维持收支相抵的状态)，但2007年之后冰盖消融的主导地位日益显著，尤其是2017年和2020年暖季(当年11月份至次年3月份)出现较大的质量亏损.格陵兰岛冰盖的时间序列显示其冰盖消融具有较强的季节性和趋势项，且在2012年夏季和2019年夏季出现了年际尺度上的加速消融现象.

图2 (a)南极和(b)格陵兰岛冰盖质量变化(等效海平面高变化)的时间序列Fig.2 Time series of ice mass changes and equivalent sea level changes in (a)Antarctica and (b)Greenland,respectively

2.2 Swarm卫星和ARIMA-MC方法在间断期的结果

如图3所示，我们分别给出了GRACE/GRACE-FO、Swarm和ARIMA-MC在南极和格陵兰岛冰盖质量变化的时间序列.2015年1月至2019年6月期间：Swarm与GRACE/GRACE-FO+ARIMA-MC(包含了GRACE/GRACE-FO与ARIMA-MC的预测结果)冰盖质量变化时间序列在南极和格陵兰区域的相关系数分别为0.58、0.72；Swarm和GRACE/GRACE-FO+ARIMA-MC观测到南极冰盖质量的趋势变化分别为-70 Gt·a-1、-53 Gt·a-1，格陵兰岛冰盖质量的趋势变化分别为-92 Gt·a-1、-138 Gt·a-1.2017—2019年期间(主要包含间断期)：Swarm与GRACE/GRACE-FO+ARIMA-MC结果在南极吻合较好，但在格陵兰岛2018年前后吻合较差.Swarm与GRACE/GRACE-FO结果在南极也存在一定的差异，Swarm结果在2015、2017、2019年的1—3月份相对于GRACE/GRACE-FO的冰盖质量亏损的幅度更大,同样在2016—2017年冰盖质量增加的幅度较GRACE/GRACE-FO结果也更大.图3b和3d的箱线图显示Swarm结果与GRACE/GRACE-FO+ARIMA-MC在2015年1月—2019年6月期间的数据分布情况：其中Swarm数据与GRACE/GRACE-FO+ARIMA-MC结果在南极的差异较大，GRACE/GRACE-FO+ARIMA-MC主要集中在-400～-1300 Gt，而Swarm数据波动的范围较大；两者的数据分布在格陵兰岛较为接近，主要集中在-1000～-1800 Gt.通过与GRACE/GRACE-FO+ARIMA-MC结果的对比可知，Swarm时变数据在格陵兰岛的表现更好.

图3 (a)和(c)分别为南极和格陵兰岛冰盖质量变化的时间序列，(b)和(d)分别为南极和格陵兰岛在2015年1月—2019年6月期间的不同数据分布Fig.3 Time series of ice mass changes in Antarctica (a)and Greenland (c),respectively.The boxplots (b)and (d)show the data distribution from Swarm and ARIMA-MC from January 2015 to June 2019

图4a和4b显示在2015年1月至2019年6月期间，Swarm与GRACE/GRACE-FO+ARIMA-MC两种方法在南极和格陵兰岛冰盖质量变化结果的空间相关性.两者在南极吻合较好的区域主要位于西南极的大部分地区、北部的莫德女王地区和东部的威尔克斯地区.在格陵兰岛吻合较好的区域主要在整个南部地区和西北部地区.如图4c—4f所示，两者在南极和格陵兰岛质量变化趋势的空间分布吻合较好.南极冰盖的消融区域主要位于西南极的阿蒙森海域和威尔克斯地区，而南极的玛丽伯德地区、北部的莫德女王地区和东部的威尔克斯地区呈现冰川累积信号.格陵兰岛的消融区域主要表现在其西部、东南部以及东北部.Swarm和GRACE/GRACE-FO+ARIMA-MC的空间趋势分布也有一定的差异.GRACE/GRACE-FO+ARIMA-MC趋势变化空间分布结果显示，南极的玛丽伯德地区的冰川累积较Swarm的结果更多，威尔克斯地区冰川质量亏损较Swarm结果更少.两者在格陵兰岛北部以及东南部地区的空间分布也存在一定差异.因此在GRACE和GRACE-FO的观测空白期，Swarm卫星和ARIMA-MC方法的结果均可有效填补南极和格陵兰岛冰盖质量变化信息，但两者之间也存在一定差异.

图4 在2015年1月—2019年6月期间，Swarm和GRACE/GRACE-FO+ARIMA-MC结果分别在(a)南极和(b)格陵兰岛的相关性，白色区域表示相关系数在统计学上未通过95%的置信区间；(c)ARIMA-MCGRACE/GRACE-FO+ARIMA-MC和(d)Swarm 结果在南极的空间趋势项；GRACE/GRACE-FO+ARIMA-MC (e)和Swarm (f)结果在格陵兰岛的空间趋势项.图中GFO为GRACE-FO的缩写Fig.4 From January 2015 to June 2019,the correlation coefficients between the results from Swarm and GRACE/GRACE-FO+ARIMA-MC in (a)Antarctica and (b)Greenland,respectively.The white indicates that the correlation coefficients don′t statistically pass the 95% confidence interval.The Antarctic spatial trend is from GRACE/GRACE-FO+ARIMA-MC (c)and Swarm (d),respectively.The spatial trend of Greenland is from GRACE/GRACE-FO+ARIMA-MC (e)and Swarm (f),respectively.GFO is the abbreviation of GRACE-FO

2.3 Swarm和ARIMA-MC结果的不确定性

图5a和5b显示的是Swarm结果在南极和格陵兰岛的总误差分布，该误差包括了公式(3)中的Swarm观测误差、泄露误差和尺度因子误差.Swarm在南极和格陵兰岛格网点上的总误差主要集中在0～500 mm区间，主要位于冰川质量稳定区域；总误差的最大值约为2000 mm，且主要位于冰川消融和冰川累积区域.在图5c和5d中，南极和格陵兰岛格网上的尺度因子主要集中在-5～5之间，其最大值约为15且主要集中在冰川消融和冰川累积区域.总而言之，Swarm数据的尺度因子和总误差主要体现在冰川消融和冰川累积区域.

图5 Swarm结果在(a)南极和(b)格陵兰岛的总误差，用于恢复(c)南极和(d)格陵兰岛真实质量变化所采用的尺度因子Fig.5 The total error of Swarm′s results in (a)Antarctica and (b)Greenland,respectively.The scale factors are used to recover the true mass in (c)Antarctica and (d)Greenland,respectively

如图6a和6b所示，利用ARIMA在南极和格陵兰岛建立质量变化模型对应的最小化AIC的空间分布较为均匀，这表明ARIMA模型和预测结果的稳健性较强.图6c和6d表示的是利用MC方法预测间断期南极和格陵兰岛冰盖质量变化的不确定度.MC方法的不确定度范围主要集中在0～500 mm区间，其最大值约为2000 mm且主要位于冰川消融和冰川累积区域.总之，ARIMA-MC预测方法具有较强的稳健性，且预测结果的不确定度与Swarm总误差分布较为类似，均集中在冰川消融和冰川累积的区域.

图6 ARIMA-MC方法在(a)南极和(b)格陵兰岛预测中的AIC值，ARIMA-MC方法在(c)南极和(d)格陵兰岛的不确定性Fig.6 The AIC of ARIMA-MC method in (a)Antarctic and (b)Greenland,respectively.The uncertainties of ARIMA-MC method in (c)Antarctica and (d)Greenland,respectively

3 讨论

3.1 南极威尔克斯地区冰川的加速消融

近些年来，南极冰盖特别是威尔克斯地区的冰川正在加速消融.图7a和7b表示的是2002—2009年期间和2010—2020年期间，GRACE/GRACE-FO+ARIMA-MC南极冰盖质量变化趋势的空间分布.南极东部的威尔克斯地区和西南极的阿蒙森海域出现了明显的加速消融，而西南极的玛丽伯德地区和北部的莫德女王地区冰川出现了加速累积信号.图7c显示的是南极威尔克斯地区冰川质量变化的时间序列.威尔克斯冰川在2002—2009年期间的消融速率为3 Gt·a-1，但在2010—2020年期间的速率增加了10倍达到30 Gt·a-1.同样Swarm结果也显示，该地区在2015年1月—2019年6月期间的消融速率达到73 Gt·a-1.Shen等(2018)基于质量输入输出模型和Lansat影像数据发现威尔克斯地区冰川加速消融.Rignot等(2019)基于较新的冰流和冰厚度数据分析得到威尔克斯地区2009—2017年间的冰川质量损失为51 Gt·a-1.Miles等(2016)认为全球大量冰川的退缩与海冰减少可能会加速温暖海水对冰川末端的侵袭，这可能会加速威尔克斯地区的冰川消融.Zhan等(2018)发现威尔克斯地区冰川的消融对厄尔尼诺影响下的海表温度异常变化较为敏感，这可能是导致威尔克斯地区冰川加速消融的原因之一.Shen等(2018)认为威尔克斯地区附近的海洋热通量会受到南极东部海域冰川动力学的影响，从而引发近年来威尔克斯地区冰川的加速消融.

图7 (a)和(b)分别为南极在2002年4月—2009年12月和2010年1月—2020年3月期间冰川质量变化趋势的空间分布，(c)是南极威尔克斯地区冰川质量变化的时间序列Fig.7 The spatial patterns of ice mass rate in Antarctic during the period of (a)April 2002—December 2009 and (b)January 2010—March 2020,respectively.(c)Time series of ice mass changes in Wilkes Land of Antarctica

3.2 夏季NAO对格陵兰岛冰盖的影响

图8a给出了格陵兰岛冰盖质量在2002年4月至2020年3月的夏季期间的年际变化和北大西洋涛动指数(North Atlantic Oscillation,NAO).格陵兰岛质量在夏季的年际变化与北大西洋涛动夏季指数呈正相关，相关系数为0.78(通过95%置信区间检验).在2012年和2019年的北大西洋涛动夏季指数负相位期间，格陵兰岛冰盖较平常夏季出现约300 Gt的质量损失加剧.图8b显示了NAO夏季指数与格陵兰岛夏季冰盖质量年际变化相关性的空间分布，格陵兰西北部和东部与NAO夏季指数呈现出较强的正相关，而格陵兰岛中央地区与NAO夏季指数呈较强的负相关.van Angelen等(2013)和Bevis等(2019)认为格陵兰岛冰盖在2012年夏季明显的加速消融与NAO有关.当NAO夏季指数出现负相位时，格陵兰岛上空会出现高压和较好的能见度，这会增强太阳辐射的表面吸收、减少降雪和气温升高，从而导致冰川在更长时间段内加速消融(Fettweis et al.,2013;Van Angelen et al.,2014;Bevis et al.,2019).因此2019年夏季的明显加速消融与2012年夏季的消融类似，是由于同时期负相位的NAO影响所致.在经历2019年负相位的夏季NAO之后，我们推测2021年夏季格陵兰岛冰盖的年际变化会类似2013年夏季，较正常年份的夏季出现消融减弱.

图8 (a)NAO夏季指数与格陵兰岛冰盖质量年际变化；(b)NAO夏季指数与格陵兰岛质量年际变化相关性的空间分布.JJAS表示选用6月、7月、8月和9月作为本文研究的夏季Fig.8 (a)The summertime NAO index (JJAS)and interannual variations of ice mass (JJAS)in Greenland.(b)The correlation coefficients between summertime NAO index and interannual variations of Greenland′s ice mass in spatial domain.JJAS denotes the months of June,July,August,and September in summer

3.3 南极和格陵兰岛冰盖质量变化趋势的估计

表1列举了本文和不同学者关于南极、格陵兰岛冰盖质量变化趋势估计的研究结果.其中南极趋势项误差考虑了拟合误差和GIA 20%的不确定性，格陵兰岛的趋势项误差参考Ran等(2018).在2003—2012、2003—2013和2003—2014年三个时间段内，本文基于球谐系数估计的趋势项与CSRM RL06 v2、JPLM RL06 v2的结果吻合较好.通过对比Mu等(2017)结果发现，本文对于2003—2012年南极冰盖质量变化趋势的估计约有～10 Gt·a-1的差异，但在格陵兰岛变化趋势的估计吻合较好，这可能是由于使用不同GIA模型对南极趋势项的估计产生了一定差异.通过对比Schrama等(2014)、Velicogna等(2014)、Watkins等(2015)和Groh等(2019)结果发现，本文对于2003—2013年南极变化趋势的估计处于所有估计值的中值，但不同学者由于采用不同的GIA模型使得南极的变化趋势估计值具有很大的不同.不同研究对格陵兰岛冰盖质量变化趋势估计值的差异大约小于10 Gt·a-1，这说明不同方法和不同GIA模型对格陵兰岛冰盖质量变化趋势估计的影响较小.本文对2003—2014年南极、格陵兰岛冰盖质量变化趋势的估计均大于Harig和Simons(2015,2016)的结果，但与Mascon估计的结果吻合较好.基于以上不同时间段的对比结果发现，Harig和Simons(2015,2016)使用的Slepian方法可能会低估南极、格陵兰岛冰盖质量变化的趋势项，这可能是因为Slepian函数在外扩南极、格陵兰岛陆海边界信号时存在信号损失等问题(高春春等,2019).总之，使用不同GIA模型对南极冰盖质量变化趋势的估计会造成较大的差异，而不同的后处理方法和GIA模型对格陵兰岛冰盖质量变化趋势估计的影响相对较小.

表1 南极和格陵兰岛冰盖消融质量变化速率的研究结果(单位：Gt·a-1)Table 1 The ice mass loss rates in Antarctica and Greenland from recent studies (unit:Gt·a-1)

4 结论

本文结合Swarm数据和ARIMA-MC预测方法填补了GRACE/GRACE-FO间断期南极和格陵兰岛冰盖质量变化的时间序列，并详细分析了其时空变化规律.结论如下：

(1)GRACE/GRACE-FO卫星探测到2002年4月至2020年3月期间南极和格陵兰岛冰盖质量损失速率分别为-119±23 Gt·a-1和-259±20 Gt·a-1，等效于全球平均海平面高每年约上升0.33 mm和0.72 mm.南极地区的冰盖消融主要在南极半岛、西南极阿蒙森海域以及威尔克斯地区，南极莫德女王地区和玛丽伯德地区则呈现出冰川累积.格陵兰岛冰盖的主要消融区域位于西北、西南、南部和东南的沿海区域，其中央区域存在微弱的冰川累积.

(2)Swarm卫星数据和ARIMA-MC预测方法均能有效填补GRACE/GRACE-FO间隔期南极和格陵兰岛冰盖质量变化的时间序列空白.Swarm卫星和ARIMA-MC方法可探测到空窗期南极、格陵兰岛冰盖的消融信号，它们在主要的冰川消融区域呈现出较强的正相关，但在数据分布和质量变化趋势的估计中存在一定差异.ARIMA-MC预测结果的不确定度与Swarm总误差分布类似，它们的最大值分布在冰川消融和冰川累积区域.

(3)南极威尔克斯地区2010—2020年期间的冰川融化速率较2002—2009年期间增加了10倍，这可能与其附近的海洋温度变化有关.夏季格陵兰岛冰盖年际变化与NAO夏季指数相关系数为0.78，格陵兰岛2019年夏季明显的加速消融与2012年夏季消融类似，均与负相位夏季北大西洋涛动有关.不同GIA模型会显著影响南极冰盖质量变化趋势项的估计，但对格陵兰岛冰盖质量变化趋势项估计的影响则相对较小.

致谢感谢ICEGM (International Centre for Global Earth Models)网站提供的Swarm和GRACE/GRACE-FO卫星观测数据(http:∥icgem.gfz-potsdam.de/home).感谢美国俄亥俄州立大学C.K.Shum教授和南方科技大学冉将军助理教授对文章给出的建议.感谢审稿专家对本文提出的宝贵修改意见.