農(nóng)田土壤重金屬空間變異多尺度研究

劉 偉， 郜允兵， 潘瑜春

[1.中國(guó)礦業(yè)大學(xué)(北京)地球科學(xué)與測(cè)繪工程學(xué)院，北京 100083； 2.國(guó)家農(nóng)業(yè)信息化工程技術(shù)研究中心，北京 100097；3.農(nóng)業(yè)部農(nóng)業(yè)信息技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京 100097； 4.北京市農(nóng)業(yè)物聯(lián)網(wǎng)工程技術(shù)研究中心，北京 100097]

土壤重金屬是表征土壤環(huán)境質(zhì)量的重要因素，如何掌握重金屬含量的空間異質(zhì)性及其空間結(jié)構(gòu)，對(duì)區(qū)域農(nóng)田土壤環(huán)境的監(jiān)控具有重要的生產(chǎn)意義。受錯(cuò)綜復(fù)雜的自然因素和長(zhǎng)期小規(guī)模、破碎化的土地利用生產(chǎn)方式的影響，我國(guó)土壤重金屬空間變異大，空間結(jié)構(gòu)復(fù)雜，當(dāng)前的采樣調(diào)查結(jié)果通常不足以準(zhǔn)確地反映實(shí)際情況，更難用于對(duì)某個(gè)具體地塊的污染范圍劃定和防治指導(dǎo)[1]。這些特征意味著：采樣尺度過大時(shí)，小尺度的空間變異信息常常被減弱或丟失[2]；采樣尺度過小時(shí)，大尺度上的空間變異由于比較微弱而被作為“隨機(jī)成分”處理[3]，同時(shí)也會(huì)導(dǎo)致采樣工作成本高、效率低下。因此，采樣尺度設(shè)計(jì)的合理性是能否精準(zhǔn)刻畫土壤重金屬空間異質(zhì)性的關(guān)鍵[4]。與此同時(shí)，影響農(nóng)田土壤重金屬空間分布的環(huán)境因子的尺度作用范圍存在差異，土壤重金屬亦表現(xiàn)出多尺度效應(yīng)[5]。而采樣問題的復(fù)雜性根本上源于農(nóng)田土壤重金屬含量空間變異的多尺度特征[2]，單純某一尺度上的研究并無法真實(shí)地揭示土壤重金屬的空間結(jié)構(gòu)特征[5-6]，有必要將土壤重金屬含量空間變異的多尺度分析作為重點(diǎn)研究方向，其要點(diǎn)是確定土壤重金屬含量空間異質(zhì)性的特征尺度。國(guó)內(nèi)外研究學(xué)者通過多尺度的采樣和對(duì)比分析，對(duì)土壤重金屬含量空間變異的多尺度結(jié)構(gòu)作諸多的探討[7-10]。然而，這些研究多集中在土壤重金屬在不同采樣尺度上的空間變異規(guī)律，將重金屬的空間異質(zhì)性定量地刻畫到相應(yīng)的尺度水平上的研究則涉及較少，尺度定量方面研究略顯不足，這些不足對(duì)重金屬污染區(qū)域范圍的確定與風(fēng)險(xiǎn)評(píng)估顯然是不利的[11]。

特征尺度是一定區(qū)域內(nèi)能夠表征土壤重金屬空間異質(zhì)性的最優(yōu)空間單元[12]，一旦確定了特征尺度，就可以此確定適宜的采樣間距[13]。本研究以北京市順義區(qū)農(nóng)田土壤重金屬砷(As)采樣點(diǎn)為數(shù)據(jù)源，采用尺度方差法研究該區(qū)域土壤As含量的空間結(jié)構(gòu)特征，確定土壤As含量異質(zhì)性的特征尺度，以期較全面地揭示農(nóng)田土壤As含量的空間多尺度特征，并為特定的采樣需求確定適宜的采樣尺度提供參考依據(jù)。

1 材料與方法

1.1 研究區(qū)概況

順義區(qū)地處北京市城區(qū)東北方向，地理位置為 40°00′～40°18′N、116°28′～116°58′E，總面積 1 019.89 km2。該地區(qū)屬于華北平原北端，北接燕山南麓，全境地勢(shì)北高南低，海拔24～637 m。境內(nèi)平原為潮白河沖積扇下段，為河流洪水?dāng)y帶沉積物質(zhì)造成，表面堆積物主要是沙、亞沙土；土壤類型主要為潮土和褐土，且潮土分布于全區(qū)大部分地區(qū)、褐土主要分布在東北部。順義區(qū)是京郊重要的農(nóng)業(yè)生產(chǎn)區(qū)，據(jù)統(tǒng)計(jì)該區(qū)2004年農(nóng)業(yè)用地面積6.36萬hm2，主要的種植作物為小麥、玉米、蔬菜//水果，冬小麥—夏玉米為主要的糧食種植模式，糧食生產(chǎn)用地占據(jù)大部分農(nóng)業(yè)用地，其中南部是主要的菜地集中區(qū)域，東北部則是果園的聚集地。

1.2 數(shù)據(jù)來源與測(cè)定

土壤樣點(diǎn)數(shù)據(jù)于2007年春季采集，采集耕層(0～20 cm)土壤樣品294個(gè)，采用GPS定位記錄樣點(diǎn)中心位置，采樣點(diǎn)主要分布于糧田、菜地、果園、苗圃等農(nóng)業(yè)用地，采用ArcMap 10.1將原始采樣點(diǎn)數(shù)據(jù)編輯生成樣點(diǎn)分布圖(圖1)。所有土樣在室內(nèi)自然風(fēng)干，碾壓磨碎后，過100目尼龍網(wǎng)篩，采用原子吸收法測(cè)定土壤重金屬As的含量。

1.3 研究方法

1.3.1 尺度方差分析 尺度方差是一種空間等級(jí)分析方法，對(duì)空間變量的多尺度結(jié)構(gòu)比較敏感[14]。該方法是將研究對(duì)象的方差按尺度等級(jí)或尺度嵌套的水平逐步分解， 觀察空間變量的尺度方差隨尺度增加是否會(huì)發(fā)生突變。一般地說，尺度方差發(fā)生突變的尺度也是空間變異性突出的尺度，表征了該等級(jí)水平上的特征尺度，與此同時(shí)尺度方差突變的相對(duì)大小可以反映不同尺度水平空間變異對(duì)系統(tǒng)總體變異的貢獻(xiàn)程度[15]。尺度方差的計(jì)算公式為：

Xijk…z=μ+α1+βij+γijk+…+ωijk…z。

(1)

式中：Xijk…z為等級(jí)系統(tǒng)最低層次上組成單元的值；μ表示當(dāng)前層次上等級(jí)系統(tǒng)基本組成單元的總體平均值；αi、βij、γijk、ωijk…z為不同尺度層次的方差。選擇與變量空間梯度變化相似的尺度劃分方式能比較合理地刻畫空間變量的尺度特征[5]。對(duì)實(shí)測(cè)樣點(diǎn)進(jìn)行趨勢(shì)分析，結(jié)果表明，土壤As含量在東西和南北方向上均有梯度變化，因此采用W-E和N-S等2種尺度劃分方式，設(shè)置最小尺度單元為0.36 km2，共12個(gè)尺度水平(表1)。

表1 尺度方差分析的劃區(qū)方案

注：a=0.6 km。

1.3.2 多環(huán)緩沖區(qū)分析 分別以河流、畜禽養(yǎng)殖基地為中心按照不同緩沖半徑進(jìn)行多環(huán)分析，發(fā)現(xiàn)以河流為中心時(shí)，0.5 km 的緩沖半徑可以使各緩沖區(qū)內(nèi)有適宜數(shù)量的樣點(diǎn)，以畜禽養(yǎng)殖基地為中心時(shí)0.3 km比較合適。因此，分別以 0.5、0.3 km的緩沖半徑對(duì)河流、畜禽養(yǎng)殖基地進(jìn)行多環(huán)緩沖區(qū)分析，將落入同一緩沖區(qū)的土壤樣點(diǎn)歸為一類，并進(jìn)行統(tǒng)計(jì)分析。

1.3.3 數(shù)據(jù)處理 本研究數(shù)據(jù)處理方法及步驟如下：(1)運(yùn)用SPSS 18.0對(duì)土壤重金屬As含量樣點(diǎn)集進(jìn)行基礎(chǔ)統(tǒng)計(jì)分析和正態(tài)分布檢驗(yàn)。(2)使用GS+7.0軟件分別對(duì)As含量樣點(diǎn)進(jìn)行半變異函數(shù)分析和高斯序貫?zāi)M，其中序貫高斯模擬設(shè)置不同的種子(seed)進(jìn)行8次單獨(dú)模擬(G1～G8表示第1次至第8次模擬)。(3)將土壤As含量的模擬數(shù)據(jù)導(dǎo)入ArcGIS生成點(diǎn)位圖后，與農(nóng)田空間分布圖疊置分析，篩選落入農(nóng)田區(qū)域內(nèi)的有效數(shù)據(jù)點(diǎn)，然后采用尺度方差法對(duì)農(nóng)田區(qū)域內(nèi)有效數(shù)據(jù)點(diǎn)進(jìn)行多尺度分析。尺度方差分析通過R語言實(shí)現(xiàn)。(4)采用SPSS進(jìn)行單因素方差分析，研究土壤類型、土壤質(zhì)地、土地利用方式等因子對(duì)研究區(qū)土壤As含量的影響；應(yīng)用ArcGIS 10.1進(jìn)行多環(huán)緩沖區(qū)分析，研究河流和畜禽養(yǎng)殖基地對(duì)研究區(qū)土壤As含量的影響。

2 結(jié)果與分析

2.1 土壤砷實(shí)測(cè)數(shù)據(jù)及條件模擬的統(tǒng)計(jì)特征

由表2可知，順義區(qū)土壤As含量變化范圍為3.85～17.34 mg/kg，變異系數(shù)為27.25%，屬于中等強(qiáng)度變異，與現(xiàn)有研究成果[16]相近。As含量的平均值為7.79 mg/kg，高于北京市土壤背景值7.09 mg/kg，說明順義區(qū)As的含量受到人類活動(dòng)的影響，有一定的累積。自然背景下土壤重金屬含量通常符合正態(tài)分布[17]，由于受外源的影響，As含量的偏態(tài)系數(shù)大于0.5，呈現(xiàn)出正偏態(tài)分布，峰度系數(shù)亦達(dá)到1.97。參照Kolmogorov-Smirnov(K-S)檢驗(yàn)結(jié)果，K-S值為1.22，雙側(cè)顯著性大于0.05，符合正態(tài)分布要求。

由圖2可知，模擬結(jié)果和實(shí)測(cè)數(shù)據(jù)的分布相近，變化范圍相同，較充分地表征了土壤As含量的極值數(shù)據(jù)。模擬數(shù)據(jù)的均值在實(shí)測(cè)均值附近波動(dòng)，這體現(xiàn)了條件模擬對(duì)空間數(shù)據(jù)的概率模擬，同時(shí)也表征出模擬數(shù)據(jù)和實(shí)測(cè)數(shù)據(jù)的分布總體是相同的。

表2 土壤As含量的統(tǒng)計(jì)結(jié)果

注：顯著性水平為0.05，K-S項(xiàng)括號(hào)中為雙側(cè)顯著性。

2.2 土壤砷空間結(jié)構(gòu)分析

在土壤屬性的多尺度分析中，要充分考慮變量的空間結(jié)構(gòu)特征，即采樣點(diǎn)的空間自相關(guān)性和空間變異等級(jí)結(jié)構(gòu)[2]。土壤重金屬As含量半方差函數(shù)參數(shù)如表3所示。采用GS+軟件多次擬合半方差模型，得到Spherical為最佳模型，塊金值C0為 1.853，基臺(tái)值(C0+C)為4.705，基底比C0/(C0+C)為 39.4%，變程為7 km。塊金值表示空間變量受隨機(jī)因素引起的變異，基臺(tái)值是空間變量總的變異?；妆菴0/(C0+C) 表示變量隨機(jī)部分引起的變異與變量總變異的比例。從結(jié)果可知，采樣點(diǎn)基底比為0.394，具有中等空間自相關(guān)性，可進(jìn)行空間變異分析。

表3 As含量半方差函數(shù)參數(shù)

概率累積曲線圖不僅可以判別變量是否符合正態(tài)分布特征(正態(tài)分布情況下概率累積曲線表現(xiàn)為一條直線)，而且通過曲線的突變拐點(diǎn)，某種程度上能夠揭示土壤重金屬含量是否來自不同的樣本總體，從而定性地判斷土壤重金屬含量是否具有等級(jí)結(jié)構(gòu)[2，17]。由土壤重金屬As含量的對(duì)數(shù)概率累計(jì)分布圖(圖3)可知，As含量累計(jì)概率曲線存在明顯的拐點(diǎn)，初步判斷土壤中As含量來自2個(gè)特征明顯不同的總體：A總體含量水平較低，初步定位于自然背景來源；B總體含量水平較高，可能來源于人類活動(dòng)的排放。其次，半方差函數(shù)圖也可以識(shí)別出土壤As含量的空間等級(jí)結(jié)構(gòu)特征，Robertson等通過半方差分析識(shí)別了農(nóng)田土壤中pH值的空間變異等級(jí)結(jié)構(gòu)[18]。由于空間變量的變異可能存在巢式等級(jí)結(jié)構(gòu)，因此該變量的半方差值隨著距離的增加表現(xiàn)為臺(tái)階式上升的趨勢(shì)，而半方差值突變的拐點(diǎn)則刻畫了不同水平上的特征尺度[14]。由圖4可知，2.5 km處為半方差圖突變轉(zhuǎn)折點(diǎn)，說明土壤As含量的空間變異在2.5 km處發(fā)生了突變，表明土壤As含量空間變異存在等級(jí)結(jié)構(gòu)。

2.3 多尺度空間變異識(shí)別

如上所述，研究區(qū)存在空間多尺度結(jié)構(gòu)，因此，本研究采用尺度方差法計(jì)算土壤As含量在不同等級(jí)水平上的方差，并識(shí)別土壤As含量的特征尺度。在Excel導(dǎo)入8次條件模擬在不同尺度水平上的方差(散點(diǎn))，并計(jì)算出8次條件模擬的方差在同一水平上的平均值(折線)，結(jié)果如圖5所示。

由圖5可知，尺度方差隨著尺度的增大而表征出不同的特征，即不同尺度上具有不同的空間異質(zhì)性。整體上尺度方差曲線呈現(xiàn)倒U形，在小尺度上隨尺度增加空間變異性增大，大尺度上則相反。先看W-E劃區(qū)(圖5-A)的尺度方差圖，隨著尺度的增大，尺度方差在等級(jí)5(4a×4a)和等級(jí)9(16a×16a)處有明顯的峰值，說明土壤As含量在等級(jí)5和等級(jí)9的空間異質(zhì)性發(fā)生了突變，因此識(shí)別出特征尺度為 2.4、9.6 km；再看N-S劃區(qū)(圖5-B)的尺度方差圖，尺度方差在等級(jí)3(2a×2a)、等級(jí)5(4a×4a)和等級(jí)9(16a×16a)處表現(xiàn)為峰值，但是等級(jí)3突變程度相對(duì)較弱，特征尺度不明顯，因此只識(shí)別出特征尺度2.4、9.6 km。綜合得到小尺度上特征尺度為2.4 km，與半方差函數(shù)2.5 km拐點(diǎn)相吻合；大尺度上特征尺度為9.6 km，與As含量半方差函數(shù)的變程 7 km 相接近。

2.4 土壤As含量多尺度效應(yīng)成因分析

土壤重金屬空間異質(zhì)性大小的尺度效應(yīng)受控于不同尺度下控制土壤重金屬含量變異的各種生態(tài)過程的重要程度，即區(qū)域內(nèi)不同影響因素對(duì)重金屬As含量變異的作用范圍不同，實(shí)際采樣工作中可以根據(jù)取樣的目的和關(guān)注的影響因素，選擇接近影響因素作用范圍的尺度作為采樣尺度。例如，交通運(yùn)輸污染源對(duì)土壤重金屬Pb含量的影響區(qū)域一般在交通線兩側(cè)幾十米的范圍內(nèi)，影響尺度比較小[19-20]。也就是說，如果某一地區(qū)的主要污染來源于交通運(yùn)輸，那么道路兩側(cè)土壤重金屬含量的空間變異范圍(變程)是在污染源的影響距之內(nèi)的，為了通過采樣來全面地刻畫重金屬的空間變異特征，采樣間距須小于影響距，這一點(diǎn)在柳云龍等的研究[21]中已得到證實(shí)。

本研究前述的統(tǒng)計(jì)分析及半方差分析都顯示As含量存在中等強(qiáng)度空間變異，且As含量存在等級(jí)結(jié)構(gòu)，空間變異受到了土壤自然背景等結(jié)構(gòu)因子和人類活動(dòng)等隨機(jī)因子的共同影響。因此，本研究以土壤類型、土壤質(zhì)地、土地利用、離河邊距、離畜禽養(yǎng)殖基地距離等探討各影響因素對(duì)土壤As含量空間變異的影響范圍。其中，土壤類型、土壤質(zhì)地、土地利用采用單因素方差分析(ANOVA)，結(jié)果顯示土壤質(zhì)地、土地利用對(duì)土壤As含量的空間變異存在顯著影響(P<0.05)；離河邊距、離畜禽養(yǎng)殖基地距離則采用多環(huán)緩沖區(qū)分析。

2.4.1 土壤質(zhì)地 研究區(qū)土壤質(zhì)地類型主要包括輕壤質(zhì)、中壤質(zhì)、沙壤質(zhì)、沙質(zhì)，294個(gè)樣點(diǎn)中有2個(gè)樣點(diǎn)位于其他幾種類型中，為了更好地進(jìn)行方差分析，本研究剔除了這2個(gè)樣點(diǎn)。由表4可知，沙質(zhì)土壤顯著低于其他幾類，且經(jīng)過多重比較分析得出輕壤質(zhì)、中壤質(zhì)差異不顯著，沙壤質(zhì)、沙質(zhì)差異也無顯著差異。因此，將輕壤質(zhì)、中壤質(zhì)歸為一類，沙壤質(zhì)、沙質(zhì)歸一類，并將圖斑合并(圖6)。采用ArcGIS計(jì)算矢量斑塊邊長(zhǎng)，結(jié)果顯示土壤質(zhì)地主要斑塊邊長(zhǎng)范圍4～19 km，斑塊平均邊長(zhǎng)為9.1 km，與As含量空間變異在9.6 km的特征尺度相近，說明在大尺度上土壤質(zhì)地是導(dǎo)致土壤As含量空間變異的主要因素。

表4 不同土壤質(zhì)地中As的含量

注：同列數(shù)據(jù)后不同字母表示差異顯著(P<0.05)。表5同。

2.4.2 土地利用類型 由表5可知，不同土地利用類型的土壤As含量依次為糧田>果園>苗圃>菜地，果園和糧田的As含量顯著高于菜地。因此，將果園、糧田歸為一類，并將相鄰圖斑合并。采用ArcGIS計(jì)算矢量斑塊邊長(zhǎng)，結(jié)果顯示，土地利用類型主要斑塊邊長(zhǎng)范圍為300～4 860 m，而土壤As含量小尺度上的特征尺度為2.5 km，說明在小尺度上土地利用類型是導(dǎo)致土壤As含量空間變異的重要因素之一。

表5 不同土地利用土壤中As的含量

2.4.3 河流 河流沿岸土地由于開發(fā)歷史較長(zhǎng)，土壤環(huán)境受人為的干擾較劇烈，河岸帶表層土壤受重金屬的污染[22]。本研究采用多環(huán)緩沖區(qū)分析法對(duì)樣點(diǎn)進(jìn)行分組統(tǒng)計(jì)以描述采樣點(diǎn)與河流距離的變化關(guān)系(圖7)，整體上隨著樣點(diǎn)與河流距離的增加，As的含量逐漸遞減，擬合曲線在0～2.6 km范圍內(nèi)迅速下降，2.6 km之后逐漸趨于平穩(wěn)?；诳臻g分析的結(jié)果顯示，河流對(duì)土壤As含量的空間分布產(chǎn)生了一定的影響，且控制范圍約為離河流2.6 km的區(qū)域內(nèi)。

2.4.4 畜禽養(yǎng)殖 規(guī)?；?、集約化的畜禽養(yǎng)殖會(huì)產(chǎn)生大量的排泄物(主要包括畜禽糞便)，由于運(yùn)輸距離有限，大量未經(jīng)處理的畜禽糞便主要以有機(jī)肥的形式施用在養(yǎng)殖基地附近的農(nóng)田[23]。有研究表明，華北地區(qū)畜禽糞便中As含量存在超標(biāo)現(xiàn)象[24]。圖8為采樣點(diǎn)與畜禽養(yǎng)殖不同距離范圍內(nèi)土壤As含量變化曲線，總體上隨著距離的增加，重金屬含量逐漸下降，擬合曲線在0～2.2 km范圍內(nèi)下降比較快，2.2 km之后趨于平緩。由此推斷，畜禽糞便主要投放在以畜禽養(yǎng)殖場(chǎng)為中心，方圓2.2 km的范圍內(nèi)，對(duì)土壤As含量空間變異的影響范圍在 2.2 km 左右。

3 結(jié)論與討論

本研究經(jīng)過多尺度分析得到以下結(jié)論：(1)順義區(qū)土壤As含量變化范圍為3.85～17.34 mg/kg，平均值為 7.79 mg/kg，變異系數(shù)為27.25%，屬于中等強(qiáng)度變異。(2)順義區(qū)農(nóng)田土壤重金屬As含量具有較強(qiáng)的空間自相關(guān)性，空間異質(zhì)性存在多尺度結(jié)構(gòu)。(3)順義區(qū)農(nóng)田土壤重金屬As含量空間異質(zhì)性存在2個(gè)特征尺度，分別是2.4、9.6 km，并且不同的特征尺度受控于不同的影響因素，在小尺度上受土地利用、畜禽養(yǎng)殖和河流的影響，大尺度上主要受土壤質(zhì)地的影響。

雖然本研究采用尺度方差研究了順義區(qū)農(nóng)業(yè)土壤中As含量空間變異的尺度效應(yīng)，但是也存在以下不足：(1)本研究采用的基本劃區(qū)尺度為0.6 km，每1級(jí)尺度上推的間距劃分都是下1級(jí)的2倍，這就造成了尺度等級(jí)的不連續(xù)性，這可能會(huì)使部分特征尺度沒有被挖掘出來或者特征尺度所在的尺度等級(jí)無法精確地表達(dá)，這種情況在大尺度上時(shí)會(huì)對(duì)特征尺度識(shí)別十分不利。因此，在接下來的研究可以通過改變基本劃區(qū)尺度，進(jìn)一步多尺度分析，以提高特征尺度識(shí)別的精度。(2)尺度涵蓋粒度和幅度2個(gè)方面，本研究所述尺度僅考慮了空間粒度。為了更加全面地掌握采樣尺度與土壤重金屬空間異質(zhì)性間的相互關(guān)系，后續(xù)的研究將就不同采樣幅度下的情況進(jìn)行探討，當(dāng)幅度充足時(shí)，有利于全面揭示農(nóng)田土壤重金屬含量的空間結(jié)構(gòu)信息。(3)由于相關(guān)數(shù)據(jù)的獲取難度比較大，如大氣沉降對(duì)土壤As含量的影響，本研究對(duì)順義區(qū)As含量影響因素的分析不夠充分，接下來的研究中將就此展開深入的探討。