基于光譜指數和機器學習的土壤有機質含量反演

朱傳梅 王宏衛 謝霞 馬利剛 仝雁軍 古力孜熱·買買提

摘要：土壤有機質（SOM）是衡量礦區生態環境質量的重要指標之一。為了提高SOM含量的估算精度，在已有二波段指數的基礎上加入第3個波段，構建新的三波段指數，利用極限學習機（ELM）和隨機森林（RF）分別建立SOM含量的預測模型。在新疆準東煤田采集168個土壤樣點，在室內進行SOM含量、光譜的測定。對光譜數據進行平滑和預處理后，構建多個兩波段、三波段光譜指數，隨后分析不同維度光譜數據與SOM含量的敏感程度和敏感區域。ELM和RF被用于對每個維度最優光譜參數建立預測模型。研究結果顯示，無論采用哪種方式建模，每個維度的光譜數據與SOM含量的敏感程度和建模精度均隨信息維度的增加而增加，即三波段指數（TBI）>二波段指數>一維光譜數據。在三波段指數中，ELM的預測效果要優于RF，其中（TBI-4）-ELM的預測效果最好，決定系數（r2）=0.87，均方根誤差（RMSEP）=4.07，相對分析誤差（RPD）=2.63。三波段指數與ELM的結合，可以很好地減弱土壤信息噪聲，提高SOM含量的預測精度。

關鍵詞：遙感;光譜分析;土壤有機質;光譜指數;機器學習

中圖分類號：S153.6+21文獻標志碼：A

文章編號：1002-1302（2020）22-0233-09

作者簡介：朱傳梅（1993—），女，新疆昌吉人，碩士研究生，主要從事綠洲聚落景觀與土壤有機質研究。E-mail：aspiration818@163.com。

通信作者：王宏衛，博士，教授，碩士生導師，主要從事干旱區綠洲聚落生態安全格局與土壤重金屬研究。E-mail：wanghw_777@163.com。

礦產資源的開采和加工可以帶來更多的經濟效益，但是開采煤礦會干擾土層，破壞植被，使土壤失去利用價值，這對土地資源的可持續發展和生態環境構成嚴重威脅[1]。在我國，大型的露天煤礦多集中分布于干旱、半干旱生態脆弱的地帶，該地區土壤自身修復能力極弱，生態敏感性極強，再加上礦產資源長期大量被開采和加工，導致當地環境問題和生態修復問題日益突出[2]。土壤是陸地生態系統中許多生態過程（例如養分循環、水平衡、凋落物分解等）的基礎，土壤有機質狀態是衡量退化生態系統中生態功能恢復和維持的關鍵指標[3]。因此，無損地監測土壤有機質（SOM）含量有助于礦區環境管理和生態恢復。一般的SOM含量測定多基于大量的野外土壤采樣和繁瑣的室內化學分析方法，較費時、費力、耗資，無法滿足現代精細農業發展的需要[3]。與傳統的方法相比，遙感技術是一種很有前景的土壤屬性定量評估方法，具有快速響應、成本低、采集快等特點，可以很好地用來描述、評估不同尺度下表層土壤的各種特征[4]。因此，基于不同的光譜反射和吸收特性，可見光-近紅外光譜（VIS-NIR）分析技術可以作為一種替代方法，保證SOM含量的準確估算。

土壤是由多種物質組成的混合物，其反射光譜常包括背景噪聲、基線漂移、傾斜等干擾信息，直接用來進行SOM含量的估算得到的結果并不理想[5]。相關研究表明，光譜預處理方法能夠較好地移除噪聲、突出光譜曲線特征、去除或減弱其他因素的影響，為建立具有較高精度的估測模型提供可能[5]。光譜微分是增加信噪比的主要技術手段之一，其中一階微分（FD）可去除不同的背景噪聲和基線漂移，擴大樣本間光譜的差異，反映被測物體的本質特征。連續統去除（CR）法可減小散射對目標光譜的影響，同時可放大微弱光譜的吸收特性[5]。這2種預處理技術被廣泛地應用于VIS-NIR分析中，對建立預測能力強、穩健性好的分析模型至關重要。

以往對于SOM的研究是在一維層面（全波段反射率或對應的數學變換）上選取單個敏感波段或多個敏感波段進行建模，該方法僅考慮了SOM與光譜間的關系，并沒有考慮光譜間的重疊吸收或相互影響[3-4]。光譜指數是由幾個窄波段或寬波段組合而成，可通過分析特定波段間的相互作用，提高對待測屬性的敏感程度。Wang等采用最優兩波段指數對土壤鹽分含量進行分析和預測，取得了較顯著的成果[3]。然而，Tian等在水稻葉片氮濃度的定量估測中，對比了兩波段指數和三波段指數的估算能力[4]，這些研究和提出的指標表明，通過兩波段指數評估某些參數存在不足。對于土壤這種組成極為復雜的物質，兩波段指數能否很好地消除或減弱土壤中其他物質產生光的散射和分子的非特征吸收的干擾有待于進一步研究。

機器學習算法在解析非線性問題時表現能力較好，常用于土壤屬性的定量化估測，其中，極限學習機（ELM）和隨機森林（RF）更是研究的焦點。Douglas等在估算土壤中總石油烴（TPH）含量時，發現與線性偏最小二乘回歸（PLSR）法相比，RF模型能更好地反映土壤光譜的非線性響應，從而提供更高的預測精度[6]。然而，ELM和RF能否在較多的土壤信息噪聲中（如嚴重的人為影響）和較低的SOM水平（如干旱區嚴重的荒漠化影響）下，建立SOM含量和光譜參數的聯系，并達到一定的預測精度，有待進一步研究。

本研究的目的：（1）利用波段優化算法，構建新的三波段光譜指數;（2）量化不同維度的光譜參數對SOM的響應;（3）通過比較SOM的預測精度，尋找最有效的建模方法。

1材料與方法

1.1研究區介紹和土壤樣本的制備

研究區為準東煤田，位于我國新疆準噶爾盆地東南緣（43°45′～45°00′N，88°45′～91°10′E）為13000km2的露天煤田，煤炭儲量預估可達到3900億t[2]。它是世界上最大的綜合煤田，被譽為“中國工業糧倉”。該地為極端干燥的大陸性氣候，年平均降水量、溫度分別為140～183mm、5.3～7.3℃。土壤母質為第四紀沖積沉積物，地表植被稀疏。準東煤田的主要土地利用和土地覆被類型為荒地、草地、農田等。自2006年準東煤礦啟動以來，大量的工礦活動已造成了生態失衡和嚴重的環境污染，土壤性質可能正在發生變化。

2014年6月，筆者所在研究團隊收集了168個土壤樣本，期間未發生極端天氣（如大雨或強風）。研究區主要地形為丘陵，因此該試驗的主要設計路線是沿著道路進行。在確保車輛和人員安全的前提下，選擇距離道路較平坦的區域（>300m）作為采樣地點。在每個采樣點，用木鏟收集直徑在10m以內的5個子樣本的表土樣本（土壤深度為0～20cm），并混合成復合樣本（質量約1.5kg）。將復合樣品立即裝入帶標簽的防水自封袋中，并使用手持全球定位系統（GPS）和植被覆蓋率儀記錄樣品的坐標、高程、植被覆蓋率等信息。將樣品送回實驗室進行室內空氣干燥2周（室溫26～28℃），小心去除非土壤物質（礫石、植物根、其他物質），然后用瑪瑙研缽研磨并通過1.5mm篩，以減小粒徑的影響。SOM含量通過重鉻酸鉀外部加熱法[5]進行測定。

1.2光譜數據的采集和預處理

利用ASD地物光譜儀進行土壤樣本的光譜采集，該儀器在350～1000nm與1000～2500nm區間的采樣間隔分別為1.4、2.0nm，重采樣間隔為1.0nm，輸出波段數為2151個。土壤光譜的測定應避免外界光源的干擾，因此試驗在暗室中進行，選用50W的鹵素燈，探頭的視場角為5°，探頭距離土壤表面30cm，每測定10次進行1次白板定標，每個土壤樣本重復測定10次，剔除異常光譜后取其算術平均值作為該樣品的反射率。

位于兩端（350～399nm和2401～2500nm）處的光譜數據通常存在高頻噪聲，因此對其進行剔除。對168個土壤樣本的400～2400nm反射率光譜進行多項式階數為2、平滑窗口為9點的Savitzky-Golay平滑。高光譜數據的維數較高，存在的多余重復信息和相鄰波段間高度相關。已有文獻表明，重采樣高光譜數據能夠降低信息冗余，提高計算的準確度、速度、反演精度[4]。因此，對光譜數據重采樣到10nm，每個光譜由此獲得201個波段，降采樣后的光譜被定義為原始光譜（R）。

本研究引入了2種光譜預處理方法：用于去除地物背景噪聲、修復基線漂移的一階微分和用于消除散射、突出光譜特性的連續統去除法。預處理工作均在MATLAB2018b中完成。

1.3二維、三維相關系數

相關性分析有助于揭示SOM含量與光譜間線性相關關系的強弱，這種關系通常以一維線性數據的形式展示。近年來，許多研究以可視化的形式展示了常見光譜指數與待測屬性的相關關系，該方法不僅將光譜特性從一維擴展到了二維，而且充分考慮了光譜間的相互關系[3-5]，公式如下：

選取特定敏感區域的第三波段疊加在兩波段光譜指數上，往往能夠增加指數的精確性，增強抗干擾性，清除常用兩波段指數存在的易飽和現象等。本研究在所選3種指數的形式基礎上增加第3個波段λ3，使用MATLAB2018b做400～2400nm的全波段循環，并選擇最優三波段指數構建估測模型。本研究具體涉及如下：

1.4建模預測

為了保證建模集和驗證集數據盡可能地涵蓋SOM含量的全部范圍，本研究使用Kennard-Stone（K-S）算法進行樣本集合的劃分。該算法根據已選擇對象的重要主成分數量來將歐氏距離（Euclideandistance）最大化。168個土壤樣本被分成2個部分，其中112個樣本用于建模，剩余56個樣本用于驗證建模結果。

1.4.1極限學習機

極限學習機（ELM）是由Huang等提出的基于單隱含層前饋神經網絡的機器算法[7]。與傳統的神經網絡模型例如反向傳播前饋（BP）神經網絡模型相比，其學習方式不同，BP神經網絡模型利用梯度下降算法，通過反向傳播方式進行學習，在學習過程中須要不斷地進行迭代來更新權值和閾值。而ELM僅須要添加隱含層節點的個數，利用正則化計算網絡的輸出權值來達到學習的目的，網絡中輸入層和輸出層的權值和閾值通過隨機初始化得到，且不影響網絡的收斂能力。因此，它的學習速度比BP等提高了數千倍，且不會出現過擬合現象。

ELM中，連接層的激活函數采用更加接近生物學激活模型的softplus函數代替傳統的sigmoid函數，已有在不同數據庫的大量實驗證明了softplus激活函數具有更好的泛化性能。本研究通過將步長從3增加到150，每次增加3個步長來確定隱含層的最優節點數，每個模型重復50次以減少隨機誤差。本研究利用MATLAB2018b對Huang等編寫的ELM源碼[7]進行了部分修改和調試。

1.4.2隨機森林

隨機森林（RF）是用于分類和回歸的一種機器學習算法。RF以決策樹學習和簡單平均算法為基礎，根據每個二叉樹上的節點數（m）和自助法（Bootstrap）抽換選取n個樣本構建決策樹，利用未被選取的樣本對每顆樹進行預測，由于RF隨機選擇特征和變量，使模型不容易陷入過度擬合。

為了建立RF模型，對模型中二叉樹節點處的變量個數（mtry）和決策樹數量（decisiontree）進行調整。選擇均方差（MSE）的倒數作為適應度函數值，即最優模型的適應度函數值越小，它與期望值間的誤差就越小。mtry以1～9為間隔，每次移動的步長為1，決策樹數量以100～2000為間隔，每次移動的步長為100。

1.5模型的精度檢驗

本研究從擬合程度、估算能力、準確性等3個方面對模型進行驗證。模型的擬合度用決定系數（r2）評價，r2越接近1說明模型的擬合程度越好。模型的估算能力用驗證集均方根誤差（RMSEP）評價，RMSEP越小說明模型的估算能力越好。估算模型的準確性用相對分析誤差（RPD）來評價。

2結果與分析

2.1不同的SOM含量和反射光譜曲線

表1顯示，168個樣點的SOM含量為0.255～45.708g/kg，平均值為7.461g/kg，標準差為8.747g/kg，變異系數為117.23%，為強變異，這表明準東煤礦SOM含量空間變異較大，這可能由不同土地利用方式和人為影響造成的。建模集和驗證集間的F-test得出P=0.72>0.05，表明集合的劃分比較相似，但差異不顯著。

由原始光譜曲線（圖1-a）可知，研究區不同SOM含量的光譜曲線形態基本一致，該地區土壤光譜在可見光波段陡峭上升，反射率與SOM含量呈明顯的負相關關系，即SOM含量越高反射率越低，土壤光譜曲線在1400、1900、2200nm附近都存在水分吸收谷，在1900nm處較為明顯[5]。如圖1-b所示，在580、1350、1900nm附近分別出現了不同的正負峰值，可見一階微分有利于顯示原始光譜的部分肩峰。如圖1-c所示，反射率經過連續統去除吸收特征得到了明顯的放大，其中500、1900nm處的為弱吸收帶在連續統去除曲線中可以觀察出來，而在反射率曲線中則不明顯。

2.2一維相關性分析

由圖2可知，R與SOM含量的相關性曲線較為平滑，且沒有波段通過0.01水平上的顯著性檢驗，說明R與SOM含量的敏感性較低。然而R經過FD處理后，顯著性明顯提升，尤其是在640～1110nm處，最大相關系數在840nm處，為0.52，因為該波段附近存在C—H吸收帶，與SOM含量直接相關[5]。R經過CR處理后，在1110～1380nm出現了1個相關系數高峰，相關系數為0.3～0.4，使原來弱的吸收特征得到了增強。600、840、2250nm分別為R、FD、CR相關系數絕對值最大的波段，本研究提取了不同處理中相關性最強的波段數據，用于后續的對比分析和估算建模（表2）。

2.3二維相關性分析

圖3為基于不同光譜處理得到的兩波段光譜指數（NDI、RI、DI）的二維相關系數圖，圖中橫軸代表λ1，[CM（20*2]波長范圍為400～2400nm，豎軸代表λ2，波長范圍為400～2400nm，右邊的顏色軸代表相關系數值到顏色圖的映射，顏色軸的上限和下限為最大正相關系數和最大負相關系數。與圖2相比，二維相關系數值的閾值范圍均大于一維相關系數。其中以SOM含量與CR的差值指數（DI）的相關性為例（圖3-i），最大相關系數與（圖2）中CR的一維最大相關系數相比提高了0.12。R和CR構建的兩波段指數其相關性較好的波段主要集中在可見光和短波近紅外區域（圖3-a、圖3-b），這與Zhang等的研究結果[5]一致。在每個兩波段指數中，FD對SOM含量的敏感性最強，其中FD-NDI、FD-DI（圖3-b、圖3-h）得到的相關系數較大，約為0.55。提取每個指數的最優波段組合并計算其與SOM含量的相關系數，結果列于表2中。

2.4三維相關性分析

為了充分挖掘光譜數據，根據前人的研究[5]進展，本研究在兩波段（λ1、λ2）指數的基礎上增加第3個波段λ3，構建新的三波段指數。三波段指數在兩波段光譜指數的基礎上放大了與SOM含量相關系數的閾值范圍（圖4至圖7）。同時，本研究發現在R所構建的三波段指數中（圖5至圖7），最優指數的波長基本都在可見光和近紅外短波范圍內，820～850nm為主要敏感區，已有研究發現820nm附近存在C—H吸收帶，這與SOM含量直接相關[5];在FD中（圖4至圖7），最優指數集中在長波近紅外范圍，每個指數中均有1個或多個波靠近1450、1950nm的水分吸收帶。同時，除TBI-4外，其余3個三波段指數中均有1個波段靠近2300nm，已有研究證明2300nm附近存在1個與SOM含量相關的C—H特征峰[8-9]，在每個三波段指數中，FD的敏感性均最強，但從切片圖（圖4至圖7中的d、e、f）中來看，敏感區域較為分散，不集中;在CR中（圖4至圖7），780～1100nm的短波近紅外為敏感區域，880、950、1020nm為主要的敏感波段，在TBI-2和TBI-3中均有出現。在所有三波段指數中，TBI-4（FD2020，FD2360，FD1150）的效果最優，最大相關系數為0.65，相比于FD的兩波段指數（圖3-b、圖3-e、圖3-h）提升了0.1。因此，優化后的三波段指數是估測SOM含量較為有前景的指標，同時也表明兩波段指數存在一定的不足。提取每個三波段指數的最優波段組合并計算其與SOM含量的相關系數，結果列于表2中。

2.5建模與預測

利用ELM和RF2種機器算法分別建立SOM含量預測模型，不同維度的最優光譜參數（表2）和建模集中的SOM含量作為機器算法的訓練樣本（表1），驗證集中的SOM含量作為機器算法的預測樣本（表1）。

結果表明，無論采用哪一種模型，一維光譜數據的預測效果明顯不如多維指數的建模預測結果，r2（0.55，0.51）、RPD（1.41，1.36）較低，RMSE（7.57，7.84）相對較高（表3），雖然ELM的估算精度略高于RF，但該模型僅具有區分樣本數據含量高低的能力，量化能力較弱。

兩波段指數的建模效果相較一維光譜數據有所提升，r2為0.59～0.72，RMSE為5.89～7.09，RPD為1.51～1.82。三波段指數的估算效果均優于前兩者，r2為0.78～0.87，RMSE為4.07～5.2，RPD為2.03～2.63。此外，在三波段指數的建模結果中，ELM對SOM含量的估算能力普遍優于RF。其中使用三波段指數4-ELM的預測效果最好（r2=0.87、RMSE=4.07、RPD=2.63），其擬合程度、估算能力和準確性均優于本研究中的所有模型。

3結論與討論

光譜指數根據地物的光譜特性，對波段進行簡單組合，以達到對地表參量簡單有效的度量[10]。前人研究發現由于研究人員自身客觀條件的局限性，研究所用的試驗數據往往來自某一特定地區，使光譜指數往往具有一定的局限性[11-14]。本研究參考植被光譜分析方法，采用波段優化算法分析了在400～2400nm范圍內任意2個波段的NDI、RI、DI與SOM含量間的關系，發現其對SOM含量的響應大大提升。

本研究發現，光譜指數法可在一維線性數據的基礎上更為深層次挖掘光譜信息，充分考慮波段間的相互關系和提高建模預測的精度。通過在已有的兩波段指數中添加第3個敏感波段，構建新的三波段指數，發現三波段光譜指數放大了與SOM含量相關系數的閾值范圍，其中FD的三波段指數4效果最優，最大相關系數為0.65。利用ELM和RF對不同維度的最優光譜參數分別建立預測模型，發現無論采用哪種模型，預測精度均跟據光譜參數維度的增加而增加，并且在三波段指數中，ELM的建模效果均優于RF。其中，三波段光譜指數4-ELM的預測能力最優，r2=0.87，RPD=2.63。該研究結果對減少土壤中其他信息噪聲的干擾和土壤屬性的定量估算具有一定的參考價值，為估測土壤中其他屬性提供了新的思路。

參考文獻：

[1]BaoNS，WuLX，YeBY，etal.Assessingsoilorganicmatterofreclaimedsoilfromalargesurfacecoalmineusingafieldspectroradiometerinlaboratory[J].Geoderma，2017，288：47-55.

[2]彭麗萍，戴岳，師慶東.新疆準東荒漠區5種典型植物水分來源[J].干旱區研究，2018，35（5）：1146-1152.

[3]WangXP，ZhangF，DingJL，etal.Estimationofsoilsaltcontent（SSC）intheEbinurLakeWetlandNationalNatureReserve（ELWNNR），NorthwestChina，basedonaBootstrap-BPneuralnetworkmodelandoptimalspectralindices[J].ScienceoftheTotalEnvironment，2018，615：918-930.

[4]TianYC，YaoX，YangJ，etal.Assessingnewlydevelopedandpublishedvegetationindicesforestimatingriceleafnitrogenconcentrationwithground-andspace-basedhyperspectralreflectance[J].FieldCropsResearch，2011，120（2）：299-310.

[5]ZhangZP，DingJL，WangJZ，etal.PredictionofsoilorganicmatterinnorthwesternChinausingfractional-orderderivativespectroscopyandmodifiednormalizeddifferenceindices[J].CATENA，2020，187：104257.

[6]DouglasRK，NawarS，CipulloS，etal.Evaluationofvis-NIRreflectancespectroscopysensitivitytoweatheringforenhancedassessmentofoilcontaminatedsoils[J].ScienceoftheTotalEnvironment，2018，626：1108-1120.

[7]HuangGB，DianHW，YuanL.Extremelearningmachines：asurvey[J].InternationalJournalofMachineLearning&Cybernetics，2011，2：107-122.

[8]史舟，王乾龍，彭杰，等.中國主要土壤高光譜反射特性分類與有機質光譜預測模型[J].中國科學：地球科學，2014，44（5）：978-988.

[9]鄭曼迪，熊黑鋼，喬娟峰，等.基于綜合光譜指數的不同程度人類干擾下土壤有機質含量預測[J].江蘇農業學報，2018，34（5）：1048-1056.

[10]張子鵬，丁建麗，王敬哲.基于諧波分析算法的干旱區綠洲土壤光譜特性研究[J].光學學報，2019，39（2）：391-401.

[11]張賢龍，張飛，張海威，等.基于光譜變換的高光譜指數土壤鹽分反演模型優選[J].農業工程學報，2018，34（1）：110-117.

[12]尼加提·卡斯木，茹克亞·薩吾提，師慶東，等.基于優化光譜指數的土壤有機質含量估算[J].農業機械報，2018，49（11）：155-163.

[13]章濤，于雷，易軍，等.高光譜小波能量特征估測土壤有機質含量[J].光譜學與光譜分析，2019，39（10）：3217-3222.

[14]徐麗華，謝德體.土壤有機質含量預測精度對光譜預處理和特征波段的響應[J].江蘇農業學報，2019，35（6）：1340-1345.