應用GPS垂向位移定量分析2011—2020年云南省極端干旱時空特征

楊興海， 袁林果*， 姜中山， 馮顯杰， 湯苗， 李長海

1 西南交通大學地球科學與環境工程學院， 成都 611756 2 云南大學國際河流與生態安全研究院， 昆明 650091

0 引言

干旱是一種漸發性的自然災害，根據影響干旱的自然因素將其分為氣象干旱、農業干旱、水文干旱與社會經濟干旱，持續的氣象干旱是引發農業干旱、水文干旱與社會經濟干旱的主要因素.干旱嚴重制約我國社會經濟的發展，近年來極端干旱事件的爆發給我國造成了巨大的經濟損失.分析干旱的嚴重程度與時空分布特征，對政府制定抗旱與災后生態環境修復決策具有重要意義.干旱指數是評估干旱特征的重要指標(Zhao et al.，2017a；Palmer，1965)，傳統的氣象干旱指數通常基于地面觀測站和遙感衛星所觀測的氣候變量、水文變量，如標準化降水指數SPI(Standardized Precipitation Index，McKee et al.，1993)、標準化徑流指數(Shukla and Wood，2008)、標準化地下水位指數(Bloomfield and Marchant，2013)、帕爾默干旱指數PDSI(Palmer Drought Severity Index，Palmer，1965)、自校正帕默爾干旱指數SCPDSI(Self-Calibrating Palmer Drought Severity Index，Wells et al.，2004；van der Schrier et al.，2013)、標準化降水蒸散指數SPEI(Standardized Precipitation Evapotranspiration Index，Vicente-Serrano et al.，2010)與植被條件指數(Kogan，1995)等.但傳統地面觀測站存在空間分布不均、區域覆蓋率低等制約，難以反映大尺度干旱特征，遙感衛星能觀測到地表淺層(土壤濕度、植被覆蓋等)的水文成分，但無法監測更深層次的水文成分.隨著大地測量技術的發展，GRACE系列重力衛星(GRACE / GRACE-FO)提供了全球重力場變化的觀測數據(Zhou et al.，2018；Duan et al.，2022)，該數據能有效記錄全球范圍內陸地水儲量的變化(Syed et al.，2008；Zhao et al.，2018；Rao and Sun.，2022).大量研究表明，GRACE數據反演的陸地水儲量與降水量、土壤水及地下水的時空相關性較好(李瓊等，2013；Han et al.，2019)，因此GRACE陸地水儲量可作為研究干旱時空分布的補充數據(Ma et al.，2017).為更準確定量評估干旱嚴重程度，Zhao等(2017b)提出GRACE干旱指數GRACE-DSI(GRACE Drought Severity Index).GRACE-DSI能反映大尺度陸地水儲量虧損情況，有更廣泛的適用性，補充了傳統氣象干旱指數的不足(Liu et al.，2020)，同時與氣象干旱指數相關性好(冉艷紅等，2021；Satish Kumar et al.，2021).但GRACE的空間分辨率為250～300 km，因此GRACE-DSI更適用于評估全球范圍的干旱事件.此外，GRACE時間分辨率為1個月，還存在較多數據缺失，這限制了GRACE-DSI監測短期干旱事件.

陸地水的質量變遷會引起地球的彈性響應.中國陸態網應用GPS大地測量技術近實時地高精度監測地表位移(Jin et al.，2019；Zhang et al.，2019)，為中國區域陸地水儲量變化的相關研究提供了新的數據來源.大量研究表明，陸態網提供的GPS垂向位移時間序列包含線性趨勢信號、階躍信號、季節性信號及白噪聲.其中，季節性信號包含大氣、海洋與陸地水等環境負荷的變化、基巖的熱膨脹效應(Yan et al.，2009)及共模誤差等.Blewitt和Lavallée(2002)與姜衛平等(2018)對GNSS時間序列的線性項、季節項與階躍項進行建模，并利用最小二乘方法求解時間序列模型的各項參數，使季節性信號得以分離.此外，Yuan等(2018a，b)與Ming等(2017)分別使用主成分分析與獨立成分分析(ICA, Independent Component Analysis)時空濾波，來降低陸態網GPS時間序列中的共模誤差等有色噪聲信號，閆昊明等(2010)對基巖熱膨脹效應引起的位移進行建模，使得季節性信號中的環境負荷位移得以分離.在此基礎上扣除環境負荷位移中的大氣負荷位移與海洋負荷位移后，陸地水負荷位移信號得以恢復(Fu and Freymueller，2012；Fu et al.，2012；Chen et al.，2013；Chew and Small，2014；胡順強等，2021).由于陸地水負荷與地表形變的關系可通過格林函數(Farrell，1972)或球諧函數(沈迎春等，2017)建立，因此國內外學者基于負荷位移理論、使用GPS垂向位移反演陸地水儲量(Argus et al.，2014；Jiang et al.，2021a；何思源等，2018；成帥等，2021).GPS陸地水儲量為單天解，相比GRACE陸地水儲量，GPS陸地水儲量的時間分辨率更高，同時GPS垂向位移對小尺度的水文負荷變化很敏感，因此更適用于研究小尺度的水文干旱事件.有學者基于GPS陸地水儲量計算了水儲量異常(Argus et al.，2017；Jiang et al.，2021b)，以反映陸地水儲量的虧損情況，進而評估干旱事件的時空分布特征.Jiang等(2021a)參考GRACE-DSI提出GPS干旱指數GPS-DSI(GPS Drought Severity Index)，并將其作為分析干旱嚴重程度的指標.

云南省位于中國版圖的西南部(97.31°—106.11°E，21.80°—29.15°N)，地處青藏高原東南側，屬于低緯度地區，海拔落差明顯，地勢自西北向東南呈階梯狀遞減，各階梯的溫差較大.區域內分布著眾多水系，以金沙江、瀾滄江與怒江為主要代表.同時，該地區屬于典型的亞熱帶高原季風氣候，受到南亞季風與東亞季風的影響，具有干季、濕季分明的特點.全區地表水資源豐富、年平均降水量充沛，但降水的時間和空間分布不均，且常出現季風活動和其他氣候條件的異常變化，導致旱災時有發生.例如，在全球變暖、厄爾尼諾與異常反氣旋的影響下，2019—2020年云南省發生了嚴重干旱，期間降水量不足歷史平均值的一半，氣溫達到歷史新高(Ding and Gao，2020).此次干旱是近十年來最嚴重的一次，其干旱覆蓋范圍廣、干旱持續時間長，期間主要河流湖泊水位急劇下降、多數水源干涸，使得飲用水短缺、農產品大幅減產、經濟損失嚴重，也是導致森林火災頻繁發生的重要原因(Du et al.，2021)，因此研究該地區干旱的時空分布特征具有重要意義.此前，已有學者基于GRACE時變重力場數據，分析西南地區干旱事件中陸地水儲量的變化(李瓊等，2013；Long et al.，2014；Wu et al.，2020)，同時西南地區密集的陸態網GPS站也為分析干旱事件提供了有力的數據基礎.Jiang等(2017)基于GPS水汽反演技術獲得可降水量，結合GPS垂向位移研究了云南省的干旱事件，結果表明云南省的GPS可降水量、GPS垂向位移與GRACE等效水高GRACE-EWH(GRACE Equivalent Water Height)的時間變化趨勢一致，驗證了利用GPS研究干旱的可行性.也有學者結合GPS垂向位移與多種氣象干旱指數，分析西南地區干旱事件(姚朝龍等，2019)，研究表明該地區的GPS垂向位移對干旱事件有顯著的響應、與氣象干旱指數一致性較高.Ma等(2021)利用GPS可降水量計算標準化降水換算指數，結合SPEI及綜合氣象干旱指數，定量評估云南省的干旱事件.然而，GPS可降水量無法反映地面實際的水儲量變化，因此在研究水文干旱時具有不確定性.何思源等(2018)與成帥等(2021)分別使用格林函數法與Slepian球諧系數法，利用GPS垂向位移反演云南地區的陸地水儲量，但沒有評估干旱事件對陸地水儲量的影響.Jiang等(2021b)基于FastICA(Fast Independent Component Analysis)反演模式，利用GPS垂向位移分析了云南地區2011—2019年的陸地水儲量虧損情況，但沒有研究2020年極端干旱事件及定量分析水文干旱特征(嚴重程度、時空分布).

本文旨在利用GPS垂向位移反演云南省時變陸地水儲量，結合GRACE產品和SCPDSI產品，研究云南省2011—2020年干旱的時空分布趨勢，并定量分析2019—2020年極端干旱的嚴重程度及時空分布.首先，利用FastICA反演模式，反演每日的GPS-EWH(GPS Equivalent Water Height).其次，將云南省按照地形地貌特征劃分為五個區域，根據GPS-EWH、GRACE-EWH和SCPDSI計算各區域的GPS-DSI、GRACE-DSI與歸一化SCPDSI(SCPDSI-Z)，并比較其相關性.最后，根據GPS-EWH計算五個區域的GPS-EWH異常時間序列，評估干旱事件中陸地水的虧損情況，同時結合GPS-DSI、GRACE-DSI與SCPDSI-Z，分析云南省2019—2020年極端干旱事件的時空分布.

1 數據來源及預處理

1.1 GPS時間序列數據

本文采用國家地震科學數據中心(https:∥data.earthquake.cn/)提供的中國陸態網GPS時間序列，其結果由GAMIT/GLOBK軟件解算得到，解算時扣除了固體潮、海潮與極潮的影響(Davis，2004).為保證觀測數據能夠覆蓋整個研究區域，采用2011年1月1日至2020年12月31日云南省內及其周邊共44個觀測站(如圖1所示)的位移時間序列.GPS時間序列包含構造運動引起的長期線性趨勢項、地球物理因素引起的季節項及儀器更換或地震引起的階躍項(姜衛平等，2018)，其中季節性信號包含由大氣、海洋與陸地水等物質的質量負荷變化引起的位移.由于陸地水負荷引起的垂向位移比水平位移更顯著(Argus et al.，2014)，因此本文使用GPS垂向位移時間序列進行研究，后文中的位移時間序列均指垂向位移的時間序列.

為消除GPS位移時間序列中由其他地球物理現象引起的位移信號、改正地震與儀器更換引起的階躍信號，利用德國地學研究中心GFZ(GeoForschungsZentrum，https:∥www.gfz-potsdam.de/en/home/)提供的非潮汐大氣負荷與非潮汐海洋負荷產品(Dill and Dobslaw，2013)扣除GPS位移時間序列中的大氣負荷與海洋負荷位移；采用歐洲中期氣象預報中心ECMWF(European Centre for Medium-Range Weather Forecasts)提供的第五代全球大氣再分析日平均地表溫度產品(https:∥cds.climate.copernicus.eu/cdsapp#!/home)，使用Yan等(2009)與閆昊明等(2010)提出的GPS基巖熱膨脹模型，改正GPS時間序列中基巖熱膨脹效應引起的位移.最后，使用最小二乘法與GPS時間序列模型(姜衛平等，2018)對各項參數進行建模，并扣除趨勢項及階躍項.采用的擬合模型如下：

(1)

其中，y(ti)為ti時刻的位移觀測值；y0為測站的標準位置；a0為線性項速度；Ck與Sk為季節項的振幅，k=1時對應周年項，k=2時對應半周年項；gj為階躍項，H為海維西特階梯函數，階躍前為0，階躍后為1，Tgj階躍的時間；e(ti)為隨機噪聲.

數據處理后得到垂向水文負荷位移時間序列，為保證反演精度還需篩選站點數據.移除數據缺失時間超過所研究時間域30%的GPS站數據，同時比較鄰近站點的時間序列信號特征，去除含有異常信號的站點數據.篩選過程中發現MMMZ站缺失數據超過5年、YNGM站含有異常信號，因此移除這兩個站的數據，得到42個GPS站的位移時間序列.

為保證反演的時變陸地水儲量的完整性，期望選擇的GPS站在每個歷元都有數據，但實際部分歷元的數據會存在缺失.因此，本文使用Liu等(2018)開發的GNSS缺失數據插值軟件(GMIS)，對數據缺失較少的GPS站進行插值.該軟件基于克里金-卡爾曼插值方法(Mardia et al.，1998)，將克里金空間域插值方法與卡爾曼濾波時間域最優估計結合.如圖2所示，云南地區所有的GPS站點都有明顯的周年振幅，表明該地區的水文負荷季節性變化顯著.GPS站的周年振幅平均值為6.9 mm，最大可達10.1 mm，最小為3.3 mm，周年振幅較大的站主要集中在云南西南部.GPS站最大沉降量對應的時間主要集中在9月下旬與10月初，該時間比最大降水量對應的月份(7月與8月)滯后約2個月.

圖2 GPS的垂向水文負荷位移周年振幅.箭頭方向為最大沉降量對應的時間，顏色圖為振幅內插結果Fig.2 Annual amplitude of GPS vertical hydrological loading displacements. The direction of arrows represents the time corresponding to the maximum settlement, and the color plot represents the amplitude interpolation result

1.2 GRACE數據

GRACE地球時變重力場數據可以反映地表質量負荷的時空變化，在陸地水儲量監測中得到廣泛的應用.美國德克薩斯大學空間研究中心CSR(Center for Space Research，http:∥www2.csr.utexas.edu/grace/)提供了時變重力場數據，其中2017年6月前的數據由GRACE重力衛星觀測得到，2018年6月之后的數據由GRACE-FO重力衛星觀測所得.為驗證GPS陸地水儲量結果的可靠性，本文采用2011年1月至2020年12月的GRACE時變重力場Mascon解產品(RL06 v02版本)與之對比.Mascon產品解算時將原始的C20系數替換為衛星激光測距的C20系數(Loomis et al.，2019)，同時還將GRACE-FO的C30系數替換為衛星激光測距的C30系數，基于Richard Peltier等(2018)提出的ICE6G-D模型改正冰后回彈.相較于球諧系數產品而言，Mascon產品不需要進行濾波、平滑或者縮放處理，已被廣泛應用于陸地水儲量變化相關研究(Save et al.，2016).Mascon解被估計在1°×1°的六邊形格網中，以0.25°×0.25°的格網形式提供，但其實際的分辨率由衛星軌道高度決定，約為250～300 km(Save et al.，2016).

1.3 降水及干旱指數

降水是陸地水補給的一個重要途徑，直接關系著氣象干旱的發生與否.本文采用國家氣象科學數據中心(http:∥data.cma.cn/)提供的月平均降水量格網數據，其時間范圍與GPS位移時間序列保持一致.該格網數據由氣象觀測站數據空間內插所得，云南省內共均勻分布109個氣象站點，所得格網空間分辨率為0.5°×0.5°.

干旱指數可直觀地反映復雜的干旱情況，有利于研究干旱嚴重程度及時空分布情況.SCPDSI由Wells等(2004)提出，是PDSI的改進方案.在SCPDSI計算過程中將參考作物替換成實際的植被覆蓋，并考慮了潛在蒸散量與積雪的動態變化，使干旱的空間特征更突出，因此SCPDSI被廣泛應用于氣象干旱的相關研究.為研究云南省水文干旱與氣象干旱的相關性，本文采用東英吉利大學氣候研究組CRU(Climatic Research Unit，https:∥lr1.uea.ac.uk/cru/data)提供的SCPDSI，時間分辨率為1個月，空間分辨率為0.5°×0.5°.

2 研究方法

2.1 時變陸地水儲量反演方法

固體地球是彈性體，陸地水的質量變化會引起地球表面的彈性形變，Farrell(1972)使用格林函數確定地表質量變化與形變之間的關系.Argus等(2014)基于圓盤負荷理論，提出利用GPS垂向位移反演陸地水儲量，其最小二乘反演模型如下：

‖Ax-b‖2+β2‖L(x)‖2→min,

(2)

其中，A為格林函數系數矩陣，基于Wang等(2012)提出的負荷勒夫數(PREM地球模型)計算得到；x為待估計的EWH向量；b為GPS垂向位移觀測向量；β為平滑因子；L為拉普拉斯平滑矩陣.式(2)理論上應考慮反演區域內與區域外所有質量負荷的影響，何思源等(2018)研究表明，基于格林函數反演方法、利用GPS垂向位移反演云南地區陸地水儲量變化時，將邊界向外擴充3°，場外負荷源的影響可以忽略不計.

Jiang等(2021b)采用FastICA反演模式，反演了云南省的時變陸地水儲量.該方法不同于單個歷元反演模式，無論觀測數據的歷元多少，其進行反演的次數都等于獨立分量IC(Independent Component)個數，因此能有效地提高反演的時間效率.并且，GPS垂向位移時間序列含有許多噪聲信號，FastICA方法可以有效過濾噪聲，具有收斂速度快、穩定性強的優點(Hsu et al.，2020；Jiang et al.，2021b；Hyv?rinen and Oja，2000).其反演步驟如下：

(1)利用GPS垂向位移時間序列構建觀測矩陣，行號為歷元的序號，列號為GPS站的序號.

(2)使用FastICA算法將觀測矩陣分解為空間函數與時間函數IC，并根據每個IC對濾波后數據的方差貢獻進行排序，以確定IC的個數，分解公式為

Xr×p=Vr×nWn×p,

(3)

其中，X為GPS垂向位移觀測矩陣；V為時間函數，垂向位移與反演的EWH共用一套時間函數；W為垂向位移空間函數；r為歷元個數；n為IC個數；p為站點個數.

(3)根據式(2)將每個垂向位移空間函數IC分別反演為EWH空間函數IC，反演過程中將邊界向外擴充3°，并采用交叉驗證方法(Jiang et al.，2021b；Matthews and Segall，1993)確定最優的平滑因子為0.012.

(4)各時間函數IC與對應的EWH空間函數IC乘積之和為等效水高總量：

Yr×m=Vr×nUn×m,

(4)

其中，Y為EWH向量；V為時間函數；U為EWH空間函數；m為等效水高格網個數.

本文沿用上述反演模式，使用2個IC來分析，共減少90%的數據方差，IC1與IC2對濾波后數據的方差貢獻率(S1,S2)分別為63.4%與36.6%.如圖3b—c所示，85.7%的GPS站IC1垂向位移空間函數絕對值大于10 mm，IC2為61.9%.其中，IC2的空間函數在滇西南較大，最大可達29 mm；而IC1則相對分布均勻，穩定在10～20 mm之間.此外，IC1的時間函數在每年第228天附近達到最小值，對應了最大降水月份(7月)；IC2的時間函數則在每年第318天附近達到最小值，滯后于最大降水月份約兩個月.兩個IC的時間函數都有短期下降或抬升的趨勢，可能與水文極端事件(極端干旱、極端降水)有關.

圖3 IC1與IC2的時間函數 (a) 與垂向位移空間函數 (b—c)Fig.3 Temporal functions (a) and vertical displacements spatial functions (b—c) of IC1 and IC2

2.2 干旱指數

用干旱指數反映水儲量的虧損情況，能直觀地評估干旱的發展過程及嚴重程度(Van Loon，2015).常用的PDSI、SPI等氣象干旱指數通常僅考慮了少數水文或氣象變量，不能完整地反映陸地水儲量的虧損情況，同時各種干旱指數考慮的水文與氣象變量不完全相同，因此各種干旱指數存在差異(Shi et al.，2021；Wang et al.，2015).Zhao等(2017a，b)計算的干旱指數GRACE-DSI，可以反映陸地水儲量的虧損情況，但受限于時空分辨率低與數據缺失.Jiang等(2021a)基于GPS-EWH計算的干旱指數GPS-DSI，同樣能反映陸地水儲量的虧損情況，且不存在數據缺失，可以作為評估干旱的一種有效補充手段.GPS-DSI計算方式如下：

(5)

為更準確地評估干旱，本文采用與式(5)相同的標準化方法，將GRACE-EWH和SCPDSI計算為GRACE-DSI與SCPDSI-Z，結合GPS-DSI、GRACE-DSI和SCPDSI-Z分析干旱特征.參考Zhao等(2017b)的GRACE-DSI分類標準，將文本中使用的干旱指數劃分為表1中的11個類別.

表1 干旱程度分類Table 1 Classification of drought

3 結果分析

3.1 陸地水儲量變化時空分布

本文估計了研究區2011—2020年的GPS-EWH日變化量，圖4展示了GPS-EWH與GRACE-EWH的周年振幅及多年平均降水量，三種數據周年變化的空間分布趨勢一致，即滇西南的周年變化最大，逐漸向滇東北方向減小.其中滇西南的GPS-EWH周年振幅最大值約為300 mm，GRACE-EWH約為220 mm；滇東北的GPS-EWH周年振幅最小值約為20 mm，GRACE-EWH約為40 mm.GPS-EWH與GRACE-EWH存在差異，其主要原因可能是GRACE的空間分辨率低，對局部區域的陸地水儲量變化不夠敏感(成帥等，2021；何思源等，2018；Argus et al.，2014)，及GPS時間序列中還含有未建模的非水文負荷位移信號(Fu et al.，2015).圖4c中多年平均降水量在滇西南達到最大值，約1750 mm；在滇東北減小到約900 mm；在滇西北達到最小值，約600 mm.多年平均降水量的總體空間分布趨勢與GPS-EWH和GRAEC-EWH周年振幅一致，表明該地區陸地水儲量變化受降水驅動.此外，滇西北降水量最小，但GPS-EWH的周年振幅較大，說明該地區陸地水儲量可能受其他因素(徑流、蒸散)影響較大.

本文為研究云南地區的陸地水儲量時間變化特征，比較了研究區域內GPS-EWH、GRACE-EWH與降水量的月度時間序列.如圖5a所示，GPS-EWH與GRACE-EWH的相關系數為0.82，兩者的時間變化趨勢具有高度一致性.但GPS-EWH與GRACE-EWH的時間序列周年振幅分別約為127 mm與64 mm，GRACE-EWH比GPS-EWH低估了50.4%的陸地水儲量.

GPS-EWH與GRACE-EWH都是陸地水的狀態量，而降水屬于通量，兩者是不同的物理量，陸地水儲量與降水的關系可由下式簡單表示：

(6)

其中，t為時間，TWSA為陸地水儲量，P為降水，R

為地表徑流量，E為蒸發量.由式(6)可知陸地水儲量由降水、徑流與蒸發量累積得到，因此圖5中GPS-EWH與GRACE-EWH的相位滯后于降水2個月，在其他相關研究中有相同結論(Jiang et al.，2017；何思源等，2018；成帥等，2021).

本文將水儲量等效水高進行微分，以此統一陸地水與降水的物理量，如圖5b所示，GPS-EWH微分、GRACE-EWH微分與降水相關性較高、時間變化趨勢一致.但三者也存在差異：(1)GRACE-EWH微分比GPS-EWH微分的變化幅度小，其原因可能是GRACE空間分辨率低；(2)GPS-EWH微分與降水的波谷不一致，結合公式(6)分析，可能是干季時陸地水儲量主要受蒸發與地表徑流影響.

圖5 陸地水儲量的時間變化.EWH時間序列(a)，EWH微分時間序列(b)，EWH異常時間序列(c).R、R1、R2與R3分別為GPS-EWH與GRACE-EWH、GPS-EWH微分與降水量、GPS-EWH微分與GRACE-EWH微分及GRACE-EWH微分與降水量的相關系數Fig.5 Temporal variations of terrestrial water storage. EWH (a), EWH differentials (b), and EWH anomalies (c). R, R1, R2, and R3 are the correlation coefficients between GPS-EWH and GRACE-EWH, GPS-EWH differential and precipitation, GPS-EWH differential and GRACE-EWH differential, as well as GRACE-EWH differential and precipitation, respectively

3.2 等效水高與干旱指數及其相關性分析

本文按照地形地貌特征將研究區域分為圖1中的五個區域(滇東北、滇東南、滇中、滇西北與滇西南)，并探究各區域的GPS-DSI、GRACE-DSI與SCPDSI-Z之間的相關性.如圖6a所示，GPS-DSI與GRACE-DSI在滇西南與滇西北地區相關性最高，兩者在滇東的相關性最低.對比圖6a與圖6b可知，在滇中與滇西地區GPS-DSI與SCPDSI-Z的相關性更高.如圖6b所示，GPS-DSI與SCPDSI-Z在滇中和滇西南地區相關性最高，相關系數約為0.85，表明在該地區GPS-DSI能作為評估干旱的補充數據.但在滇東地區GPS-DSI與SCPDSI-Z呈現出負相關，對比圖2中GPS站點分布，可能是滇東地區GPS站空間分布稀疏導致GPS-DSI不確定性較大.

圖7a展示了五個區域GPS-EWH與GRACE-EWH的時間序列及其相關系數.GPS-EWH與GRACE-EWH表現出高度的一致性，兩者在滇西北、滇西南與滇中地區的相關系數最高，分別為0.87、0.86與0.79.在滇東南與滇東北GPS-EWH與GRACE-EWH的相關系數分別為0.56與0.42，呈現出弱相關性，同時GPS-EWH與GRACE-EWH的時間變化趨勢略有區別.結合圖2分析，可能是GPS站空間分布稀疏導致了滇東南與滇東北地區的GPS-EWH具有較大不確定性.

圖7b展示了云南省五個區域GPS-DSI、GRACE-DSI與SCPDSI-Z的時間序列及它們間的相關系數R1、R2與R3.滇中、滇西南與滇西北的R2分別為0.72、0.71與0.56，大于R3和R1，表明GPS-DSI對小區域的干旱事件更敏感.此外，五個區域的GRACE-DSI與SCPDSI-Z相關性較好，R3均約為0.5，表明水文干旱與氣候干旱的一致性較好.但滇東地區的R2明顯小于R3，表明該地區的GPS-DSI可能具有較大不確定性，這與圖6b的結果相同.除滇東地區的GPS-DSI外，三種干旱指數在2019年出現大幅度下降，一直持續到2020年，表明期間有嚴重的干旱事件發生.

圖6 GPS-DSI與GRACE-DSI (a) 及GPS-DSI與SCPDSI-Z (b) 的相關系數Fig.6 Correlation coefficients between GPS-DSI and GRACE-DSI (a), as well as GPS-DSI and SCPDSI-Z (b)

圖7 云南省五個區域的EWH (a) 與干旱指數 (b) 時間序列.R、R1、R2和R3分別為GPS-EWH與GRACE-EWH、GPS-DSI與GRACE-DSI、GPS-DSI與SCPDSI-Z及GRACE-DSI與SCPDSI-Z的相關系數Fig.7 Time series of EWH (a) and drought severity indices (b) in the five regions of Yunnan Province. R, R1, R2, and R3 are the correlation coefficients between GPS-EWH and GRACE-EWH, GPS-DSI and GRACE-DSI, GPS-DSI and SCPDSI-Z, as well as GRACE-DSI and SCPDSI-Z, respectively

綜上可知，三種干旱指數存在一定的差異，其主要原因可能是：(1)滇東地區GPS站點空間分布稀疏，導致估計的GPS-EWH與GPS-DSI具有較大不確定性；(2)GPS時間序列中含有未建模的非水文負荷位移信號；(3)GRACE的空間分辨率較低，因此在小區域內GPS-DSI比GRACE-DSI更敏感；(4)GRACE衛星在任務最后幾個月數據質量不佳，導致該時期GRACE-DSI與GPS-DSI、SCPDSI-Z的差異較大；(5)SCPDSI-Z屬于氣象干旱指數，與水文干旱指數考慮的因素不同，使得它們表現出時間差異、干旱程度差異.因此，結合多種干旱指數，使其優勢互補，能更好地評估干旱事件.

3.3 干旱的時空分布特征

圖8展示了干旱指數分別在DSI≤-0.8(中度干旱及以上)、DSI≤-1.3(重度干旱及以上)與DSI≤-1.6(極端干旱及以上)干旱區間的頻率.如圖8a—c所示，云南省大部分地區發生中度及以上干旱事件的頻率超過0.2.GPS-DSI、GRACE-DSI與SCPDSI-Z在中度及以上干旱區間且頻率超過0.2的格網分別有87%、43%與74%.如圖8d—f所示，GPS-DSI與SCPDSI-Z屬重度及以上干旱區間且頻率大于0.1的格網分別有53%與38%，該類型干旱主要發生在滇中與滇南.圖8g—i表明GPS-DSI與SCPDSI-Z在極端干旱區間的頻率較低，且極端干旱主要集中在滇西南地區，GRACE-DSI在該種干旱區間的頻率最低.

圖8 DSI≤-0.8 (a—c)、DSI≤-1.3 (d—f)與DSI≤-1.6 (g—i)的頻率空間分布Fig.8 Spatial distributions of frequency of DSI≤-0.8 (a—c), DSI≤-1.3 (d—f) and DSI≤-1.6 (g—i)

圖9 云南省五個區域GPS-EWH (a)、GPS-DSI異常與降水異常時間序列(b)Fig.9 Time series of GPS-EWH (a), GPS-DSI anomalies, and precipitation anomalies (b) in the five regions of Yunnan Province

圖9b展示了五個區域的降水異常與GPS-DSI時間序列，降水異常用于評價降水量的虧損程度.滇西南、滇西北、滇中的GPS-DSI與降水異常整體變化趨勢相似，降水異常顯著上升時GPS-DSI會出現上升趨勢，降水異常顯著下降后GPS-DSI也會下降.但滇東南2018年的GPS-DSI存在異常的波谷，滇東北GPS-DSI與降水異常的變化趨勢一致性較差，因此分析這兩個地區的水文干旱時需參考GRACE-DSI.

根據圖9與圖7b中的GPS-DSI與GRACE-DSI時間序列，云南省在2011—2020年間發生了三次明顯的水文干旱事件.第一次水文干旱開始于2011年8月份左右，該次水文干旱持續時間短、干旱程度較輕，在2012年4月左右恢復為無旱.第二次開始于2012年10月左右，此次水文干旱程度與第一次相似，并于2013年12月左右恢復為無旱.第三次水文干旱開始于2019年4月左右，持續到2020年10月左右出現緩解，而后又加劇.Ma等(2017)與Jiang等(2017)分別利用GRACE陸地水儲量和GPS可降水量分析云南省干旱，均發現第一次與第二次水文干旱.Shi等(2021)通過分析滇中、滇西南與滇西北的GRACE陸地水儲量異常，發現了這三次水文干旱事件.在這三次干旱事件中，第三次水文干旱事件在10年內持續時間最長、干旱程度最嚴重.此次干旱事件中降水異常大多數月份為虧損，雖然有少許月份的降水異常為正值，但降水持續時間短，并沒有緩解水文干旱，同時在滇西南降水量虧損的月份最多，其虧損量在2019年8月達到最大值.持續的降水虧損導致GPS-DSI與GRACE-DSI大幅度下降，2020年5月大部分地區都達到異常干旱.

綜上所述，云南省時常有中度及以上的干旱發生，偶爾有重度干旱與極端干旱發生.2011—2020年共有3次顯著的干旱事件發生，其中極端干旱開始于2019年4月，持續到2020年底.

3.4 2019—2020年極端干旱的空間分布

為進一步評估2019—2020年云南省極端干旱事件的空間分布，圖10展示了此次干旱最嚴重的時期(2019年6—11月、2020年3—8月)的GPS-DSI、GRACE-DSI與SCPDSI-Z的空間分布.如圖10所示，三種干旱指數整體的干旱程度相似，2019年云南大部分地區屬于中度干旱，少部分地區達到重度干旱，2020年云南大部分地區達到重度干旱與異常干旱，直到2020年8月才出現消退的跡象.

圖10 2019年(左3列)與2020年(右3列)GPS-DSI、GRACE-DSI與SCPDSI-Z空間分布Fig.10 Spatial distributions of GPS-DSI, GRACE-DSI and SCPDSI-Z in 2019 (left 3 columns) and 2020 (right 3 columns)

具體分析圖10中GPS-DSI空間分布可知，2019年6—7月滇西北與滇中地區的干旱程度最嚴重，可達到重度干旱與極端干旱.2019年雨季后，8—10月極端干旱消退，11月重度干旱也消退，但大部分地區仍屬于中度干旱.2019年冬季與次年春季降水量偏小，陸地水沒有得到足夠的補充，導致2020年春季干旱再次加重.2020年3月滇中、滇西北和滇西南絕大部分地區發展為極端干旱，4—7月干旱范圍繼續擴大，使云南大部分地區發展為極端干旱與異常干旱.2020年8月，得到降水補充后，干旱程度有所減輕，滇西北大部分地區降低為中度干旱，僅滇西南存在部分地區屬于異常干旱.

由于滇東地區的GPS-DSI不確定性較大，因此研究此區域的干旱還需參考GRACE-DSI與SCPDSI-Z.如圖10中的GRACE-DSI所示，2019年6—8月干旱的覆蓋面積最廣，大部分地區為中度干旱，9—11月滇西北與滇東北逐漸減輕為無旱，而滇西南出現干旱加重趨勢，此結論與冉艷紅等(2021)的研究結果一致.2020年3月，滇南的GRACE-DSI仍為干旱，4—7月干旱范圍擴大、干旱程度加重，大部分地區發展為極端干旱與異常干旱，8月干旱由南向北逐漸減弱.此外，根據SCPDSI-Z可知，2019年6—11月大部分地區都為中度干旱，6—10月小部分地區出現重度干旱或極端干旱.除滇東的GPS-DSI外，2020年3—7月SCPDSI-Z的空間分布與GPS-DSI、GRACE-DSI相似，滇西北地區為無旱，而其他地區干旱覆蓋范圍廣、干旱程度高，直到8月干旱有所緩解.

綜上可知，2019—2020年云南省經歷了嚴重的干旱，其覆蓋范圍廣、持續時間長.2019年6—11月云南全區以中度干旱為主，2020年3—8月云南大部分地區以極端干旱為主，其中滇西南2020年可達異常干旱.此外，GPS-DSI、GRACE-DSI與SCPDSI-Z存在差異：(1)滇東地區的GPS-DSI比GRACE-DSI和SCPDSI-Z更濕潤；(2)2019年8—11月滇西南地區的GRACE-DSI有極端干旱和異常干旱，但GPS-DSI與SCPDSI-Z都沒有出現此類干旱；(3)2020年SCPDSI-Z的干旱程度相比GPS-DSI和GRACE-DSI低.本文3.2節已詳細討論上述差異的原因.

4 結論

本文反演了云南省2011—2020年的GPS-EWH，并結合GPS-EWH、GRACE Mascon產品和SCPDSI產品計算了GPS-DSI、GRACE-DSI和SCPDSI-Z，從而定量分析云南省干旱的時空分布特征.為分析干旱事件中水儲量虧損情況，本文基于GPS-EWH和降水產品計算了云南省五個區域的GPS-EWH異常和降水異常，并分析其時間變化趨勢.此外，展示了云南省2019年6—11月及2020年3—8月GPS-DSI、GRACE-DSI以及SCPDSI-Z的空間分布，評估了此次干旱事件的時空變化.本文的主要研究結論如下：

(1)GPS-EWH與GRACE-EWH周年振幅的空間分布一致，均在滇西南最大，并逐漸向滇東北減少.GPS-EWH與GRACE-EWH時間序列的周期性變化趨勢一致，均在每年8—10月最大、4—6月最小，但GRACE-EWH時間序列的周年振幅比GPS-EWH小50.4%.此外，GPS-EWH微分、GRACE-EWH微分與降水量相關性較高，表明陸地水儲量受降水的驅動.

(2)滇西北、滇西南與滇中GPS-DSI與SCPDSI-Z的相關系數分別為0.56、0.71與0.72，高于GPS-DSI與GRACE-DSI、GRACE-DSI與SCPDSI-Z的相關系數，表明在局部地區GPS-DSI比GRACE-DSI對干旱更敏感.但滇東南與滇東北地區GPS-DSI與SCPDSI-Z的相關性較低，對比GPS的站點空間分布，其原因可能與滇東南與滇東北地區GPS站點分布稀疏有關.

(3)分析GPS-DSI、GRACE-DSI與SCPDSI-Z的頻率結果可知，云南省時常發生中度干旱，偶爾發生重度與極端干旱.同時三種干旱指數、GPS-EWH異常與降水異常結果表明，2011—2020年共有3次顯著的干旱事件發生，其中極端干旱事件發生在2019—2020年.極端干旱期間滇中、滇西南與滇西北的GPS-EWH虧損值達10年內最大值，同時GPS-DSI、GRACE-DSI與SCPDSI-Z均大幅下降.

(4)結合GPS-DSI、GRACE-DSI和SCPDSI-Z分析2019—2020年極端干旱，2019年6—11月滇中及滇西大部分地區為中度干旱，少部分地區為重度干旱，2020年3—8月干旱加劇，云南大部分地區達到重度干旱與異常干旱.此次干旱在這10年內持續時間最長、干旱程度最嚴重.

實驗結果表明，GPS-DSI具有較高的時空分辨率，可以作為傳統干旱監測手段的補充.但GPS-DSI受限于GPS站的空間分布密度，因此結合GPS-DSI、GRACE-DSI與SCPDSI-Z可以更有效地評估干旱的嚴重程度與時空分布.此外，隨著我國北斗導航定位系統的全面建成及其他GNSS導航衛星系統的不斷發展，將來有望聯合GNSS多星座觀測數據，進行陸地水儲量及水文干旱研究.

致謝感謝國家地震科學數據中心提供的GPS坐標時間序列，感謝美國德克薩斯大學空間研究中心提供的時變重力場Mascon解產品，感謝德國地學研究中心提供的負荷位移產品，感謝中國氣象局提供的降水數據，感謝東英吉利大學氣候研究組提供的SCPDSI數據.文中部分圖使用GMT6.1開源軟件繪制 (Wessel et al., 2019).