胡旭東 ,聶晶晶 ,梅紅波 ,李遠遠,任曉杰,李夢迪
(1.中國地質大學(武漢)資源學院,武漢 430074; 2.云南省環(huán)境監(jiān)測中心站,昆明 650034)
地下水脆弱性一詞最早是在1968年由法國的Albinet和Margat首次提出[1],目前國際上比較公認的定義是1993年美國國家科學研究委員會給定的“地下水脆弱性是污染物到達最上層含水層之上某特定位置的傾向性與可能性”[2]。國內對地下水脆弱性的研究始于20世紀90年代中期[3],雖然起步較晚,但是近十幾年發(fā)展較快。地下水資源對于人們生活和社會發(fā)展起著不可替代的作用,已經(jīng)越來越成為人們的關注的焦點。
地下水脆弱性評價方法主要有三類:迭置指數(shù)法、過程數(shù)學模擬法和統(tǒng)計方法[4]。三類方法各有優(yōu)缺點,其中迭置指數(shù)法以其低成本、操作簡單、數(shù)據(jù)易獲取、表達直觀等優(yōu)點得到廣泛的運用。迭置指數(shù)法中最常用也是最經(jīng)典的評價模型是美國環(huán)保局在1987年提出的DRASTIC法[5,6]。該方法選取了7種與地下水脆弱性相關的指標,根據(jù)指標自身的屬性建立評分系統(tǒng),然后按照不同指標對脆弱性影響程度的大小對其分配權重,將評價區(qū)劃分為若干個單元,對各單元的各評價指標進行加權求和,即得到脆弱性評價值,如公式(1)所示。
(1)
式中:Rij表示第i個評價單元第j個指標的評分;ωj表示第j個指標相應的權重。
DRASTIC模型的評價指標在傳統(tǒng)權重分配方法上是綜合了美國多位水文地質學專家的經(jīng)驗給出的。而在實際運用中,照搬這一套評權重值顯然不合理,因為各個評價區(qū)的地理、地質以及水文地質條件都不盡相同,需要根據(jù)評價區(qū)的具體情況,結合對評價區(qū)屬性的了解,合理選擇評價指標以及判斷出各個指標的相對重要程度。指標權重賦值直接決定評價結果的精度。地下水脆弱性研究中賦權方法主要包括主觀賦權法和客觀賦權法。主觀賦權法對經(jīng)驗具有較強的依賴性,權值大小取決于分析者對實際情況的了解程度,決策結果具有較強的主觀隨意性,客觀性較差;客觀賦權法具有較強的數(shù)學理論依據(jù),克服了主觀意愿,然而賦權結果依賴于數(shù)據(jù)質量和數(shù)量,往往導致評價結果與實際情況不符,比如最重要的評價指標并未賦予最大的權重值。綜合來看,主觀方法和客觀方法各利弊,而基于DRASTIC模型的地下水脆弱性評價對評價指標的賦權過程具有主觀特點和客觀特點[7],因此本文結合兩者進行綜合賦權,其中主觀方法選取的是改進的層次分析法(Analytic Hierarchy Process,AHP),客觀方法選取的是因子分析法(Factor Analysis,F(xiàn)A)。將層次分析法與因子分析法綜合運用于地下水脆弱性評價之中,旨在確定相對合理的指標權重,最終提高評價結果的可信度。
層次分析法是1973年美國運籌學家Saaty創(chuàng)建的一種非常重要的決策方法[8]。作為主觀賦權法[9],層次分析法通過主觀驅動,使用數(shù)學轉換將性質明確但難以量化的問題定量化表達,具有操作簡單、實用性和適應性強等優(yōu)點[10,11]。層次分析法計算過程包括建立層次結構模型、構建判斷矩陣、計算各指標的權重、一致性檢驗等步驟,其中首要的步驟是構建判斷矩陣。傳統(tǒng)的層次分析法以九尺度原則確定判斷矩陣,然而指標間對應的重要程度不好把握,難以確定合理的判斷矩陣。本文采用三標度法簡化模型[12,13],有利于比較相對關系,經(jīng)過數(shù)學運算可由標度矩陣間接得到判斷矩陣。此外,判斷矩陣的一致性檢驗也是層次分析法中的核心問題,一致性檢驗通過加入平均隨機一致性指標的檢驗系數(shù)CR衡量[14],CR小于0.1即通過一致性檢驗。反之若CR大于0.1,表示偏離了一致性,則此時求得的權重矩陣W不可靠。
然而三標度法并未逃離層次分析法固有的缺陷,即一致性檢驗需要多次調試,調試過程中比較盲目,并不能確切地掌握是哪些因素導致的不一致結果,若一直不能通過一致性檢驗則要耗費大量的工作量[15]。因此本文在三標度的基礎上引入最優(yōu)傳遞矩陣的方法[16],采用自調節(jié)方式將構建的標度矩陣轉化為一致性矩陣,從而計算各個指標的權重值而不需要再進行一致性檢驗。具體計算過程如下:
利用三標度判斷指標間重要程度關系,構建標度矩陣A=(aij)n×n,矩陣A滿足aij=-aji,為反對稱矩陣。三標度準則如表1所示。
表1 三標度判斷值Tab.1 Standard based on 3 scale
設矩陣T=(tij)n×n為矩陣A的最優(yōu)傳遞矩陣,其中:
(2)
若矩陣F中的元素滿足fij=1/fji(互反矩陣)且fij=fik/fjk,即為一致性矩陣,利用傳遞矩陣確定一致性矩陣F=(fij)n×n,其中:
fij=exp(tij)
(3)
根據(jù)一致性判斷矩陣F計算出各指標權重。F矩陣的最大特征值對應的特征向量作為各個指標的相對權重向量,即FX=λmaxX。其中,X=[X1,X2,…,Xn]T為特征向量,其近似值作為評價指標的相對權重。本文采用根法計算特征向量的近似值:
(4)
歸一化得各個指標權重:
(5)
客觀方法選用因子分析法,因子分析法是一種將高維數(shù)據(jù)降維的方法[17]。在地下水脆弱性評價中,一種被稱為主成分-因子分析的方法被廣泛運用于評價指標賦權方面[18,19],其本質上還是因子分析法,只不過在因子荷載矩陣的計算上使用了主成分法。除了主成分法,還有主軸因子法、最小二乘法、極大似然法等方法可以計算出因子荷載矩陣[20]。
因子分析法是在主成分分析法基礎上發(fā)展而來的,其基本思想是,在保證原始信息量丟失最少的情況下,使用降維方法,分析原始變量之間的公共關系,將原始變量進行歸類,歸類后的因子稱為公共因子,使得原本需要許多信息重疊的變量描述的問題,通過因子分析法后變成少數(shù)幾個帶有主要信息的公共因子即可描述。
因子分析法的步驟是:對原始變量進行標準化,計算樣本相關系數(shù)矩陣R、求解R的特征根λi(λ1≥λ2≥…≥λp>0)和相應的標準正交的特征向量li,確定公共因子數(shù)m(按前m個特征值之和占特征值總和的百分比來確定),計算主因子載荷矩陣A=[aij];根據(jù)需要對載荷矩陣進行旋轉,得到能被更好解釋的公共因子,最后計算因子得分矩陣。整個計算過程可以通過軟件SPSS快速實現(xiàn)。因子分析的主要目的是計算可解釋的公共因子、各個樣品在公共因子上的得分及各變量的權重。由得分矩陣可得出各主因子的得分,如公式(6)所示。
(6)
式中:fi為第i個主因子的得分;βij為相應的得分系數(shù)。
由于主因子可以代替原來的變量研究問題,而每個主因子fi對事件的貢獻度可以用相關系數(shù)矩陣的特征根λi衡量,則原始變量的權重值是由每個主因子的特征值做權數(shù)對每個主因子得分系數(shù)進行加權求和得到,如公式(7)所示。
(7)
由公式(6)和(7)可以計算各變量的權重,如公式(8)所示。
(8)
歸一化后得最終的權重值,如公式(9)所示。
(9)
需要指出的是,由主成分分析求出的因子荷載矩陣只是一個很基礎的矩陣,理論上因子荷載矩陣有很多種,一般求解過程會進行因子旋轉,目的是使因子荷載值兩極分化,趨向于0或者±1,相關系數(shù)矩陣的特征值則會發(fā)生變化,即方差貢獻度會變化,但是主因子的總貢獻度保持不變。之后由因子荷載矩陣求得的因子得分系數(shù)以及后續(xù)得到的權重值發(fā)生改變。上述內容表明理論上得出的權重值不唯一,介于因子旋轉之后更符合實際意義和解釋,本文在求解權重時使用因子旋轉之后的結果。
綜合賦權結合了主觀賦權和客觀賦權,綜合兩者信息,既尊重主觀又不失客觀[21]。計算方法一般包括積法以及和法。和法的公式添加了一個人為估計的參數(shù),用來衡量主觀權值和客觀權值的重要程度,導致該方法操作較為復雜,而且對于參數(shù)尺度也不好把握,因此本文選擇積法求取綜合賦權值,如公式(10)所示[22]。
(10)
昆明安寧工業(yè)園草鋪片區(qū)位于云南省安寧市草鋪鎮(zhèn)中部,占地面積約為62 km2。研究區(qū)屬中亞熱帶低緯度高海拔地區(qū),平均海拔1 800 m,年平均氣溫14.8 ℃,年平均降雨量1 087.2 mm。區(qū)內地表水系發(fā)育,主要為鳴矣河和螳螂川,鳴矣河是螳螂川的一級支流,在連然鎮(zhèn)通仙橋匯入螳螂川。地下水類型主要包括碎屑巖裂隙水、火成巖裂隙水以及孔隙水,富水性為中等。區(qū)內地下水主要接受大氣降水補給,局部受地表水補給,多以分散的隙流或片狀散流形式向當?shù)刈畹团判够鶞拭骟氪ǜ伞⒅Я髋判埂?/p>
經(jīng)過近15年的發(fā)展,草鋪片區(qū)已成為云南省主要石油及磷化工基地,由于工業(yè)的快速發(fā)展,以氟化物、煙(粉)塵以及二氧化硫等為主的企業(yè)污染物對當?shù)卮髿猸h(huán)境、土壤環(huán)境、植物農(nóng)作物等造成不同程度的影響[23]。研究區(qū)基本處于地下水的排泄區(qū),且大部分區(qū)域淺層地下水埋深較淺,為易受污染區(qū)。當?shù)氐叵滤Y源正受到嚴重的潛在威脅。
(1)評價指標選取。本文基于DRASTIC模型,結合研究區(qū)地質、地理、水文等條件選取適合的指標??紤]到評價區(qū)的實際情況以及指標數(shù)據(jù)的可獲得性,選取地下水埋深、河網(wǎng)密度、降雨量、地形坡度、含水層厚度、包氣帶介質以及含水層滲透系數(shù)等7個指標來描述評價區(qū)地下水脆弱性。傳統(tǒng)的DRASTIC模型指標還應包括含土壤類型,研究區(qū)內土壤類型與包氣帶介質分布比較一致,該指標對研究區(qū)地下水脆弱性評價并無太大意義,故舍棄。區(qū)內水系發(fā)達,河網(wǎng)密集分布,河網(wǎng)切割侵蝕易影響包氣帶介質,進而影響地下水脆弱性[24],故選用河網(wǎng)密度代替土壤介質。各個指標評分如表2所示。
表2 指標等級劃分以及評分Tab.2 Ranges and ratings for indexes
(2)權重計算。根據(jù)評價區(qū)水文地質條件特點,結合專家意見,確定該研究區(qū)指標對地下水脆弱性影響程度關系為:地下水埋深>降雨量=包氣帶介質>含水層厚度>含水層滲透系數(shù)>河網(wǎng)密度=地形坡度。依據(jù)三標度原則確定標度矩陣(表3),構造傳遞矩陣進而間接求出層次分析法的一致性判斷矩陣(表4),計算出的AHP權重值如表7所示。本次脆弱性評價借助GIS平臺[25],將評價區(qū)劃分為24 801個單元(單元大小為50 m×50 m),把每個單元內包含的7種指標屬性導入到SPSS軟件,包氣帶介質為定性數(shù)據(jù),沒有具體的數(shù)值即選用評分值代替,算得的因子旋轉后的特征值與方差貢獻度如表5所示,選擇特征值大于1且累加貢獻度大于80%的為前3個主因子,得分系數(shù)矩陣如表6所示,則由因子分析法得到的權重如表7所示。將兩種權重值帶入公式(10),得到最終的權重值(表7)。
表3 標度矩陣Tab.3 Scale matrix
表5 特征值與方差貢獻度Tab.5 Eigenvalues and variance contribution rate
表6 得分系數(shù)矩陣Tab.6 Factor score matrix
表7 評價指標權重分配統(tǒng)計表Tab.7 The table of index weight distribution
(3)地下水脆弱性評價。利用得到的綜合權重值,按照公式(1)將各個指標的評分進行加權求和,得到研究區(qū)地下水脆弱性評價值。理論上脆弱性評價值范圍為1~10,按照等間距分為5個級別,即低級別、較低級別、中等級別、較高級別、高級別,評價結果如圖1所示。級別為低的區(qū)域占研究區(qū)面積的18.45%,主要分布在研究區(qū)東側;級別為較低的區(qū)域分布最多,占研究區(qū)的50.64%;中等級別的區(qū)域分布也較多,約占28.29%,主要位于研究區(qū)中部以及少量分布在西側;較高級別的區(qū)域分布較少,約占1.71%,主要在研究區(qū)中心地帶,此區(qū)域地下水埋深淺且含水層厚度??;高級別的區(qū)域分布最少,約占0.91%。研究區(qū)地下水脆弱性整體情況比較良好。此外,為驗證綜合權重法的合理性,本文分別使用層次分析法以及因子分析法得到各指標的權重值并計算出研究區(qū)地下水脆弱性分級進行比較(圖2)。
圖1 研究區(qū)地下水脆弱性評價結果圖(綜合賦權)Fig.1 Groundwater vulnerability assessment of study area(combination)
圖2 研究區(qū)地下水脆弱性評價結果圖Fig.2 Groundwater vulnerability assessment of study area
(4)評價結果驗證與對比。采用單點特征污染物濃度和該點地下水脆弱性指數(shù)之間聯(lián)系的強弱的斯皮爾曼相關系數(shù)ρ進行評價結果驗證[26],如公式(11)所示。
(11)
式中:ui和vi分別為兩變量按等級大小或優(yōu)劣排序后的秩;n為樣本數(shù)量;ρ為斯皮爾曼相關系數(shù)。
評價區(qū)共布置了14個地下水監(jiān)測點(圖1),監(jiān)測點信息見表8。云南省環(huán)境監(jiān)測中心站在2017年4月25日-4月26日進行了地下水采樣,采集了pH、氯化物、硫化鹽、硝酸鹽氮、亞硝酸鹽、氨氮、石油類、揮發(fā)酚、鉛等污染物數(shù)據(jù)信息。氨氮濃度是衡量水體受污染程度的重要指標之一[27],選取氨氮為特征污染物,參照《地下水質量標準》(GBT-14848-2017)中的規(guī)定分級,結合脆弱性評價結果,得到各監(jiān)測點氨氮含量等級及脆弱性等級分類對照表(表9)。根據(jù)公式(11)分別算出3種賦權方法的斯皮爾曼相關系數(shù),結果顯示綜合權重法的相關系數(shù)最高,在95%的置信度下為0.604;因子分析法的相關系數(shù)次之,為0.373;而層次分析法最低,僅為0.308。上述結果表明相比主觀權重法以及客觀權重法而言,綜合權重法的評價結果精度得到提高。并且經(jīng)過綜合權重法后研究區(qū)脆弱性評價等級與氨氮含量呈中度相關性,評價結果具有可靠性。
表8 研究區(qū)地下水監(jiān)測點氨氮數(shù)據(jù)信息Tab.8 Information of surveillance site in study area
本文基于DRASTIC模型就脆弱性評價中評價指標的權重問題展開討論。指標權重大小反映該指標對地下水脆弱性的影響程度,其分配結果直接影響評價結果的合理性,賦權過程需要考慮主觀判斷和客觀信息,分別采用層次分析法和因子分析法進行賦權,并結合兩者得出最終的權重值。通過在云南省安寧市工業(yè)園草鋪片區(qū)的實例應用,驗證了綜合賦權方法的可靠性。本文主要結論如下:
(1)傳統(tǒng)層次分析法不便于使用,故引入三標度準則以及傳遞矩陣,有利于判斷指標間重要程度關系和簡化一致性檢驗,推廣了層次分析法的使用。
表9 各監(jiān)測點氨氮含量等級及脆弱性等級分類表Tab.9 Ratings for ammonia nitrogen and groundwater vulnerability at surveillance site
(2)綜合賦權結果顯示,地下水埋深和降雨量占有較大比重,而河網(wǎng)密度和地形坡度比重較小,符合研究區(qū)實際情況。主觀賦權結果與客觀賦權結果差異最大的3個指標分別是地下水埋深、河網(wǎng)密度以及地形坡度。根據(jù)主觀判斷,地下水埋深應具有最大的權重值而河網(wǎng)密度和地形坡度的權值應最??;然而客觀賦權結果顯示地下水埋深對脆弱性的影響并不大,反而河網(wǎng)密度以及地形坡度的作用得到大幅增長。就脆弱性評價結果而言,主觀賦權的結果略顯平緩,基本處于較低級別,而客觀賦權的結果顯示研究區(qū)地下水脆弱性基本處于較高級別,造成這一現(xiàn)象可能是由于兩種賦權方法給予河網(wǎng)密度以及地形坡度的權重差異造成的??陀^賦權法中脆弱性評價級別較高的區(qū)域里,河網(wǎng)密度值以及坡度值都比較大,加上較高權重的影響因此最終的評價值較高,反觀主觀賦權由于權值的稀釋而造成這些區(qū)域的評價值不高。
(3)以云南省安寧市工業(yè)園草鋪片區(qū)地下水脆弱性評價為例,經(jīng)過斯皮爾曼系數(shù)的驗證,綜合權重法的相關系數(shù)最高,氨氮含量等級與脆弱性等級呈現(xiàn)中度相關性。相比之下,主觀賦權及客觀賦權下的氨氮含量等級與脆弱性等級僅為弱相關,表明在地下水脆弱性研究中綜合權重法具有一定的合理性,有必要結合主觀判斷以及客觀信息綜合考慮。而主觀賦權的相關系數(shù)最低,可見決策時決不能僅憑主觀臆斷,應充分尊重事物的客觀性。
安寧市工業(yè)園草鋪片區(qū)地下水整體情況比較良好,基本處于中等級別及以下,然而情況也不容樂觀。排污量較大的石油煉化產(chǎn)業(yè)和磷鹽化工產(chǎn)業(yè)基本位于中等級別區(qū)域,這些區(qū)域地下水埋深較淺,包氣帶介質吸附凈化能力較差,地下水資源受污染的潛在可能性仍然較大,相關政府部門及企業(yè)應做好防護處理,盡量避免在工業(yè)正常生產(chǎn)條件下,發(fā)生廢料下滲污染地下水的情況。