基于光譜FOD與優(yōu)化指數(shù)的銀川平原土壤有機質(zhì)含量反演

張俊華 尚天浩 陳睿華 王怡婧 丁啟東 李小林

(1.寧夏大學(xué)生態(tài)環(huán)境學(xué)院, 銀川 750021；2.西北土地退化與生態(tài)恢復(fù)國家重點實驗室培育基地, 銀川 750021；3.西安煤航遙感信息有限公司，西安 710199；4.寧夏大學(xué)地理科學(xué)與規(guī)劃學(xué)院，銀川 750021)

0 引言

土壤有機質(zhì)(Soil organic matter，SOM)是植物營養(yǎng)物質(zhì)的主要來源之一，在土壤肥力質(zhì)量評價中具有重要地位[1]。SOM成分復(fù)雜，光譜特性差異大，但總體呈現(xiàn)出SOM含量升高，整個譜線反射率逐漸降低的趨勢[2]。但當(dāng)SOM含量較低，光譜波段被吸收部分的能量過低時，光譜特性分析效果不佳[3]。ABERGAZ等[4]指出，SOM含量(質(zhì)量比)小于20 g/kg時，光譜曲線易受其他成分的影響，不能準確反映由官能團引起的吸收特征。因此，如何提高低含量SOM的反演精度是相關(guān)領(lǐng)域的研究熱點之一。

利用不同的反射率變換方式篩選SOM敏感波段和采用不同建模方法都是提高低含量SOM反演模型精度的有效手段。XIE等[5]研究指出，與反射率的一階微分、倒數(shù)的一階微分、平方根的一階微分和對數(shù)的一階微分變換方式相比，特征波段選擇一階微分可以顯著提高低含量SOM的反演精度。LIU等[6]在一階微分、二階微分和吸光度對數(shù)、吸光度對數(shù)的一階和二階微分及連續(xù)統(tǒng)去除等變換方式的基礎(chǔ)上，采用競爭自適應(yīng)加權(quán)算法篩選出15～40個不等的特征波段，發(fā)現(xiàn)采用連續(xù)統(tǒng)去除-隨機森林建立反演陜西省靖邊縣SOM模型精度達0.96，相對分析誤差(Residual predictive derivation，RPD)為3.02。在建模方法中，最小二乘法(Least square，LS)[7-8]、支持向量機(Support vector machine，SVM)[7-8]、支持向量機分類(Support vector machine classification，SVMC)[9]、遺傳算法(Genetic algorithm，GA)[10]、最小二乘回歸(Least square regression，LSR)[11]、反向傳播神經(jīng)網(wǎng)絡(luò)(Back propagation neural network，BPNN)[12-13]、隨機森林(Random forest，RF)[14-15]等一種或多種方法被用于巴西、以色列荒漠區(qū)以及我國河北、新疆、西藏、寧夏等地區(qū)低含量SOM的反演，都取得了較好的效果。此外，引入其他變量也可以提高SOM反演精度[13，16-17]。雖然變量的引入可以顯著提高SOM含量的反演精度，但這些變量的獲取技術(shù)難度較大，所以其應(yīng)用和推廣有一定困難。

為提高模型運算效率，減少建模變量輸入個數(shù)，有學(xué)者在敏感波段篩選的基礎(chǔ)上，用兩個或兩個以上特征波段分別進行波段耦合運算，構(gòu)建比值指數(shù)(Ratio index，RI)、差值指數(shù)(Difference index，DI)、土壤沙化指數(shù)(Soil desertification index，SDI1)、土壤退化指數(shù)(Soil degradation index，SDI2)、修正歸一化差分指數(shù)(Corrected normalized difference index，CNDI)等[18-20]，也有學(xué)者嘗試在二維空間建立光譜數(shù)據(jù)與SOM間的全波段相關(guān)性運算，以解決單一敏感波段在一維空間構(gòu)建光譜指數(shù)所造成的光譜信息缺失問題[21-22]。二維空間所建模型精度較一維空間所建模型精度得到顯著提升，但因研究區(qū)域不同，導(dǎo)致所選最佳建模光譜指數(shù)存在較大差異。綜合來看，已有研究或基于不同二維光譜指數(shù)分別建模，或以分數(shù)階微分(Fractional-order derivarives，F(xiàn)OD)、整數(shù)階微分聯(lián)合單一光譜指數(shù)建模，將FOD聯(lián)合優(yōu)化指數(shù)進行較低SOM含量反演的報道相對較少。

1 材料與方法

1.1 研究區(qū)概況

銀川平原(37°50′～39°23′N，104°17′～107°39′E)地處寧夏北部，是我國糧食主產(chǎn)區(qū)之一，屬溫帶干旱氣候，年均降水量約190 mm，全年蒸發(fā)量約 1 750 mm。土壤類型多以灌淤土、灰鈣土和鹽堿土為主[24]，SOM含量普遍偏低。

1.2 土壤樣本采集

根據(jù)銀川平原氣候、水文、地質(zhì)地貌和土壤等基礎(chǔ)資料數(shù)據(jù)，主要以農(nóng)田土壤為研究對象，利用ArcGIS 10.4預(yù)設(shè)5 km×5 km 網(wǎng)格。于2019年5月以梅花五點法采集裸露地表0～20 cm的土樣(以光譜測定點為中心，半徑為1 m距離4個方向處各選1個點，共5個點)，混合均勻后保留0.5 kg作為該點土樣。采樣過程中受建筑、道路和水體的影響，適當(dāng)調(diào)整采樣點(圖1)。土樣風(fēng)干后采用外加熱-重鉻酸鉀法測定SOM含量，共獲得有效土壤樣本187個。

圖1 研究區(qū)及樣點分布

1.3 野外光譜采集及預(yù)處理

采用ASD FieldSpc4型地物光譜儀進行野外原位土壤的高光譜測定，波段為350～2 500 nm，每次測定前均進行白板校正，光譜探頭始終保持在距地面約80 cm，每一樣點重復(fù)測定5次，取平均值作為該樣點的光譜反射率。測定時間為10:00—14:00，天氣狀況良好，晴朗無風(fēng)。

為消除儀器噪聲和環(huán)境背景干擾，去除噪聲過大的邊緣波段(350～399 nm和2 401～2 500 nm)。在Matlab 2018b軟件中對400～2 400 nm范圍內(nèi)的高光譜數(shù)據(jù)作S-G平滑去噪處理。FOD作為一種光譜信號處理技術(shù)[22]，通過細化階數(shù)改變光譜曲線中的斜率和曲率，使不同分數(shù)階下光譜曲線間的特征波段處微弱差異得到明顯放大[25]。本研究以整數(shù)階微分定義推廣出的G-L算法進行分數(shù)階微分計算，計算式為

(1)

式中n——上下限微分階數(shù)之差

f(λ)——以波長λ處反射率為自變量的函數(shù)

v——階數(shù)

本研究分數(shù)階設(shè)置為0～2階，以0.20階為間隔，共11階。

1.4 優(yōu)化光譜指數(shù)

基于高光譜遙感信息進行表層土壤屬性快速反演，其中光譜信息多以連續(xù)窄波段為主。本文在全波段400～2 400 nm范圍，通過對差值指數(shù)、比值指數(shù)、歸一化指數(shù)(Normalized difference index，NDI)、再歸一化差值指數(shù)(Renormalized difference index，RDI)和廣義差值指數(shù)(Generalized difference index，GDI)進行兩兩算法優(yōu)化(表1)，建立全波段(400～2 400 nm)與SOM含量間的二維相關(guān)性圖。

表1 優(yōu)化光譜指數(shù)及公式

1.5 建模方法與精度評價

SVM模型基于結(jié)構(gòu)風(fēng)險最小化處理，能夠很好地解決模型中存在的局部最優(yōu)問題[26]。本研究選用決定系數(shù)R2、均方根誤差(Root mean square error，RMSE)和RPD作為所建模型的精度和預(yù)測效果評價。一般將RPD分為3類來評價模型的可靠性：RPD小于等于1.4，表明模型預(yù)測能力極差；RPD大于1.4、小于等于2.0，表明模型能夠具有較為可靠預(yù)測能力；RPD大于2.0，表明模型預(yù)測能力極好[23]。

2 結(jié)果與分析

2.1 銀川平原土壤有機質(zhì)含量特征

根據(jù)全國第二次土壤普查推薦的土壤肥力分級標準，將銀川平原187個土樣按SOM含量共分為6個等級(表2)：SOM含量為2.81～56.14 g/kg，均值12.48 g/kg、變異系數(shù)57.26%，其中處于四級～六級(SOM含量小于20 g/kg)的土樣有174個，占總數(shù)的93.05%，該范圍SOM含量平均值為11.06 g/kg，變異系數(shù)為39.30%，本文選取這174個土樣為研究對象。

表2 銀川平原土壤SOM統(tǒng)計特征

利用K-S算法[27]中的歐氏距離劃分總樣本集，其中建模和驗證樣本分別為117和57(表3)[23]，所劃分的建模數(shù)據(jù)和驗證數(shù)據(jù)較總樣本SOM含量區(qū)間值和標準差較為一致，變異系數(shù)說明SOM均屬中度變異。

表3 建模與驗證土壤樣本有機質(zhì)含量統(tǒng)計特征

2.2 基于原始光譜反射率的分數(shù)階微分光譜特征

由圖2可知，原始光譜反射率吸收帶整體較寬，吸收特征差異明顯，其中1 400、1 900 nm處存在明顯吸收峰。在400～1 800 nm范圍內(nèi)，光譜反射率隨波長的增加而增大；在1 900～2 200 nm范圍內(nèi)，光譜重疊嚴重。隨著分數(shù)階的不斷增加，重疊峰和基線漂移逐漸消除，光譜強度的幅度下降，光譜反射率不斷趨近于0，其中1 900 nm處的吸收谷逐漸由正峰變?yōu)樨摲澹唤t外波段玫紅色區(qū)域的光譜標準差差異明顯大于可見光區(qū)域。當(dāng)FOD從0階增加到1.6階時，土壤光譜反射率逐漸降低并最終穩(wěn)定在-0.07～0.07；當(dāng)從1.6階增加到2.0階時，土壤光譜反射率在-0.06～0.06穩(wěn)定不變，且高分數(shù)階微分的吸收谷明顯小于低分數(shù)階微分。

圖2 土壤反射率分數(shù)階微分曲線

2.3 不同分數(shù)階下的光譜指數(shù)與土壤有機質(zhì)含量相關(guān)性及敏感波段優(yōu)選

比較土壤NDI/RDI、DI/NDI、DI/RDI、RDI/NDI、DI/GDI和RI/GDI在不同分數(shù)階微分下的最大相關(guān)系數(shù)絕對值(Maximum absolute correlation coefficient，MACC)可以看出(圖3)，優(yōu)化指數(shù)DI/RDI與SOM的MACC顯著高于其他指數(shù)，平均值高達0.997 6，其次為NDI/RDI和RDI/NDI，均值分別為0.740 6和0.730 2，DI/NDI最低，MACC均值為0.309 3。同時，DI/RDI、NDI/RDI和RDI/NDI與SOM含量間的MACC在0.8階時達到峰值后逐漸下降，此時MACC分別為0.998 6、0.762 1和0.799 3。

圖3 不同分數(shù)階微分下光譜優(yōu)化指數(shù)最大相關(guān)系數(shù)絕對值

為篩選最佳建模變量，本研究以MACC大于等于0.800 0為閾值進行優(yōu)選。土壤DI/NDI、DI/GDI、RI/GDI、NDI/RDI和RDI/NDI在0～2階范圍內(nèi)并無MACC大于等于0.800 0的敏感波段組合存在，故在后續(xù)建模中不作考慮。DI/RDI在0.2～2.0階范圍內(nèi)的MACC為0.996 5～0.998 6，全部高于0.800 0，故在后續(xù)建模中采用DI/RDI指數(shù)。

2.4 優(yōu)化光譜指數(shù)與土壤有機質(zhì)含量二維相關(guān)性分析

圖4為SOM含量與優(yōu)化光譜指數(shù)DI/RDI的二維相關(guān)系數(shù)分布圖，該圖能夠清晰展示相關(guān)系數(shù)較高特征指數(shù)的波段分布區(qū)間。與優(yōu)化指數(shù)NDI/RDI和RDI/NDI的最大MACC相比(分布圖未展示)，DI/RDI最大相關(guān)系數(shù)絕對值整體提高0.25，其敏感波段主要集中在1 450～1 750 nm和 2 100～2 400 nm之間(圖4)，這也進一步說明分數(shù)階微分聯(lián)合優(yōu)化指數(shù)在較低含量SOM的敏感波段篩選中，能夠有效消除土壤水分敏感波段在1 400 nm 和1 900 nm處的干擾[9]。

圖4 最優(yōu)DI/RDI指數(shù)與SOM含量的二維相關(guān)系數(shù)

2.5 基于支持向量機的SOM含量反演模型構(gòu)建

表4 基于DI/RDI-SVM的SOM建模與預(yù)測結(jié)果

2.6 克里格模型與銀川平原SOM含量分布反演

在ArcGIS 10.2中，利用模型殘差結(jié)合普通克里格對實測和0.2階下DI/RDI-SVM模型的結(jié)果分別作反演插值(圖5)，可以看出，銀川平原南部的吳忠市利通區(qū)、青銅峽市和巴音陶亥鎮(zhèn)黃河以西地區(qū)SOM含量明顯高于其他地區(qū)。賀蘭縣整體及惠農(nóng)區(qū)東北部SOM含量較高，而銀川市西夏區(qū)西南部、永寧縣東西兩側(cè)及石嘴山市大武口區(qū)SOM含量整體較低，這是由于該區(qū)域地勢低洼、蒸發(fā)強烈，加之多年引黃灌溉，導(dǎo)致該區(qū)域土壤鹽漬化較重，土壤結(jié)構(gòu)差，保水、保肥性能下降，故SOM含量低下。

圖5 銀川平原土壤有機質(zhì)含量實測值和0.2階下DI/RDI-SVM反演值插值圖

與實測值相比，四級水平(SOM含量 10～20 g/kg)SOM的反演值所占面積比實測值高2.80%(表5)，而五級(SOM含量6～10 g/kg)和六級(SOM含量小于6 g/kg)水平面積低了2.07%和0.73%，說明該方法反演的銀川平原SOM在相對較高水平反演值略高于實測值，而低水平反演值則略低于實測值。

表5 銀川平原不同等級有機質(zhì)實測值與反演值所占面積

3 討論

光譜指數(shù)能夠從二維光譜空間凸顯SOM的響應(yīng)特征，充分利用高光譜信息，減少其他土壤信息對SOM的影響，降低模型的復(fù)雜性、去除冗余信息變量[29]。郭燕等[21]分析了NDI、DI、RI與SOM的二維相關(guān)性，最大相關(guān)系數(shù)在0.5～0.7之間；HONG等[22]基于不同光譜變換的光譜指數(shù)與SOM的二維相關(guān)性，得到最大相關(guān)系數(shù)為0.810 0，比一維光譜與SOM的相關(guān)性有了顯著提升。尚天浩等[23]采用單一光譜指數(shù)對銀川平原SOM估算，發(fā)現(xiàn)RDI的MACC最高，MACC可達0.801 0，RDI-SVM模型估測精度最高，RPD為2.32。目前大多應(yīng)用RI、DI和NDI等單一光譜指數(shù)來反演土壤屬性，而兩兩組合很少使用[15]。本研究基于5種光譜指數(shù)的兩兩算法優(yōu)化，發(fā)現(xiàn)NDI/RDI和RDI/NDI的MACC最大值分別為0.762 1和0.769 3，說明光譜指數(shù)NDI和RDI無論作何算法，其比值優(yōu)化波段指數(shù)間的MACC范圍都差異不大；與此相比，DI/RDI的MACC顯著提升20%以上，最高達0.998 6，表明在土壤低有機質(zhì)含量條件下，DI/RDI指數(shù)作算法優(yōu)化時的作用優(yōu)于單一NDI或RDI指數(shù)。這也說明SOC對優(yōu)化光譜指數(shù)的靈敏度明顯高于使用單一光譜指數(shù)時的靈敏度，從而突出了使用光譜指數(shù)組合的優(yōu)勢[7]。在SOM估測建模方法中，SVM作為一種有效監(jiān)督機器學(xué)習(xí)的方法，通過對數(shù)據(jù)的智能化處理和數(shù)據(jù)價值的充分挖掘，在一定程度上解決了數(shù)據(jù)處理過程中的“過學(xué)習(xí)”和“離散值多”問題[23]。本研究利用分數(shù)階微分聯(lián)合不同優(yōu)化指數(shù)所建SVM模型，預(yù)測RPD普遍高于2.00，其中DI/RDI-SVM模型RPD高達4.31，這是因為SVM為機器算法，其模型在數(shù)據(jù)驅(qū)動中具有極強的自學(xué)習(xí)能力，能夠準確獲取土壤光譜中極為復(fù)雜的非線性特征[33]，SVM模型依靠核函數(shù)，經(jīng)過多次訓(xùn)練后，將輸入數(shù)據(jù)繪制到新的超空間中，在該超空間中執(zhí)行分離，最終使用ε不敏感損失函數(shù)獲得用于數(shù)據(jù)擬合和預(yù)測的最優(yōu)超空間[17]，該超空間可以容忍小于常數(shù)ε集作為閾值的誤差[34]，使得模型擬合精度不斷接近100%[35]。ZHANG等[20]還基于波段優(yōu)化算法和兩波段指數(shù)形式構(gòu)建了一個新的三波段指數(shù)——修正歸一化差分指數(shù)，并指出該三維指數(shù)在利用可見-近紅外光譜估計土壤其他生化參數(shù)方面具有很強的應(yīng)用潛力。本研究后期也可以嘗試采用三維光譜指數(shù)建立SOM估算模型。此外，還需考慮銀川平原土壤屬性中的水分、鹽分、pH值和養(yǎng)分等因素對SOM光譜反演的影響。

4 結(jié)論

(1)銀川平原土壤樣本SOM含量平均值為12.48 g/kg，屬于高度變異，其中93.05%的土樣SOM處于四級～六級水平。

(2)土壤野外原始光譜反射率吸收帶整體較寬，吸收特征差異明顯，在1 400 nm和1 900 nm處有明顯吸收峰。隨著分數(shù)階的不斷增加，光譜反射率不斷趨近于0；高分數(shù)階微分的吸收谷明顯小于低分數(shù)階微分。

(3)優(yōu)化光譜指標MACC從大到小依次為DI/RDI、NDI/RDI、RDI/NDI、DI/GDI、DI/NDI和RI/GDI，敏感波段主要集中在1 450～1 750 nm和2 100～2 400 nm。