本發(fā)明涉及土壤重金屬含量高光譜反演領(lǐng)域，特別涉及基于室內(nèi)標(biāo)樣高光譜特征的土壤重金屬cd含量反演建模及其光譜響應(yīng)特征波段識別的方法。

背景技術(shù)：

隨著城市化、工業(yè)化和農(nóng)業(yè)集約化的發(fā)展，農(nóng)用地土壤重金屬污染物因殘留時間長、隱蔽性強、遷移性小、毒性大、不易被微生物降解、可以通過食物鏈富集等特點，逐漸受到國內(nèi)外研究者的普遍關(guān)注。然而，因大范圍高密采樣成本耗費巨大、常態(tài)化快速高效動態(tài)監(jiān)測難以實施，大范圍開展重金屬污染治理至今仍然存在“重金屬污染安全隱患點不清、污染空間格局不明”兩大突出問題，已成為擺在土壤重金屬污染全面與深度治理面前的一道“理論與方法”鴻溝，無法滿足國家在廣域地理空間快速高效查清重金屬污染空間格局、發(fā)現(xiàn)流域重金屬污染時空演化規(guī)律與形成機制、控制污染范圍擴大和轉(zhuǎn)移、合理規(guī)劃農(nóng)業(yè)生產(chǎn)、減少污染產(chǎn)品損害國民健康的重大需求。

與此同時，新興的遙感探測技術(shù)具有能通過非接觸的方式對目標(biāo)物體進行快速、動態(tài)、大范圍監(jiān)測的特征，高光譜更是因其多且連續(xù)的光譜波段優(yōu)勢逐漸被廣泛應(yīng)用于土壤重金屬污染監(jiān)測研究中，為定量反演土壤重金屬含量提供了新的契機，使得大范圍、高時空分辨率、高精度快速探測土壤重金屬污染靶區(qū)的含量空間分布成為可能。

近年來，國內(nèi)外圍繞土壤重金屬含量高光譜反演開展了大量的研究工作。截止目前，大量研究已證實：當(dāng)土壤重金屬含量達(dá)到一定數(shù)值時，在地面可見光、近紅外和熱紅外波段光譜測試中可發(fā)現(xiàn)明顯的重金屬光譜響應(yīng)特征；在土壤重金屬含量偏低時，仍可通過地面土壤有機質(zhì)、鐵猛氧化物和粘土礦物的光譜特征間接推算重金屬含量。然而，土壤成分復(fù)雜，重金屬含量甚微，有效排除土壤水分、有機質(zhì)、鐵錳氧化物等成分光譜干擾，建立重金屬含量反演模型，準(zhǔn)確識別微量重金屬元素光譜響應(yīng)特征波段，是亟待解決的關(guān)鍵技術(shù)難題。

技術(shù)實現(xiàn)要素：

為了有效排除土壤水分、有機質(zhì)、鐵錳氧化物等成分光譜干擾，準(zhǔn)確識別微量重金屬元素光譜響應(yīng)特征波段，本發(fā)明提供一種基于室內(nèi)標(biāo)樣高光譜特征的土壤重金屬cd含量反演建模及其光譜響應(yīng)特征波段識別的方法。

為了實現(xiàn)上述技術(shù)目的，本發(fā)明的技術(shù)方案是，一種基于室內(nèi)標(biāo)樣高光譜特征的土壤重金屬cd含量反演建模及其光譜響應(yīng)特征波段識別的方法，包括以下步驟：

步驟1：土壤標(biāo)樣制作及其準(zhǔn)確重金屬cd含量與土壤光譜數(shù)據(jù)獲??；

步驟2：土壤重金屬cd含量分析與光譜數(shù)據(jù)預(yù)處理；

步驟3：基于步驟2轉(zhuǎn)換得到的光譜數(shù)據(jù)的土壤重金屬cd含量反演建模，實現(xiàn)光譜響應(yīng)特征波段識別。

所述的方法，步驟1中所述的土壤標(biāo)樣的制作及準(zhǔn)確其重金屬cd含量的獲取，包括以下步驟：

步驟1)：充分考慮土壤重金屬cd污染含量實際變化范圍及其光譜響應(yīng)效果，人為設(shè)計不同cd污染含量(>50等級)(單位：mg/kg)下的標(biāo)準(zhǔn)土壤樣本；

步驟2)：采用格網(wǎng)采集方式以1m×1m分辨率采集研究區(qū)未經(jīng)重金屬污染的背景土壤樣本若干(>50份)，留作標(biāo)樣添加土壤底樣；

步驟3)：將步驟2)采集到的背景土壤經(jīng)實驗室烘干、研磨處理，準(zhǔn)確稱量得到標(biāo)樣添加土壤底樣，每份100g；

步驟4)：取步驟3)處理后的底樣一份，采用三酸消化——原子吸收分光光度法實驗室化學(xué)分析測定背景土壤重金屬cd含量(單位：mg/kg)；

步驟5)：計算達(dá)到步驟1)中設(shè)定的cd污染含量值需向100g底樣中添加的1000mg/kg的cd標(biāo)準(zhǔn)溶液，計算公式如下：

c＝100(a-s)/ρ

式中，c代表加入100g底樣土壤中的cd標(biāo)準(zhǔn)溶液量(單位：ml)，a代表步驟1)中設(shè)計的土壤標(biāo)樣的重金屬cd含量(單位：mg/kg)，s代表步驟4)中測得的背景土壤樣本的重金屬cd含量，ρ代表cd標(biāo)準(zhǔn)溶液的密度(單位：kg/m³)。

特別地，向100g土壤底樣中人為添加cd標(biāo)準(zhǔn)溶液后，需經(jīng)風(fēng)干處理后測定土壤重金屬cd含量，標(biāo)樣cd含量最高不超過0.1mg，為底樣的10^-4數(shù)量級，可忽略不計，因此土壤總量按100g計。

步驟6)：實驗室按步驟5)中計算得到的cd標(biāo)準(zhǔn)溶液量依次添加到步驟3)獲取的土壤底樣中，制得土壤標(biāo)準(zhǔn)樣本；

步驟7)：采用三酸消化一一一原子吸收分光光度法實驗室化學(xué)分析測定步驟6)中制得的土壤標(biāo)樣，準(zhǔn)確獲取土壤標(biāo)樣重金屬cd含量。

所述的方法，步驟2中所述的土壤污染標(biāo)樣光譜與背景樣本光譜作比值的光譜預(yù)處理方法，包括以下步驟：

步驟1)：將實驗室人為添加cd污染標(biāo)樣制作的土壤標(biāo)樣光譜與未經(jīng)人為添加重金屬cd標(biāo)準(zhǔn)溶液的土壤背景樣本光譜同時做10nm間隔重采樣—標(biāo)準(zhǔn)正態(tài)變換(snv)—一階/二階微分處理，得到轉(zhuǎn)換后的光譜數(shù)據(jù)；

步驟2)：將經(jīng)步驟1)光譜變換后的土壤污染標(biāo)樣光譜與背景樣本光譜作比值，得到光譜特征數(shù)據(jù)，計算公式如下：

ri＝ai/bi(i＝1,2,…,n)

式中，n代表光譜波段數(shù)，i代表第i個波段，ai代表第i個波段的步驟1)光譜轉(zhuǎn)換后的土壤污染樣本光譜反射率，bi代表第i個波段的步驟1)光譜轉(zhuǎn)換后的土壤光譜反射率，ri代表第i個波段經(jīng)比值變換后得到的光譜反射率，作為土壤重金屬cd含量高光譜反演建模輸入變量光譜特征數(shù)據(jù)。

附圖說明

圖1為實驗室標(biāo)樣制作方法流程；

圖2為光譜基于各類光譜預(yù)處理方法的光譜變換效果實例；

圖3為偏最小二乘法建模因子提取過程的基本原理；

圖4為偏最小二乘法建模最小平方誤差理論原理。

具體實施方式

下面是對本發(fā)明一個優(yōu)選實施例，結(jié)合附圖進行的詳細(xì)說明。

1、如圖1流程圖所示，土壤標(biāo)樣制作及其準(zhǔn)確重金屬cd含量與土壤光譜數(shù)據(jù)獲取，具體原理方法如下：

步驟1)：充分考慮土壤重金屬cd污染含量實際變化范圍及其光譜響應(yīng)效果，人為設(shè)計不同cd污染含量(>50等級)(單位：mg/kg)下的標(biāo)準(zhǔn)土壤樣本；

步驟2)：采用格網(wǎng)采集方式以1m×1m分辨率采集研究區(qū)未經(jīng)重金屬污染的背景土壤樣本若干(>50份)，留作標(biāo)樣添加土壤底樣；

步驟3)：將步驟2)采集到的背景土壤經(jīng)實驗室烘干、研磨處理，準(zhǔn)確稱量得到標(biāo)樣添加土壤底樣，每份100g；

步驟4)：取步驟3)處理后的底樣一份，采用三酸消化——原子吸收分光光度法實驗室化學(xué)分析測定背景土壤重金屬cd含量(單位：mg/kg)；

步驟5)：計算達(dá)到步驟1)中設(shè)定的cd污染含量值需向100g底樣中添加的1000mg/kg的cd標(biāo)準(zhǔn)溶液，計算公式如下：

c＝100(a-s)/ρ

式中，c代表加入100g底樣土壤中的cd標(biāo)準(zhǔn)溶液量(單位：ml)，a代表步驟1)中設(shè)計的土壤標(biāo)樣的重金屬cd含量(單位：mg/kg)，s代表步驟4)中測得的背景土壤樣本的重金屬cd含量，ρ代表cd標(biāo)準(zhǔn)溶液的密度(單位：kg/m³)。

特別地，向100g土壤底樣中人為添加cd標(biāo)準(zhǔn)溶液后，需經(jīng)風(fēng)干處理后測定土壤重金屬cd含量，標(biāo)樣cd含量最高不超過0.1mg，為底樣的10^-4數(shù)量級，可忽略不計，因此土壤總量按100g計。

步驟6)：實驗室按步驟5)中計算得到的cd標(biāo)準(zhǔn)溶液量依次添加到步驟3)獲取的土壤底樣中，制得土壤標(biāo)準(zhǔn)樣本；

步驟7)：采用三酸消化一一一原子吸收分光光度法實驗室化學(xué)分析測定步驟6)中制得的土壤標(biāo)樣，準(zhǔn)確獲取土壤標(biāo)樣重金屬cd含量。

2、采用光譜重采樣—標(biāo)準(zhǔn)正態(tài)變換(snv)—一階/二階微分—背景樣比值系列光譜變換方法進行光譜數(shù)據(jù)預(yù)處理，具體原理方法如下：

步驟1)：為消除數(shù)據(jù)冗余，降低噪聲，平滑波形，得到平穩(wěn)的變化波譜，以10nm為間隔進行光譜數(shù)據(jù)的重采樣，計算公式如下：

r＝(r1+r2+…+rn)/n

式中，r代表重采樣后的光譜反射率；n代表該10nm間隔內(nèi)波段數(shù)，rn代表第n個波段的光譜反射率。

步驟2)：對重采樣后的光譜采用加權(quán)平均化方法進行標(biāo)準(zhǔn)正態(tài)變量變換(snv)。對于給定的樣本，snv計算了所有變量的標(biāo)準(zhǔn)差，整個樣本再通過該值進行歸一化，計算公式如下：

式中，n是變量的個數(shù)；xi,j是第i個樣本的第j個變量的值；δ是自定義的偏移量。

步驟3)：采用一階/二階微分處理去除部分線性或接近線性的背景值干擾，即噪聲光譜對目標(biāo)光譜的影響，增強光譜特征差異，提取光譜特征吸收帶，光譜微分計算公式為：

r′(λi)＝[r(λi)-r(λi-1)]/2δλ

r″(λi)＝[r(λi)-2r(λi-1)+r(λi-2)]/(2δλ)²

式中，λi為每個波段的光譜反射率，r′(λi)和r″(λi)分別為波長λi的一階和二階微分光譜，δλ是波長λi-1到λi的間隔，本研究間隔經(jīng)重采樣變換均為10nm。

步驟4)：將經(jīng)重采樣—snv—一階/二階微分光譜變換后的土壤污染標(biāo)樣光譜與背景樣本光譜作比值，得到光譜特征數(shù)據(jù)，計算公式如下：

ri＝ai/bi(i＝1,2,…,n)

式中，n代表光譜波段數(shù)，i代表第i個波段，ai代表第i個波段的步驟1)光譜轉(zhuǎn)換后的土壤污染樣本光譜反射率，bi代表第i個波段的步驟1)光譜轉(zhuǎn)換后的土壤光譜反射率，ri代表第i個波段經(jīng)比值變換后得到的光譜反射率，作為土壤重金屬cd含量高光譜反演建模輸入變量光譜特征數(shù)據(jù)。

3、基于數(shù)據(jù)預(yù)處理光譜轉(zhuǎn)換得到的光譜數(shù)據(jù)采用偏最小二乘回歸(plsr)方法構(gòu)建土壤重金屬含量反演模型，初步識別出土壤重金屬含量反演的特征波段土壤重金屬cd含量高光譜反演建模的步驟包括：

步驟1)：設(shè)光譜變換后得到的波段數(shù)為m，土壤樣本數(shù)量為n，t，u分別為從自變量與因變量中提取的因子，稱為偏最小二乘因子。從變換后得到的各波段光譜反射率(x1,x2,…,xm)中提取相互獨立的成分th(h＝1,2,…)，從各土壤樣本重金屬cd含量值(y1,y2,…,yn)中提取相獨立的成分uj(j＝1,2,…)，要求th和uj之間的相關(guān)程度達(dá)到最大。從原始變量集中提取第一對因子t1和u1的線性組合公式為：

式中，w1′＝(w1,1,w1,2,…,w1,m)′為模型效應(yīng)權(quán)重；v1′＝(v1,1,v1,2,…,v1,m)′為因變量權(quán)重。

步驟2)：建立初始變量對t1的方程。

其中，α′1＝(α1,1,α1,2,…,α1,m)′；β1′＝(β1,1,β1,2,…,β1,m)′為僅有一個自變量t1時的參數(shù)向量，其中α′1為模型效應(yīng)載荷量。e1和f1分別為殘差陣。根據(jù)普通最小二乘法可求得參數(shù)向量α1，β1。

步驟3)：利用多元回歸方法建立提取的成分與因變量土壤重金屬含量值之間的回歸方程。若提取的第一因子不能達(dá)到回歸模型的精度要求，則利用殘差陣e1和f1相應(yīng)取代x0和y0，重復(fù)步驟2)，步驟3)繼續(xù)提取因子，以此類推。假設(shè)最終提取了r個因子，x0和y0對r個因子的回歸方程為：

把第一步分析所得的自變量中提取因子th(h＝1,2,…r)的線性組合帶入因變量對r個因子建立的回歸方程，即得因變量的回歸方程，即為土壤重金屬含量反演模型，具體公式為：

y＝α1x1+α2x2+···+αmxm

步驟3)：提取上述模型得到的與土壤重金屬含量顯著相關(guān)的波段(xs1,xs2,…,xsm)，即為光譜響應(yīng)特征波段。