基于多源數(shù)據(jù)的華北平原地下水儲量時空變化分析

張明剛 魏長壽

1 山東科技大學測繪與空間信息學院，青島市前灣港路579號，266590 2 內蒙古科技大學礦業(yè)與煤炭學院，內蒙古自治區(qū)包頭市阿爾丁大街7號，014010

近20 a來，華北平原農(nóng)業(yè)迅猛發(fā)展，地下水超量開采，造成地下水資源稀缺。因此，獲取華北平原地區(qū)的地下水儲量變化，探究影響其變化的因素對當?shù)剞r(nóng)業(yè)的良好運營非常重要。

對華北平原地下水儲量進行監(jiān)測的研究較多[1-3]，但一般側重于獲取地下水水儲量的變化趨勢，減少后處理誤差，很少關注農(nóng)業(yè)用水量變化、降水量變化的影響，缺乏對地下水減少影響因素的具體分析。本文使用GRACE衛(wèi)星數(shù)據(jù)質量較好的2003～2015年的數(shù)據(jù)反演華北平原地下水儲量的時空變化，并利用經(jīng)驗正交函數(shù)(empirical orthogonal function,EOF)具體分析影響華北平原地下水儲量時空變化的主要因素。

1 研究區(qū)概況

華北平原位于31°36′～40°29′N、112°13′～120°53′E，屬半干旱、半濕潤季風性氣候，春季干旱少雨、蒸發(fā)強，夏季高溫多雨，冬季寒冷干燥，地區(qū)年均降水量在空間上分布不均。冬小麥和夏玉米的一年兩熟制是該地區(qū)的典型農(nóng)業(yè)種植模式[4]。其中，玉米的灌漿期在7月份左右，此時華北平原氣溫較高，雨水供應充足，地下水開采對地下水儲量的變化影響較??；華北平原冬小麥的產(chǎn)量占全國的60%～80%，高密度的冬小麥種植需要消耗大量的水資源，且其灌漿期為3～5月，恰逢干旱少雨的春季，導致水資源嚴重短缺，此時超過60%的用水來自地下水[5]。

圖1 研究區(qū)位置Fig.1 Location of the study area

2 數(shù)據(jù)與方法

2.1 監(jiān)測井數(shù)據(jù)

從中國地質環(huán)境監(jiān)測研究院獲取2003～2013年華北平原地區(qū)61口監(jiān)測井水位數(shù)據(jù)。為了將水位數(shù)據(jù)轉化為與GRACE等效的水儲量變化，參考文獻[6]，取華北平原的平均給水度為0.06，將水位數(shù)據(jù)乘以區(qū)域給水度求得水儲量變化的等效水柱高。

2.2 GRACE重力衛(wèi)星數(shù)據(jù)獲取陸地水儲量變化

使用CSR提供的GRACE level-2 RL06數(shù)據(jù)產(chǎn)品，時間跨度為2003-01～2015-12，分辨率為1個月，使用三次樣條差值方法補齊缺失數(shù)據(jù)[7]。對GRACE數(shù)據(jù)進行如下預處理：1)用SLR對應項替代原始的1階項和C20項；2)對生成的全球格網(wǎng)數(shù)據(jù)進行去平均處理；3)采用P4M15去相關濾波去除“南-北條帶誤差”；4)采用300 km的高斯濾波對去條帶之后的球諧系數(shù)進行平滑，降低高階項球諧系數(shù)的噪聲；5)采用正向建模法改正信號的泄露誤差[8]。最后反演獲得陸地水儲量(terrestrial water storage，TWS)變化。

2.3 GLDAS水文模型數(shù)據(jù)

選取2003～2015年共156個月GLDAS Noah陸地模型輸出的0～200 cm土壤含水量數(shù)據(jù)和雪水當量月解數(shù)據(jù)來估計地表水儲量的變化，其空間分辨率為1°×1°，并采用與GRACE數(shù)據(jù)相同的處理方式獲取土壤水儲量(soil moisture storage，SMS)變化與雪水當量(snow water equivation,SWE)變化。

2.4 地下水儲量變化ΔGWS的獲取

簡化后的區(qū)域水量平衡方程為[9]：

ΔGWS=ΔTWS-ΔSMS-ΔSWE

利用GRACE獲取的ΔTWS，減去GLDAS水文模型數(shù)據(jù)獲取的ΔSWS和ΔSWE，即可得地下水儲量(ground water storage，GWS)的變化ΔGWS。

2.5 TRMM衛(wèi)星降雨數(shù)據(jù)

降雨對地下水儲量研究具有重要意義。本文使用熱帶降雨測量任務(tropical rainfall measuring mission,TRMM)衛(wèi)星數(shù)據(jù)(https:∥disc.gsfc.nasa.gov/datasets/TRMM_3B43_7/summary)估算華北平原的年均降水量。

2.6 小麥藍水足跡

小麥藍水足跡WFblue指在小麥生長過程中用于灌溉消耗的地表與地下水的總量[10-11]。由于華北平原降水量在季節(jié)上分布不均，使得冬小麥在生長期間主要依靠地下水補充灌溉。因此，本文利用冬小麥WFblue數(shù)據(jù)來說明華北平原年均灌溉用水量的空間變化及其對地下水儲量的影響。冬小麥WFblue數(shù)據(jù)可以從水足跡網(wǎng)站(https://waterfootprint.org/en/resources/waterstat/wf-crop-production-and-consumption-china/)獲取，因網(wǎng)站數(shù)據(jù)僅截止到2009年，故本文使用的冬小麥WFblue數(shù)據(jù)時間跨度為2003～2009年。

2.7 EOF的原理

EOF的基本原理是對包含p個空間點的變量隨時間進行分解。假設有n個點的資料，則任一空間點i和任一時間點j的變量值(xij)可以看成由p個空間函數(shù)EOFik和時間函數(shù)tkj(k=1,2,…,p)的線性組合。具體原理參見文獻[12]，本文不再贅述。分解后使用North檢驗對各模態(tài)進行檢核，以保證各模態(tài)相互獨立[13]。

3 結果與分析

3.1 反演結果驗證

圖2、3分別為2003～2013年月尺度和季尺度GRACE獲取的ΔGWS與監(jiān)測井實測ΔGWS的時變序列，二者振幅和趨勢大致相同。從圖2可知，月尺度下，2個時間序列之間的相關系數(shù)為0.725 3；從圖3可知，季尺度下，2個時變序列之間的相關系數(shù)為0.793 3。由此可得，基于GRACE計算得到的地下水儲量變化與實測數(shù)據(jù)之間在時間上有較強的相關性，用GRACE研究華北平原的地下水儲量變化是可行的，且時間尺度越大，結果越準確。

圖2 月尺度GRACE反演和實測ΔGWS時間序列Fig.2 Monthly time series of ΔGWS inversed by GRACE and measured ΔGWS

圖3 季尺度GRACE反演和實測ΔGWS時間序列Fig.3 Seasonal time series of ΔGWS inversed by GRACE and measured ΔGWS

在2003～2005年初，GRACE獲取的地下水儲量變化相比監(jiān)測井實測數(shù)據(jù)整體偏低，這可能是由于現(xiàn)有的監(jiān)測井數(shù)量有限、有些月份缺失、空間分辨率不足所致。

3.2 地下水儲量變化影響因素分析

圖4為2003～2013年監(jiān)測井實測ΔGWS、2003～2015年GRACE獲取的ΔGWS、降水以及2003～2009年冬小麥WFblue隨時間的變化。其中，地下水對降水的響應延遲2～3個月，故將降水時序后移3個月。

圖4 監(jiān)測井實測ΔGWS、GRACE反演的ΔGWS、冬小麥WFblue和降水隨時間的變化Fig.4 ΔGWS obtained from monitoring wells,ΔGWS inversed by GRACE,winter wheat WFblue and precipitation vary with time

可以看出，2003～2015年華北平原地下水儲量以2.62±0.40 cm/a的速率銳減，可分為2個階段：1)2003～2008年地下水儲量處于上升狀態(tài)，速率為3.50±1.47 cm/a；2)2009年開始，地下水儲量由增轉減，2009～2015年地下水儲量的下降速率為4.44±0.54 cm/a，2014～2015年下降最快，達到8.96±4.52 cm/a。

圖4還表明，2003～2005年降水較多，冬小麥WFblue處于較低值，冬小麥灌溉對地下水的依賴程度低，因此地下水儲量在此期間處于上升階段；2006年降水較少，不能滿足冬小麥灌溉的需要，地下水被大量用于灌溉生產(chǎn)，冬小麥WFblue處于2003～2009年的最大值；到了干旱少雨的2009年，降水峰值未超過200 mm，灌溉超采情況進一步加劇，地下水儲量開始急劇下降；2012年和2013年降雨量有所上升，但因農(nóng)業(yè)灌溉常年超采地下水，導致地下水資源補給、排泄的平衡難以恢復，地下水儲量并未回升；2014年和2015年降水量最少，使得地下水儲量進一步減少，地下水儲量下降速率達到峰值。由以上分析可知，冬小麥的灌溉活動和降水與華北平原地下水儲量下降關系密切。

3.3 年際ΔGWS的EOF分解

為清晰地展現(xiàn)華北平原地下水儲量的年際變化，進一步分析影響地下水儲量變化的主要因素，將GRACE獲取的年際ΔGWS進行EOF分解，表1為EOF分解后前5個模態(tài)的方差貢獻率?？梢钥闯觯?個模態(tài)的累計方差貢獻率達到93.09%，基本可以解釋總體的地下水儲量變化特征。

表1 華北平原年際ΔGWS EOF分解后前5個模態(tài)方差貢獻率Tab.1 Accumulated variance contribution rate of the frist five EOF modes of ΔGWS in North China plain

表2為進行North檢驗的前5個模態(tài)的特征值λ、特征值的誤差范圍ej和檢驗值λj-λj+1。從表2可以看出，EOF分析結果均符合North檢驗標準，所對應的EOF是有價值的。綜上考慮，本文選取前2個整體模態(tài)，結合2003～2009年冬小麥年均WFblue數(shù)據(jù)與2003～2015年年均降水數(shù)據(jù)進一步分析華北平原ΔGWS的時空變化特征。

表2 North檢驗結果Tab.2 Results of North test

圖5為第1特征向量與冬小麥年均WFblue對比。從圖5(a)中可以看出，2003～2009年第1特征向量場與研究區(qū)冬小麥年均WFblue空間分布相關系數(shù)為-0.69，呈現(xiàn)良好的負相關性。第1特征向量場在黃河流域與海河流域的大部分地區(qū)為負值，在邯鄲、邢臺、鶴壁、濮陽、新鄉(xiāng)一帶ΔGWS長期處于低值區(qū)，地下水虧損嚴重，與此區(qū)域較高的冬小麥WFblue相對應；只有在淮河流域的南部第1特征向量場為正值，地下水儲量增加，與農(nóng)業(yè)灌溉開采地下水較少有關。不同區(qū)域地下水儲量的變化幅度不同，黃河流域地下水儲量下降最快，以黃河為界，向南北兩側下降速率減小。從圖5(b)中可以看出，2003～2009年第1時間序列與冬小麥年均WFblue時變序列呈現(xiàn)良好的一致性。這種與冬小麥年均WFblue一致的特性占總體方差的80.04%，基本表達了地下水儲量變化場的主要結構，表明冬小麥的灌溉用水是影響華北平原地下水儲量變化的主要因素。

圖5 第1特征向量與冬小麥年均WFblue對比Fig.5 Comparison of the first feature vector and annual average WFblue of the winter wheat

圖6(a)為2003～2015年第2特征向量場與年均降水的空間分布?？梢钥闯?，華北平原的年均降水量由北向南隨緯度遞減，在南部的淮河流域降雨量較大，而中部的黃河下游平原和北部的京津冀地區(qū)則較少?？傮w而言，降水量空間分布不均，以黃河為界，呈現(xiàn)出南多北少的格局，在黃河以北地區(qū)，太行山山前平原整體上降水量相對較少。這與第2特征向量場有著強一致性，二者之間的相關系數(shù)為0.93。結合圖6(b)可知，二者的時變序列趨勢保持一致，表明降水是影響華北平原地下水儲量年際變化的次要因素。

圖6 第2特征向量與年均降水對比Fig.6 Comparison of the second fecture vector and annual average precipitation

4 結 語

本文利用GRACE重力衛(wèi)星數(shù)據(jù)反演了華北平原2003～2015年地下水儲量的年際變化，使用EOF分解方法進行時空分解，并將第1特征向量和第2特征向量結合冬小麥年均WFblue和降水進行對比分析，得到如下結論：

1)GRACE反演結果與監(jiān)測井實測數(shù)據(jù)相關性好，證明了使用GRACE與GLDAS水文模型結合反演地下水儲量時空變化的可行性。

2)2003～2015年華北平原地下水儲量以2.62±0.40 cm/a的速率下降，使用EOF將華北平原地下水儲量的年際變化分解為2個主要模態(tài)，第1模態(tài)占總體變化的80.04%，第2模態(tài)占總體變化的13.05%，二者對總體變化的解釋率可以達到93.09%。

3)第1特征向量場呈現(xiàn)以黃河流域為中心、向南北兩側逐漸變緩的趨勢，與2003～2009年冬小麥年均WFblue空間分布之間有良好的負相關性，相關系數(shù)為-0.69，可見冬小麥的灌溉用水是華北平原地下水儲量減少的主要因素。

4)第2特征向量場方差貢獻率為13.05%，在時間上沒有明顯的趨勢性變化，由北向南隨緯度遞減，與年均降水的空間分布之間的相關系數(shù)為0.93，可知降水是影響華北平原地下水儲量變化的第二大因素。