GNSS/GRACE/GRACE-FO/气象数据结合反演干旱指数

姚朝龙,陈涌鑫,罗志才,李 琼,叶雪淼,温进杰

1.华南农业大学资源环境学院,广东 广州 510642; 2.自然资源部建设用地再开发重点实验室,广东 广州 510642; 3.广东省土地信息工程技术研究中心,广东 广州 510642; 4.华中科技大学物理学院地球物理研究所,湖北 武汉 430074; 5.西南石油大学土木工程与测绘学院,成都 四川 610500

频发的严重干旱事件导致区域乃至全球水资源和粮食短缺,引发经济和社会危机。干旱指数是监测、评估与预报干旱的基础。目前,国内外建立的干旱指数已达上百种,广泛应用于气象、农业和水文干旱监测[1-5]。然而,由于影响干旱的因素众多,干旱发生发展过程复杂,基于单一数据的干旱指数难以综合反映干旱事件的多尺度特征及其综合影响,容易出现干旱情势的误报和漏报[6]。因此,利用多种水循环要素[7]发展融合多源数据的综合干旱指数,以全面刻画干旱的时空演变过程及其影响成为当前研究的热点问题[4-5,8]。

由于学术界对综合干旱指数的概念未达成共识,至今没有统一明确的定义[5]。世界气象组织将通过加权或建模方法合并不同干旱指数得到的结果统称为综合干旱指数。综合干旱指数的构建主要是通过权重组合、多变量联合分布和机器学习等方法融合降雨、蒸散发、径流、土壤水、地下水等实测或模拟水文气象多源观测信息[1-2,4-5,8]。然而,由于传统地面监测、水文模型、多光谱和高光谱卫星遥感等技术手段难以准确获取大范围深层土壤水和地下水信息,一定程度上制约了现有综合干旱指数监测干旱的能力[5]。

近20年来,GNSS和GRACE及GRACE Follow-On(GRACE-FO)卫星重力测量成为监测水资源整体变化(包括地表水、土壤水、地下水、积雪等)的新兴技术手段[9],为构建理想的水文干旱指数提供了新的数据。全球大部分GNSS测站垂向形变主要是由水文负荷变化引起的。因此,扣除大气、潮汐、板块运动等非陆地水文负荷形变效应后,GNSS垂向形变可用于反演区域水储量变化[10-11]、监测干旱过程[12-13]和构建干旱指数[14-15],在测站密集的区域,空间分辨率可达数十千米。自2002年4月以来,GRACE/GRACE-FO卫星重力测定的水储量变化具有全球精度一致的优势,被广泛应用于监测大范围干旱[16-17]和构建干旱指数[18-20]。

然而,受数据处理策略[21-22]及测站局地复杂环境因素[23-26]的影响,相邻GNSS站点在反映大尺度水文负荷形变时可能出现较大差异[14,23],特别是构造运动活跃的区域,从而影响GNSS监测干旱的精度。卫星重力数据的空间分辨率仅为300 km左右,限制了其监测小范围干旱事件的能力[16]。因此,为了充分发挥密集GNSS站网较高的空间分辨率、GRACE/GRACE-FO数据较好的空间一致性的优势,同时顾及气象干旱与水文干旱的综合影响,本文基于我国西南地区陆态网GNSS站点的垂向形变、GRACE/GRACE-FO数据和我国使用的综合干旱指数(composite index,CI)数据,利用干旱指数融合常用的主成分分析(PCA)方法[27-29],构建能从水文气象角度全面描述干旱的综合干旱指数(combined drought index,CDI),并与改进的帕默尔干旱指数(scPDSI)进行对比验证。

1 GNSS垂向形变干旱指数构建

本文利用中国地震局GNSS数据产品服务平台提供的30个GNSS站点(图1)每日的垂向形变时间序列产品对中国西南地区进行分析。研究区域(99°E—104°E,22°N—30°N)包含云南和四川两个省份。两个站点之间的最小距离约为20 km,平均距离约为100 km。

图1 研究区GNSS与气象站点空间分布

GNSS垂向形变产品已扣除地球固体潮、极潮和海潮的影响[11]。在此基础上,采用二阶Butterworth滤波削弱随机噪声和非连续粗差的影响[13]。针对出现连续粗差的GNSS台站,进一步剔除各站点垂向形变绝对值大于3倍标准差的数据。各GNSS站点的大气负荷形变和非潮汐海洋负荷形变均由德国地学研究中心(GFZ)提供的格网数据产品插值得到,并将其从GNSS垂向形变中扣除,得到水文负荷形变。最后扣除包含板块运动和冰川均衡调整信号的线性趋势,并把每日的GNSS数据转换为月间隔。

在GNSS水文负荷形变的基础上,通过扣除气候平均值,并进行标准化处理得到GNSS垂向位移干旱指数(GNSS vertical deformation-derived hydrological drought index,GNSSVD-HDI)[14]

(1)

2 基于PCA的综合干旱指数构建方法

参与CDI指数构建的数据包括GRACE干旱强度指数(GRACE drought severity index,GRACE-DSI)和中国气象局国家气候中心提供的综合干旱指数(CI)数据。GRACE-DSI是基于JPL提供的GRACE/GRACE-FOMascon格网数据产品,利用文献[19]的方法计算得到,反映了所有水文要素的变化信息。作为国家气候中心用于气象干旱监测的指标,CI站点位置如图1所示,利用近30天(月尺度)和近90天(季尺度)标准化降水指数,以及近30天湿润指数计算而得。首先通过最近点插值得到GNSS站点上的GRACE-DSI和CI。然后,在各个站点上利用PCA将二者与GNSSVD-HDI进行融合得到CDI。需要指出的是,本文仅对GNSSVD-HDI、GRACE-DSI和CI 3种指数均有数据的情况下进行融合,对缺失的数据未做插补处理。

PCA方法通过将高维数据经线性投影到低维空间,同时尽可能多地保留信息量(方差最大方向),达到数据降维的目的。计算步骤[25]如下。

(1) 考虑到不同干旱指数量级可能存在不一致的现象,首先对3种干旱指数的原始数据分别进行标准化处理

(2)

(2) 计算协方差阵

(3)

(3) 求解主成分。计算协方差矩阵Sij的p个特征值(设为λ1>λ2>…>λp)和特征向量v1,v2,…,vp,则样本的p个主成分为

(4)

式中,i=1,2,…,p;j=1,2,…,n。

(4) 计算综合干旱指数CDI为

(5)

式中,第i个主成分Di和主成分权重系数wi分别为

Di=v1GHSSVD-HDI+v2GRACE-DSI+

v3CI

(6)

(7)

由于第一主成分保留原数据集绝大部分的信息[28-29],因此本文通过PCA方法计算第一主成分的权重来构建CDI。为了验证CDI的精度和可靠性,本文计算其与scPDSI的相关系数和均方根误差(RMSE)。scPDSI包含了8种水循环要素的信息,被认为是一种干旱分析的多变量指标[27]。该指数根据当地条件重新校准,能更好地根据当地气候进行评估干旱,已用于验证GRACE-DSI[17]和基于GNSS的干旱指数[28]。由于该数据更新较慢,scPDSI与其他干旱指数相比,缺少2020年的数据,数据详细信息见表1。

表1 采用的数据信息

3 结果与分析

3.1 单一干旱指数结果分析

图2为30个站点GNSSVD-HDI、GRACE-DSI、CI 3种单一干旱指数两两之间的相关系数。由图2可知,3种单一干旱指数在西南地区的适用性差异较为明显。其中,GNSSVD-HDI与GRACE-DSI、GNSSVD-HDI与CI、GRACE-DSI与CI相关性最弱的站点分别位于研究区东南部(图2(a))、北部(图2(b))、南部和北部的SCJL站点(图2(c)),且在个别站点出现了负相关的情况。结合表2可以看出,整体上GRACE-DSI与CI的相关性最强,30个站点的相关系数平均值为0.31,最大值达到0.69(YNYL站),且在17个站点的相关系数大于GNSSVD-HDI与GRACE-DSI/CI的相关系数,但相关系数的波动范围也最大,在-0.18～0.69之间。30个站点中,GNSSVD-HDI与GRACE-DSI/CI的相关系数高于GRACE-DSI与CI相关系数的站点分别为5个和8个。

表2 GNSSVD-HDI、GRACE-DSI、CI在30个站点两两之间相关系数统计

图3给出了不同情况下3个站点上两种干旱指数与第3种干旱指数差异较为明显的时间序列。在YNGM站点,GNSSVD-HDI与GRACE-DSI符合较好,相关系数为0.52,而与CI差异较大,CI与GNSSVD-HDI/GRACE-DSI的相关系数分别为0.04和0.10。在YNMZ站点,GNSSVD-HDI与CI符合较好,相关系数为0.43,而与GRACE-DSI差异较为明显,GRACE-DSI与GNSSVD-HDI/CI的相关系数分别为-0.17和-0.18。在YNYL站点,虽然3种单一干旱指数的相关性均相对较强,但GRACE-DSI与CI符合更好,相关系数为0.69,而与GNSSVD-HDI存在一定的差异,GNSSVD-HDI与GRACE-DSI/CI的相关系数分别为0.41和0.37。

图3 3个站点上两种干旱指数与第3种干旱指数差异较为明显的时间序列

图4给出了3种干旱指数在30个站点的时间序列和均值序列。由图4可知,3种干旱指数均能很好地反映2011—2013年和2019—2020年的长期偏干旱的状态,以及2016年至2017年6月长期偏湿润的状态,但在反映干旱的严重程度方面存在差异。其中,GNSSVD-HDI与GRACE-DSI显示2011年末至2012年初的干旱比2012年末的干旱严重、2020年的干旱比2019年的干旱严重,而CI的结果则相反。这主要是由于CI反映的是气象干旱,而GNSSVD-HDI与GRACE-DSI反映的是水文干旱。相比于GNSSVD-HDI与CI干旱指数,由于卫星重力观测具有空间精度一致的优势,30个站的GRACE-DSI更趋一致,但数据缺失严重,这将影响其评估干旱的准确性和可用性。大部分GNSS站点的GNSSVD-HDI具有较好的一致性,但受局地因素的影响,个别测站偏离较为明显。CI作为一种气象干旱指数主要反映降雨的变化,由于西南地区地形变化大造成显著的降雨空间分布差异,使得各站CI指数存在较大的波动。但相比于GNSSVD-HDI和GRACE-DSI水文干旱指数,CI无法很好地反映干旱对水资源的整体影响。例如,CI未能很好地探测到2020年的干旱事件。综合以上分析表明,3种干旱指数在西南地区的应用各有优缺点,使用单一指数均难以在所有测站获得满意的干旱监测效果,而构建融合3种干旱指数的综合指数将可能实现优势互补,提升区域干旱监测的能力。

图4 30个站GNSSVD-HDI、GRACE-DSI、CI干旱指数与均值时间序列

3.2 综合干旱指数(CDI)结果分析

为了分析单一干旱指数与CDI指数的差异,本文计算了30个站4种干旱指数与scPDSI的相关系数和RMSE,以及CDI与scPDSI相关系数与其他3种单一干旱指数与scPDSI相关系数的差值(图5)。由图5(a)可以看出,4种干旱指数与scPDSI在各站的相关系数差异较为明显,这主要是由于各指数所反映的水循环要素以及时空分辨率不同导致的。例如,虽然GNSSVD-HDI和GRACE-DSI均反映了水资源的整体变化信息,但二者的空间分辨率以及数据所受影响因素不同。GNSS数据的空间分辨率高于GRACE/GRACE-FO数据,但其受到局地地壳运动、人类活动等非水文负荷的影响更多,使其与scPDSI的相关系数较低。CI与scPDSI在各站的相关系数差异较大的原因除了二者所反映的水循环要素不完全相同之外,较大地形起伏也使得单站CI指数与格网scPDSI存在一定的空间差异。

图5 CDI、GNSSVD-HDI、GRACE-DSI、CI与scPDSI的相关系数与均方根误差

由表3可知,相比于单一干旱指数,CDI与scPDSI的相关系数整体上均有提升,最大相关系数达到0.84,30个站点的平均相关系数从0.34(GNSSVD-HDI)、0.57(GRACE-DSI)、0.48(CI)提升到0.64。由于scPDSI比CI指数包含了更多水循环要素的信息,相关系数的提高说明融合3种干旱指数能提升CI指数在西南地区的干旱监测能力,但对于不同干旱指数,融合后的结果提升程度略有差别。相比于单一的CI指数,CDI与scPDSI的相关系数在所有站点均有提高,相关系数最大提高0.36(图5(b)中的SCMB站),平均提高0.16,说明经数据融合后CI指数监测干旱综合影响的水平全面提升。相比于单一的GNSSVD-HDI指数,除了YNLA、YNMZ两个站点(图5(b)),融合后的相关系数在28个站明显提高,最大提高达0.64(图5(b)中的YNTH站),相关系数平均提高0.29,说明数据融合后显著改善了GNSSVD-HDI在绝大部分测站的干旱监测能力。对于GRACE-DSI干旱指数,融合后67%的站点(20个站)的相关系数有所提高,最大提升0.31(图5(b)中的YNLA站),相关系数平均提高0.07。从RMSE结果来看,3种单一干旱指数、CDI与scPDSI的偏差差异较小,均小于0.3(图5(c)和表2)。

表3 CDI、GNSSVD-HDI、GRACE-DSI、CI与scPDSI的相关性和RMSE统计

进一步分析图5(b)中相关系数下降最大的YNMZ站可以发现,相比于GNSSVD-HDI与scPDSI的相关系数(0.59),GRACE-DSI、CI与scPDSI的相关系数相对较低,约为0.37(图5(a))。图6(c)也显示GRACE-DSI与CI的一致性较低(相关系数为0.10),且二者与GNSSVD-HDI的相关性最大仅为0.28(图6(a)和图6(b)),说明当3种干旱指数中两种指数一致性较差,且与第3种指数的相关性也较低时,PCA提取的主要信息受一致性较差的两种干旱指数的影响更大,数据降维可能会引起有用信息的损失,进而降低融合后的效果。因此,文献[30]采用核熵成分分析方法构建综合干旱指数,该方法能以最小化的特征集最大限度地保留输入数据集的信息量。而在相关系数提高最大的YNTH站,图6(c)显示GRACE-DSI与CI具有较强的相关性(相关系数为0.48),而GNSSVD-HDI与二者相关系数最大值仅为0.08(图6(c)和图6(d)),但是融合后3种指数与scPDSI的相关性均提高,且GNSSVD-HDI提高最大。说明3种干旱指数中两种指数一致性较好时,即使与第3种干旱指数的相关性较低,PCA方法提取的主要信息受一致性较好的两种干旱指数的影响更大,融合后可以显著提升第3种干旱指数的监测能力。

4 结 论

本文在构建基于陆态网30个站点GNSS垂向形变和GRACE/GRACE-FO数据的单一干旱指数的基础上,利用PCA方法将二者与国内常用的CI指数进行融合,在西南地区建立了顾及气象干旱和水文干旱的综合干旱指数(CDI)。研究结果表明,融合GNSS、GRACE/GRACE-FO和气象数据的CDI指数提升了西南地区干旱监测的整体水平。相比于3种单一干旱指数,CDI指数与scPDSI的相关性分别在30、28和20个站点得到提升,平均相关系数分别从0.34、0.48、0.34、0.57提高到0.64。研究结果为融合大地测量观测与其他数据构建高精度新型综合干旱指数提供了新的途径,促进大地测量学与水文学等学科的交叉融合。

对于融合效果欠佳的部分测站,今后的研究需要进一步分析不同数据包含的干旱特征信息,创新数据融合方法,达到更佳的融合效果,进一步提升干旱的监测能力。此外,本文仅对同时存在GNSS、GRACE/GRACE-FO和CI 3种数据的情况进行数据融合,对缺失数据进行有效补充再进行融合,将能获得时间连续性更好的综合干旱指数,有利于从干旱发生发展过程和机制的角度全面评价综合干旱指数的优势。

致谢：感谢中国地震局GNSS数据产品服务平台(http:∥www.cgps.ac.cn/)提供数据支撑。