侯沁言，張世熔*，馬小杰，潘小梅，王貴胤，鐘欽梅

（1.四川農(nóng)業(yè)大學(xué)環(huán)境學(xué)院，成都 611130；2.四川省土壤環(huán)境保護(hù)重點實驗室，成都 611130）

近年來，由于人類活動強度的提高，農(nóng)田土壤重金屬污染已成為全球關(guān)注的環(huán)境問題之一[1-2]。在污染土壤的重金屬元素中，Cd、Pb、Hg、As、Cr、Zn、Cu、Ni等8種元素的研究報道較多[3-6]。值得關(guān)注的是，即使在較低含量水平條件下，土壤Cd、Pb、Hg和As也不僅會影響農(nóng)產(chǎn)品的產(chǎn)量和品質(zhì)，還會通過食物鏈傳遞，對人體健康造成潛在威脅[5]。作為我國糧食、油料等作物主要生產(chǎn)區(qū)之一的四川盆地丘陵地區(qū)，其農(nóng)田土壤重金屬的研究已有少量報道[7-9]。例如，在尚二萍等[7]的研究結(jié)果中，四川盆地耕地點位超標(biāo)率為43.55%，而嚴(yán)明書等[8]和劉潔等[9]對重慶都市經(jīng)濟圈土壤及沱江中游土壤質(zhì)量的研究表明，其點位重金屬超標(biāo)率分別為13.59%和0.6%。盡管樣點超標(biāo)率可以描述土壤重金屬的污染狀況，但不同研究結(jié)果較大的空間異質(zhì)性也增大了準(zhǔn)確評估區(qū)域尺度重金屬污染分布特征的難度。因此，進(jìn)一步開展該區(qū)域土壤重金屬空間分布及污染評價研究，為區(qū)域土壤環(huán)境質(zhì)量管理提供參考依據(jù)是十分必要的。

關(guān)于區(qū)域土壤重金屬污染評價方法，傳統(tǒng)上多采用基于樣點的單因子污染指數(shù)法（Pi）、內(nèi)梅羅綜合指數(shù)法（PN）、地累積指數(shù)法（Igeo）和潛在生態(tài)風(fēng)險指數(shù)法（RI）等方法進(jìn)行評價[10-12]。然而，重金屬在土壤中的移動性較差，基于樣點的指數(shù)法難于提供其污染空間分布的真實信息[13]。特別是在丘陵地區(qū)，重金屬在同一農(nóng)田地塊內(nèi)及不同地塊間的遷移易受地形的限制。隨著丘陵地區(qū)土壤重金屬污染防治進(jìn)入規(guī)劃與工程實施階段，提供更準(zhǔn)確的污染空間分布信息成為必然要求。因此，近年來一些研究者基于GIS、GPS和地統(tǒng)計學(xué)等空間分析方法研究了區(qū)域土壤重金屬的空間分布特征[9，14-16]。同時，也有研究者將這些空間分析法與指數(shù)法結(jié)合對其進(jìn)行評價[17]。然而，迄今基于點位的污染超標(biāo)率與基于空間分布的超標(biāo)率之間的關(guān)系研究仍少見報道。此外，在相同的研究尺度（丘陵區(qū)）、采樣方案及軟件平臺下，系統(tǒng)性地對比不同空間插值方法對農(nóng)田土壤重金屬超標(biāo)率的影響仍鮮見報道。

凱江流域為四川盆地丘陵區(qū)涪江流域的一級支流，是四川盆地典型的糧油及蔬菜生產(chǎn)基地之一。流域內(nèi)種植業(yè)和畜禽養(yǎng)殖業(yè)發(fā)達(dá)，人口密度高達(dá)500～700人·km-2。近年來，隨著農(nóng)業(yè)集約化程度的提高和鄉(xiāng)鎮(zhèn)企業(yè)的發(fā)展，其農(nóng)田土壤重金屬的含量水平及污染狀況尚不清楚。因此，本文根據(jù)我國《土壤環(huán)境質(zhì)量農(nóng)用地土壤污染風(fēng)險管控標(biāo)準(zhǔn)（試行）》（GB 15618—2018）[18]，采用單因子指數(shù)法，在GIS和GPS軟件平臺上，系統(tǒng)對比研究反比距離權(quán)重法（Inverse distance weighted，IDW）、徑向基函數(shù)法（Radial basis function，RBF）、普通克里格法（Ordinary kriging，OK）、析取克里格法（Disjunctive kriging，DK）和概率克里格法（Probability kriging，PK）5種空間插值方法對該區(qū)域農(nóng)田土壤重金屬污染評價結(jié)果的影響，旨在為區(qū)域土壤污染治理提供科學(xué)參考。

1 材料與方法

1.1 研究區(qū)域概況

研究區(qū)域地處東經(jīng) 104°40′～105°6′和北緯 30°56′～31°4′，面積約 471 km2，屬亞熱帶濕潤性季風(fēng)性氣候，年平均氣溫15℃，年均降水量800 mm，具有氣候溫和、四季分明、雨水充沛等特點。研究河段處于四川盆地丘陵低山地區(qū)，其土地利用方式主要為耕地，以種植水稻、小麥、棉花和油菜等為主。成土母質(zhì)以侏羅系蓬萊組紫色砂泥巖、白堊系城墻巖群黃色砂巖、厚砂薄頁巖風(fēng)化坡殘積物為主，凱江主流及支流沿岸有部分地區(qū)為第四系新沖積物。

1.2 土壤樣品采集與分析方法

根據(jù)該區(qū)地形地貌、成土母質(zhì)、土壤類型、土地利用等條件，按隨機網(wǎng)格法和分層抽樣法進(jìn)行樣點的布設(shè)，共187個采樣點及20個檢查點（圖1）。在每個樣點按五點采樣法采集同一地塊0～20 cm表層土樣，混合均勻后用四分法保留土樣約1 kg；樣點位置采用GPS進(jìn)行定位，并調(diào)查記錄地形、母質(zhì)、土壤類型等自然環(huán)境條件，以及作物、灌溉、施肥等社會經(jīng)濟信息。

圖1 樣點空間分布圖Figure 1 Distribution of soil sampling points

土壤樣品于實驗室自然風(fēng)干后，去掉礫石、動植物殘骸等，進(jìn)行研磨，再分別過20目、100目篩，待測。土壤pH值在水土比為2.5∶1條件下采用復(fù)合電極法進(jìn)行測定；土壤Cd、Pb的總量用HNO3-HF-HClO4體系（V∶V∶V=1∶2∶2）消解，消解液冷卻后，加入5 mL 2 mol·L-1HCl溶解，并轉(zhuǎn)移至25 mL容量瓶，定容、過濾。待測液采用石墨爐原子吸收分光光度法測定；土壤Hg和As總量采用HCl與HNO3混合液體系（V∶V=1∶1）消解，原子熒光法測定[19]。在分析過程中所用試劑均為分析純，所有土壤樣品中每個重金屬全量測定時均設(shè)置3個重復(fù)，并設(shè)置空白樣，檢測過程中采用國家標(biāo)準(zhǔn)土壤樣品（GSS-3）進(jìn)行質(zhì)量控制，Cd、Hg、As和Pb的平均回收率分別為106.8%、92.5%、96.8%和102.3%。

1.3 單因子污染指數(shù)法

單因子污染指數(shù)法作為最簡單的評價方法之一，通過實際含量和評價標(biāo)準(zhǔn)的比值直接表現(xiàn)出重金屬的污染級別及危害程度。其計算公式[7]為：

式中：Pi為第i種土壤重金屬的環(huán)境質(zhì)量指數(shù)；Ci為土樣中第i種重金屬的實測含量，mg·kg-1；Si為第i種重金屬的農(nóng)用地土壤污染風(fēng)險篩選值（試行）（GB 15618—2018）[18]，mg·kg-1。本文評價分為兩個等級：當(dāng) Pi＜1時，土樣重金屬含量低于篩選值，樣點基本未受到污染；當(dāng)Pi≥1時，土壤樣點受到污染，指數(shù)值越高，污染的程度越高。

1.4空間插值方法

本研究基于樣點GPS坐標(biāo)和ArcGIS 10.2地理信息系統(tǒng)軟件平臺，采用地統(tǒng)計學(xué)的普通克里格法進(jìn)行空間插值獲得區(qū)域土壤4種重金屬含量的空間分布圖。同時，為系統(tǒng)對比研究不同插值方法對空間分布的影響，在相同的樣點GPS坐標(biāo)和ArcGIS 10.2平臺上，分別采用反比距離權(quán)重法和徑向基函數(shù)法兩種確定性方法，以及普通克里格法和析取克里格法兩種地統(tǒng)計學(xué)方法對單因子指數(shù)Pi進(jìn)行空間插值[20]，以研究不同空間插值方法對該區(qū)域農(nóng)田土壤重金屬污染評價結(jié)果的影響。其中，基于地統(tǒng)計學(xué)的普通克里格法和析取克里格法均要求數(shù)據(jù)為正態(tài)分布。此外，本研究還采用概率克里格法[20]對單因子指數(shù)Pi進(jìn)行了不同污染水平的風(fēng)險分析。最后，根據(jù)以上5種空間插值方法獲得的信息，綜合討論該區(qū)域農(nóng)田土壤重金屬污染超標(biāo)狀況。

1.5 數(shù)據(jù)處理方法

本研究所有數(shù)據(jù)統(tǒng)計均在SPSS 19.0平臺上進(jìn)行；樣點分布圖和各類空間插值圖利用ArcGIS 10.2軟件繪制。

2 結(jié)果與討論

2.1 土壤中重金屬含量統(tǒng)計

研究區(qū)土壤樣點pH范圍為4.90～8.25。土壤中Cd、Hg、As和Pb 4種重金屬均呈正態(tài)分布（表1），平均含量分別為 0.38、0.10、7.77 mg·kg-1和 22.54 mg·kg-1。其中，土壤Cd、Hg、As和Pb的平均含量分別為當(dāng)?shù)赝寥拉h(huán)境背景值的1.58、1.43、1.27倍和1.11倍，表明該區(qū)域存在不同程度的重金屬累積情況。根據(jù)GB 15618—2018[18]，4種重金屬元素的平均值均低于其風(fēng)險篩選值，但有部分土壤樣點的Cd和As含量超過其污染風(fēng)險篩選值。變異系數(shù)（CV）相較于標(biāo)準(zhǔn)差更能反映出數(shù)據(jù)的離散性[21]，土壤Cd、Hg、As、Pb的變異系數(shù)分別為50%、40%、29%和13%，均處于中度變異程度（10%＜CV＜100%）[13，22]，表明這幾種元素數(shù)據(jù)的離散性較小，分布相對較為均勻。

表1 土壤重金屬含量統(tǒng)計分析Table 1 Statistical analysis of soil heavy metals

2.2 土壤重金屬含量空間特征分析

由于該區(qū)域土壤Cd、Hg、As和Pb含量均符合正態(tài)分布，因此本研究在ArcGIS 10.2平臺上選擇地統(tǒng)計學(xué)模塊中的普通克里格法進(jìn)行空間插值，獲得它們含量的空間分布圖（圖2）。研究區(qū)土壤Cd（圖2a）含量在西北角和凱河一帶有兩個高值中心（＞0.8 mg·kg-1），同時呈現(xiàn)向四周遞減的趨勢。研究區(qū)西北角靠近縣城，工業(yè)活動頻繁；同時，在凱河附近曾有一養(yǎng)豬場，其養(yǎng)殖廢棄物可能是造成污染的主要原因。土壤Hg的空間分布（圖2b）呈現(xiàn)以南華北部為高值區(qū)（＞0.17 mg·kg-1），向東、向北逐漸減少到 0.06～0.10 mg·kg-1的低值區(qū)，并在東南部最低（＜0.06 mg·kg-1）。由于研究區(qū)土壤Hg平均含量較低，因此該分布可能與土壤母質(zhì)因素有較大關(guān)系。土壤As的空間分布（圖2c）呈現(xiàn)以南華和萬安的兩個高值中心（＞11 mg·kg-1），向四周遞減的趨勢。研究區(qū)內(nèi)土壤Pb（圖2d）的高值中心（＞25 mg·kg-1）位于南華、凱河處，并向周圍呈現(xiàn)遞減趨勢，其東北及西北部屬于低值區(qū)（＜20 mg·kg-1）。有研究表明，土壤中As和Pb主要源于交通污染和大氣沉降[23-24]。研究區(qū)土壤As和Pb的高值區(qū)均處于交通較好的區(qū)域，同時在南華附近存在一個垃圾焚燒廠，這也可能是造成高值的原因。

圖2 凱江流域土壤重金屬空間分布圖Figure 2 Spatial distribution of soil heavy metal in the Kaijiang Valley

2.3 土壤重金屬污染評價

2.3.1 基于樣點的單因子污染指數(shù)評價

研究區(qū)內(nèi)土壤Cd、Hg、As和Pb的單項污染指數(shù)的平均值分別為0.487、0.111、0.381和0.102（表2），其中僅土壤Cd和As存在Pi≥1的超過風(fēng)險篩選值的樣點，其占比分別為7.5%和0.5%。因此，Cd是研究區(qū)農(nóng)田土壤中主要的重金屬污染物。

2.3.2 基于空間的單因子污染指數(shù)評價

為進(jìn)一步探究研究區(qū)重金屬污染指數(shù)的空間分布特征，分別采用反比距離權(quán)重法、徑向基函數(shù)法、普通克里格法和析取克里格法對4種重金屬的Pi進(jìn)行空間插值，并分級統(tǒng)計其面積占總面積的百分比（表3）。

總體上4種插值方法的平均誤差（ME）、標(biāo)準(zhǔn)化平均誤差（MSE）和20個校驗點的檢驗誤差（MEP）均較?。ū?），表明它們用于插值是可行的，但普通克里格法的誤差相對較小。從土壤重金屬種類分析，Hg和Pb的4種空間插值分析方法的統(tǒng)計結(jié)果完全相同，全區(qū)域均為Pi＜0.4等級的面積分布。因此，這兩種土壤重金屬在研究區(qū)內(nèi)無污染風(fēng)險。然而，不同空間插值方法所得的土壤Cd和As的面積統(tǒng)計結(jié)果卻呈現(xiàn)一定差異（表3）。

土壤 As僅存在 Pi＜0.4、0.4≤Pi＜0.6 和 0.6≤Pi＜0.8這3個污染指數(shù)等級的空間分布。因此，盡管這4種空間插值方法所得的分布等級相同，但在前3個等級各自分布的面積卻呈現(xiàn)一定差異。例如，在0.6≤Pi＜0.8等級的分布中，占區(qū)域面積百分?jǐn)?shù)由小到大依次為析取克里格法面積（僅0.07%）＜徑向基函數(shù)法（0.10%）＜普通克里格法（0.67%）＜反比距離權(quán)重法（0.89%），與基于樣點污染指數(shù)的評價結(jié)果（2.7%）不同。此外，4種空間插值方法皆因僅有1個超標(biāo)樣點而無法形成超標(biāo)等級的面積分布。

表2 基于樣點的單因子指數(shù)評價Table 2 Assessment of single factor indexes based on sampling points

表3 基于不同插值方法的單因子指數(shù)分級面積比Table 3 Area ratios of single factor index based on the different interpolation methods

土壤Cd污染指數(shù)的空間分布則因插值方法不同呈現(xiàn)明顯差異。但無論哪一種空間插值方法，其超標(biāo)區(qū)域占研究區(qū)總面積的百分比（0.14%～1.06%）均明顯低于基于樣點統(tǒng)計的百分比（7.5%）。

為進(jìn)一步觀察4種空間插值方法對其污染指數(shù)，特別是超標(biāo)指數(shù)空間分布的影響，本文選取在土壤質(zhì)量環(huán)境評價研究中采用最頻繁的反比距離權(quán)重法[25]和普通克里格法[9]作圖，進(jìn)行輔助分析（圖3和圖4）。對比圖3、圖4和表3可知，反比距離權(quán)重法和徑向基函數(shù)法中土壤Cd的Pi≥1超標(biāo)等級面積（1.06%和0.67%）明顯高于普通克里格法（0.42%）和析取克里格法（0.14%），這是由于算法公式上的差異所致[26-27]。前兩種方法屬確定性插值方法，均突出高值點對其鄰近空間區(qū)域影響，其權(quán)重更大，所以插值土壤Cd的超標(biāo)面積偏大（圖3a）；反之，基于地統(tǒng)計學(xué)的普通克里格法和析取克里格法并不強調(diào)高值點對其鄰近空間區(qū)域的影響權(quán)重，且不同含量樣點取等權(quán)重，因而插值結(jié)果具有一定的平滑效應(yīng)[20，28]。特別是析取克里格法，其平滑效應(yīng)明顯，導(dǎo)致其超標(biāo)面積僅為反比距離權(quán)重法和普通克里格法的13.2%和33.3%。

圖3 基于反比距離權(quán)重法的單因子指數(shù)空間分布Figure 3 Spatial distribution of single factor index for soil heavy metals by inverse distance weighting method

2.4 單因子指數(shù)的風(fēng)險概率評價

與前面4種研究單因子污染指數(shù)值空間分布特征的插值方法不同，概率克里格法[20]能基于污染指數(shù)的某一閾值水平來觀察研究區(qū)超過該閾值水平概率的空間分布特征[29]。毛竹等[30]利用該方法分析了漢源富泉鉛鋅礦區(qū)土壤重金屬的概率空間分布，預(yù)測出研究區(qū)各區(qū)域重金屬的超標(biāo)概率。本研究基于土壤Cd、Hg、As、Pb的含量水平及單因子污染指數(shù)統(tǒng)計結(jié)果（表2），通過預(yù)測概率分析，將Cd和As的閾值設(shè)定為1.0，而Hg和Pb的最大閾值分別為0.3和0.2，以超過這些閾值為條件，用概率克里格法依次進(jìn)行空間插值獲得它們的風(fēng)險概率分布圖（圖5）。其中，Hg在Pi＞0.3及Pb在Pi＞0.2條件下均無高風(fēng)險概率，因此整個研究區(qū)農(nóng)田土壤受Hg和Pb污染的概率很低，該高值區(qū)域評價結(jié)果與前面基于樣點或空間的單因子指數(shù)評價結(jié)果均一致。

土壤Cd與As污染超標(biāo)（Pi＞1）的風(fēng)險概率區(qū)域占整個研究區(qū)面積的比例則與前面幾種評價不完全一致。其中，在萬安附近，有一小塊農(nóng)田土壤As的超標(biāo)概率達(dá)到0.8～1.0的水平（圖5c，占研究區(qū)面積的1.05%），即該區(qū)域土壤有80%的概率受到As污染。其評價結(jié)果與基于樣點的單因子污染指數(shù)評價相似，而未像基于空間的前4種插值方法一樣漏掉有1個樣點土壤As超標(biāo)的信息。同時，土壤Cd存在兩個概率為1.0的高值區(qū)域，其空間分布與Pi值的超標(biāo)區(qū)域空間分布特征相似，但面積約為區(qū)域面積的0.16%，更接近基于空間的析取克里格法的評價結(jié)果。

3 結(jié)論

（1）研究區(qū)土壤Cd、Hg、As和Pb均呈正態(tài)分布，平均含量依次是當(dāng)?shù)赝寥拉h(huán)境背景值的1.58、1.43、1.27倍和1.11倍，存在不同程度的重金屬累積情況，其中有部分土壤點位Cd和As含量超過GB 15618—2018中的風(fēng)險篩選值。

（2）采用普通克里格法獲得區(qū)域土壤中4種重金屬含量的空間分布圖。研究區(qū)內(nèi)土壤Cd、As和Pb的含量總體上呈現(xiàn)以兩個高值中心向四周遞減的趨勢；土壤Hg則總體呈現(xiàn)西部高、東部低的趨勢。

（3）基于樣點的單因子指數(shù)評價得出土壤污染程度為 Cd＞As＞Hg＞Pb?；诳臻g的面積統(tǒng)計結(jié)果表明，研究區(qū)僅土壤Cd存在超標(biāo)面積。4種方法中，普通克里格插值方法的誤差最小，更接近于真實情況。

（4）概率克里格空間插值表明研究區(qū)僅有土壤Cd和As存在高污染風(fēng)險概率區(qū)域，其空間分布特征與Pi的高值區(qū)域空間分布特征基本一致。

圖4 基于普通克里格法的土壤重金屬單因子指數(shù)空間分布Figure 4 Spatial distribution of single factor index for soil heavy metals by ordinary kriging method

圖5 基于概率克里格法的土壤重金屬單因子指數(shù)概率空間分布Figure 5 Probability distribution of single factor index for soil heavy metals by probability