基于地統(tǒng)計學的寶雞市區(qū)地下水位空間變異特征研究

許義和，魏曉妹，蔡明科，胡國杰，徐萬林

（西北農(nóng)林科技大學 水利與建筑工程學院，陜西 楊凌 712100）

地下水是我國北方城市重要的供水水源，研究地下水的時空變異特征，對合理調(diào)控與管理地下水資源、保障城市用水安全具有重要的意義［1］。目前，國內(nèi)有關地下水時空變異特征的研究較多，但大多數(shù)集中在地下水水質(zhì)方面［2－6］，采用地統(tǒng)計學方法對地下水水位進行時空變異特征的研究還尚不多見。寶雞市位于關中平原的西部，是陜西省第二大工業(yè)城市，市區(qū)地下水開采歷史悠久。2000年初，市區(qū)普查登記水源井503眼，地下水開采量達6984 萬m3，占城市總供水量的66%，超采1106 萬m3。長期超采使市區(qū)地下水位大幅度下降、形成區(qū)域性降落漏斗，在集中開采地段導致淺層承壓含水層被局部疏干。為了合理開發(fā)利用地下水資源，保護生態(tài)環(huán)境，寶雞市政府決定從2001年起逐年關閉一定數(shù)量的單位自備井，以控制市區(qū)地下水開采量。到2004年，自備井關閉數(shù)量已達80眼。同年10月，利用渭河市區(qū)段修筑的攔河閘開始蓄水，形成了總面積約140萬m2市區(qū)湖泊，既創(chuàng)造了優(yōu)美的市區(qū)水域景觀，又涵養(yǎng)了地下水源。本文基于寶雞市DEM數(shù)據(jù)，利用ArcGIS軟件，確定研究區(qū)范圍，通過野外調(diào)查、資料收集，運用地統(tǒng)計學方法對寶雞市區(qū)自備井關閉前后的地下水位時空變異規(guī)律及其特征進行研究，旨在為寶雞市區(qū)地下水資源的保護及可持續(xù)利用提供參考。

1 材料與方法

1.1 研究區(qū)概況

本次研究區(qū)范圍包括寶雞市區(qū)建成區(qū)，面積約100km2。市區(qū)北臨黃土臺塬，南依秦嶺山脈，西為群山所阻，東去進入寬闊的關中平原，呈西高東低、長條狀分布在渭河兩岸的高漫灘和一級階地之上。市區(qū)主要的地貌單元是河谷階地，另外還分布有沖洪積扇和滑坡堆積單元。市區(qū)地下水比較豐富，按地下水埋藏條件和含水介質(zhì)特征分為第四系松散巖類孔隙潛水和第三系碎屑巖類孔隙裂隙承壓水。市區(qū)多年平均降水量為673.1mm，其50%～60%集中在7－9月份。區(qū)內(nèi)水系密布，渭河橫貫全區(qū)，近些年由于降水量偏小兼寶雞峽引渭河水用于灌溉，渭河經(jīng)常發(fā)生斷流。區(qū)內(nèi)渭河兩側(cè)有大小支流9條，除清姜河、金陵河外，其他河流降水后多呈短暫洪水排入渭河，對地下水補給有限。

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

1.2.1 分離研究區(qū)DEM數(shù)據(jù) 本文的地理信息資料來源于互聯(lián)網(wǎng)——國際科學數(shù)據(jù)服務平臺，選用ASTER GDEM30m分辨率數(shù)字高程數(shù)據(jù)。軟件工具采用ArcGIS 9.3以及Excel 2003。由于原始數(shù)據(jù)覆蓋區(qū)域遠大于研究區(qū)，需借助ArcGIS平臺，剪切出包含研究區(qū)的DEM數(shù)據(jù)（見圖1）。

圖1 研究區(qū)DEM

1.2.2 生成矢量研究圖層 通過提取等高線操作，形成研究區(qū)矢量圖層。利用水文分析模塊，提取并分離研究區(qū)范圍內(nèi)渭河及其主要支流圖層。地下水監(jiān)測井圖層依據(jù)寶雞市區(qū)24眼監(jiān)測井1990年、1995年、2000年以及2004年的逐年地下水位埋深和各井的測點高程以及經(jīng)緯度等數(shù)據(jù)建立。建立該點圖層時，將監(jiān)測井的多個屬性數(shù)據(jù)首先保存為Excel文件，然后在ArcGIS軟件中利用其經(jīng)緯度將Excel文件轉(zhuǎn)換為點圖層。導入研究區(qū)和渭河圖層，通過投影變換，使其與監(jiān)測井的投影相一致。

依據(jù)市區(qū)地下水開采現(xiàn)狀調(diào)查，已形成包括姜潭水源地、福臨堡水源地、市區(qū)水源地、石壩河水源地、八里橋水源地、十里鋪水源地、下馬營水源地和臥龍寺等8個水源地在內(nèi)的地下水集中開采區(qū)。市區(qū)水源地及地下水監(jiān)測井空間分布如圖2所示。

圖2 市區(qū)水源地及監(jiān)測井分布圖

1.3 研究方法

由于研究區(qū)范圍較大，地下水動態(tài)監(jiān)測井數(shù)目又相對較少，且分布不均。在這種情況下，要得到整個研究區(qū)的地下水位空間分布，必須要進行插值［7］。本文選用克里金法（Kriging）對寶雞市區(qū)地下水位進行插值。Kriging方法用下面的方程描述一個區(qū)域化變量z（x）［8］，即

式中:u（x）——確定性的全局趨勢分量。ε（x）——自相關隨機的分量，Krining方法假定該分量具有均值為0，滿足二階平穩(wěn)假定的性質(zhì)，它是通過變異函數(shù)進行描述。

變異函數(shù)是地統(tǒng)計學的基本工具，為區(qū)域變量z（xi）和z（xi＋h）增量平方的數(shù)學期望，即區(qū)域化變量的方差，其通式為

式中:γ（h）——變異函數(shù)；h——步長，即為了減少各樣點組合對的空間距離個數(shù)而對其進行分類的樣點空間間隔距離；N（h）——間隔距離為h時的樣點對數(shù)；z（xi）——變量z在空間位置xi的取值；z（xi＋h）——變量z在空間位置xi＋h上的取值。

在本研究中步長的劃分以800m為間距，即800，1600 ，…，8000 m，共10組組距；8000 m 小于樣本距離（17km）的一半，符合地統(tǒng)計分析要求。在后續(xù)模型模擬過程中，綜合考慮曲線擬合情況、誤差大小等因素，對步長予以修正。進行Kriging插值的最關鍵一步在于確定變量的變異函數(shù)。通過對不同模型的運行比較表明，寶雞市區(qū)地下水位的實際變異函數(shù)值與球狀模型（spherical model，公式3）的擬合效果較好，該模型的一般公式為

式中:γ（h）、h含義同式2；C0——塊金值，反映數(shù)據(jù)點的測量誤差；C0＋C——基臺值；a——相關距離，當數(shù)據(jù)點之間的距離大于a時，則認為其是不相關的。

式（3）為各向同性時球狀模型表達式，由于區(qū)域化變量往往具有各向異性，既存在幾何各向異性又存在區(qū)域各向異性［8］，因此，還需確定異性比ρ和方向角θ，前者定義為變程橢圓短軸與長軸之比，后者定義為長軸與y軸（或x軸）的夾角，順時針方向為正。

本文運用上述地質(zhì)統(tǒng)計學方法，借助ArcGIS的地統(tǒng)計分析模塊（Geostatistical Analyst），通過檢驗數(shù)據(jù)正態(tài)性、剔除數(shù)據(jù)趨勢、擬合球狀變異曲線以及交叉驗證等步驟，得到寶雞市區(qū)各年地下水位動態(tài)插值結果，據(jù)此分析市區(qū)地下水位的空間變異特征。

2 結果與分析

2.1 地下水位的趨勢性分析

變異函數(shù)的計算，首先要求數(shù)據(jù)符合正態(tài)分布，否則可能存在比例效應［9］。單樣本柯爾莫哥洛夫·斯米諾夫（One－sample Kolmogorov－Smimov（KS））檢驗表明，1990年、1995年、2000年和2004年的地下水位數(shù)據(jù)符合正態(tài)分布；其次，要求區(qū)域變量滿足二階平穩(wěn)假設［10］，然而地下水位多是非平穩(wěn)的，因此需分離出其趨勢性分量，使數(shù)據(jù)具有較好的正態(tài)性，更好地滿足二階平穩(wěn)假定。在ArcGIS地統(tǒng)計分析模塊中，將4個年份的數(shù)據(jù)分別向2個正交方向投影，以此來分析地下水位整體變化趨勢，如圖3所示。

圖3 地下水位變化趨勢

從圖3中可看出，在各個年份中，市區(qū)地下水位在y方向（南北方向）上均呈現(xiàn)明顯的“U”型趨勢，這與研究區(qū)的地形地貌較為一致，研究區(qū)北邊是渭北黃土臺塬，南邊是秦嶺山脈，中間是地勢較低的渭河漫灘及一級階地，因此呈現(xiàn)南北高、中間低的“U”型態(tài)勢。在x方向上（西東方向），即沿著渭河水流方向，相對當年y方向來說，地下水位的趨勢變化稍顯平穩(wěn)。整體上可以選擇一個二階多項式對地下水位變化趨勢進行較好的模擬。

1990－2004年，南北方向上的地下水位變化趨勢基本沒有發(fā)生改變，而在西東方向上則發(fā)生了變化，1990年地下水位西高東低，2000年轉(zhuǎn)變?yōu)闁|高西低，2004年形成了東西兩頭低，中間高的趨勢。這種變化趨勢與寶雞市區(qū)地形地貌、不同時期地下水開采強度及人工湖建成蓄水等因素有關。1990－2000年，寶雞市區(qū)地下水開采主要集中在西邊的姜潭水源地、石壩河水源地、市區(qū)水源地以及中部的下馬營水源地，東邊的水源地開采強度相對較小。10a開采使得研究區(qū)地下水流向發(fā)生了改變，這種改變主要表現(xiàn)在東西方向上，地下水由南北匯流、單一向東的天然流向轉(zhuǎn)變?yōu)橐韵埋R營水源地東部為界，東西部流向相反的局面，在中西部地下水匯流至各個水源地由于超采形成的漏斗區(qū)，在東部則順渭河而下，與地表水流向一致。2000－2004年，在寶雞市區(qū)中部，由于關閉自備井措施得力，加上人工湖蓄水對地下水產(chǎn)生補給，使地下水位上升速度較東西兩側(cè)快。

2.2 空間變異特征分析

利用二次漂移假定將趨勢性去掉后，采用各向異性球狀模型對研究區(qū)4個年份的地下水位變異函數(shù)進行擬合，擬合后的模型參數(shù)見表1。

表1 各年份球狀模型參數(shù)表

塊金值代表隨機變異量，較大的塊金值表明較小尺度上的某種過程不可忽視?；_值代表變量空間變異的結構性方差，其值越大表示總的空間異質(zhì)性程度越高。當不同的區(qū)域化變量相比較時，塊金值和基臺值并不十分有效，用塊金值與基臺值之比（塊金系數(shù)）來反映塊金方差占總空間異質(zhì)性變異的大小則非常有意義［9，11］。塊金值與基臺值的比值小于0.25時，空間相關性強；在0.25～0.75之間時，空間相關性中等；大于0.75時，空間相關性弱［12］。表1中地下水位的塊金值與基臺值的比值都在0.25～0.75之間，表明具有中等的空間相關性；各年的塊金值均為正值，且值較大，說明存在由采樣誤差、短距離變異、隨機變異等因素引起的塊金效應，而且在800m以下尺度上的短距離變異所占比重較大。

變程反映了在該搜索方向上的空間相關程度［12］，它與觀測尺度以及在取樣尺度上影響地下水特征的各種生態(tài)過程相互作用有關［9］。在變程之內(nèi)，變量具有空間自相關特性，反之不存在。表1中各年份模型長軸變程在5～6km，小于研究區(qū)域的尺度，說明研究區(qū)各年份地下水的空間連續(xù)性差，地下水位的空間自相關距離較短；異性比（短長軸變程之比）從1990年的0.713減小到2004年的0.171，表明市區(qū)地下水位空間各向異性逐漸增強。長軸所在角度與渭河走向基本一致，反映出西東方向上地下水位呈現(xiàn)較長的相關尺度；而短軸所在角度與渭河大約呈垂直方向，表明地下水沿南北方向向渭河排泄。

2.3 地下水位空間插值結果及模型交叉驗證

研究區(qū)地下水位空間插值最終結果如圖4所示。插值過程中充分考慮地下水位的趨勢性和各向異性，并利用模型交叉驗證對插值結果進行評價。模型交叉驗證統(tǒng)計值見表2。一般來說，模型擬合良好的判別原則為［12］:實測樣本數(shù)據(jù)和估計值之間平均誤差絕對值接近0，平均標準誤差接近0，平均預測標準差與均方根誤差最接近，以及標準化均方根預測誤差接近1。因此，各年份的模型擬合結果較為理想。

從插值結果（圖4）來看，南北方向上地下水位是兩頭高中間低，并且地下水流都以近似垂直的角度流向中部，與地形吻合；東西方向上因地形地貌、開采強度等的空間差異，地下水位由西高東低逐漸轉(zhuǎn)變?yōu)闁|西兩頭低，中間高，水流方向亦隨之改變。從整體上來看，市區(qū)地下水位的空間插值結果與地下水位變化趨勢分析是一致的。圖4顯示，1990年研究區(qū)最低地下水位是522m，1995年下降到492m，到2000年上升為502m，2004年則繼續(xù)上升到533m，最低水位分布區(qū)域由最東部逐漸轉(zhuǎn)移到最西部，說明14a間市區(qū)不同地方地下水位變幅較大，且空間分布變化明顯。地下水位的變化受開采強度的影響較大，地下水資源的采補難以保持平衡。根據(jù)統(tǒng)計數(shù)據(jù)分析，1990－1995年，開采系數(shù)持續(xù)增大且大于1的水源地有福臨堡水源地、市區(qū)水源地和石壩河水源地，超強度開采使得三水源地部分低水位區(qū)域相連，并形成了降落漏斗。姜潭水源地和十里鋪水源地持續(xù)強開采則形成了次一級低水位區(qū)域，從整體上看，市區(qū)地下水處于嚴重超采狀態(tài)；到2000年，雖然姜潭水源地、石壩河水源地以及臥龍寺水源地仍處在強開采或超強開采的狀態(tài)，但由于其他水源地開采強度已大幅度消減，地下水補給大于開采，市區(qū)地下水處于強開采狀態(tài)。2000年后，隨著市區(qū)自備井逐年關閉，地下水開采量也在逐年減小，到2004年，地下水開采量僅為4712 .4萬m3，較1990年減少了36.2%，這使得市區(qū)整體地下水位迅速回升，地下水位最低值由2000年的502m上升到2004年的533m，上升幅度達30m，地下水補給繼續(xù)大于開采，市區(qū)地下水處于中等開采強度。

圖4 插值結果

綜上所述，寶雞市區(qū)地下水位在短時期內(nèi)呈現(xiàn)出較為劇烈的波動既受制于區(qū)域特殊的地形及地質(zhì)條件，同時又與市區(qū)地下水的補排條件變化密切相關。

表2 模型交叉驗證統(tǒng)計值

3 結論

（1）各年份球狀模型擬合的塊金值與基臺值的比值都在0.25～0.75之間，因此寶雞市區(qū)地下水位具有中等的空間相關性。各年份的塊金值均為正值，且值較大，說明在800m以下尺度上的短距離變異所占比重較大。各年份模型長軸變程小于研究區(qū)域的尺度，說明研究區(qū)各年份地下水的空間連續(xù)性差，地下水位的空間自相關距離較短。異性比（短長軸變程之比）從1990年的0.713減小到2004年的0.171，說明研究區(qū)地下水位空間各向異性逐漸增強。

（2）交叉驗證顯示各年份的模型擬合結果較為理想。插值結果表明，市區(qū)不同地方的地下水位變幅較大，且空間分布變化明顯，說明研究區(qū)局部地下水位變化受開采強度變化的影響較大。

［1］李新波，郝晉珉，胡克林，等.集約化農(nóng)業(yè)生產(chǎn)區(qū)淺層地下水埋深的時空變異規(guī)律［J］.農(nóng)業(yè)工程學報，2008，24（4）:95－98.

［2］米日姑·買買提，海米提·依米提，古麗娜爾·托合提，等.伊犁河流域地下水埋深與TDS時空變異及水化學特征分析［J］.水土保持學報，2010，24（2）:174－178.

［3］陳淑峰，李帷，胡克林，等.基于GIS的華北高產(chǎn)糧區(qū)地下水硝態(tài)氮含量時空變異特征［J］.環(huán)境科學，2009，30（12）:3541－3547.

［4］賈艷紅，趙傳燕，南忠仁.黑河下游地下水波動帶土壤鹽份空間變異性研究［J］.土壤學報，2008，45（3）:420－430.

［5］姚榮江，楊勁松，劉廣明.土壤鹽分和含水量的空間變異性及其CoKriging估值:以黃河三角洲地區(qū)典型地塊為例［J］.水土保持學報，2006，20（5）:133－138.

［6］蘇里坦，宋郁東，張展翔.新疆三工河流域地下水礦化度的時空變異及其分形特征［J］.地質(zhì)科技情報，2005，24（1）:85－90.

［7］李新，程國棟，盧玲.空間內(nèi)插方法比較［J］.地球科學進展，2000，15（3）:260－295.

［8］曹洪，駱冠勇，廖建三，等.廣州城區(qū)地下空間開發(fā)對地下水環(huán)境的影響研究［J］.巖石力學與工程學報，2006，25（2）:3347－3356.

［9］肖篤寧，李小玉，宋冬梅，等.民勤綠洲地下水開采時空動態(tài)模擬［J］.中國科學，2006，36（6）:567－578.

［10］王仁鐸，胡光道.線性地質(zhì)統(tǒng)計學［M］.北京:地質(zhì)出版社，1988.

［11］韓業(yè)珍，魏曉妹，李立.基于地統(tǒng)計學的地下水位時空變異特征研究［J］.人民黃河，2010，32（5）:52－53.

［12］阮本清，許鳳冉，蔣任飛.基于球狀模型參數(shù)的地下水水位空間變異特性及其演化規(guī)律分析［J］.水利學報，2008，39（5）:573－579.