基于GRACE 衛星的黃淮海地區地下水儲量變化與歸因分析

王 坤，李玉龍，李琛曦，賴燕玲，楊傳國,3

(1.河海大學 水文水資源學院, 江蘇 南京 210098； 2.河北省地質資源環境監測與保護重點實驗室, 河北 石家莊050021； 3.河海大學 水災害防御全國重點實驗室, 江蘇 南京 210098)

地下水是水資源的重要組成部分，由于具有分布廣、水質好、可持續開發利用等許多優點，已成為最重要的和最佳的供水選擇[1]，特別是在水資源相對匱乏的北方地區。黃淮海地區是中國重要的糧食生產基地[2]，是人均水資源量較少的地區之一，也是人類活動對水資源影響最嚴重的區域[3]。大規模開發地下水資源導致地下水位下降嚴重，產生了水質惡化、地面沉降、海水入侵等一系列水安全問題[4]，已成為制約地區發展的主要因素。全面監測地下水并掌握其變化規律是當務之急。傳統人工監測地下水的方法耗費人力物力，且由于監測站點布設的困難，難以實現地下水的大面積整體監測。隨著重力衛星探測技術的發展，通過測量地球重力場的變化反演全球陸地水儲量變化，為實現地下水的區域整體監測提供了新途徑[5]。國內外諸多學者利用GRACE （Gravity Recovery and Climate ）重力衛星在不同區域開展了地下水研究，證實了在大尺度范圍內反演地下水資源量變化的可靠性[6]。Famiglietti 等[7]和Scanlon 等[8]利用GRACE 數據反演加利福尼亞州中央山谷地下水儲量變化，與研究區實測地下水位驗證結果大體一致。國內學者驗證了GRACE 與GLDAS （Global Land Data Assimilation System）數據反演黃河流域[9]、海河流域[10]、三江源[11]等流域地下水儲量變化具有較好適用性。

黃淮海地區地下水儲量變化引起越來越多的研究者關注，對其影響的主導因素開展了相關研究。Guo 等[12]研究發現人類活動在海河流域地下水儲量變化中占主導作用，農業灌溉是地下水儲量消耗的主要因素，而降水和蒸發量相對穩定；林雪敏等[13]研究表明農業用水量對京津冀地區地下水儲量變化的影響程度明顯大于降水量與蒸散量；秦歡歡等[14]利用模型設定不同情景模擬華北平原水資源利用狀況，結果表明農業節水措施和南水北調工程可以有效恢復地下水儲量。南水北調通水后黃淮海地區供水結構發生了改變，相關研究較少，亟需采用多源數據深入分析區域地下水儲量變化的特征與成因。因此，本文利用最新的GRACE 數據和GLDAS 陸面同化數據，定量計算地下水儲量變化，評估黃淮海地區地下水儲量變化的時空特征，結合降水、蒸散發、地下水供水、跨流域調水等自然和人為因素，探討影響地下水儲量變化的原因，從而為科學評估與合理利用地下水提供依據。

1 研究區概況

本文選擇黃淮海地區，包括海河流域、淮河流域及黃河流域花園口以下地區為研究區域（圖1）。地理位置位于北緯31°～43°，東經112°～122°，北起燕山南麓，南抵桐柏山、大別山，西起太行山，東臨渤海、黃海，總面積為67.8 萬km2。黃淮海地區屬溫帶季風氣候，季節變化明顯，冬季干燥寒冷，夏季高溫多雨[15]；年平均氣溫為13.1 ℃；年平均降水量為500～800 mm；人均水資源量為462 m3，為全國平均水平的五分之一，水資源承載能力與經濟社會發展矛盾突出，是南水北調東線和中線工程的受水區[16]。黃淮海地區土層深厚，土質肥沃，農業規模大，主要種植小麥和玉米。地下水開采量的持續增加，導致一些地區的地下水位下降嚴重。

圖1 研究區位置Fig.1 The location of the study area

2 數據與方法

2.1 數據來源

（1）GRACE 重力衛星數據。GRACE 計劃是美國國家航空航天局和德國航空中心的合作項目，通過衛星監測地球重力場，估算陸地水儲量的變化。本文選用美國德克薩斯大學研究中心（CSR）最新發布的RL06 Mascon 版本（https://www2.csr.utexas.edu/grace/）提供的重力衛星陸地水儲量變化數據。相比于傳統的球諧系數，Mascon 產品的最大特點是無需做任何后處理，可以直接提供格網產品。該數據以2004—2009 年的陸地水儲量等效水深為基準，空間分辨率為0.25°×0.25°，時間分辨率為1 個月，時間序列從2003 年2 月至2021 年12 月，個別缺測時段采用奇異譜分析法[17]插補缺失值。

（2）GLDAS 陸面同化數據。全球陸面數據同化系統（GLDAS）包含4 個陸面模型（Noah、CLM、Mosaic、VIC），本文采用的是GLDAS-2.1 Noah 模型陸面同化數據，可解決內陸水域無效數據的問題。從中選取總冠層水量（Canoplnt）、雪水當量（SWE）、土壤含水量（SMS）和蒸散發（Evap）數據，其中SMS 為4 層土壤含水量數據之和（0～10 cm、10～40 cm、40～100 cm、100～200 cm）。時間長度及空間分辨率與GRACE 數據一致。

（2）降水和供水量數據。數據來源于水利部和流域部委水資源公報，時間序列為2003—2021 年。選取各流域降水、地表水供水量、地下水供水量、跨流域調水和其他水源等數據，其他水源包括污水處理回水、集雨工程供水、海水淡化供水等。

2.2 研究方法

2.2.1 地下水儲量變化估算 陸地水儲量包括地表水、地下水、土壤水、生物水、冰雪、植被冠層水等，其中地表水與生物水相對于其他成分，其多年變化量可以忽略不計[18-19]。因此，陸地水儲量變化可表示為：

由此可得到地下水儲量變化量的計算：

式中：ΔGWS為地下水儲量變化（cm）；ΔTWS為陸地水儲量變化（cm）；ΔSMS為土壤含水量變化（cm）；ΔSWE為雪水當量變化（cm）；ΔCanoplnt為冠層水量變化（cm）。陸地水儲量變化由GRACE 重力衛星反演得到，土壤含水量、雪水當量和植被冠層水的變化由GLDAS 數據計算。

2.2.2 Theil-Sen Median 和Mann-Kendall 法 Theil-Sen Median 方法又稱Sen 斜率估計，是一種穩健的非參數統計的趨勢計算方法[20]，用于估計地下水儲量的變化速率 β，計算公式為：

式中：x為 各時段的ΔGWS；i和j為時段月份數；Median 表示取中值。若β大于0，則表示地下水儲量為增長趨勢，反之則為下降趨勢。

Mann-Kendall 檢驗是一種非參數的時間序列趨勢性檢驗方法，廣泛應用于水文氣象數據分析，定義檢驗統計量S：

使用檢驗統計量Z進行趨勢檢驗，Z值計算方法為：

式中：Svar為統計量S的方差。取顯著性水平α=0.05，當|Z|≥1.96 時，表示該序列通過95%置信度檢驗，地下水儲量變化趨勢顯著，反之變化趨勢不顯著。

2.2.3 隨機森林模型 隨機森林（Random Forest）是一種基于統計學習理論的組合分類智能算法?；締卧菦Q策樹，使用自抽樣（Bootstrap）法從原始數據集中隨機重復采樣，得到多個決策樹模型，并將所有結果的平均值作為最終回歸結果。根據基尼指數判斷多個自變量對因變量的所做貢獻，可以有效進行特征重要性排序，隨機森林算法具有訓練速度快、精度高、不易過擬合等優點，在水文地質學領域應用越來越廣泛[21]。本文采用python 語言中sklearn 庫建立隨機森林模型，以地下水儲量變化作為因變量，將影響地下水儲量變化的因素作為自變量，計算影響因素的貢獻率。

基尼指數代表了模型的不純度，基尼指數越小，則不純度越低，特征越好。其計算公式為：

式中：K為數據集分類總數；pk為某個子集屬于第k類的概率。

特征A在內部某個節點m的重要性Vim,m為分枝前后基尼指數的變化量：

式中：Gini,m為節點m處的基尼指數，Gini,l、Gini,r為分枝后出現新節點處的基尼指數；l、r為節點m處分枝后的新節點。

則特征A在隨機森林模型中基尼指數評分為：

式中：Vim,i為第i棵決策樹各節點處分支前后基尼指數變化量，n為節點數。

3 結果與分析

3.1 時間變化特征

根據GRACE 和GLDAS 數據，采用式（2）計算研究區各網格地下水儲量的逐月變化。整個黃淮海地區相對于基準期，地下水儲量在2003—2021 年總體波動范圍為?24.94～5.99 cm（見圖2），整體呈減少趨勢，近20 年內共減少23.18 cm，平均每年減少1.22 cm。其中在2003—2010 年平均減少0.45 cm/a，在2011—2019 年間減少速度加快為1.68 cm/a，其原因可能與經濟快速發展及人口增長有關，在2020 年以后以1.05 cm/a 的速度開始增加。

圖2 黃淮海地區地下水儲量時間變化趨勢Fig.2 Temporal fluctuations in groundwater storage within the Huang-Huai-Hai region

逐一分析比較黃淮海地區內各二級子流域地下水儲量變化（見圖3），具體變化速度見表1。從表1 和圖3 可見研究期內各子流域的地下水儲量均處于持續減少狀態，海河流域平均下降速率為1.75 cm/a，淮河流域平均下降速率為0.57 cm/a。在2003—2010 年間，地下水儲量下降最快的是海河南系，達1.94 cm/a，而淮河流域整體處于緩慢上升階段。2011—2019 年研究區地下水儲量虧損程度擴大，虧損地區增多，除淮河下游變化不明顯，各子流域呈現不同程度的下降趨勢，花園口以下地區以2.92 cm/a 的速率下降，海河南系下降速率為2.86 cm/a。2020 年后，部分子流域地下水儲量開始上升，花園口以下地區上升速率為7.22 cm/a。

表1 黃淮海地區子流域地下水儲量變化速度Tab.1 The rates of change in groundwater storage within the Huang-Huai-Hai region

圖3 黃淮海地區各子流域地下水儲量變化趨勢Fig.3 Temporal fluctuations in groundwater storage within the sub-basins of the Huang-Huai-Hai region

各子流域地下水儲量變化在年內變化表現出差異性（圖4）。灤河及冀東沿海、海河北系地區地下水儲量逐月變化較為穩定，山東半島沿海諸島、淮河下游地區地下水儲量在4—8 月份處于高位。海河南系、徒駭馬頰河、花園口以下、淮河上游等主要農作物耕種地區的地下水儲量下降明顯，且在6—7 月份明顯偏少。黃淮海地區水儲量持續減少主要是因農業灌溉進行的大規模開采地下水[22]，用于農作物灌溉的地下水供水量占地下水開采總量的80%。黃淮海地區的農作物以冬小麥和夏玉米為主，5—8 月為作物需水關鍵期，需水量占整個生長發育過程的79%[23]，地下水下降最快，地下水儲量處于最低狀態。冬小麥開始大量播種時間為10 月中旬，越冬期不需要灌溉，降水基本能滿足需要，此時地下水儲量得到恢復。提高農業用水效率和減少地下水開采量是緩解區域地下水資源短缺的重要措施。

圖4 黃淮海地區地下水儲量逐月變化Fig.4 The monthly average of groundwater storage in the Huang-Huai-Hai region

3.2 空間變化趨勢

從空間分布看，研究區地下水儲量變化存在明顯的空間差異性，整體上由東至西呈現下降幅度增大趨勢，變化范圍為?21.9～2.27 cm，負值為相對于基準期處于減少狀態。圖5 表明地下水儲量減幅最大的是海河南系，主要集中在邯鄲、安陽、鶴壁，淮河流域上游較為盈余；最低值出現在鶴壁市，最高值出現在淮安市。位于海河南系的白洋淀流域是華北最大的淡水湖，已成為雄安新區發展的重要生態水體，該地區地下水儲量最低值為?10.75 cm，最高值為?4.07 cm，虧損較為嚴重。

圖5 黃淮海地區多年平均地下水儲量變化分布Fig.5 The spatial pattern of groundwater storage within the Huang-Huai-Hai region

基于Sen 斜率估計和Mann-Kendall 法分析黃淮海地區2003—2021 年地下水儲量變化趨勢，計算出每個格網的β斜率和顯著性Z值，研究區88.87%區域通過置信水平0.05 的檢驗，其中94.2%的區域呈現出顯著的負向變化，即地下水儲量減少。為了進一步比較各個地區地下水儲量的變化情況，將研究區域劃分為極顯著增加、顯著增加、微顯著增加、不顯著增加、不顯著減少、微顯著減少、顯著減少、極顯著減少等8 個等級（圖6）。從圖6 可以看出：黃淮海地區地下水儲量減少速率整體表現出自南向北逐步增大的特點，海河南系大部分區域處于極顯著減少狀態，水儲量增加的地區位于淮河流域上游及中游部分區域。

圖6 Sen+M-K 計算地下水儲量變化趨勢分布Fig.6 The trend of variations in groundwater storage calculated using Sen+M-K methods

3.3 地下水儲量變化成因

3.3.1 降水、蒸散發變化影響 降水入滲補給是影響地下水儲量變化最直接的因子?；春恿饔蚪邓S沛，年均降水量為873.3 mm，充足的降水可促進地下水補給，有利于恢復地下水儲量。海河流域年均降水量為532.7 mm，降雨徑流經水庫截留等損失后流域下滲較少，對地下水儲量變化影響的并不明顯（圖7）。同時，由于海河流域農業部門耗水最多，降水量變化會影響灌溉量，枯水年灌溉量增多，地下水開采量會隨之增多。蒸散發過程會消耗土壤水分，減少地下水的補給。黃淮海地區多年平均蒸散量表現為增加趨勢，淮河流域的蒸散量由2003—2015 年的722.4 mm 增加至2016—2021 年的782.7 mm，海河流域的蒸散量由533.9 mm 增加至597.3 mm，花園口以下的蒸散量由629.7 mm 增加至678.0 mm。

圖7 降水、蒸散發與地下水儲量變化關系Fig.7 The interrelationship among precipitation, evapotranspiration, and groundwater storage

降水在黃淮海地區時間空間上的分布差異明顯?；春恿饔蚪邓觌H變化更大，直接影響地下水儲量變化速度，地下水儲量變化與降水變化趨勢具有較明顯的一致性。而海河流域降水對地下水儲量變化的影響不顯著，整個海河流域年降水量與地下水儲量之間的相關系數僅為0.02，灤河及冀東沿海、海河北系甚至出現負相關，這表明影響海河流域地下水儲量變化的主要因素不是降水量，而是人類開采的地下水量，這與之前學者研究結果一致[24]。

3.3.2 供水量變化影響 流域供水來源包括當地地表水、地下水、跨流域調水等，其中淮河流域外調水為引長江水與黃河水，海河流域外調水在2015 年之前為引黃河水，南水北調引長江水主要從2015 年開始?；春恿饔虻叵滤┧壤鹉昃徛郎p少，從2003 年占總供水量的32.3%降至2021 年的22%，主要原因是跨流域調水量呈增加趨勢；地表水多年平均供水比例為55.8%，是其最主要的供水來源，地下水次之，平均占比26.9%（圖8）。黃河花園口以下地區地下水供水比例呈減小趨勢，但總供水量呈增加趨勢。

圖8 黃淮海地區年供水量變化Fig.8 The annual fluctuations in water supply within the Huang-Huai-Hai region

海河流域地下水是主要的供水來源，2015 年以前地下水供水占比64.2%，2015 年南水北調工程通水后開始顯著減少，2019—2021 年已降到39.1%，但仍是海河流域最主要的供水來源。南水北調工程通水后，海河流域供水結構發生了改變，地下水供水量明顯減少，但由于地下水多年持續超采嚴重，目前仍維持在近150 億m3/a，地下水儲量減少速率有所減緩，但下降趨勢未得到根本遏制[10]。

3.3.3 成因分析 以降水量、蒸散量、地下水供水量和跨流域調水量為自變量建立隨機森林模型，將計算所得影響因素重要性評分進行歸一化處理，定量分析比較不同因子對黃淮海地區地下水儲量變化的相對影響程度（圖9）。結果表明，淮河流域地下水儲量變化的影響因子相對貢獻率為降水＞蒸散發＞地下水供水＞跨流域調水，而海河流域影響因子的相對貢獻率為地下水供水＞跨流域調水＞降水＞蒸散發，黃河流域花園口以下的相對貢獻率為地下水供水＞蒸散發＞降水。各個流域地下水儲量變化的主導因素不同，降水是淮河流域的主導因素，貢獻率28.3%；地下水供水是海河流域和黃河流域花園口以下地區的主導因素，貢獻率分別達52.6%和58.0%。

圖9 黃淮海地區地下水儲量變化影響因素的相對貢獻度Fig.9 The relative contributions of factors influencing variations in groundwater storage within the Huang-Huai-Hai region

降水和蒸發等自然因素對淮河流域地下水儲量變化的影響程度顯著大于海河流域，而地下水供水是影響海河流域和黃河流域花園口以下地區地下水儲量下降的主導因素。研究表明人為因素對黃淮海地區地下水儲量變化影響逐漸增大[25-26]，2003—2021 年間，地下水消耗遠遠超出了地下水系統的自然補給能力，導致整個海河流域和淮河流域北部地下水儲量大幅下降。海河流域跨流域調水的相對貢獻率大于淮河流域，隨著地下水超采治理及南水北調供水工程持續推進，地下水儲量下降速率近年來有所減緩，北京、天津等受水區城市地下水位上升明顯[27-28]，但農業區地下水位仍在下降。

4 結 語

（1）2003—2021 年間，黃淮海地區地下水儲量呈波動減少趨勢，平均變化速率為?1.22 cm/a，其中海河流域為?1.75 cm/a，淮河流域為?0.57 cm/a；海河南系虧損程度最嚴重，淮河上游和下游變化不顯著。2010 年代的下降速率高于2000 年代，隨著南水北調工程調水量的增加，地下水儲量下降速率有所減緩，沿線部分城市地下水儲量有所回升。

（2）地下水儲量變化存在明顯的空間差異性，海河流域虧損程度大于淮河流域。黃淮海地區94.2%的區域地下水儲量呈下降狀態，下降速率整體表現出自南向北逐步增加的特征；海河南系為極顯著減少趨勢，淮河流域上游及中游部分地區呈增加趨勢。

（3）黃淮海地區地下水儲量變化的主導因素不同，降水對淮河流域地下水儲量變化具有主導作用，而影響海河流域和黃河流域花園口以下地區最主要的因素是地下水供水量，多年持續超采造成地下水儲量大幅減少。南水北調工程通水后海河流域地下水開采量減少，改善了區域供水結構，地下水儲量下降得到緩解，但還需進一步提高農業用水效率，加強地下水資源管理。