姚朝龍,陳涌鑫,羅志才,李 瓊,葉雪淼,溫進(jìn)杰

1.華南農(nóng)業(yè)大學(xué)資源環(huán)境學(xué)院,廣東 廣州 510642; 2.自然資源部建設(shè)用地再開(kāi)發(fā)重點(diǎn)實(shí)驗(yàn)室,廣東 廣州 510642; 3.廣東省土地信息工程技術(shù)研究中心,廣東 廣州 510642; 4.華中科技大學(xué)物理學(xué)院地球物理研究所,湖北 武漢 430074; 5.西南石油大學(xué)土木工程與測(cè)繪學(xué)院,成都 四川 610500

頻發(fā)的嚴(yán)重干旱事件導(dǎo)致區(qū)域乃至全球水資源和糧食短缺,引發(fā)經(jīng)濟(jì)和社會(huì)危機(jī)。干旱指數(shù)是監(jiān)測(cè)、評(píng)估與預(yù)報(bào)干旱的基礎(chǔ)。目前,國(guó)內(nèi)外建立的干旱指數(shù)已達(dá)上百種,廣泛應(yīng)用于氣象、農(nóng)業(yè)和水文干旱監(jiān)測(cè)[1-5]。然而,由于影響干旱的因素眾多,干旱發(fā)生發(fā)展過(guò)程復(fù)雜,基于單一數(shù)據(jù)的干旱指數(shù)難以綜合反映干旱事件的多尺度特征及其綜合影響,容易出現(xiàn)干旱情勢(shì)的誤報(bào)和漏報(bào)[6]。因此,利用多種水循環(huán)要素[7]發(fā)展融合多源數(shù)據(jù)的綜合干旱指數(shù),以全面刻畫干旱的時(shí)空演變過(guò)程及其影響成為當(dāng)前研究的熱點(diǎn)問(wèn)題[4-5,8]。

由于學(xué)術(shù)界對(duì)綜合干旱指數(shù)的概念未達(dá)成共識(shí),至今沒(méi)有統(tǒng)一明確的定義[5]。世界氣象組織將通過(guò)加權(quán)或建模方法合并不同干旱指數(shù)得到的結(jié)果統(tǒng)稱為綜合干旱指數(shù)。綜合干旱指數(shù)的構(gòu)建主要是通過(guò)權(quán)重組合、多變量聯(lián)合分布和機(jī)器學(xué)習(xí)等方法融合降雨、蒸散發(fā)、徑流、土壤水、地下水等實(shí)測(cè)或模擬水文氣象多源觀測(cè)信息[1-2,4-5,8]。然而,由于傳統(tǒng)地面監(jiān)測(cè)、水文模型、多光譜和高光譜衛(wèi)星遙感等技術(shù)手段難以準(zhǔn)確獲取大范圍深層土壤水和地下水信息,一定程度上制約了現(xiàn)有綜合干旱指數(shù)監(jiān)測(cè)干旱的能力[5]。

近20年來(lái),GNSS和GRACE及GRACE Follow-On(GRACE-FO)衛(wèi)星重力測(cè)量成為監(jiān)測(cè)水資源整體變化(包括地表水、土壤水、地下水、積雪等)的新興技術(shù)手段[9],為構(gòu)建理想的水文干旱指數(shù)提供了新的數(shù)據(jù)。全球大部分GNSS測(cè)站垂向形變主要是由水文負(fù)荷變化引起的。因此,扣除大氣、潮汐、板塊運(yùn)動(dòng)等非陸地水文負(fù)荷形變效應(yīng)后,GNSS垂向形變可用于反演區(qū)域水儲(chǔ)量變化[10-11]、監(jiān)測(cè)干旱過(guò)程[12-13]和構(gòu)建干旱指數(shù)[14-15],在測(cè)站密集的區(qū)域,空間分辨率可達(dá)數(shù)十千米。自2002年4月以來(lái),GRACE/GRACE-FO衛(wèi)星重力測(cè)定的水儲(chǔ)量變化具有全球精度一致的優(yōu)勢(shì),被廣泛應(yīng)用于監(jiān)測(cè)大范圍干旱[16-17]和構(gòu)建干旱指數(shù)[18-20]。

然而,受數(shù)據(jù)處理策略[21-22]及測(cè)站局地復(fù)雜環(huán)境因素[23-26]的影響,相鄰GNSS站點(diǎn)在反映大尺度水文負(fù)荷形變時(shí)可能出現(xiàn)較大差異[14,23],特別是構(gòu)造運(yùn)動(dòng)活躍的區(qū)域,從而影響GNSS監(jiān)測(cè)干旱的精度。衛(wèi)星重力數(shù)據(jù)的空間分辨率僅為300 km左右,限制了其監(jiān)測(cè)小范圍干旱事件的能力[16]。因此,為了充分發(fā)揮密集GNSS站網(wǎng)較高的空間分辨率、GRACE/GRACE-FO數(shù)據(jù)較好的空間一致性的優(yōu)勢(shì),同時(shí)顧及氣象干旱與水文干旱的綜合影響,本文基于我國(guó)西南地區(qū)陸態(tài)網(wǎng)GNSS站點(diǎn)的垂向形變、GRACE/GRACE-FO數(shù)據(jù)和我國(guó)使用的綜合干旱指數(shù)(composite index,CI)數(shù)據(jù),利用干旱指數(shù)融合常用的主成分分析(PCA)方法[27-29],構(gòu)建能從水文氣象角度全面描述干旱的綜合干旱指數(shù)(combined drought index,CDI),并與改進(jìn)的帕默爾干旱指數(shù)(scPDSI)進(jìn)行對(duì)比驗(yàn)證。

1 GNSS垂向形變干旱指數(shù)構(gòu)建

本文利用中國(guó)地震局GNSS數(shù)據(jù)產(chǎn)品服務(wù)平臺(tái)提供的30個(gè)GNSS站點(diǎn)(圖1)每日的垂向形變時(shí)間序列產(chǎn)品對(duì)中國(guó)西南地區(qū)進(jìn)行分析。研究區(qū)域(99°E—104°E,22°N—30°N)包含云南和四川兩個(gè)省份。兩個(gè)站點(diǎn)之間的最小距離約為20 km,平均距離約為100 km。

圖1 研究區(qū)GNSS與氣象站點(diǎn)空間分布

GNSS垂向形變產(chǎn)品已扣除地球固體潮、極潮和海潮的影響[11]。在此基礎(chǔ)上,采用二階Butterworth濾波削弱隨機(jī)噪聲和非連續(xù)粗差的影響[13]。針對(duì)出現(xiàn)連續(xù)粗差的GNSS臺(tái)站,進(jìn)一步剔除各站點(diǎn)垂向形變絕對(duì)值大于3倍標(biāo)準(zhǔn)差的數(shù)據(jù)。各GNSS站點(diǎn)的大氣負(fù)荷形變和非潮汐海洋負(fù)荷形變均由德國(guó)地學(xué)研究中心(GFZ)提供的格網(wǎng)數(shù)據(jù)產(chǎn)品插值得到,并將其從GNSS垂向形變中扣除,得到水文負(fù)荷形變。最后扣除包含板塊運(yùn)動(dòng)和冰川均衡調(diào)整信號(hào)的線性趨勢(shì),并把每日的GNSS數(shù)據(jù)轉(zhuǎn)換為月間隔。

在GNSS水文負(fù)荷形變的基礎(chǔ)上,通過(guò)扣除氣候平均值,并進(jìn)行標(biāo)準(zhǔn)化處理得到GNSS垂向位移干旱指數(shù)(GNSS vertical deformation-derived hydrological drought index,GNSSVD-HDI)[14]

(1)

2 基于PCA的綜合干旱指數(shù)構(gòu)建方法

參與CDI指數(shù)構(gòu)建的數(shù)據(jù)包括GRACE干旱強(qiáng)度指數(shù)(GRACE drought severity index,GRACE-DSI)和中國(guó)氣象局國(guó)家氣候中心提供的綜合干旱指數(shù)(CI)數(shù)據(jù)。GRACE-DSI是基于JPL提供的GRACE/GRACE-FOMascon格網(wǎng)數(shù)據(jù)產(chǎn)品,利用文獻(xiàn)[19]的方法計(jì)算得到,反映了所有水文要素的變化信息。作為國(guó)家氣候中心用于氣象干旱監(jiān)測(cè)的指標(biāo),CI站點(diǎn)位置如圖1所示,利用近30天(月尺度)和近90天(季尺度)標(biāo)準(zhǔn)化降水指數(shù),以及近30天濕潤(rùn)指數(shù)計(jì)算而得。首先通過(guò)最近點(diǎn)插值得到GNSS站點(diǎn)上的GRACE-DSI和CI。然后,在各個(gè)站點(diǎn)上利用PCA將二者與GNSSVD-HDI進(jìn)行融合得到CDI。需要指出的是,本文僅對(duì)GNSSVD-HDI、GRACE-DSI和CI 3種指數(shù)均有數(shù)據(jù)的情況下進(jìn)行融合,對(duì)缺失的數(shù)據(jù)未做插補(bǔ)處理。

PCA方法通過(guò)將高維數(shù)據(jù)經(jīng)線性投影到低維空間,同時(shí)盡可能多地保留信息量(方差最大方向),達(dá)到數(shù)據(jù)降維的目的。計(jì)算步驟[25]如下。

(1) 考慮到不同干旱指數(shù)量級(jí)可能存在不一致的現(xiàn)象,首先對(duì)3種干旱指數(shù)的原始數(shù)據(jù)分別進(jìn)行標(biāo)準(zhǔn)化處理

(2)

(2) 計(jì)算協(xié)方差陣

(3)

(3) 求解主成分。計(jì)算協(xié)方差矩陣Sij的p個(gè)特征值(設(shè)為λ1>λ2>…>λp)和特征向量v1,v2,…,vp,則樣本的p個(gè)主成分為

(4)

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

(4) 計(jì)算綜合干旱指數(shù)CDI為

(5)

式中,第i個(gè)主成分Di和主成分權(quán)重系數(shù)wi分別為

Di=v1GHSSVD-HDI+v2GRACE-DSI+

v3CI

(6)

(7)

由于第一主成分保留原數(shù)據(jù)集絕大部分的信息[28-29],因此本文通過(guò)PCA方法計(jì)算第一主成分的權(quán)重來(lái)構(gòu)建CDI。為了驗(yàn)證CDI的精度和可靠性,本文計(jì)算其與scPDSI的相關(guān)系數(shù)和均方根誤差(RMSE)。scPDSI包含了8種水循環(huán)要素的信息,被認(rèn)為是一種干旱分析的多變量指標(biāo)[27]。該指數(shù)根據(jù)當(dāng)?shù)貤l件重新校準(zhǔn),能更好地根據(jù)當(dāng)?shù)貧夂蜻M(jìn)行評(píng)估干旱,已用于驗(yàn)證GRACE-DSI[17]和基于GNSS的干旱指數(shù)[28]。由于該數(shù)據(jù)更新較慢,scPDSI與其他干旱指數(shù)相比,缺少2020年的數(shù)據(jù),數(shù)據(jù)詳細(xì)信息見(jiàn)表1。

表1 采用的數(shù)據(jù)信息

3 結(jié)果與分析

3.1 單一干旱指數(shù)結(jié)果分析

圖2為30個(gè)站點(diǎn)GNSSVD-HDI、GRACE-DSI、CI 3種單一干旱指數(shù)兩兩之間的相關(guān)系數(shù)。由圖2可知,3種單一干旱指數(shù)在西南地區(qū)的適用性差異較為明顯。其中,GNSSVD-HDI與GRACE-DSI、GNSSVD-HDI與CI、GRACE-DSI與CI相關(guān)性最弱的站點(diǎn)分別位于研究區(qū)東南部(圖2(a))、北部(圖2(b))、南部和北部的SCJL站點(diǎn)(圖2(c)),且在個(gè)別站點(diǎn)出現(xiàn)了負(fù)相關(guān)的情況。結(jié)合表2可以看出,整體上GRACE-DSI與CI的相關(guān)性最強(qiáng),30個(gè)站點(diǎn)的相關(guān)系數(shù)平均值為0.31,最大值達(dá)到0.69(YNYL站),且在17個(gè)站點(diǎn)的相關(guān)系數(shù)大于GNSSVD-HDI與GRACE-DSI/CI的相關(guān)系數(shù),但相關(guān)系數(shù)的波動(dòng)范圍也最大,在-0.18～0.69之間。30個(gè)站點(diǎn)中,GNSSVD-HDI與GRACE-DSI/CI的相關(guān)系數(shù)高于GRACE-DSI與CI相關(guān)系數(shù)的站點(diǎn)分別為5個(gè)和8個(gè)。

表2 GNSSVD-HDI、GRACE-DSI、CI在30個(gè)站點(diǎn)兩兩之間相關(guān)系數(shù)統(tǒng)計(jì)

圖3給出了不同情況下3個(gè)站點(diǎn)上兩種干旱指數(shù)與第3種干旱指數(shù)差異較為明顯的時(shí)間序列。在YNGM站點(diǎn),GNSSVD-HDI與GRACE-DSI符合較好,相關(guān)系數(shù)為0.52,而與CI差異較大,CI與GNSSVD-HDI/GRACE-DSI的相關(guān)系數(shù)分別為0.04和0.10。在YNMZ站點(diǎn),GNSSVD-HDI與CI符合較好,相關(guān)系數(shù)為0.43,而與GRACE-DSI差異較為明顯,GRACE-DSI與GNSSVD-HDI/CI的相關(guān)系數(shù)分別為-0.17和-0.18。在YNYL站點(diǎn),雖然3種單一干旱指數(shù)的相關(guān)性均相對(duì)較強(qiáng),但GRACE-DSI與CI符合更好,相關(guān)系數(shù)為0.69,而與GNSSVD-HDI存在一定的差異,GNSSVD-HDI與GRACE-DSI/CI的相關(guān)系數(shù)分別為0.41和0.37。

圖3 3個(gè)站點(diǎn)上兩種干旱指數(shù)與第3種干旱指數(shù)差異較為明顯的時(shí)間序列

圖4給出了3種干旱指數(shù)在30個(gè)站點(diǎn)的時(shí)間序列和均值序列。由圖4可知,3種干旱指數(shù)均能很好地反映2011—2013年和2019—2020年的長(zhǎng)期偏干旱的狀態(tài),以及2016年至2017年6月長(zhǎng)期偏濕潤(rùn)的狀態(tài),但在反映干旱的嚴(yán)重程度方面存在差異。其中,GNSSVD-HDI與GRACE-DSI顯示2011年末至2012年初的干旱比2012年末的干旱嚴(yán)重、2020年的干旱比2019年的干旱嚴(yán)重,而CI的結(jié)果則相反。這主要是由于CI反映的是氣象干旱,而GNSSVD-HDI與GRACE-DSI反映的是水文干旱。相比于GNSSVD-HDI與CI干旱指數(shù),由于衛(wèi)星重力觀測(cè)具有空間精度一致的優(yōu)勢(shì),30個(gè)站的GRACE-DSI更趨一致,但數(shù)據(jù)缺失嚴(yán)重,這將影響其評(píng)估干旱的準(zhǔn)確性和可用性。大部分GNSS站點(diǎn)的GNSSVD-HDI具有較好的一致性,但受局地因素的影響,個(gè)別測(cè)站偏離較為明顯。CI作為一種氣象干旱指數(shù)主要反映降雨的變化,由于西南地區(qū)地形變化大造成顯著的降雨空間分布差異,使得各站CI指數(shù)存在較大的波動(dòng)。但相比于GNSSVD-HDI和GRACE-DSI水文干旱指數(shù),CI無(wú)法很好地反映干旱對(duì)水資源的整體影響。例如,CI未能很好地探測(cè)到2020年的干旱事件。綜合以上分析表明,3種干旱指數(shù)在西南地區(qū)的應(yīng)用各有優(yōu)缺點(diǎn),使用單一指數(shù)均難以在所有測(cè)站獲得滿意的干旱監(jiān)測(cè)效果,而構(gòu)建融合3種干旱指數(shù)的綜合指數(shù)將可能實(shí)現(xiàn)優(yōu)勢(shì)互補(bǔ),提升區(qū)域干旱監(jiān)測(cè)的能力。

圖4 30個(gè)站GNSSVD-HDI、GRACE-DSI、CI干旱指數(shù)與均值時(shí)間序列

3.2 綜合干旱指數(shù)(CDI)結(jié)果分析

為了分析單一干旱指數(shù)與CDI指數(shù)的差異,本文計(jì)算了30個(gè)站4種干旱指數(shù)與scPDSI的相關(guān)系數(shù)和RMSE,以及CDI與scPDSI相關(guān)系數(shù)與其他3種單一干旱指數(shù)與scPDSI相關(guān)系數(shù)的差值(圖5)。由圖5(a)可以看出,4種干旱指數(shù)與scPDSI在各站的相關(guān)系數(shù)差異較為明顯,這主要是由于各指數(shù)所反映的水循環(huán)要素以及時(shí)空分辨率不同導(dǎo)致的。例如,雖然GNSSVD-HDI和GRACE-DSI均反映了水資源的整體變化信息,但二者的空間分辨率以及數(shù)據(jù)所受影響因素不同。GNSS數(shù)據(jù)的空間分辨率高于GRACE/GRACE-FO數(shù)據(jù),但其受到局地地殼運(yùn)動(dòng)、人類活動(dòng)等非水文負(fù)荷的影響更多,使其與scPDSI的相關(guān)系數(shù)較低。CI與scPDSI在各站的相關(guān)系數(shù)差異較大的原因除了二者所反映的水循環(huán)要素不完全相同之外,較大地形起伏也使得單站CI指數(shù)與格網(wǎng)scPDSI存在一定的空間差異。

圖5 CDI、GNSSVD-HDI、GRACE-DSI、CI與scPDSI的相關(guān)系數(shù)與均方根誤差

由表3可知,相比于單一干旱指數(shù),CDI與scPDSI的相關(guān)系數(shù)整體上均有提升,最大相關(guān)系數(shù)達(dá)到0.84,30個(gè)站點(diǎn)的平均相關(guān)系數(shù)從0.34(GNSSVD-HDI)、0.57(GRACE-DSI)、0.48(CI)提升到0.64。由于scPDSI比CI指數(shù)包含了更多水循環(huán)要素的信息,相關(guān)系數(shù)的提高說(shuō)明融合3種干旱指數(shù)能提升CI指數(shù)在西南地區(qū)的干旱監(jiān)測(cè)能力,但對(duì)于不同干旱指數(shù),融合后的結(jié)果提升程度略有差別。相比于單一的CI指數(shù),CDI與scPDSI的相關(guān)系數(shù)在所有站點(diǎn)均有提高,相關(guān)系數(shù)最大提高0.36(圖5(b)中的SCMB站),平均提高0.16,說(shuō)明經(jīng)數(shù)據(jù)融合后CI指數(shù)監(jiān)測(cè)干旱綜合影響的水平全面提升。相比于單一的GNSSVD-HDI指數(shù),除了YNLA、YNMZ兩個(gè)站點(diǎn)(圖5(b)),融合后的相關(guān)系數(shù)在28個(gè)站明顯提高,最大提高達(dá)0.64(圖5(b)中的YNTH站),相關(guān)系數(shù)平均提高0.29,說(shuō)明數(shù)據(jù)融合后顯著改善了GNSSVD-HDI在絕大部分測(cè)站的干旱監(jiān)測(cè)能力。對(duì)于GRACE-DSI干旱指數(shù),融合后67%的站點(diǎn)(20個(gè)站)的相關(guān)系數(shù)有所提高,最大提升0.31(圖5(b)中的YNLA站),相關(guān)系數(shù)平均提高0.07。從RMSE結(jié)果來(lái)看,3種單一干旱指數(shù)、CDI與scPDSI的偏差差異較小,均小于0.3(圖5(c)和表2)。

表3 CDI、GNSSVD-HDI、GRACE-DSI、CI與scPDSI的相關(guān)性和RMSE統(tǒng)計(jì)

進(jìn)一步分析圖5(b)中相關(guān)系數(shù)下降最大的YNMZ站可以發(fā)現(xiàn),相比于GNSSVD-HDI與scPDSI的相關(guān)系數(shù)(0.59),GRACE-DSI、CI與scPDSI的相關(guān)系數(shù)相對(duì)較低,約為0.37(圖5(a))。圖6(c)也顯示GRACE-DSI與CI的一致性較低(相關(guān)系數(shù)為0.10),且二者與GNSSVD-HDI的相關(guān)性最大僅為0.28(圖6(a)和圖6(b)),說(shuō)明當(dāng)3種干旱指數(shù)中兩種指數(shù)一致性較差,且與第3種指數(shù)的相關(guān)性也較低時(shí),PCA提取的主要信息受一致性較差的兩種干旱指數(shù)的影響更大,數(shù)據(jù)降維可能會(huì)引起有用信息的損失,進(jìn)而降低融合后的效果。因此,文獻(xiàn)[30]采用核熵成分分析方法構(gòu)建綜合干旱指數(shù),該方法能以最小化的特征集最大限度地保留輸入數(shù)據(jù)集的信息量。而在相關(guān)系數(shù)提高最大的YNTH站,圖6(c)顯示GRACE-DSI與CI具有較強(qiáng)的相關(guān)性(相關(guān)系數(shù)為0.48),而GNSSVD-HDI與二者相關(guān)系數(shù)最大值僅為0.08(圖6(c)和圖6(d)),但是融合后3種指數(shù)與scPDSI的相關(guān)性均提高,且GNSSVD-HDI提高最大。說(shuō)明3種干旱指數(shù)中兩種指數(shù)一致性較好時(shí),即使與第3種干旱指數(shù)的相關(guān)性較低,PCA方法提取的主要信息受一致性較好的兩種干旱指數(shù)的影響更大,融合后可以顯著提升第3種干旱指數(shù)的監(jiān)測(cè)能力。

4 結(jié) 論

本文在構(gòu)建基于陸態(tài)網(wǎng)30個(gè)站點(diǎn)GNSS垂向形變和GRACE/GRACE-FO數(shù)據(jù)的單一干旱指數(shù)的基礎(chǔ)上,利用PCA方法將二者與國(guó)內(nèi)常用的CI指數(shù)進(jìn)行融合,在西南地區(qū)建立了顧及氣象干旱和水文干旱的綜合干旱指數(shù)(CDI)。研究結(jié)果表明,融合GNSS、GRACE/GRACE-FO和氣象數(shù)據(jù)的CDI指數(shù)提升了西南地區(qū)干旱監(jiān)測(cè)的整體水平。相比于3種單一干旱指數(shù),CDI指數(shù)與scPDSI的相關(guān)性分別在30、28和20個(gè)站點(diǎn)得到提升,平均相關(guān)系數(shù)分別從0.34、0.48、0.34、0.57提高到0.64。研究結(jié)果為融合大地測(cè)量觀測(cè)與其他數(shù)據(jù)構(gòu)建高精度新型綜合干旱指數(shù)提供了新的途徑,促進(jìn)大地測(cè)量學(xué)與水文學(xué)等學(xué)科的交叉融合。

對(duì)于融合效果欠佳的部分測(cè)站,今后的研究需要進(jìn)一步分析不同數(shù)據(jù)包含的干旱特征信息,創(chuàng)新數(shù)據(jù)融合方法,達(dá)到更佳的融合效果,進(jìn)一步提升干旱的監(jiān)測(cè)能力。此外,本文僅對(duì)同時(shí)存在GNSS、GRACE/GRACE-FO和CI 3種數(shù)據(jù)的情況進(jìn)行數(shù)據(jù)融合,對(duì)缺失數(shù)據(jù)進(jìn)行有效補(bǔ)充再進(jìn)行融合,將能獲得時(shí)間連續(xù)性更好的綜合干旱指數(shù),有利于從干旱發(fā)生發(fā)展過(guò)程和機(jī)制的角度全面評(píng)價(jià)綜合干旱指數(shù)的優(yōu)勢(shì)。

致謝：感謝中國(guó)地震局GNSS數(shù)據(jù)產(chǎn)品服務(wù)平臺(tái)(http:∥www.cgps.ac.cn/)提供數(shù)據(jù)支撐。