黃河流域地下水儲量時空演變特征及影響因素分析

關鍵詞：ＧＲＡＣＥ；地下水儲量；耦合算法；氣候變化；植被；黃河流域

中圖分類號：Ｐ６４１．７；ＴＶ８８２．１ 文獻標志碼：Ａ ｄｏｉ：１０．３９６９／ ｊ．ｉｓｓｎ．１０００－１３７９．２０２４．０９．０１０

引用格式：劉家宏，周美林，邵薇薇．黃河流域地下水儲量時空演變特征及影響因素分析［Ｊ］．人民黃河，２０２４，４６（９）：６７－７３．

地下水是重要的淡水資源，對干旱和半干旱區(qū)工農(nóng)業(yè)生產(chǎn)、生活和生態(tài)系統(tǒng)的維持具有重要的支撐［１］作用。過去幾十年，受氣候變化和人類活動影響，地下水灌溉區(qū)域地下水位下降明顯［２］ 。氣候變化對地下水系統(tǒng)具有直接和間接影響，直接影響表現(xiàn)為改變地下水補給量或消耗量，間接影響為改變地下水利用方式。這些影響同時受土地利用等因素的顯著調(diào)控［３］ 。鑒于缺乏足夠的地下水位監(jiān)測井和對復雜水文地質(zhì)信息的全面了解，量化長時間尺度地下水變化及分析氣候和人類活動對其造成的影響具有挑戰(zhàn)性。

重力恢復及氣候?qū)嶒灒ǎ牵遥粒茫牛┬l(wèi)星于２００２ 年４月成功發(fā)射以來，有效監(jiān)測地球重力場變化。其繼任衛(wèi)星ＧＲＡＣＥ Ｆｏｌｌｏｗ－Ｏｎ 于２０１８ 年發(fā)射，為地表質(zhì)量變化持續(xù)監(jiān)測提供支撐。依據(jù)ＧＲＡＣＥ／ ＧＲＡＣＥ Ｆｏｌｌｏｗ－Ｏｎ監(jiān)測數(shù)據(jù)通過對地表質(zhì)量變化的反演進而得到陸地水儲量變化，是全球和區(qū)域尺度陸地水儲量變化監(jiān)測的重要手段［４－６］ 。陸地水儲量包含地表水、土壤水、地下水、冰川積雪等，基于水量平衡原理去除其他分量后，可得到地下水儲量的變化情況，該模式在地下水研究中得到廣泛應用［７］ 。Ｆｅｎｇ 等［８］ 根據(jù)２００３—２０１０ 年ＧＲＡＣＥ 衛(wèi)星監(jiān)測數(shù)據(jù)和地面測量數(shù)據(jù)，基于水量平衡法估算了華北地區(qū)地下水儲量變化，估算結(jié)果與同一時期監(jiān)測井站估計結(jié)果一致。Ｚｈｏｎｇ 等［９］ 利用水量平衡法估算了２００５—２０１１ 年遼河西部地下水儲量的變化，評估了該地區(qū)地下水枯竭狀況。Ｘｕ 等［１０］ 利用ＧＲＡＣＥ 和全球陸面數(shù)據(jù)同化系統(tǒng)（Ｇｌｏｂａｌ Ｌａｎｄ Ｄａｔａ Ａｓｓｉｍｉｌａｔｉｏｎ Ｓｙｓｔｅｍ，ＧＬＤＡＳ）數(shù)據(jù)評估了華北平原地下水長期變化趨勢，并分析了降水和南水北調(diào)工程的影響。

黃河流域位于干旱和半干旱氣候區(qū)，水資源供需矛盾十分突出。流域地表水資源開發(fā)利用率已達８０％，超過流域水資源承載能力［１１］ 。地下水資源在流域水資源安全和生態(tài)保護方面具有至關重要的作用。氣候變化背景下黃河流域氣溫顯著升高，降水量減小。作為生態(tài)系統(tǒng)敏感區(qū)域，加之退耕還林還草工程的實施，黃河流域植被覆蓋面積顯著增加［１２］ 。許多學者對黃河流域地下水變化進行了相關研究。Ｌｉｎ 等［１３］ 依據(jù)水量平衡原理，利用時間穩(wěn)定性分析識別了２００３—２０１６ 年黃河流域地下水波動的熱點和代表性監(jiān)測區(qū)域。王希棟等［１４］ 采用ＧＲＡＣＥ／ ＧＲＡＣＥ Ｆｏｌｌｏｗ－Ｏｎ 和ＧＬＤＡＳ 結(jié)合的方式評估了２００３—２０２１ 年黃河流域地下水儲量變化，定性分析了降水和植被對流域上游、中游和下游水儲量變化的貢獻。然而，現(xiàn)有研究對地下水時空變化特征討論較多，對地下水儲量空間變化主要驅(qū)動因素及驅(qū)動過程的認識存在不足。

考慮到ＧＲＡＣＥ 和ＧＲＡＣＥ Ｆｏｌｌｏｗ－Ｏｎ 衛(wèi)星間有近１ａ時間缺測，常用的線性和樣條插值在捕捉時空變化特征方面存在不足，本研究基于時間序列分解與隨機森林算法，構建像元尺度流域地下水儲量變化耦合模型，分析近２０ ａ 黃河流域地下水儲量時空演變規(guī)律，量化影響因素的貢獻，揭示氣候和人類活動對地下水變化的驅(qū)動過程，以期為黃河流域地下水資源利用和相關政策制定提供參考。

１流域概況及數(shù)據(jù)來源

１．１研究區(qū)域

黃河流域覆蓋青藏高原、黃土高原、內(nèi)蒙古高原和華北平原，地勢復雜，貫穿東、中、西三大地勢階梯，流域面積７９．５ 萬ｋｍ^２。流域多年平均氣溫為５．８ ℃、降水量為４６６．５ ｍｍ、潛在蒸散發(fā)量為６５３．６ ｍｍ［１５］ 。流域地下水分為松散巖類平原盆地孔隙水、松散巖類黃土高原孔隙水、碳酸鹽巖裂隙巖溶水、碎屑巖孔隙裂隙水、巖漿巖與變質(zhì)巖類裂隙水、凍結(jié)層水，對流域水資源循環(huán)演化有重要調(diào)節(jié)作用。流域地下水資源劃分為４ 個二級區(qū)、１１ 個三級區(qū)、３１ 個四級區(qū)、１０９ 個五級區(qū)［１６］ 。本文以黃河流域花園口水文站以上區(qū)域為研究對象（后文統(tǒng)稱黃河流域），分析地下水儲量變化，根據(jù)水沙特性和地形地質(zhì)，將黃河流域分為７ 個子流域，分別為龍羊峽以上、龍羊峽—蘭州、蘭州—頭道拐、頭道拐—龍門、龍門—三門峽、三門峽—花園口、內(nèi)流區(qū)。

１.２數(shù)據(jù)來源

１）ＧＲＡＣＥ 和ＧＲＡＣＥ Ｆｏｌｌｏｗ－Ｏｎ 監(jiān)測數(shù)據(jù)。采用得克薩斯州大學空間研究中心（Ｕｎｉｖｅｒｓｉｔｙ ｏｆ ＴｅｘａｓＣｅｎｔｅｒ ｆｏｒ Ｓｐａｃｅ Ｒｅｓｅａｒｃｈ，ＣＳＲ）、美國噴氣推進實驗室（Ｊｅｔ Ｐｒｏｐｕｌｓｉｏｎ Ｌａｂｏｒａｔｏｒｙ，ＪＰＬ） 發(fā)布的Ｍａｓｃｏｎ ＲＬ０６Ｖ２ 等效水高數(shù)據(jù)產(chǎn)品，時間序列為２００３ 年１ 月到２０２２ 年１２ 月，空間分辨率分別為０．２５°×０．２５°和０．５°×０．５°。不同數(shù)據(jù)產(chǎn)品對流域水儲量變化研究各有優(yōu)勢，相關研究表明，融合多個數(shù)據(jù)產(chǎn)品可有效降低水儲量變化監(jiān)測的不確定性［１７］ 。本研究以ＣＳＲ 和ＪＰＬ 數(shù)據(jù)的算術平均值作為地下水儲量分析的基礎，考慮到ＣＳＲ 和ＪＰＬ 數(shù)據(jù)在空間分辨率上的差異，為了更好地表征空間變異性，對ＪＰＬ 產(chǎn)品重采樣至分辨率為０．２５°后，再與ＣＳＲ 數(shù)據(jù)進行算術平均。對于傳感器故障導致的ＧＲＡＣＥ 數(shù)據(jù)序列缺失，本研究采用三次樣條法逐像元進行插值。

２） ＧＬＤＡＳ 數(shù)據(jù)。由美國航空航天局和國家海洋與大氣管理局合作開發(fā)，采用陸地表面模型和數(shù)據(jù)同化技術，結(jié)合衛(wèi)星和地面觀測數(shù)據(jù)，生成最佳陸面通量和狀態(tài)量［１８］ 。以ＧＬＤＡＳ Ｎｏａｈ－Ｖ２．１ 模型輸出的２ ｍ 深度土壤含水量為基礎，評估黃河流域的土壤水儲量變化。該數(shù)據(jù)的空間分辨率為０．２５°，時間分辨率為月。

３）氣象數(shù)據(jù)。采用國家青藏高原科學數(shù)據(jù)中心發(fā)布的中國１ ｋｍ 分辨率逐月降水量、逐月平均氣溫和潛在蒸散發(fā)量表征流域氣候特征。上述數(shù)據(jù)集的氣溫和降水來源于ＣＲＵ 發(fā)布的全球０．５°分辨率氣候數(shù)據(jù)集和ＷｏｒｌｄＣｌｉｍ 發(fā)布的全球高分辨率氣候數(shù)據(jù)集，通過Ｄｅｌｔａ 空間降尺度方案在中國地區(qū)進行降尺度，并采用中國４９６ 個獨立氣象觀測點的數(shù)據(jù)對降尺度結(jié)果進行驗證。潛在蒸散發(fā)量以最低氣溫、最高氣溫和月均氣溫為輸入，基于Ｈａｒｇｒｅａｖｅｓ 公式計算得到。

４）ＭＯＤＩＳ 數(shù)據(jù)。采用ＭＯＤ１３Ｃ１ 歸一化植被指數(shù)（ＮＤＶＩ）和ＭＯＤ１５Ａ２Ｈ 葉面積指數(shù)（ＬＡＩ）表征地表植被覆蓋和植被生長情況。數(shù)據(jù)的時間分辨率為８ ｄ，空間分辨率為５００ ｍ。以最大值法合成月尺度和年際尺度ＮＤＶＩ 和ＬＡＩ，并使用Ｓ－Ｇ 濾波算法對ＮＤＶＩ 和ＬＡＩ進行時間濾波重建，以減少噪聲并插補缺失數(shù)據(jù)。

５）ＧＬＥＡＭ 數(shù)據(jù)。ＧＬＥＡＭ（Ｇｌｏｂａｌ Ｌａｎｄ ＥｖａｐｏｒａｔｉｏｎＡｍｓｔｅｒｄａｍ Ｍｏｄｅｌ）是一個評估陸地蒸散發(fā)的模型，利用Ｐｒｉｅｓｔｌｅｙ－Ｔａｙｌｏｒ 公式和微波觀測數(shù)據(jù)來分析潛在和實際蒸散發(fā)。該模型結(jié)合了地表凈輻射、近地表氣溫、植被光學特征和土壤水分，基于Ｇａｓｈ 模型計算降水攔截損失，從而評估水體及冰雪覆蓋區(qū)的蒸發(fā)量［１９］ 。本研究采用ＧＬＥＡＭ ３．７ ａ數(shù)據(jù)集，該數(shù)據(jù)以遙感和再分析數(shù)據(jù)為驅(qū)動，具有月尺度時間分辨率和０．２５°×０．２５°的空間分辨率。

２研究方法

陸地水儲量包含土壤水、地下水、雪水當量、冰川積雪和冠層含水等。冰川積雪及冠層含水對黃河流域年際水儲量變化影響較小， 本研究不予考慮。以ＧＲＡＣＥ 衛(wèi)星陸地水儲量異常與ＧＬＤＡＳ 模擬土壤水儲量異常的差異作為區(qū)域地下水儲量異常。由于ＧＬＤＡＳ模擬輸出土壤水為當前土壤水的狀態(tài)值，因此為了與水儲量異常保持一致，以研究時段土壤水儲量相對于２００４—２００９年多年平均土壤水儲量的差值作為土壤水儲量的異常。

ＧＲＡＣＥ 和ＧＲＡＣＥ Ｆｏｌｌｏｗ－Ｏｎ 數(shù)據(jù)序列間存在１１個月（２０１７ 年７ 月—２０１８ 年５ 月）的數(shù)據(jù)缺失，缺失數(shù)據(jù)序列較長，傳統(tǒng)的插值方法在人類活動影響顯著的干旱半干旱區(qū)域插值精度不足［２０］ 。隨機森林算法是常用的高效機器學習算法，在地球科學、遙感領域得到廣泛應用，尤其在全球和區(qū)域尺度應用效果較好［２１］ 。同時，隨機森林算法的樹結(jié)構模型通過量化變量對樹節(jié)點的貢獻評估變量重要性，為分析影響因素對地下水儲量變化的貢獻提供了新的視角。

ＳＨＡＰ（ＳＨａｐｌｅｙ Ａｄｄｉｔｉｖｅ ｅｘＰｌａｎａｔｉｏｎｓ）模型是一種解釋機器學習模型輸出的方法，該模型借鑒博弈論中的Ｓｈａｐｌｅｙ 值概念，可量化每個特征對模型預測結(jié)果的貢獻，關于該模型的詳細介紹見參考文獻［２３］。為了提高模型的可解釋性，采用ＳＨＡＰ 方法評估不同時間尺度隨機森林模型各自變量對因變量的貢獻，并將不同時間尺度貢獻之和作為變量最終貢獻。考慮到各因素對地下水儲量變化影響存在明顯的交互作用［２４］ ，模型構建時對變量進行了歸一化處理，并假定貢獻率大于０．４ 的因素為地下水儲量變化的主導因素。

３結(jié)果與分析

３．１地下水儲量時空變化特征

圖１ 為２００３—２０２２年黃河流域中上游地區(qū)地下水儲量異常值年際變化情況，近２０ ａ 地下水儲量持續(xù)顯著減小，變率為－３．３ｍｍ／ａ。各子流域地下水儲量均呈顯著減小趨勢，變率為－０．９～－５．７ ｍｍ／ａ，其中三門峽—花園口地下水儲量減小幅度最大，龍羊峽—蘭州地下水儲量減小幅度較?。ㄒ妶D２、表１）。２００８年之后，蘭州以下相鄰子流域地下水儲量差異明顯增大，地下水儲量空間變異性顯著增強。三門峽—花園口地下水儲量變化在２００４ 年、２０１２年、２０１８年和２０２２年均出現(xiàn)極大值，其他時期地下水儲量以持續(xù)減小為主，２０２０年地下水儲量最小。內(nèi)流區(qū)和龍門—三門峽，龍羊峽以上和龍羊峽—蘭州，地下水儲量時間變化特征基本一致。

黃河流域地下水儲量年際變化空間格局如圖３所示，從東到西地下水儲量年際變化幅度表現(xiàn)為先減小后增大的格局，流域東部地下水減小速度最快。

具體而言，黃土高原中心區(qū)域地下水減小幅度較大，尤其在汾河盆地和關中平原，地下水儲量減小速度為－１５～－３０ ｍｍ／ａ。流域南部邊緣（甘肅境內(nèi)）和寧夏部分區(qū)域地下水儲量呈增加趨勢，增加速度為０ ～６ ｍｍ／ａ。對年際變率進行顯著性分析，結(jié)果表明趨勢擬合決定系數(shù)（Ｒ^２）大于０．５ 的網(wǎng)格占比約為７２％，主要分布在蘭州以下區(qū)域、黃河源區(qū)中部及扎陵湖和鄂陵湖區(qū)域。

３.２地下水儲量影響因素分析

圖４為黃河流域降水年際時空變化情況（圖中①為頭道拐—龍門、②為內(nèi)流區(qū)、③為龍門—三門峽、④為龍羊峽—蘭州、⑤為蘭州—頭道拐、⑥為三門峽—花園口、⑦為龍羊峽以上，下同），２００３—２０２２年黃河流域年降水量無顯著趨勢，其中２０１５—２０２２年波動較大。相較于２０１５年，２０１９年流域平均年降水量增加了１２２ｍｍ，到２０２２年減小至接近２０１５年值。

相較于降水的變化，地下水儲量僅在２０１７年略微增大，其他時期均持續(xù)減小。這一現(xiàn)象說明地下水對降水的響應具有明顯的滯后性，并且降水增加對地下水儲量增加的補給效應不顯著。黃河流域降水具有明顯的空間變異性，其中龍門—三門峽、三門峽—花園口２ 個子流域降水年際波動較為明顯，２００４ 年、２０１１年、２０１６ 年和２０２１ 年為降水顯著增加年份，除２００４ 年外，地下水儲量均在后一年表現(xiàn)為增加，說明地下水儲量變化對降水變化的響應具有時間尺度效應。

氣溫通過蒸發(fā)影響地下水補給，潛在蒸散發(fā)表征區(qū)域最大蒸散發(fā)能力。當氣溫升高時，潛在蒸散發(fā)量的增加意味著更多的水分可能從地表和土壤中蒸發(fā)，減少了能夠滲透到地下的水量，從而影響地下水的自然補給。黃河流域氣溫和潛在蒸散發(fā)量年際變化如圖５ 和圖６ 所示。

研究時段氣溫具有明顯的升高趨勢，在２００６ 年、２０１６ 年和２０２２ 年流域平均和子流域氣溫達到時段極值。潛在蒸散發(fā)量無顯著變化趨勢，波動明顯，與氣溫一致性較好。氣溫升高促進土壤和植被蒸發(fā)，因而相應時期的地下水儲量表現(xiàn)為減小趨勢，地下水儲量對氣溫升高的響應速度大于對降水增加的。

黃河流域陸地蒸散發(fā)量、歸一化植被指數(shù)（ＮＤＶＩ）、葉面積指數(shù)（ＬＡＩ）年際時空變化如圖７～圖９所示。自１９９９ 年退耕還林工程實施后，黃河流域植被綠度顯著提高［２５］ 。本研究也表明２００３—２０２２ 年黃河流域平均ＮＤＶＩ 和ＬＡＩ 均呈顯著增大趨勢，植被覆蓋狀況和植被生長狀況持續(xù)向好。黃河流域陸地蒸散發(fā)量也表現(xiàn)出明顯的增加趨勢。空間區(qū)域上，龍門—花園口區(qū)域即黃土高原退耕還林核心區(qū)ＮＤＶＩ 和ＬＡＩ 顯著增大。相應的龍門到花園口區(qū)域陸地蒸散發(fā)量也表現(xiàn)出明顯的增加趨勢。陸地蒸散發(fā)量變化特征與前述氣溫和降水變化基本一致，２００４ 年、２００８ 年、２０１２ 年和２０２２年蒸散發(fā)量顯著減小，對應年份降水量減小、氣溫下降，植被對地下水儲量的影響同時受氣候調(diào)節(jié)。

３．３影響因素貢獻量化分析

圖１０ 為各子流域降水、氣溫、ＮＤＶＩ、ＬＡＩ、蒸散發(fā)量、潛在蒸散發(fā)量的ＳＨＡＰ 貢獻值。整體而言，各因素對黃河流域地下水儲量變化的貢獻從大到小依次為ＮＤＶＩ（０．８０）、氣溫（０．５５）、蒸散發(fā)量（０．５１）、潛在蒸散發(fā)量（０．５１）、ＬＡＩ（０．５０）和降水（０．４３），表明黃河流域地下水儲量減小主要由ＮＤＶＩ 增大、氣溫升高和蒸散發(fā)量增加等共同驅(qū)動。黃河流域不同區(qū)域地下水儲量變化主要驅(qū)動因素表現(xiàn)出顯著的空間差異性。龍羊峽以上地下水儲量變化主要驅(qū)動因素為陸地蒸散發(fā)量，說明該子流域地下水儲量變化主要依賴大氣條件和植被的相互作用。龍羊峽—蘭州地下水儲量變化主要受植被狀況和覆蓋（ＮＤＶＩ 和ＬＡＩ）、氣溫、降水、蒸散發(fā)量和潛在蒸散發(fā)量共同影響，這些因素組合反映了植被狀況、氣候變化對水資源的綜合影響。蘭州—頭道拐地下水儲量變化主導因素為ＮＤＶＩ，植被覆蓋變化是該區(qū)域地下水變化的關鍵。頭道拐—三門峽地下水儲量變化主要由ＮＤＶＩ、ＬＡＩ 和氣溫等因素共同驅(qū)動，表明該子流域植被覆蓋、植被生長狀況和氣候條件對地下水儲量變化影響顯著。內(nèi)流區(qū)地下水儲量變化主要由ＮＤＶＩ、潛在蒸散發(fā)量、氣溫和陸地蒸散發(fā)量驅(qū)動，這些因素反映了復雜地表－大氣相互作用、植被狀況對地下水的影響。三門峽—花園口ＮＤＶＩ 是地下水儲量變化主要驅(qū)動因素，揭示了植被綠度提高對該區(qū)域地下水資源具有不可忽視的影響。

３．４討論與建議

地下水儲量變化受地下水補給量和消耗量共同影響，研究發(fā)現(xiàn)近２０ ａ 黃河流域氣候呈暖干化趨勢，降水無明顯增加，氣溫上升顯著。加之陸地表面蒸發(fā)量顯著增加，地下水消耗量大于補給量，導致地下水儲量顯著減小。雖然２０１５ 年后黃河流域豐水年頻次增加，但是地表植被覆蓋面積增加、氣溫上升導致水汽壓差增大、蒸散發(fā)量增大，加之地下水對降水響應的滯后性，地下水儲量不能得到及時恢復，而存儲的土壤水在下一個干旱年被消耗。因此，黃河流域地下水儲量減小是氣候變暖和植被變綠共同作用的結(jié)果。在流域不同區(qū)域植被對地下水儲量變化的貢獻均較為顯著，其他因素貢獻不同，與植被的作用方式也存在差異，應對黃河流域地下水儲量減小應該采取與地域特點相一致的策略。

４結(jié)束語

分析了近２０ ａ 黃河流域地下水儲量時空演變規(guī)律，量化了各因素的貢獻，揭示了氣候和人類活動對地下水變化的驅(qū)動過程，主要結(jié)論如下。

１）地下水儲量時空變化特征：時間上，２００３—２０２２年黃河流域地下水儲量顯著減小，變率為－３．３ｍｍ／ａ，其中頭道拐—龍門和三門峽—花園口變率為－５．７ ｍｍ／ ａ 和－５．５ ｍｍ／ａ。空間上，地下水儲量變幅表現(xiàn)為自東向西先減小后增大的格局，在汾河盆地和關中平原減小幅度較大（－１５～－３０ ｍｍ／ ａ）。

２）因素變化特征：黃河流域氣候呈暖干化趨勢。２００３—２０２２年降水變化趨勢不顯著，降水量明顯波動。２０１５年后豐水年頻次增加，氣溫呈現(xiàn)出明顯的上升趨勢，且２０１５年之后平均氣溫顯著高于之前。相比之下，潛在蒸散發(fā)量未表現(xiàn)出一致的變化特征。地表植被覆蓋和植被生長狀況持續(xù)向好，地表蒸散發(fā)量顯著增大。

３）變化驅(qū)動機制：黃河流域地下水儲量的減小是植被增加、氣溫升高和蒸散發(fā)量變化共同驅(qū)動的結(jié)果。地下水儲量空間變化主導因素差異顯著，龍羊峽以上由陸地蒸散發(fā)過程驅(qū)動，龍羊峽—蘭州由植被覆蓋、氣溫、降水和蒸散發(fā)共同驅(qū)動，蘭州—頭道拐由植被覆蓋變化驅(qū)動，頭道拐—三門峽由植被覆蓋變化和氣溫驅(qū)動，內(nèi)流區(qū)由植被和氣候變化驅(qū)動，三門峽—花園口由ＮＤＶＩ 驅(qū)動。

本研究未考慮地下水抽取及原煤開采引起的地下水排泄等人類活動對地下水儲量變化的影響。同時，黃河流域植被變綠與生態(tài)恢復工程密切相關，需構建合適的人類活動量化指標，綜合考慮人類生產(chǎn)生活用水、提取地下水及煤炭開采和生態(tài)恢復工程，分離人類活動的正負效應。另外，本研究在子流域劃分時以地表水為依據(jù)，對水文地質(zhì)條件考慮不足，后續(xù)研究應以黃河流域地下水資源分區(qū)為標準，考慮地質(zhì)水文參數(shù)和土壤入滲率，從降水－土壤水－地下水三水轉(zhuǎn)換關系及機理層面揭示地下水儲量變化的機制。