王雪梅，玉米提·買明，黃曉宇，李 銳，劉 東

1.新疆師范大學(xué)地理科學(xué)與旅游學(xué)院, 新疆 烏魯木齊 830054 2.新疆維吾爾自治區(qū)重點(diǎn)實(shí)驗(yàn)室“新疆干旱區(qū)湖泊環(huán)境與資源實(shí)驗(yàn)室”, 新疆 烏魯木齊 830054

引 言

隨著城市化進(jìn)程的加劇，工業(yè)和生活污水的排放以及農(nóng)藥等大量使用導(dǎo)致土壤重金屬污染日益嚴(yán)重[1]。近幾十年來，新疆石油化工產(chǎn)業(yè)和綠洲農(nóng)業(yè)的快速發(fā)展，使得土壤中的鉻、汞、鉛、銅、鎘和砷等重金屬元素在土壤環(huán)境中普遍存在。渭干河—庫車河三角洲綠洲是新疆南疆重要的農(nóng)業(yè)灌溉區(qū)和石油化工園區(qū)，土壤受到重金屬污染的風(fēng)險(xiǎn)較大。在自然界中，有害元素砷(As)常以離子形態(tài)存在于土壤溶液中，對動(dòng)植物和人體健康危害極大，而以膠狀物質(zhì)存在于土壤中的砷容易被植物吸收并經(jīng)食物鏈最終到達(dá)人體，危害生命健康[2]。開展土壤重金屬As含量的檢測對防治土壤重金屬污染、保障人體健康，以及構(gòu)建健康生態(tài)環(huán)境具有十分重要的現(xiàn)實(shí)意義[3]。

檢測土壤有害元素As的傳統(tǒng)方法主要采用分光光度法、化學(xué)分析法以及原子熒光光譜法等實(shí)驗(yàn)手段進(jìn)行，雖檢測精度高，但儀器設(shè)備和人力耗費(fèi)大且不易推廣，隨著高光譜遙感技術(shù)的快速發(fā)展，為土壤有害元素砷的遙感估測提供了有力的技術(shù)手段[4]。土壤樣品的高光譜采集受眾多環(huán)境因素的影響，光譜數(shù)據(jù)中存在許多細(xì)節(jié)噪聲。研究最初嘗試對土壤原始光譜數(shù)據(jù)進(jìn)行倒數(shù)、倒數(shù)對數(shù)以及一階微分和二階微分等多種數(shù)學(xué)變換處理，從而消除土壤樣品中的環(huán)境噪聲，增強(qiáng)土壤光譜數(shù)據(jù)中的有效信號[5-6]。隨著研究的深入，發(fā)現(xiàn)單純采用傳統(tǒng)的數(shù)學(xué)變換方法對原始光譜數(shù)據(jù)進(jìn)行處理，雖在一定程度上提高了模型的估測精度，卻難以有效去除白噪聲[7]。相關(guān)研究顯示，經(jīng)過連續(xù)小波變換(continuous wavelet transformation，CWT)后的光譜反射率與土壤屬性具有更強(qiáng)的關(guān)聯(lián)性，可實(shí)現(xiàn)光譜信號的近似特征和細(xì)節(jié)特征的有效分離，構(gòu)建的模型估測精度更高[8-9]。Zhang[10]、肖艷[11]等研究認(rèn)為光譜數(shù)據(jù)的連續(xù)小波變換對有效波段的挖掘和反演模型的估測都明顯優(yōu)于傳統(tǒng)的數(shù)學(xué)變換。為了探討在連續(xù)小波變換下進(jìn)行高光譜遙感估測干旱區(qū)綠洲土壤As含量的最佳小波基函數(shù)及最優(yōu)分解尺度，本研究基于野外調(diào)查的98個(gè)土壤樣品的光譜數(shù)據(jù)以及有害元素砷含量，通過4種小波基函數(shù)對原始光譜采用10個(gè)尺度進(jìn)行連續(xù)小波變換處理，篩選出與土壤As含量密切相關(guān)的敏感小波系數(shù)。以敏感小波系數(shù)為自變量，土壤As為因變量，采用偏最小二乘回歸(partial least squares regression，PLSR)、支持向量機(jī)回歸(super vector machine regression，SVMR)、BP神經(jīng)網(wǎng)絡(luò)(back propagation neural network，BPNN)和隨機(jī)森林回歸(random forest regression，RFR)方法構(gòu)建干旱區(qū)綠洲土壤有害元素As含量的反演模型，為有害元素As的準(zhǔn)確估算提供依據(jù)。

1 實(shí)驗(yàn)部分

1.1 土樣采集與制備

渭干河—庫車河三角洲綠洲位于新疆塔里木盆地北緣，是一個(gè)典型的干旱區(qū)扇形平原綠洲，轄區(qū)包括新疆維吾爾自治區(qū)阿克蘇地區(qū)的庫車市、沙雅縣和新和縣，地理位置為北緯39°30′—42°40′，東經(jīng)81°27′—84°07′(圖1)。該區(qū)多年平均降水量為51.6 mm，多年平均蒸發(fā)量可達(dá)2 123 mm，屬干旱與極干旱地區(qū)。綠洲植被主要以小麥(TriticumaestivumL)、棉花(Gossypiumspp)和玉米(ZeamaysLinn)等農(nóng)作物，以及檉柳(Tamarixramosissima)、鹽爪爪(Kalidiumfoliatum)、鹽穗木(Halostachyscaspica)、駱駝刺(Alhagisparsifolia)、蘆葦(PhragmitesaustraliasTrin)和胡楊(Populuseuphratica)等荒漠植被為主。土壤類型包括棕漠土、灌淤土、沼澤土和鹽土等多種類型。以遙感影像為參考，在室內(nèi)提前合理布點(diǎn)，2019年7月中到下旬進(jìn)行野外調(diào)查和土壤樣品采集。土壤剖面深度為0～20 cm，共采集98個(gè)土壤樣品，并分別采用GPS系統(tǒng)對樣點(diǎn)進(jìn)行準(zhǔn)確定位。采集的樣品裝入做好標(biāo)記的土樣袋中。土壤樣品在實(shí)驗(yàn)室經(jīng)風(fēng)干處理，挑出雜物，研磨過1 mm篩，送至新疆維吾爾自治區(qū)分析測試研究院，由工作人員依據(jù)GB/T 22105.2—2008第2部分測定標(biāo)準(zhǔn)，采用原子熒光光譜法進(jìn)行土壤總砷的測定。

圖1 采樣點(diǎn)分布圖

1.2 光譜數(shù)據(jù)的測定與預(yù)處理

采用ASDFieldSpec3光譜儀進(jìn)行研究區(qū)土壤樣品的高光譜數(shù)據(jù)采集。在野外測定光譜時(shí)，為避免溫差干擾數(shù)據(jù)的采集，首先要對光譜儀進(jìn)行預(yù)熱和白板校正。由于環(huán)境因素會使光譜數(shù)據(jù)在采集時(shí)產(chǎn)生干擾噪聲，故要遠(yuǎn)離高大物體。按照標(biāo)簽順序?qū)⑼寥罉悠肪鶆驍傞_并完全覆蓋在50 cm×50 cm的牛皮紙上，探頭與采樣點(diǎn)的垂直上方距離為15 cm。為了使采集的光譜數(shù)據(jù)更為準(zhǔn)確，每個(gè)樣品需要重復(fù)觀測10次，取平均值作為該樣品的光譜曲線值。為了消除儀器兩端產(chǎn)生的較大噪聲，剔除350～399和2 451～2 500 nm的光譜數(shù)據(jù)，采集的有效光譜范圍為400～2 450 nm。進(jìn)一步對98個(gè)樣品的光譜反射率進(jìn)行9點(diǎn)Savitzky-Golay濾波平滑處理，通過觀察98條光譜曲線，將受水分吸收和大氣影響的1 341～1 400和1 811～1 950 nm波段的異常數(shù)值進(jìn)一步剔除。

1.3 基于光譜數(shù)據(jù)的連續(xù)小波變換

連續(xù)小波變換(continuous wavelet transformation，CWT)也稱為積分小波變換(integral wavelet transform，IWT)，常應(yīng)用于各類信號的處理，具有較高的分辨性和適應(yīng)性，其變換公式見式(1)和式(2)

(1)

(2)

式中：a為尺度因子；b為平移因子；Wf(a，b)為小波系數(shù)，采用二維矩陣表示；λ為波長；f(λ)為土壤光譜反射率；Ψa，b為小波基函數(shù)。CWT作為一種信號處理方法，可有效去除白噪音，然而在光譜分解過程中，小波基函數(shù)的選擇一直沒有統(tǒng)一的說法。相關(guān)研究發(fā)現(xiàn)，采用不同的小波基函數(shù)對光譜數(shù)據(jù)進(jìn)行各種尺度的分解，會直接導(dǎo)致實(shí)驗(yàn)結(jié)果不一致，目前多以經(jīng)驗(yàn)和反復(fù)實(shí)驗(yàn)比較進(jìn)行選擇[12]。研究中使用4種小波基函數(shù)bior1.3，db4，gaus4與mexh對原始光譜反射率R進(jìn)行10個(gè)尺度的光譜分解，選取敏感小波系數(shù)作為自變量構(gòu)建干旱區(qū)土壤有害物質(zhì)砷含量的高光譜估測模型。

1.4 模型的建立與檢驗(yàn)

采用偏最小二乘回歸(PLSR)、支持向量機(jī)回歸(SVMR)、BP神經(jīng)網(wǎng)絡(luò)(BPNN)和隨機(jī)森林回歸(RFR)4種建模方法對研究區(qū)土壤有害物質(zhì)砷含量進(jìn)行高光譜反演，其中PLSR模型的最佳潛在變量(Latent variables)按照小波系數(shù)的數(shù)量而定；SVMR模型以徑向基函數(shù)radial作為核函數(shù)，參數(shù)Cost、Epsilon和Gamma可通過反復(fù)訓(xùn)練和調(diào)試最終確定最佳值；BP神經(jīng)網(wǎng)絡(luò)模型的訓(xùn)練函數(shù)選擇newff，訓(xùn)練過程中最大迭代次數(shù)為500，訓(xùn)練目標(biāo)最小誤差為0.001，學(xué)習(xí)率為0.01，通過不斷改變隱含層的神經(jīng)元個(gè)數(shù)訓(xùn)練該模型以獲得最佳效果；隨機(jī)森林回歸(RFR)采用randomForest函數(shù)，通過指定決策樹的數(shù)目ntree來反復(fù)多次訓(xùn)練尋找最優(yōu)的模型參數(shù)mtry。模型的精度和預(yù)測能力通過決定系數(shù)(R2)、均方根誤差(RMSE)和相對分析誤差(RPD)進(jìn)行評價(jià)。R2和RPD值越大，RMSE值越小，則模型的精度越高，估測能力越強(qiáng)。

2 結(jié)果與討論

2.1 有害物質(zhì)砷含量與光譜特征分析

對研究區(qū)98個(gè)土壤樣本的有害物質(zhì)砷含量進(jìn)行基本統(tǒng)計(jì)(表1)，有害元素砷含量的最大值為22.8 mg·kg-1，最小值為4.35 mg·kg-1，平均水平為11.67 mg·kg-1，屬于中等空間變異。由偏度系數(shù)和峰度系數(shù)發(fā)現(xiàn)，砷含量總體偏高，且呈高聚集狀態(tài)。根據(jù)建模需要將總樣本隨機(jī)分為68個(gè)訓(xùn)練樣本和30個(gè)驗(yàn)證樣本。與新疆土壤砷的背景值和國家標(biāo)準(zhǔn)比較后發(fā)現(xiàn)，研究區(qū)As含量的平均水平已超過新疆土壤背景值(11.2 mg·kg-1)，但還未超出國家標(biāo)準(zhǔn)(40 mg·kg-1)，說明As已過量沉積形成輕度污染。分析認(rèn)為隨著研究區(qū)工業(yè)化進(jìn)程的加劇以及農(nóng)藥的過量使用，使耕層土壤中的有害物質(zhì)As含量在土壤中不斷蓄積，成為影響農(nóng)業(yè)可持續(xù)發(fā)展和危及人體健康的重要因素。

表1 土壤樣品的基本統(tǒng)計(jì)特征

為了分析砷含量與土壤光譜反射率之間的關(guān)聯(lián)性，將98個(gè)樣本按照砷含量由低到高劃分為5個(gè)等級，Ⅰ級：As<8 mg·kg-1，Ⅱ級：8 mg·kg-1≤As<10 mg·kg-1，Ⅲ級：10 mg·kg-1≤As<12 mg·kg-1，Ⅳ級：12 mg·kg-1≤As≤14 mg·kg-1，Ⅴ級：As≥14 mg·kg-1。通過分析不同砷含量等級的光譜反射率曲線(圖3)，發(fā)現(xiàn)隨著砷含量的增加，光譜曲線呈現(xiàn)出下降的趨勢，說明土壤中砷元素的含量直接影響到光譜反射率，采用高光譜數(shù)據(jù)進(jìn)行土壤砷含量估測具有一定的依據(jù)。

圖2 不同砷含量等級的光譜反射率曲線

2.2 CWT下的光譜反射率與有害物質(zhì)砷含量的相關(guān)分析

通過對原始光譜數(shù)據(jù)分別采用4種小波基函數(shù)進(jìn)行10個(gè)尺度的連續(xù)小波變換，變換后的小波系數(shù)與As含量進(jìn)行相關(guān)分析并繪制相關(guān)譜圖。圖3(a—d)是土壤As含量與不同小波基函數(shù)變換下10個(gè)尺度的相關(guān)系數(shù)二維圖，橫坐標(biāo)為波長(400～2 450 nm)，縱坐標(biāo)為不同分解尺度(2n，n=1, 2, 3, …, 10)，顏色表示相關(guān)系數(shù)的大小。從圖3可以看出，有害物質(zhì)As含量與CWT處理后的光譜反射率存在一定程度的相關(guān)性。其中，在可見光400～700nm以及近紅外的1 100～1 700和2 200～2 400 nm附近具有較強(qiáng)的相關(guān)關(guān)系，且表現(xiàn)出3～8尺度可極大程度提高光譜與土壤屬性之間的相關(guān)性，其中4～6尺度具有最佳的光譜分解能力。

進(jìn)一步分析圖3發(fā)現(xiàn)采用bior1.3小波基函數(shù)進(jìn)行第5尺度分解的小波系數(shù)與As含量的相關(guān)性最強(qiáng)，能夠通過p<0.01顯著性檢驗(yàn)的光譜波段有507個(gè)，且主要集中分布在可見光400～600 nm與分散分布在近紅外1 100～2 400 nm。為了得到建模所需的敏感小波系數(shù)，進(jìn)一步對第5尺度分解結(jié)果通過p<0.001相關(guān)顯著性檢驗(yàn)，篩選出相關(guān)系數(shù)絕對值的敏感波段共12個(gè)，波長范圍主要集中在2 343～2 354 nm，相關(guān)系數(shù)r最高可達(dá)0.687。通過對db4小波基函數(shù)的各尺度分解結(jié)果分析，發(fā)現(xiàn)第6尺度與As的相關(guān)性最強(qiáng)，通過p<0.01顯著性檢驗(yàn)的光譜波段有203個(gè)，主要集中分布在可見光的500～700 nm和近紅外的2 300～2 400 nm附近；通過進(jìn)一步篩選的敏感波段共有28個(gè)(p<0.001)，相關(guān)系數(shù)r的絕對值最高可達(dá)0.686。采用gaus4小波基函數(shù)進(jìn)行10個(gè)尺度的光譜分解結(jié)果顯示，通過p<0.01顯著性檢驗(yàn)的光譜波段在第5尺度最多，達(dá)到280個(gè)，且主要集中分布在400～600 nm以及2 300～2 400 nm；與As具有極顯著相關(guān)性且的敏感波段有16個(gè)(p<0.001)，相關(guān)系數(shù)絕對值最高可達(dá)0.669。mexh小波基函數(shù)的第4尺度的光譜分解能力最強(qiáng)，通過p<0.01顯著性檢驗(yàn)的光譜波段達(dá)到216個(gè)，主要集中在400～600，1 100～1 300以及2 300～2 400 nm，通過p<0.001相關(guān)顯著性檢驗(yàn)且的敏感波段有24個(gè)，相關(guān)系數(shù)絕對值最高可達(dá)0.689。

圖3 砷含量與小波系數(shù)間的相關(guān)性

通過上述分析可以看出，采用4種小波基函數(shù)對光譜數(shù)據(jù)進(jìn)行連續(xù)小波變換后，小波系數(shù)與As含量具有較顯著的相關(guān)性，表明基于不同小波基函數(shù)的CWT能較好地放大各個(gè)波段范圍內(nèi)微弱的土壤特征信息并更有效地挖掘信息。對比4種小波基函數(shù)的分解結(jié)果發(fā)現(xiàn)，經(jīng)bior1.3小波處理后與As相關(guān)性顯著的小波系數(shù)數(shù)量遠(yuǎn)高于其他3種小波函數(shù)(p<0.01)，說明bior1.3小波變換可更大程度上實(shí)現(xiàn)光譜信號的近似特征和細(xì)節(jié)特征的有效分離，提高光譜與As之間的相關(guān)性。進(jìn)一步進(jìn)行p<0.001相關(guān)極顯著性檢驗(yàn)發(fā)現(xiàn)，小波基函數(shù)db4和mexh的分解效果要優(yōu)于gaus4。通過相關(guān)系數(shù)絕對值|r|>0.6可篩選出4種小波基函數(shù)下的敏感小波系數(shù)作為As含量估算的建模變量(見表2)。

表2 篩選的敏感小波系數(shù)

2.3 反演模型的建立

為進(jìn)一步確定敏感小波系數(shù)與土壤As含量之間的定量關(guān)系，以各小波基函數(shù)CWT處理后篩選出的敏感小波系數(shù)為自變量，采用PLSR，SVMR，BPNN和RFR方法構(gòu)建土壤As含量的高光譜估測模型，并通過決定系數(shù)(R2)、均方根誤差(RMSE)和相對分析誤差(RPD)對模型進(jìn)行精度評價(jià)和可靠性檢驗(yàn)。分析表3反演結(jié)果認(rèn)為，4種模型整體估測效果良好，與其他3種機(jī)器學(xué)習(xí)算法相比，PLSR模型的估測精度相對較低；對比SVMR，BPNN與RFR模型發(fā)現(xiàn)，SVMR和BPNN模型訓(xùn)練集的R2和RPD值較為接近，驗(yàn)證集則表現(xiàn)為BPNN模型的估測能力要優(yōu)于SVMR模型，而RFR模型的訓(xùn)練集和驗(yàn)證集估測精度均較高，模型更穩(wěn)定。

表3 不同反演模型的估算精度

對比4種小波基函數(shù)的估測效果發(fā)現(xiàn)，PLSR模型中mexh函數(shù)的估測精度要高于其他3種小波基函數(shù)；在SVMR模型中，4種小波基函數(shù)的光譜分解能力一般，模型估測效果依次為mexh-24>db4-26>bior1.3-25>gaus4-25；而在BPNN和RFR模型中，以bior1.3小波基函數(shù)構(gòu)建的模型具有更強(qiáng)的估測能力，其建模集和驗(yàn)證集的R2均在0.6以上，RMSE在1.9 mg·kg-1以下，RPD值均大于1.6，說明bior1.3小波基函數(shù)對光譜數(shù)據(jù)的分解效果較佳，構(gòu)建的模型具有較好的估測能力。通過對比4種小波基函數(shù)對原始光譜數(shù)據(jù)中有效信號的辨識能力，可看出bior1.3和mexh的變換效果明顯優(yōu)于db4和gaus4函數(shù)。通過綜合比較分析認(rèn)為，bior1.3-25-RFR模型對As的估測精度最高，具有最強(qiáng)的穩(wěn)定性。該模型參數(shù)ntree和mtry的最佳值分別為500和8，其建模集和驗(yàn)證集的R2分別為0.893和0.639，RMSE為1.075和1.651 mg·kg-1，RPD值分別為2.89和1.64，均大于1.6，模型具有較好的估測能力。

分別對4種建模方法下有害物質(zhì)As的最佳反演結(jié)果繪制散點(diǎn)圖[圖4(a—d)]，發(fā)現(xiàn)各模型的估測值與實(shí)測值構(gòu)成的散點(diǎn)絕大多數(shù)都分布在1∶1線附近，相對較為集中，說明這4種建模方法對有害物質(zhì)As均具有良好的解釋能力。對比分析認(rèn)為bior1.3-25-RFR模型的散點(diǎn)在1∶1線附近分布的更為集中，總體估測誤差最小，說明該模型能更好地反演土壤有害物質(zhì)As的含量。

圖4 實(shí)測值與估測值的散點(diǎn)圖

3 結(jié) 論

適當(dāng)?shù)墓庾V變換方法有利于土壤信息的提取和光譜數(shù)據(jù)中有效信息的識別。作為一種光譜數(shù)據(jù)處理方法，CWT可以有效去除土壤高光譜數(shù)據(jù)中的噪音，檢測有效信號的效果明顯，可更好地挖掘出光譜數(shù)據(jù)中的微弱卻有效的隱藏信息，對土壤高光譜反演起著重要的作用[13]。Li等[14]的研究表明，土壤原始光譜反射率經(jīng)CWT處理在一定程度上提高了與土壤有機(jī)質(zhì)含量間的相關(guān)性，說明光譜的CWT變換可消除不確定因素對光譜信息的干擾，從土壤光譜數(shù)據(jù)中分解出穩(wěn)定的高頻信息。Zhang等[10]對4種元素(Cr，As，Ni和Cd)進(jìn)行高光譜反演結(jié)果顯示，與傳統(tǒng)數(shù)學(xué)變換方法相比，CWT變換能更好地提取光譜中的有效信息，具有良好的光譜分解能力。Hong等[15]研究顯示，采用CWT變換的隨機(jī)森林模型對土壤有機(jī)碳具有更好的估算能力。由于土壤性質(zhì)的復(fù)雜性以及土壤類型的多樣性，不同研究區(qū)域的土壤屬性存在較大空間差異，如何采用CWT方法選擇合適的小波基函數(shù)和最優(yōu)的分解尺度，以及小波基函數(shù)是否穩(wěn)定和具有普遍適用性仍需做進(jìn)一步探討。

研究基于bior1.3，db4，gaus4和mexh 4種小波基函數(shù)進(jìn)行連續(xù)小波變換處理，結(jié)合土壤As含量進(jìn)行相關(guān)分析，通過篩選的敏感小波系數(shù)采用PLSR，SVMR，BPNN和RFR模型對干旱區(qū)土壤As含量進(jìn)行估算，結(jié)果表明：(1)在CWT過程中，4種小波基函數(shù)在4～6尺度上的光譜分解效果明顯優(yōu)于其他尺度，與土壤As具有顯著相關(guān)的小波系數(shù)數(shù)量更多(p<0.01)；(2)對比4種小波基函數(shù)的光譜變換效果，bior1.3和mexh對光譜數(shù)據(jù)中有效信息的提取能力要強(qiáng)于db4和gaus4，說明這兩種小波基函數(shù)更利于挖掘光譜數(shù)據(jù)中的隱藏信息，可在一定程度上放大土壤中的有效信息；(3)與PLSR，SVMR和BPNN模型相比，RFR模型具有更高的估測精度和較好的穩(wěn)定性；土壤中As的最佳估算模型為bior1.3-25-RFR模型，其訓(xùn)練集和驗(yàn)證集的R2均大于0.6，均方根誤差RMSE小于1.7 mg·kg-1，RPD值大于1.6，估測能力較好。