楊 靜,肖天昀,李海波,王全榮 (中國地質大學(武漢)環境學院,湖北 武漢 430074)
地下水硝酸鹽污染已經成為世界性的環境問題,尤其是集約化種植區氮肥投入的逐年加大,污染更為嚴重[1-4].硝酸鹽污染直接威脅到人體健康,飲用水中硝酸鹽含量過高會影響血液中氧的傳輸能力,引發高鐵血紅蛋白癥[5].江漢平原作為我國重要的糧食生產基地,采取輪作耕作制度,屬典型集約農區,農業面源污染嚴重.該區地下水水位埋深較淺,地下水水質受地表環境和人類活動的影響較大[6].已有研究表明,江漢平原地下水中硝酸鹽含量季節變化明顯,地表水體氮素含量變化與雨季、農業施肥有一定的同步性[7-8].其主要原因是硝酸鹽具有良好的溶解性和移動性,農業生產中施用的氮肥可以通過土壤進入含水層中,造成地下水中硝酸鹽含量增加.另外,生活污水、垃圾與糞便的下滲水也是地下水硝酸鹽污染的主要來源[9].
目前針對地下水中硝酸鹽污染特征的研究主要是利用地質統計學方法.地質統計學是一門研究具有一定結構和區域化變量的數學地質方法.通過地質統計對未采樣區域變量取值進行線性無偏最優估計,分析其空間分布特征[10].部分學者結合GIS和數學模型對影響地下水中硝酸鹽含量的因素進行了辨識,認為地下水位埋深、含水層巖性、土地利用方式等是影響地下水中硝酸鹽含量的主要因素[11-12].此外,地表氮負荷量、降雨入滲補給量也是影響地下水中硝酸鹽含量的重要因素[13-14].然而上述研究大多集中在影響因子的定性識別上,對于構建地下硝酸鹽污染預測模型,定量分析不同因子對硝酸鹽污染的影響關系還鮮見報道.
本文基于2011~2014年地下水采樣測試數據,利用半方差函數模型和普通克里格插值方法,分析江漢平原地下水中硝酸鹽含量的空間變異和分布特征,并以證據權重法模型為基礎,構建江漢平原地下水硝酸鹽污染預測模型,并對模型中的因子進行統計分析,定量評價不同因子對地下水硝酸鹽污染的影響關系,以期為江漢平原地下硝酸鹽污染的治理提供一定的科學依據和參考.
江漢平原地處長江中游,湖北省中南部,位于東經111.3°~114.6°,北緯29.4°~31.4°之間,區域面積約3.87×104km2,主要由長江和漢江及其支流沖積物形成(見圖1).屬亞熱帶濕潤季風氣候區,年均氣溫15~17℃,年均降雨量約1269mm.含水層在垂向上分為淺層全新統孔隙潛水含水巖組、上部上、中更新統孔隙承壓含水巖組及下部下更新統、第三系裂隙孔隙承壓含水巖組[15].第四系全新統松散孔隙潛水是江漢平原潛水的主要類型,含水層厚度一般小于20m,巖性以粘土、粉質粘土、粉砂為主,局部地段有砂礫石層,滲透系數一般小于1.5m/d.地下水位埋深較淺,平原腹地一般0.5~2.0m.大氣降雨、地表水體的滲流及田間灌溉水的入滲是其主要補給來源,該層是農村分散住戶人畜用水的主要供水水源.上部孔隙承壓含水巖組在區內分布最為廣泛,巖性以粉砂、砂、砂礫石為主,普遍含有淤泥.含水層厚度變化較大,平原腹地最大厚度90~100m,邊緣區最薄的只有10~15m.下部裂隙孔隙承壓含水巖組主要分布在江漢平原腹部,多不連續且呈透鏡狀.巖性以粉砂、泥質粉砂、細砂為主,自上而下普遍含有淤泥,局部含少量礫石[16].研究區土壤類型以人為填土和沖積土為主,分別占全區面積的47.0%和30.4%.人為填土主要分布在研究區北部靠山地的平原、丘陵地帶,沖積土則主要集中在長江、漢江河流兩岸沖洪積階地.

圖1 江漢平原地理位置和采樣點分布Fig.1 Location of Jianghan Plain and distribution of the sampling points
2011~2014年對江漢平原淺層地下水進行隨機布點取樣,采樣點的設計主要考慮上覆土壤類型、水文、土地利用方式等,用GPS定位采樣點坐標,最終取樣825件(見圖1).取樣點類型多為機井和壓把井,取樣層位為淺層全新統孔隙潛水含水巖組.采集前先用潛水泵或壓把井將井中殘留水排盡,待井水清澈,再取新鮮水樣.現場采用雷磁DZB-718型便攜式多參數水質分析儀測定水溫、pH值、電導率、溶解氧、總硬度等.取樣后用注射器通過0.45μm微孔濾膜對水樣進行抽濾,抽濾后的水樣分別裝入2個50mL聚乙烯塑料瓶中,1瓶用優級純硝酸酸化至pH<2,用于常規陽離子分析,1瓶用于常規陰離子分析.另取一瓶500mL未過濾的水樣于24h內用滴定法進行堿度分析.取樣瓶先在實驗室內用去離子水沖洗3次,在野外取樣時再用樣品水充分潤洗3次.樣品于一周內送至中國地質大學(武漢)生物地質與環境地質國家重點實驗室進行檢測分析.Ca2+、Mg2+、K+、Na+等陽離子采用電感耦合等離子體光譜儀(ICP-AES, Thermo IRIS Intrepid II XSP)測定,F-、Cl-、、等陰離子采用離子色譜儀(Dionex ICS-1000)分析.采用標準樣、平行樣和空白樣多種監控手段進行質量控制,保證數據的準確度和精度.
用PHREEQC軟件對825個樣品檢測結果進行陰陽離子平衡分析,33個樣品的離子平衡誤差大于5%,不符合天然水中電荷平衡的一般要求,后續分析過程中予以剔除.地下水硝酸鹽檢出點589個,檢出率74.37%,并將樣品含量測定結果換算為與其對應含量.江漢平原地下水中濃度范圍為0.05~23.11mg/L,平均值5.39mg/L,標準差5.94mg/L,17個采樣點地下水含量超過 20mg/L. Box-Cox變換是一種廣泛應用的將數據變換到正態的方法,可以根據數據自動尋找“最佳”變換函數[17].采用Box-Cox變換(λ=0.21)對數據進行預處理后,使用SPSS 20.0對變換后的數據進行單樣本柯爾莫哥洛夫-斯米諾夫檢驗,結果表明,變換后濃度數據符合正態分布,顯著性水平α=0.76,滿足地學統計要求.

表1 不同半方差函數擬合主要參數及交叉驗證結果Table 1 Critical parameters of semi-variogram function and cross validation results

圖2 江漢平原地下水含量分布Fig.2 Distribution of the concentration in groundwater in Jianghan Plain

表2 證據因子分類及權重計算結果Table 2 The classes and weights of evidence factors
證據權重法是一種基于概率不確定性與貝葉斯定理的空間決策支持模型[23].該方法通過研究已知類型礦床(訓練點或響應因子)和二元圖層(證據因子或預測因子)之間的空間關系,計算證據因子的權重值和對比度等參數,對優選后的證據因子加權疊加,得到響應因子發生的后驗概率.證據權重法實質是一種空間統計模型,該方法不僅可以根據后驗概率的高低對整個區域進行遠景預測,而且能夠分析響應因子和證據因子之間的空間關系.結合GIS技術,基于證據權重法的地下水資源評價均取得了較好的研究成果[24-26].
3.2.1 先驗概率計算 假設研究區面積為A{ T},將研究區劃分為N{ T}個等大小的網格,保證每個網格至多有一個響應因子,每個網格面積為μ,定義D為響應因子,統計整個研究區內發生D的網格數N{ D},則區域內任意一個網格發生D的先驗概率P{ D}可以通過式(1)求得:

由ArcSDM模塊,定義每個網格面積為25km2,計算其先驗概率P{D}=6.1305%.
3.2.2 證據權計算 對任意證據二值(或多值)因子,假設Bj表示該因子的某種模式存在,表示該因子某種模式不存在,根據貝葉斯定理,響應因子D發生的后驗概率可以由式(2)、(3)計算:



證據因子Bj的一對正負權重定義為:當先驗概率和后驗概率相等時,==0,表示對該模式而言響應因子是隨機分布的,即該模式對響應因子的分布具有不可預見性.各證據因子權重計算結果見表2.部分因子分類借鑒DRASTIC方法提供的分類標準[28].
3.2.3 證據因子優選 對任意證據因子Bj,該證據因子和響應因子之間的空間關聯性強弱可以用對比度(正負權重的差值)Cj衡量.Cj越大,表示證據因子Bj對響應因子D的預見性越好.引入學生化對比度Stud(C)確定最優切值,可以克服權重的不確定性和數據缺失出現的差異,其統計學含義是對C的顯著性檢驗.采用綜合考慮Stud(C)>1.96和N{ B∩D}最大原則對證據因子進行優選,結果見表2.
貝葉斯定理以條件獨立性假設為前提,因此在計算后驗概率之前,必需研究各證據因子之間的獨立性.綜合檢驗法以實際統計出的響應因子點數N{ D}與預測出來的響應因子點數Np{D}的比值作為衡量標準,若N{ D}Np{D}>0.85,則認為各證據因子滿足條件獨立性要求[29].經計算,本次預測模型N{ D}Np{D}=0.87,符合檢驗要求.
3.2.4 響應因子后驗概率分布 將優選后的證據因子轉化為二值圖層后疊加,根據式(6)~(9)就可以得到響應因子D的后驗概率分布:

將所得后驗概率分為4個等級:低概率(<P{D})、較低概率(P{D}~1.5P{D})、較高概率(1.5P{D}~2.0P{D})及高概率(>2.0P{D}),最終結果如圖3.

圖3 江漢平原地下水后驗概率分布Fig.3 The posteriori probability distribution map of the concentration in groundwater in Jianghan plain
成功率曲線可以評價模型的執行性能和預測精度[30].利用不同等級的脆弱性區域累積面積比例(按后驗概率從高到低)和響應因子(超過濃度閾值的硝酸鹽采樣點)落入該等級的累積比例關系,繪制的成功率曲線見圖4. 20%的高后驗概率區域具有84%的響應因子,后驗概率大于先驗概率的區域中響應因子達90%,曲線線下面積為0.91,說明本次預測結果較好.

圖4 成功率曲線Fig.4 The success rate curve
對比度可以衡量證據因子對響應因子的影響情況.研究各證據因子不同分類類別下對比度,見圖5(a)~5(g),分析其響應關系和原因.
降雨入滲補給量:降雨入滲補給量和其對比度呈明顯的正響應關系.當降雨入滲補給量范圍為101.6~177.8mm/a時,對比度接近0,指示當前范圍降雨入滲補給量對地下水受硝酸鹽污染的可能性(含量大于10mg/L)影響較小.降雨入滲補給不但在包氣帶垂向傳輸硝酸鹽,還控制著硝酸鹽的彌散和稀釋作用.補給量越大,地下水受硝酸鹽污染的可能性也越大,但補給量大到一定程度以致硝酸鹽被稀釋時,地下水受污染的可能性將變小.本次研究結果表明,對比度隨著降雨入滲補給量的增加而增加,并沒有呈現減小的趨勢,說明降雨入滲補給量可能并未超過“稀釋閾值”.
地下水位埋深:地下水位埋深和其對比度呈明顯的負響應關系,對比度隨地下水位埋深的增加而減小.地下水位埋深小于1.5m時,對比度為正,地下水位埋深大于1.5m時,對比度為負.地下水位埋深決定地表污染物到達含水層的運移路徑.地下水位埋深越大,污染物的運移路徑越長,污染物與氧氣接觸的時間越長,被稀釋的機會越大,到達含水層的可能性越小.
滲透系數:滲透系數和其對比度呈明顯正響應關系.滲透系數小于0.5m/d時對比度為負,大于0.5m/d時對比度為正.給定水力梯度條件下,滲透系數控制地下水的流動速率,進而決定硝酸鹽進入含水層之后的遷移速率.負的對比度指示了當滲透系數小于0.5m/d時,地下水中硝酸鹽遷移速率較慢,含水層相對不易受硝酸鹽的污染.
包氣帶巖性:將研究區包氣帶巖性分為粘土、亞粘土、粉質粘土、粉土、砂礫石5類.粘土和亞粘土的對比度為負值,粉質粘土的對比度接近0,粉土和砂礫石的對比度為正值.表明包氣帶巖性為砂礫石時,地下水受硝酸鹽污染可能性較大,包氣帶巖性為粘土和亞粘土時,地下水不易受硝酸鹽污染.

圖5 證據因子分析Fig.5 The analysis of the evidence factors
土地利用類型:利用ENVI遙感影像處理軟件,基于2013年TM遙感數據,采用非監督分類方法將研究區土地利用類型分為耕地、建設用地、水體、草地及林地5類.水體覆蓋范圍內沒有采樣點,對比度為0.林地和草地的對比度為負值,耕地和建設用地對比度為正值.分析認為研究區居民區生活污水、垃圾、糞便的下滲濾液和農業氮肥等是地下水中硝酸鹽的重要污染來源.此外,林地可以有效“狙擊”地下水硝酸鹽污染.
土壤類型:沖積土的對比度為正值,潛育土、人為填土、飽和粘盤土和部分薄層土對比度為負值.區內沖積土多由長江、漢江近期沖積物演化而來,土質疏松,排水通氣良好,有機質含量較多,硝化作用強烈,地下水更易受硝酸鹽污染.
土壤全氮含量:隨著土壤全氮含量的增加,對比度由負值變為正值,二者呈明顯正響應關系.在1.0~1.5g/kg時其對比度接近0,說明當前范圍土壤全氮含量對地下水受硝酸鹽污染的可能性基本無指示意義.土壤中全氮的累積,特別是氮素投入超過作物需求和吸收時,未被吸收利用的硝態氮極易向下淋洗,造成地下水硝酸鹽污染.
證據因子正權重最大值可以反映其對響應因子的影響程度,即“貢獻度”.從圖5(h)可以看出,響應因子后驗概率達到最大值時,證據因子“貢獻度”從大到小依次是:包氣帶巖性、土壤全氮含量、地下水位埋深、土地利用類型、滲透系數、降雨入滲補給量和土壤類型.
5.1 經Box-Cox變換后, 研究區地下水NO3-N含量服從正態分布,半方差函數模型分析結果確定其最優擬合模型為球狀模型.NO3-N含量在一定范圍存在空間自相關性,相關距離為68.02km,空間相關性為87.93% ,指示人類活動(生活污水污染、農業化肥污染及污水回灌等)對地下水中硝酸鹽含量分布影響較大.
5.2 普通克里格插值結果顯示研究區地下水NO3-N含量具有南部高、北部及西部山前平原低的特點.研究區南部長江沿岸NO3-N含量普遍大于10.0mg/L.
5.3 成功率曲線顯示模型預測精度為0.91.后驗概率結果表明研究區最容易發生硝酸鹽污染的區域主要分布在研究區南部長江沿岸,約占全區面積12.65%,與實測NO3-N含量高值區分布一致.
5.4 包氣帶巖性對江漢平原地下水硝酸鹽脆弱性的“貢獻度”最大,其次是土壤全氮含量,“貢獻度”最小的是土壤類型.地下水受硝酸鹽污染的可能性隨著降雨入滲補給量、滲透系數和土壤全氮含量的增加而增加,隨著地下水埋深的增大而減小.當包氣帶巖性為砂礫石,土地利用類型為建設用地和土壤類型為沖積土覆蓋條件時,其對比度均達到最大,表明此時地下水極易受硝酸鹽污染.
[1] Masetti M, Poli S, Sterlacchini S, et al. Spatial and statistical assessment of factors influencing nitrate contamination in groundwater [J]. Journal of Environmental Management, 2008,86(1):272-281.
[2] 孟志龍,楊永剛,秦作棟,等.汾河下游流域水體硝酸鹽污染過程同位素示蹤 [J]. 中國環境科學, 2017,37(3):1066-1072.
[3] Hu K, Huang Y, Li H, et al. Spatial variability of shallow groundwater level, electrical conductivity and nitrate concentration, and risk assessment of nitrate contamination in North China Plain [J]. Environment International, 2005,31(6):896-903.
[4] 劉興權,許晶玉,江麗華,等.山東省種植區地下水硝酸鹽污染空間變異及分布規律研究 [J]. 農業環境科學學報, 2010,29(6):1172-1179.
[5] Greer F, Shannon M. Infant methemoglobinemia: the role of dietary nitrate in food and water [J]. Pediatrics, 2005,116(3):784-786.
[6] 王建偉,張彩香,潘真真,等.江漢平原地下水中有機磷農藥的分布特征及影響因素 [J]. 中國環境科學, 2016,36(10):3089-3098.
[7] 張 婷,陳世儉,傅嬌鳳.四湖地區地下水“三氮”含量及時空分布特征分析 [J]. 長江流域資源與環境, 2014,23(9):1295-1330.
[8] 田月明,盧碧林,鄭林章.江漢平原集約農區水體氮通量分析 [J].湖北農業科學, 2013,52(19):4612-4618.
[9] Wakida F, Lerner D. Non-agricultural sources of groundwater nitrate: a review and case study [J]. Water Research, 2005,39(1):3-16.
[10] 鄧 林,王文科,楊曉婷,等.關中盆地地下水硝酸鹽含量的空間變異特征 [J]. 干旱區資源與環境, 2008,22(10):152-155.
[11] Nolan B, Hitt K, Ruddy B. Probability of nitrate contamination of recently recharged groundwaters in the conterminous United States [J]. Environmental Science & Technology, 2002,36(10):2138-2145.
[12] Gardner K, Vogel R. Predicting ground water nitrate concentration from land use [J]. Ground Water, 2005,43(3):343-352.
[13] Wang M, Liu G, Wu W, et al. Prediction of agriculture derived groundwater nitrate distribution in North China Plain with GIS-based BPNN [J]. Environmental Geology, 2006,50(5):637-644.
[14] Baker D. Groundwater quality assessment through cooperative private well testing - an Ohio Example [J]. Journal of Soil and Water Conservation, 1990,45(2):230-235.
[15] 趙德君.江漢平原地下水系統三維數值模擬 [D]. 武漢:中國地質大學, 2005.
[16] Zhou Y, Wang Y, Li Y, et al. Hydrogeochemical characteristics of central Jianghan Plain, China [J]. Environmental Earth Sciences,2012,68(3):765-778.
[17] Box G, Cox D. An analysis of transformations [J]. The Royal Statistical Society. Series B (Methodological), 1994,26(2):211-252.
[18] Zhu K, Cui Z, Jiang B, et al. A DEM-based residual Kriging model for estimating groundwater levels within a large-scale domain: a study for the Fuyang River Basin [J]. Clean Technologies and Environmental Policy, 2012,15(4):687-698.
[19] 熊秋林,趙文吉,宮兆寧,等.北京城區2007~2012年細顆粒物數濃度時空演化 [J]. 中國環境科學, 2013,33(12):2123-2130.
[20] 王幼奇,白一茹,王建宇.基于GIS的銀川市不同功能區土壤重金屬污染評價及分布特征 [J]. 環境科學, 2016,37(2):710-716.
[21] 呂君偉,劉湘南.香港近岸海域懸浮固體濃度空間變異特征的地統計分析 [J]. 中國環境科學, 2014,34(3):734-741.
[22] Cambardella C, Moorman T, Novak J, et al. Field-scale variability of soil properties in central Iowa soils [J]. Soil Science Society of America Journal, 1994,58(5):1501-1511.
[23] Bonhamcarter G, Agterberg F, Wright D. Integration of geological datasets for gold exploration in Nova-Scotia [J]. Photogrammetric Engineering and Remote Sensing, 1988,54(11):1585-1592.
[24] Sorichetta A, Ballabio C, Masetti M, et al. A comparison of data-driven groundwater vulnerability assessment methods [J].Ground Water, 2013,51(6):866-879.
[25] 雷 靜.地下水環境脆弱性的研究 [D]. 北京:清華大學, 2002.
[26] 孫才志,王言鑫.基于WOE法的下遼河平原地下水硝酸鹽氮特殊脆弱性研究 [J]. 水土保持研究, 2009,16(4):80-84.
[27] 金銀龍.GB5749-2006《生活飲用水衛生標準》釋義 [M]. 北京:中國標準出版社, 2007.
[28] Yang J, Tang Z, Jiao T, et al. Combining AHP and genetic algorithms approaches to modify DRASTIC model to assess groundwater vulnerability: a case study from Jianghan Plain,China [J]. Environmental Earth Sciences, 2017,76(12):426.
[29] Carranza E. Weights of Evidence modeling of mineral potential: a case study using small number of prospects, Abra, Philippines [J].Natural Resources Research, 2004,13(3):173-187.
[30] Chung C, Fabbri A. Probabilistic prediction models for landslide hazard zonation [J]. Photogrammetric Engineering & Remote Sensing, 1999,65(12):1389-1399.