楊鑫春 游 為 萬祥禹 宋夢(mèng)芝

1 西南交通大學(xué)地球科學(xué)與環(huán)境工程學(xué)院，成都市犀安路999號(hào)，611756 2 西南科技大學(xué)環(huán)境與資源學(xué)院，四川省綿陽市青龍大道中段59號(hào)，621010

GRACE衛(wèi)星于2017-10退役，其提供的數(shù)據(jù)只更新至2017-06。但其繼任者GRACE-FO衛(wèi)星直到2018-05才發(fā)射升空，并于同年6月開始發(fā)布全球陸地水儲(chǔ)量變化(terrestrial water storage changes, TWSC)信息。兩代GRACE衛(wèi)星數(shù)據(jù)之間存在11個(gè)月的空缺期，打破了GRACE/GRACE-FO監(jiān)測TWSC的連續(xù)性，也會(huì)限制其進(jìn)一步應(yīng)用[1]。因此，填補(bǔ)11個(gè)月的數(shù)據(jù)空缺或重構(gòu)連續(xù)的TWSC尤為重要。

部分學(xué)者采用機(jī)器學(xué)習(xí)方法來重構(gòu)或預(yù)測連續(xù)的TWSC[1-2]。如Sun等[1]采用深度神經(jīng)網(wǎng)絡(luò)(deep neural network, DNN)、多變量季節(jié)性自回歸整合移動(dòng)平均模型(seasonal autoregressive integrated moving average with external variables, SARIMAX)和多元線性回歸(multiple linear regression, MLR)3種機(jī)器學(xué)習(xí)方法，重構(gòu)多種全球TWSC數(shù)據(jù)產(chǎn)品。結(jié)果表明，DNN性能略優(yōu)于SARIMAX，明顯優(yōu)于MLR。Li等[2]聯(lián)合多種統(tǒng)計(jì)分解、時(shí)間序列分解和機(jī)器學(xué)習(xí)方法重構(gòu)全球TWSC信號(hào)，與采用整體信號(hào)重構(gòu)不同，其首先分別重構(gòu)TWSC時(shí)變模態(tài)上的季節(jié)項(xiàng)、年際變化項(xiàng)和殘余項(xiàng)信號(hào)，然后將GRACE/GRACE-FO原有時(shí)變模態(tài)上的趨勢項(xiàng)信號(hào)加回到重構(gòu)的非趨勢項(xiàng)信號(hào)，以生成重構(gòu)的全頻段時(shí)變模態(tài)信號(hào)，最后將GRACE/GRACE-FO原有空間模態(tài)乘以重構(gòu)時(shí)變模型得到總的TWSC信號(hào)。同時(shí)，Li等[2]認(rèn)為在聯(lián)合使用統(tǒng)計(jì)分解方法時(shí)，相比于其他復(fù)雜的機(jī)器學(xué)習(xí)方法，雖然MLR方法在訓(xùn)練階段表現(xiàn)并非最佳，但在測試階段表現(xiàn)較好，因此MLR模型是一種可靠和穩(wěn)健的重構(gòu)方法。基于以上研究，本文將采用MLR方法，以降水、氣溫和水文模型模擬的陸地水儲(chǔ)量數(shù)據(jù)作為預(yù)測參數(shù)，重構(gòu)中國九大流域(邊界數(shù)據(jù)由資源環(huán)境科學(xué)與數(shù)據(jù)中心(http:∥www.resdc.cn)提供)連續(xù)的TWSC；同時(shí)，對(duì)比分析整體信號(hào)重構(gòu)、去趨勢項(xiàng)信號(hào)重構(gòu)和去趨勢項(xiàng)去季節(jié)項(xiàng)信號(hào)重構(gòu)等3種策略對(duì)重構(gòu)數(shù)據(jù)質(zhì)量的影響，以期提供最適合九大流域的重構(gòu)策略。

1 數(shù)據(jù)處理

1.1 GRACE/GRACE-FO數(shù)據(jù)

GRACE/GRACE-FO觀測數(shù)據(jù)產(chǎn)品具有2種不同的形式，即傳統(tǒng)的球諧系數(shù)產(chǎn)品和最新的Mascon產(chǎn)品[1]。相較于傳統(tǒng)的球諧系數(shù)產(chǎn)品，Mascon產(chǎn)品主要具有以下優(yōu)勢[3]：1)減少從陸地到海洋的泄漏誤差，增加數(shù)據(jù)的信噪比；2)在處理過程中應(yīng)用地球物理數(shù)據(jù)約束，幾乎無經(jīng)驗(yàn)后處理的平滑濾波要求，更方便非大地測量用戶的使用。本實(shí)驗(yàn)采用CSR最新發(fā)布的RL06 v02 Mascon產(chǎn)品[4]，下文簡稱CSR-M。CSR-M空間分辨率為0.25°×0.25°，所選擇的時(shí)間跨度為2002-04～2021-12，即共計(jì)183個(gè)月的GRACE數(shù)據(jù)、43個(gè)月的GRACE-FO數(shù)據(jù)以及11個(gè)月的數(shù)據(jù)空缺。需要指出的是，本文重點(diǎn)關(guān)注GRACE/GRACE-FO系統(tǒng)間11個(gè)月空缺數(shù)據(jù)的重構(gòu)問題，因此對(duì)于因GRACE/GRACE-FO衛(wèi)星儀器問題、校準(zhǔn)活動(dòng)、電池管理等原因造成的系統(tǒng)1～2個(gè)月的數(shù)據(jù)空缺采用簡單線性插值方法進(jìn)行補(bǔ)全[5]。

1.2 預(yù)測參數(shù)

降水和氣溫是陸地水儲(chǔ)量變化的重要驅(qū)動(dòng)因子[6]，故常作為重構(gòu)TWSC的主要預(yù)測參數(shù)[1-2]。此外，Sun等[1]指出雖然水文模型模擬的TWSC通常不包括地下水和表面水，但其與GRACE TWSC存在強(qiáng)相關(guān)性，增加該變量作為預(yù)測參數(shù)所重構(gòu)的TWSC性能優(yōu)于只考慮降水和氣溫作為驅(qū)動(dòng)參數(shù)的重構(gòu)數(shù)據(jù)?；诖?，本文采用ECMWF發(fā)布的月度氣候再分析數(shù)據(jù)集ERA5-Land中降水、氣溫和TWSC(包括土壤水和雪水)數(shù)據(jù)作為預(yù)測參數(shù)[7]?？紤]到3個(gè)預(yù)測參數(shù)與GRACE/GRACE-FO TWSC存在相位差[6]，選擇時(shí)間窗口前3～0個(gè)月內(nèi)的所有變量作為預(yù)測參數(shù)，即每個(gè)預(yù)測參數(shù)在2002-01～2021-12范圍內(nèi)將具有4個(gè)預(yù)測變量。同時(shí)，為了與CSR-M數(shù)據(jù)保持一致性，3個(gè)驅(qū)動(dòng)參數(shù)的空間分辨率由0.1°×0.1°重新采樣為0.25°×0.25°，并刪除2004～2009年均值。

2 重構(gòu)方法與策略

2.1 MLR模型

MLR模型是一種用于估計(jì)2個(gè)或多個(gè)自變量與1個(gè)因變量之間線性關(guān)系的統(tǒng)計(jì)方法。本文使用3個(gè)預(yù)測參數(shù)(降水、氣溫和水文模型模擬的TWSC)作為自變量，GRACE/GRACE-FO TWSC 作為因變量，其對(duì)應(yīng)的MLR模型計(jì)算公式為：

yi=β0+β1xi1+β2xi2+β3xi3+ε

(1)

式中，yi為某個(gè)流域的GRACE/GRACE-FO TWSC 時(shí)間序列；xi1、xi2、xi3為對(duì)應(yīng)流域內(nèi)3個(gè)預(yù)測參數(shù)的時(shí)間序列；β0、β1、β2、β3為待估參數(shù)；ε為模型誤差。本文使用183個(gè)月的GRACE數(shù)據(jù)，采用最小二乘方法估算MLR模型的待估參數(shù)，并采用43個(gè)月的GRACE-FO數(shù)據(jù)測試所確定MLR模型的可靠性。根據(jù)所確定的MLR模型，通過預(yù)測參數(shù)來重構(gòu)2002-04～2021-12連續(xù)的TWSC。

2.2 重構(gòu)策略

本文使用3種不同的重構(gòu)策略。重構(gòu)策略1為整體信號(hào)重構(gòu)，是對(duì)GRACE/GRACE-FO TWSC整體信號(hào)使用式(1)直接重構(gòu)生成連續(xù)的TWSC；重構(gòu)策略2為去趨勢項(xiàng)信號(hào)重構(gòu)，是在使用MLR模型重構(gòu)GRACE/GRACE-FO TWSC去長趨勢項(xiàng)信號(hào)基礎(chǔ)上加回GRACE/GRACE-FO TWSC原有的長趨勢項(xiàng)信號(hào)以生成連續(xù)的TWSC；重構(gòu)策略3為去趨勢項(xiàng)去季節(jié)項(xiàng)信號(hào)重構(gòu)，是在采用MLR模型重構(gòu)GRACE/GRACE-FO TWSC去長趨勢項(xiàng)去季節(jié)項(xiàng)信號(hào)后加回GRACE/GRACE-FO TWSC原有的長趨勢項(xiàng)和季節(jié)項(xiàng)信號(hào)以生成連續(xù)的TWSC。

為獲得信號(hào)的趨勢項(xiàng)和季節(jié)項(xiàng)，采用最小二乘分解方法對(duì)GRACE/GRACE-FO TWSC 和預(yù)測參數(shù)的時(shí)間序列進(jìn)行分解[3,8]：

Stotal=Strend+Sseasonal+Sresidual

(2)

式中，Stotal為總的原有信號(hào)；Strend為信號(hào)的長期線性趨勢項(xiàng)；Sseasonal為信號(hào)的季節(jié)項(xiàng)；Sresidual為信號(hào)的余項(xiàng)，主要包含年際變化項(xiàng)、次季節(jié)項(xiàng)和噪聲。需要說明的是，雖然11個(gè)月數(shù)據(jù)空缺會(huì)導(dǎo)致GRACE/GRACE-FO TWSC趨勢項(xiàng)值存在一定偏差，但這種偏差對(duì)重構(gòu)數(shù)據(jù)質(zhì)量的影響可以忽略[2,8]。此外，GRACE/GRACE-FO TWSC季節(jié)項(xiàng)信號(hào)在重構(gòu)時(shí)間范圍內(nèi)可認(rèn)為并未發(fā)生實(shí)質(zhì)性變化[2]，因此GRACE/GRACE-FO TWSC季節(jié)項(xiàng)信號(hào)計(jì)算也不受數(shù)據(jù)空缺的影響。

2.3 精度評(píng)定指標(biāo)

本文采用皮爾遜相關(guān)系數(shù)(Pearson correlation coefficient，CC)、歸一化均方根誤差(normalized root mean square error，NRMSE)和納什效率系數(shù)(Nash-Sutcliffe efficiency coefficient，NSE)等3種常用的精度指標(biāo)來評(píng)定重構(gòu)TWSC的性能，各指標(biāo)詳細(xì)計(jì)算公式見文獻(xiàn)[1]。CC、NRMSE和NSE的取值范圍分別為[-1,1]、[0,+∞)和(-∞,1]，CC和NSE數(shù)值越大、NRMSE數(shù)值越小，重構(gòu)數(shù)據(jù)的性能越好。此外，采用重構(gòu)的TWSC與GRACE/GRACE-FO TWSC位于測試階段(2018-06～2021-12)的RMSE值作為重構(gòu)數(shù)據(jù)的不確定性值[2]。

3 結(jié)果與討論

圖1為中國九大流域的CSR-M TWSC以及基于3種不同策略重構(gòu)的TWSC，其相應(yīng)的性能指標(biāo)和不確定性值如圖2和表1所示。基于策略1重構(gòu)的TWSC在松花江遼河、黃河、淮河、東南諸河和珠江五個(gè)流域表現(xiàn)出良好的性能(CC/NRMSE/NSE值分別優(yōu)于0.80/0.10/0.70)，但在海灤河、長江、西南諸河和內(nèi)陸河四個(gè)流域表現(xiàn)出較差的性能。這種性能差異主要是因?yàn)楸疚闹徊捎门c氣候變化有關(guān)的降水、氣溫和水文模型模擬的TWSC數(shù)據(jù)作為驅(qū)動(dòng)參數(shù)，未考慮人類活動(dòng)和冰川融化等其他影響因子。松花江遼河、黃河、淮河、東南諸河和珠江流域TWSC主要受氣候變化影響，而海灤河、長江、西南諸河和內(nèi)陸河流域TWSC除考慮氣候變化影響外，還應(yīng)考慮人類活動(dòng)或冰川融化的影響[5]。此外，重構(gòu)的TWSC在海灤河和西南諸河兩個(gè)流域存在異常的不確定性，其對(duì)應(yīng)數(shù)值分別高達(dá)13.55 cm和8.67 cm。因此，基于整體信號(hào)重構(gòu)策略重構(gòu)的TWSC不宜作為九大流域的最后重構(gòu)結(jié)果。

圖1 基于不同策略重構(gòu)的中國九大流域TWSCFig.2 The reconstructed TWSC of nine major river basins in China based on different strategies

相較于策略1，策略2可顯著提高九大流域重構(gòu)數(shù)據(jù)的性能，特別是在海灤河、長江、西南諸河以及內(nèi)陸河流域。例如在海灤河流域，重構(gòu)數(shù)據(jù)的性能指標(biāo)CC/NRMSE/NSE分別提高0.25/0.09/0.48，不確定性值降低9.71 cm。這種性能的提高得益于人類活動(dòng)和冰川融化等因素主要影響GRACE/GRACE-FO TWSC趨勢項(xiàng)信號(hào)的變化[2]，而去趨勢項(xiàng)信號(hào)重構(gòu)策略完全保留該部分信號(hào)。

相較于策略2，策略3進(jìn)一步提升重構(gòu)數(shù)據(jù)的質(zhì)量，但提升幅度并不明顯。重構(gòu)數(shù)據(jù)質(zhì)量的進(jìn)一步提升是因?yàn)椴呗?完全包含原有的GRACE/GRACE-FO TWSC 趨勢項(xiàng)和季節(jié)項(xiàng)信號(hào)，進(jìn)而能最大限度地保留人類活動(dòng)和冰川融化對(duì)TWSC的影響。此外還可以發(fā)現(xiàn)，基于策略3重構(gòu)的TWSC在海灤河、淮河和東南諸河流域存在較大的不確定性，其值分別為3.82 cm、4.01 cm和3.66 cm，這主要受GRACE/GRACE-FO 空間分辨率(從低緯度約 200 000 km2到兩極附近約 90 000 km2范圍[5])限制，較小流域的TWSC可能會(huì)存在明顯的信號(hào)泄露[9]。

圖2 中國九大流域重構(gòu)的TWSC性能指標(biāo)對(duì)比Fig.2 Comparison of performance indexes of reconstructed TWSC of nine major river basins in China

表1 不同策略重構(gòu)的中國九大流域TWSC的不確定性

根據(jù)以上分析，將基于策略3重構(gòu)的TWSC作為九大流域最后的重構(gòu)結(jié)果，表2為重構(gòu)數(shù)據(jù)的性能指標(biāo)。從表2可以看出，重構(gòu)數(shù)據(jù)的性能在不同流域存在一定差異，這主要與各自流域的GRACE/GRACE-FO 數(shù)據(jù)信噪比以及預(yù)測參數(shù)與GRACE/GRACE-FO TWSC的相關(guān)性有關(guān)。例如在信噪比以及相關(guān)性高的西南諸河流域(圖1(h)和表3)，重構(gòu)數(shù)據(jù)的性能指標(biāo)CC/NRMSE/NSE高達(dá)0.97/0.05/0.94，而在信噪比和相關(guān)性低的內(nèi)陸河流域(圖1(i)和表3)，重構(gòu)數(shù)據(jù)的性能指標(biāo)CC/NRMSE/NSE值僅有0.77/0.12/0.60。因此在對(duì)比不同流域的重構(gòu)數(shù)據(jù)性能指標(biāo)時(shí)，應(yīng)綜合考慮各流域的GRACE/GRACE-FO 數(shù)據(jù)信噪比以及預(yù)測參數(shù)與GRACE/GRACE-FO TWSC的相關(guān)性。

表2 基于策略3重構(gòu)的中國九大流域TWSC性能

表3 預(yù)測參數(shù)與GRACE TWSC在中國九大流域的最大相關(guān)系數(shù)(時(shí)間窗口為前3～0個(gè)月)

4 結(jié) 語

本文利用MLR模型，采用3種不同的重構(gòu)策略，重構(gòu)中國九大流域的TWSC，得到以下結(jié)論：

1)基于整體信號(hào)重構(gòu)策略的重構(gòu)結(jié)果僅在TWSC主要受氣候變化影響的流域(如松花江遼河、黃河、淮河、東南諸河和珠江流域)表現(xiàn)出良好的性能，而在人類活動(dòng)或冰川融化(如海灤河、長江、西南諸河和內(nèi)陸河流域)頻繁的流域表現(xiàn)較差。

2)去趨勢項(xiàng)信號(hào)重構(gòu)策略的重構(gòu)結(jié)果在九大流域均優(yōu)于整體信號(hào)重構(gòu)策略，在人類活動(dòng)和冰川融化劇烈的流域優(yōu)勢更為明顯。

3)去趨勢項(xiàng)去季節(jié)項(xiàng)信號(hào)重構(gòu)策略的重構(gòu)結(jié)果在九大流域均稍優(yōu)于去趨勢項(xiàng)信號(hào)重構(gòu)策略。

4)不同流域重構(gòu)結(jié)果的性能還與該流域的GRACE/GRACE-FO數(shù)據(jù)信噪比以及預(yù)測參數(shù)(降水、氣溫、模型模擬的TWSC)與GRACE/GRACE-FO TWSC的相關(guān)性有關(guān)，信噪比和相關(guān)性高的流域(如西南諸河流域)其重構(gòu)結(jié)果的性能優(yōu)于信噪比和相關(guān)性低的流域(如內(nèi)陸河流域)。