礦山土壤重金屬含量估算研究

饒振興，馬超振，邵鵬遠(yuǎn)，翟會君，顧展飛

(1.河南省地質(zhì)礦產(chǎn)勘查開發(fā)局 第四地質(zhì)勘查院，河南 鄭州 450001； 2.河南省地質(zhì)礦產(chǎn)勘查開發(fā)局 第二地質(zhì)環(huán)境調(diào)查院，河南 鄭州 450001； 3.鄭州航空工業(yè)管理學(xué)院，河南 鄭州 450001)

土壤資源關(guān)于農(nóng)業(yè)生產(chǎn)安全，并直接影響人類的身體健康。近年來，因煤礦開采而導(dǎo)致的土壤污染問題逐步受到人們的關(guān)注。煤礦開采雖能促進(jìn)社會和經(jīng)濟(jì)的發(fā)展[1]，但不合理的開發(fā)將會造成一系列的環(huán)境問題，包括水資源污染、地面塌陷、溫室氣體排放污染以及噪聲污染等[2-5]，嚴(yán)重影響生態(tài)環(huán)境。其中，煤矸石等廢料的長期堆放會導(dǎo)致土壤中的重金屬異常富集[6]，再經(jīng)由農(nóng)作物的吸收間接進(jìn)入人體內(nèi)，危及人類健康[7]。因此，高效地估算土壤中的重金屬含量，尤其是在礦山環(huán)境條件下，意義重大。

傳統(tǒng)的估算方法一般是通過實驗室的精密分析，但成本高昂且操作效率低[8-9]。隨著科技水平的發(fā)展，遙感技術(shù)在高光譜領(lǐng)域發(fā)展迅速，為土壤污染物的高效估算提供了可能[10-12]，但基于高光譜的重金屬含量反演推斷技術(shù)仍大多以實驗室的光譜數(shù)據(jù)分析為基礎(chǔ)，一定程度上限制了其使用效率。為此，已有不少專家針對重金屬含量和野外土壤高光譜的關(guān)系，嘗試建立直接的映射函數(shù)，并取得了部分成果[13-14]，如陶超等[15]采用支持向量機(jī)(SVM)作為分類器，將一定數(shù)量的室內(nèi)樣本數(shù)據(jù)用于分類器訓(xùn)練，同時對野外樣本數(shù)據(jù)進(jìn)行分類測試，通過光譜變換后進(jìn)行了重金屬污染定性分類模型遷移試驗，結(jié)果表明不同試驗區(qū)域的砷(As)、鉛(Pb)和鋅(Zn)等重金屬污染狀況可直接利用室外光譜采樣進(jìn)行研究。鄒濱等[16]依托湖南衡陽市某礦區(qū)，針對土壤鎘(Cd)含量測定，開展了實驗室與野外原位DS(direct standardization，直接校正法)光譜聯(lián)合反演模型的有效性驗證。試驗結(jié)果表明，將野外原位光譜DS轉(zhuǎn)換校正結(jié)合融合實驗室光譜先驗知識的模型精度最高，其決定系數(shù)R2達(dá)到了0.72，對野外原位光譜模型的可靠性提升具有重要作用。張霞等[17]針對野外高光譜實時監(jiān)測土壤重金屬含量過程中存在的光譜數(shù)據(jù)冗余問題，基于雄安新區(qū)的野外農(nóng)田土壤樣本，分別從反演機(jī)理和波段選擇算法2個角度出發(fā)，提出了相應(yīng)的反演方法，實現(xiàn)了土壤樣本中鉛含量估算精度的有效提升。綜上所述，在野外土壤高光譜與重金屬含量估算的模型關(guān)系方面，較多研究集中于可行性或有效性的驗證以及存在問題的完善方面，缺乏基于偏最小二乘模型的高光譜估算模型研究及其具體應(yīng)用方法的介紹。因此，本文以某城市煤礦環(huán)境下的周邊土地為研究基礎(chǔ)，為將土壤中隱蔽性較強(qiáng)的重金屬污染物含量進(jìn)行高效估算與表征，建立了偏最小二乘模型的高光譜估算模型，其中偏最小二乘模型的建立基礎(chǔ)為土壤反射光譜數(shù)據(jù)與光譜變換指標(biāo)，從而突破了上述基于實驗室數(shù)據(jù)進(jìn)行重金屬含量估算的局限性。

1 礦區(qū)概況

1.1 基本情況

某城市北部礦區(qū)內(nèi)土壤成土類型較多，土壤pH值在6.5～7.5。第四系松散沉積物在該城市北部礦區(qū)內(nèi)具有普遍的分布，厚度達(dá)50 m左右，部分區(qū)域超過100 m。砂壤土和砂土為區(qū)域內(nèi)主要的表層土質(zhì)，具有良好的耕性。土壤中鉀、磷、氮等含量較為豐富，整體評價為肥沃。該礦區(qū)經(jīng)過長期的煤礦開采，對周邊土地環(huán)境影響巨大，原有土地利用結(jié)構(gòu)遭到破壞，可耕地面積減少，同時大量的采空區(qū)形成改變了原有地貌，使農(nóng)田出現(xiàn)下沉。與此同時，該礦區(qū)的周邊土壤不可避免地遭到了不同程度的污染，土壤重金屬富集污染較嚴(yán)重，造成周邊居民身體健康的潛在隱患。

1.2 主要治理對策及原理

為更有針對性地采取土壤污染治理對策，實現(xiàn)可持續(xù)發(fā)展，需對土壤重金屬含量進(jìn)行估算和評價。為此，提出并綜合應(yīng)用了一種根據(jù)現(xiàn)場高光譜數(shù)據(jù)實現(xiàn)反演模型建立的方案，從而準(zhǔn)確地估算出土壤中的重金屬含量。其基本原理為：土壤中的重金屬含量較低，其對于光譜的響應(yīng)來說是非常微弱的，即無法通過光譜表面直接得出土壤重金屬的光譜反應(yīng)。但重金屬與土壤中的黏土或鐵氧化物等較易聚集在一起，且光譜響應(yīng)對此理化性質(zhì)比較敏感，故重金屬含量可先與土壤中的鐵氧化物建立一定的函數(shù)關(guān)系，再通過上述理化性質(zhì)建立高光譜模型，從而間接進(jìn)行重金屬含量的估算[18-20]。同時，在研究方法方面，為解決高光譜波段數(shù)據(jù)量較大、不同波段間相關(guān)性強(qiáng)等問題，擯棄傳統(tǒng)的逐步回歸和主成分方法，采用了偏最小二乘方法，進(jìn)一步完善了高光譜估算模型。

2 數(shù)據(jù)采集及預(yù)處理

2.1 土壤樣本制備

基于格網(wǎng)GIS技術(shù)，結(jié)合礦區(qū)實地情況，按1.5 km2范圍采集1個樣點的原則進(jìn)行控制，根據(jù)布點方案采集土壤樣本36個，如圖1所示。樣本采集時，盡量選擇土壤粒度及含水量相近的區(qū)域，以避免土壤環(huán)境對光譜的影響。

圖1 礦區(qū)土壤采樣點示意Fig.1 Schematic diagram of soil sampling points in mining area

2.2 野外光譜測定

采用手持地物光譜儀進(jìn)行本次野外光譜測定。測定時，選取太陽輻射強(qiáng)烈、晴朗無云的中午時段，風(fēng)力不大于3級且空氣能見度在10 km以上。操作人員宜穿深色衣服，操作時保持與測量的樣本相距約500 mm，同時，光譜儀跟水平面的法線的夾角需控制在一定范圍內(nèi)，該值一般為10°。采集時，每10 min定標(biāo)1次。為便于后期取平均數(shù)，在每一個采樣點處，均需要進(jìn)行10條光譜數(shù)據(jù)的測量。

2.3 土壤理化性質(zhì)測試

采集的土壤樣本作通風(fēng)陰干處理，經(jīng)壓碎研磨后，過篩備用。在試驗過程中，擬通過電感耦合等離子體原子發(fā)射光譜法對土壤的理化性質(zhì)進(jìn)行檢測，主要測定土壤中Cr、Mn、Cu、Fe、Ni、Cd、Zn、As、Pb等重金屬的含量。該方法在測定時，先進(jìn)行微波消解，冷卻后加超純水搖勻，再按提前確定好的參數(shù)檢測樣本里元素的發(fā)射強(qiáng)度，并基于校準(zhǔn)函數(shù)計算出重金屬元素的含量。該方法在操作時，不僅快速精準(zhǔn)，而且可同時對多種元素進(jìn)行測定，另外還具有線性范圍寬的優(yōu)勢，應(yīng)用廣泛。

2.4 其他相關(guān)數(shù)據(jù)采集

在土壤的野外高光譜數(shù)據(jù)采集前，采用影像融合方法，將2.5 m全色波段數(shù)據(jù)與Spot-5衛(wèi)星中10 m多光譜數(shù)據(jù)進(jìn)行融合，由此獲得礦區(qū)2.5 m高清彩色影像。同時，由于人為耕種與施肥等因素對土壤的理化性質(zhì)存在較大影響。因此，搜集了礦區(qū)不同采樣點的土地利用類型。另外，此次試驗還收集了該地區(qū)土壤類別空間分布數(shù)據(jù)，以更好地考慮重金屬元素臨界含量與容量在不同土壤類別的差異。

2.5 高光譜數(shù)據(jù)預(yù)處理

土壤光譜數(shù)據(jù)由于是野外采集的，故需要在模型擬合以及光譜指標(biāo)變換前，對相關(guān)數(shù)據(jù)進(jìn)行平滑降噪、剔除不合格光譜等預(yù)處理，并進(jìn)行均值化。

3 土壤重金屬污染估算

3.1 高光譜估算模型建立方案

在模型建立時引入偏最小二乘回歸模型，相對于傳統(tǒng)的多元線性等回歸模型，其不僅吸收了主成分分析方法中提取主成分的優(yōu)勢，在進(jìn)行回歸建模時綜合了自變量信息，而且解決了變量中多重相關(guān)的問題。在完成建模后，為更好地檢驗實測值與預(yù)測值的關(guān)系，評估方法的準(zhǔn)確性，提出了分析方法中常用的相關(guān)系數(shù)、均方根誤差等指標(biāo)參數(shù)。

3.2 提取高光譜指標(biāo)信息

在獲取土壤光譜后，其一般含有繁雜的各種信息，難以分清光譜波段的相關(guān)關(guān)系。為更高效地實現(xiàn)有效數(shù)據(jù)獲取，對土壤高光譜特征進(jìn)行更準(zhǔn)確的表征，基于如圖2所示的高光譜反射率趨勢，提取出4個關(guān)鍵的指標(biāo)，具體涵蓋有去除包絡(luò)線、反射率倒數(shù)對數(shù)以及一階、二階微分指標(biāo)。

圖2 土壤的反射率光譜走勢示意Fig.2 Schematic diagram of reflectance spectrum trend of soil

(1)去除包絡(luò)線。每個土樣的高光譜曲線整體發(fā)展趨勢類似，為有效體現(xiàn)土樣波譜的反射峰和吸收峰，使之歸一化到統(tǒng)一的背景值中(0～1)，放大微小的波譜變化，可進(jìn)行去除包絡(luò)線處理，如圖3所示。

圖3 土壤反射率光譜去除包絡(luò)線曲線 Fig.3 Removal envelope curve of soil reflectance spectrum

(2)為增強(qiáng)光譜曲線在坡度上的響應(yīng)，消除空氣中氣體對光譜的影響，并在一定程度上突出對可見光波段的光譜波段的差異性，可進(jìn)行反射光譜的倒數(shù)對數(shù)變換，如圖4所示。

圖4 土壤反射率光譜倒數(shù)對數(shù)曲線Fig.4 Reciprocal logarithm curve of soil reflectance spectrum

(3)受儀器自身的影響，光譜儀在使用中所采集的光譜基線會隨著時間而出現(xiàn)漂移，不利于后期偏最小二乘法的建模。為此，可采取光譜微分處理，不僅可以去除基線漂移現(xiàn)象，還能提取到光譜曲線中較為隱蔽的信息。同時，更高的波譜分辨率可以在一階和二階微分處理后的反射光譜中實現(xiàn)，如圖5所示。

圖5 土壤反射率光譜一階微分和二階微分曲線Fig.5 Soil reflectance spectrum first and second order differential curves

3.3 對重金屬的相關(guān)性分析

根據(jù)高光譜估算理論原理，為滿足高光譜建模條件，需對土壤中各重金屬元素含量與Fe、Mn含量的相關(guān)性進(jìn)行判定。對于重金屬元素的相關(guān)系數(shù)，一般來說，當(dāng)其值大于0.5，同時呈現(xiàn)為顯著相關(guān)的，則表明兩者呈強(qiáng)顯著相關(guān)性，據(jù)此可以實現(xiàn)高光譜的建模分析。本次土壤中幾種重金屬元素含量的相關(guān)系數(shù)矩陣見表1，可知，F(xiàn)e、Mn的相關(guān)系數(shù)為0.94，體現(xiàn)了兩者的共生關(guān)系。Ni、As、Cr、Pb、Zn等與Fe的相關(guān)系數(shù)均在0.5以上，說明均呈強(qiáng)顯著相關(guān)。而Cu與Fe、Mn的雖呈顯著相關(guān)，但相關(guān)系數(shù)未超過0.5；Cd與Fe既未呈顯著相關(guān)，相關(guān)系數(shù)也較小。因此，在上述重金屬元素中，不滿足高光譜建模要求的只有Cu與Cd元素。

表1 重金屬相關(guān)系數(shù)矩陣Tab.1 Correlation coefficient matrix of heavy metal

3.4 提取光譜特征波段

野外現(xiàn)場測定試驗中，光譜分辨率設(shè)定為1 nm，通過測定獲得的光譜波長區(qū)間在350～1 050 nm。高分辨率在提高測量數(shù)據(jù)精度的同時，也帶來了非常多的無效信息，導(dǎo)致后續(xù)建模困難。因此，可將光譜分辨率增大十倍，即從原先的1 nm變更至10 nm，達(dá)到簡化模型、提升工效的目的。在這之后，接著計算樣本重金屬含量與上述4個光譜指標(biāo)的相關(guān)函數(shù)，并按照所得關(guān)系，優(yōu)先出較為敏感的波段，以此建立相關(guān)性函數(shù)模型。其中可將波段作為高光譜建模的相關(guān)性標(biāo)準(zhǔn)如下：不論是正相關(guān)還是負(fù)相關(guān)，均需滿足0.05的顯著性水平，同時相關(guān)系數(shù)需大于0.5。

3.5 偏最小二乘回歸模型的建立

根據(jù)重金屬含量與不同光譜指標(biāo)的相關(guān)性分析結(jié)果，結(jié)合特征波段和范圍，將36個土壤樣本分為驗證樣本與建模樣本兩部分，前者設(shè)定為10個，后者建模分析用取26個。在建模完成后，對于所得結(jié)果，分別采用決定系數(shù)(R2)和均方根誤差(RMSE)進(jìn)行精度對比，見表2(考慮篇幅，僅列出元素Ni的建模結(jié)果)。其中，為防止出現(xiàn)比較離譜的數(shù)據(jù)結(jié)論，采用了決定系數(shù)的指標(biāo)，一般其值大于0.3時，表示具有一定的參考意義。而對于均方根誤差來說，其值與精度呈正相關(guān)，與誤差為負(fù)相關(guān)。在表2所示Ni元素的回歸模型中，建模精度最高的是二階微分，在建模樣本中，決定系數(shù)與均方根誤差分別為0.869 3、16.980 0；在驗證樣本中，決定系數(shù)與均方根誤差分別為0.798 3、17.990 0，對于土壤中重金屬含量的估算要求基本是可以滿足的。而對于二階微分外的其余光譜指標(biāo)，精度均不是很高，不適用于重金屬含量估算。

表2 重金屬元素Ni的建模結(jié)果Tab.2 Modeling results of Ni heavy metal elements

3.6 估算結(jié)果精度評價

為更加好地比較各光譜指標(biāo)的估算精度，直觀展示建模樣本與驗證樣本的實測值、預(yù)測值，并對建模結(jié)果進(jìn)行評價，制作了相關(guān)數(shù)據(jù)的散點圖。本文僅以元素Ni和Mn為例，如圖6所示。

圖6 元素Ni、Mn的數(shù)據(jù)散點Fig.6 Data scatter diagram of elements Ni and Mn

在圖6中，三角形為建模樣本，圓形為驗證樣本，斜線為y=x的趨勢線，模型擬合效果主要由散點與趨勢線的距離體現(xiàn)，距離越近，則表示精度越高。由圖6(a)可知，除少數(shù)樣本距離趨勢線較遠(yuǎn)外，大部分的散點與趨勢線均較近。因此，Ni元素的二階微分回歸模型擬合效果較好。由圖6(b)可看出，其大部分散點亦分布在趨勢線周圍，但上下位置分別有部分建模樣本偏離趨勢線較遠(yuǎn)，故Mn元素一階微分回歸模型的精度低于Ni。據(jù)此可同樣分析出其余重金屬元素的模型精度，本文不再贅述。

由模型精度分析可知，基于偏最小二乘法的應(yīng)用模型能有效估算重金屬含量，但最優(yōu)估算模型隨重金屬元素種類的不同而變化，對于重金屬元素Cr、Zn、As、Pb、Mn、Ni，其效果依次逐漸變好。其中，Ni元素對應(yīng)估算模型最適合的光譜指標(biāo)為二階微分，Mn元素對應(yīng)為一階微分，Pb元素對應(yīng)為反射率，As元素對應(yīng)為反射率，Zn元素對應(yīng)為二階微分，Cr元素對應(yīng)為去除包絡(luò)線。由此可知，為提高建模精度，不同的重金屬元素應(yīng)嘗試運用不同的光譜指標(biāo)；另外，元素Ni、Mn、Zn最優(yōu)模型對應(yīng)的光譜指標(biāo)均為微分變換值，說明在提高建模精度的措施方面，利用微分方法處理光譜反射率是一種有效的途徑。

4 結(jié)語

(1)土壤中的重金屬較易與黏土、鐵氧化物以及有機(jī)質(zhì)等聚集，土壤的這些理化性質(zhì)在高光譜數(shù)據(jù)中均有明顯的光譜響應(yīng)，這為利用高光譜模型間接估算重金屬含量提供了基礎(chǔ)原理。

(2)土壤的性質(zhì)與光譜的關(guān)系為非線性，因此應(yīng)用非線性的偏最小二乘等方法，可極大提高預(yù)測的精確性；采取去除包絡(luò)線、反射光譜倒數(shù)對數(shù)、微分?jǐn)?shù)據(jù)處理等方法，既可消除高光譜噪聲，也能提高建模精度。

(3)土壤重金屬估算的建模結(jié)果顯示，所述估算模型基本能滿足估算要求，同時也有良好的精度保證。另外，對光譜反射率進(jìn)行微分處理是提高建模精度的有效方法。