張清文，楊曉蕪，楊 睿，尹 軒，寇財(cái)垚

（1.華北理工大學(xué) 礦業(yè)工程學(xué)院，河北 唐山 063210；2.國(guó)家知識(shí)產(chǎn)權(quán)局專利局專利審查協(xié)作廣東中心，廣東 廣州 510700）

引言

土壤有機(jī)碳（Soil organic carbon，SOC）和土壤全氮（Soil total nitrogen，STN）作為土壤質(zhì)量和健康的主要指標(biāo)，在保持土壤質(zhì)量方面發(fā)揮著重要作用[1]。SOC 與STN 的比值（C/N）也是土壤質(zhì)量和肥力的重要指標(biāo)，反映了SOC 與STN 之間的相互作用和耦合關(guān)系[2]。此外，C/N 是影響土壤微生物群落的主要因素，因此在陸地碳氮循環(huán)中起著關(guān)鍵作用[3]。研究土壤C/N 對(duì)于土壤管理、生態(tài)環(huán)境監(jiān)測(cè)和揭示生態(tài)系統(tǒng)碳氮耦合關(guān)系有重要意義。

傳統(tǒng)的化學(xué)分析方法測(cè)定土壤C/N 需要分別測(cè)定土樣的SOC 含量和STN 含量，再進(jìn)行計(jì)算C/N，雖然精度較高，但是土壤樣品采集勞動(dòng)強(qiáng)度大，實(shí)驗(yàn)室測(cè)定人力財(cái)力消耗多[4]。近年來(lái)，高光譜技術(shù)的發(fā)展對(duì)于快速，低成本且高精度測(cè)定土壤各類(lèi)理化性質(zhì)提供了新的思路。目前，很多學(xué)者已經(jīng)使用實(shí)測(cè)高光譜對(duì)SOC 含量和STN 含量進(jìn)行精確的估測(cè)[5,6]，陳秋實(shí)利用實(shí)測(cè)土壤光譜反射率結(jié)合互花米草不同入侵年限探究了互花米草入侵下的濱海濕地SOC 含量，模型決定系數(shù)（Coefficient of determination，R2）最高達(dá)到0.81[7]，盧艷麗對(duì)農(nóng)田耕層土壤的STN 含量進(jìn)行建模分析，模型R2最高達(dá)到0.863[8]。但是將土壤光譜反射率作為自變量，探究土壤的SOC 含量和STN 含量的比值關(guān)系相關(guān)研究較少。

高光譜數(shù)據(jù)波段信息豐富，每個(gè)波段可以看作是土壤的一種特征，這些特征可以捕捉不同波段下土壤的化學(xué)和物理性質(zhì)，對(duì)土壤理化性質(zhì)的估測(cè)有重要意義[9]，但波段間反射率信息具有強(qiáng)共線性，即部分波段信息存在冗余。前人多使用Pearson 相關(guān)系數(shù)確定敏感波段[10]，李焱等選取了相關(guān)系數(shù)最大的10 個(gè)波段作為自變量，成功預(yù)測(cè)了新疆地區(qū)不同地表利用類(lèi)型的STN 含量[11]，于雷等比較了全波段與敏感波段建立的預(yù)測(cè)模型精度，發(fā)現(xiàn)敏感波段建立的模型精度略低于全波段建立模型的精度，但全波段模型的復(fù)雜度要遠(yuǎn)高于敏感波段預(yù)測(cè)模型的復(fù)雜度[10]。目前，降維和特征篩選是處理波段間的多重共線性有效方法[12]，基于這兩種思想，機(jī)器學(xué)習(xí)中的偏最小二乘回歸（Partial least square regression，PLSR）和基于決策樹(shù)的方法得到研究人員廣泛重視。

濱海濕地土壤含水量大，可達(dá)性差，野外采樣難度大。本文基于有限野外采樣點(diǎn)的SOC 和STN含量數(shù)據(jù)和實(shí)測(cè)高光譜數(shù)據(jù)，并對(duì)原始光譜進(jìn)行光譜變換以挖掘土壤信息，通過(guò)降維和特征篩選的方法確定與土壤C/N 相關(guān)的敏感波段，最后基于PLSR 和隨機(jī)森林回歸 （Random forest regression，RFR）建立濱海濕地土壤C/N 的光譜預(yù)測(cè)模型。以期為濕地碳氮循環(huán)提供幫助。

1 材料與方法

1.1 土壤數(shù)據(jù)采集分析

研究區(qū)位于江蘇省鹽城市大豐區(qū)東北部沿海，該地區(qū)年降水量為980～1070mm，多集中于夏季，年平均氣溫在13.7～14.6℃之間[13]。研究區(qū)蘆葦與堿蓬等本土植物逐漸被外來(lái)種互花米草替代，已形成面積廣泛的互花米草單一優(yōu)勢(shì)植物群落，受互花米草入侵與擴(kuò)散的影響，原有泥灘逐漸演變?yōu)辂}沼濕地[14]。本研究考慮到該區(qū)范圍內(nèi)互花米草海向入侵特征，垂直于海岸線設(shè)計(jì)采樣方案，每隔40 米進(jìn)行土壤樣本采集，記錄每個(gè)采樣點(diǎn)的經(jīng)緯度坐標(biāo)，土壤類(lèi)型。每個(gè)樣點(diǎn)處采集3 個(gè)表層（0～30cm）土壤樣品，充分混合，形成復(fù)合樣品，共獲取33 個(gè)土樣，樣點(diǎn)分布如圖1 所示。

圖1 研究區(qū)概況圖

為減少外界因素對(duì)土樣的影響，土樣室內(nèi)自然風(fēng)干，剔除異物，研磨過(guò)2mm 篩，分為3 份，分別用于凱氏蒸餾法測(cè)定STN 含量[15]，稀釋熱重鉻酸鉀容量法測(cè)定SOC 含量[16]，高光譜數(shù)據(jù)采集。土壤C/N統(tǒng)計(jì)特征如表1 所示，數(shù)據(jù)集偏度為-18.312%，為小偏度數(shù)據(jù)集且數(shù)據(jù)分布左偏。

表1 土樣C/N 描述性統(tǒng)計(jì)

1.2 高光譜測(cè)定及預(yù)處理

利用ASD FieldSpec 4 便攜式地物光譜儀測(cè)定土壤光譜反射率，其測(cè)量范圍是350～2500nm，采樣間隔為1.4nm（350～1000nm）和1.1nm（1001～2500 nm），重采樣間隔為1nm，共輸出2151 個(gè)波段。將過(guò)篩的土壤樣品放置于深2cm，長(zhǎng)寬為5cm 的樣品盒中，均勻平鋪，接觸式探頭獲取該土樣的光譜反射率。在每次測(cè)量前進(jìn)行白板校正，以減小光線漫反射影響，每個(gè)土樣采集10 條光譜曲線，ViewSpecPro計(jì)算10 條曲線平均值得到該土樣的光譜反射率。

為消除外界環(huán)境因素造成的噪聲干擾，本研究采用Savitzky-Golay 算法（多項(xiàng)式階數(shù)為2，窗口點(diǎn)數(shù)為20）對(duì)光譜反射率曲線進(jìn)行平滑處理。為充分挖掘光譜信息，對(duì)平滑后的原始光譜（R）進(jìn)行數(shù)學(xué)變換，包括倒數(shù)（1/R），對(duì)數(shù)（LGR），一階微分（R′）和二階微分（R"）[17,18]。

1.3 數(shù)據(jù)分析方法

1.3.1 偏最小二乘回歸

PLSR 是一種多元回歸分析方法，與多元線性回歸類(lèi)似，但不同之處在于PLSR 在建立模型時(shí)同時(shí)考慮了自變量和因變量之間的關(guān)系，即通過(guò)構(gòu)建主成分來(lái)最大化自變量和因變量之間的協(xié)方差[10]。對(duì)原始數(shù)據(jù)進(jìn)行正交變換構(gòu)建的主成分是原始自變量的線性組合，且主成分之間是線性無(wú)關(guān)的，所以PLSR 適用于處理自變量間多重共線性的數(shù)據(jù)。PLSR 以提取的主成分代替原始變量參與回歸建模，以達(dá)到對(duì)高維度自變量數(shù)據(jù)降維的目的[19]。綜上，PLSR 可以解決波段間的共線性問(wèn)題和樣本數(shù)少于光譜波段數(shù)的問(wèn)題，適用處理高光譜數(shù)據(jù)。本研究依據(jù)留一法（Leave-one-out，LOO）交叉驗(yàn)證計(jì)算的預(yù)測(cè)均方根誤差 （Root mean square error of prediction，RMSEP） 確定主成分?jǐn)?shù)量，RMSEP 值隨主成分個(gè)數(shù)的增加先遞減，達(dá)到最低點(diǎn)后隨主成分個(gè)數(shù)的增加出現(xiàn)微小上升或波動(dòng)，為防止模型過(guò)擬合，提取最小的RMSEP 對(duì)應(yīng)的主成分作為最終建模變量。PLSR 模型使用R 語(yǔ)言中的“pls”包構(gòu)建。

1.3.2 Boruta 特征篩選

Boruta 特征篩選的內(nèi)核是隨機(jī)森林算法。首先Boruta 算法會(huì)先根據(jù)原始數(shù)據(jù)集創(chuàng)建一個(gè)影子數(shù)據(jù)集，影子數(shù)據(jù)集與原始數(shù)據(jù)集有相同因變量，但每一個(gè)自變量被隨機(jī)打亂，使用兩個(gè)數(shù)據(jù)集同時(shí)進(jìn)行隨機(jī)森林建模，比較模型中的變量重要性，如果同一自變量的重要性在兩個(gè)模型中有顯著差異，則說(shuō)明此變量為預(yù)測(cè)土壤C/N 的敏感變量，反之，則為非敏感變量，最后輸出所有敏感變量，所有的敏感變量被用于RFR 建模[2]。Boruta 特征篩選使用R語(yǔ)言中的“Boruta”包完成。

1.3.3 隨機(jī)森林回歸

RFR 采用了集成學(xué)習(xí)的思想，由多個(gè)決策樹(shù)組成，每個(gè)決策樹(shù)是一個(gè)弱學(xué)習(xí)器，通過(guò)對(duì)這些決策樹(shù)的集成，隨機(jī)森林對(duì)于高維數(shù)據(jù)和復(fù)雜關(guān)系的建模效果良好，能夠捕捉非線性、交互作用和噪聲[20]。首先隨機(jī)森林從訓(xùn)練數(shù)據(jù)集中進(jìn)行有放回的隨機(jī)抽樣，生成多個(gè)不同的訓(xùn)練數(shù)據(jù)子集。隨機(jī)森林在每個(gè)訓(xùn)練數(shù)據(jù)子集上隨機(jī)選擇一個(gè)特征子集，隨機(jī)森林對(duì)每個(gè)訓(xùn)練數(shù)據(jù)子集和特征子集的組合構(gòu)建一個(gè)決策樹(shù)。在得到所有決策樹(shù)之后，對(duì)于新的輸入數(shù)據(jù)，每個(gè)決策樹(shù)產(chǎn)生一個(gè)預(yù)測(cè)結(jié)果，取所有決策樹(shù)預(yù)測(cè)結(jié)果平均值，作為最終的模型預(yù)測(cè)結(jié)果。所以RFR 對(duì)數(shù)據(jù)異常值不敏感，具有較高的預(yù)測(cè)準(zhǔn)確性和穩(wěn)定性。隨機(jī)森林建模中需要定義兩個(gè)參數(shù)： 每次拆分時(shí)隨機(jī)抽樣作為候選特征的數(shù)量（mtry）、森林中生成的決策樹(shù)數(shù)量（ntree），RFR 模型基于R 語(yǔ)言的“caret”包使用網(wǎng)格搜索法調(diào)整最優(yōu)參數(shù)，采用均方根誤差（Root mean square error，RMSE）最小的參數(shù)組合構(gòu)建最終模型。

1.4 精度評(píng)價(jià)指標(biāo)

本研究建立模型使用LOO 交叉驗(yàn)證以防止過(guò)擬合問(wèn)題，LOO 交叉驗(yàn)證使用n-1 個(gè)樣本作為訓(xùn)練集（n 為樣本數(shù)），沒(méi)有參與訓(xùn)練的1 個(gè)樣本作為測(cè)試集，幾乎用到了所有樣本信息，共循環(huán)n 次，保證每個(gè)樣本都被用作測(cè)試集，其充分地利用了每個(gè)樣本，更適用于小樣本數(shù)據(jù)集，驗(yàn)證結(jié)果更具代表性。模型的驗(yàn)證指標(biāo)采用預(yù)測(cè)值和實(shí)測(cè)值的R2，RMSE 和平均絕對(duì)誤差（Mean absolute error，MAE）。R2越接近于1，RMSE 和MAE 越小說(shuō)明模型精度越高。計(jì)算公式為：

2 結(jié)果與分析

2.1 不同光譜變換的曲線特征和敏感變量

對(duì)R 反射率進(jìn)行數(shù)學(xué)變換，各變換后的光譜反射率如圖2 所示。各樣本R 反射率曲線波形基本相似，在可見(jiàn)光波段（350～780nm）光譜反射率迅速上升，在近紅外波段（780～1900nm）上升趨勢(shì)減緩，在2100nm 處反射率達(dá)到峰值，之后平緩下降。在1410nm 和1920nm 處出現(xiàn)明顯的水分吸收谷。經(jīng)過(guò)對(duì)數(shù)變換后，曲線形態(tài)與變換前基本一致，但相比于R 反射率曲線，各樣本之間的曲線更加聚攏。經(jīng)倒數(shù)變換后的光譜反射率變化趨勢(shì)與R 反射率曲線完全相反且更加聚集，但曲線趨勢(shì)發(fā)生改變的節(jié)點(diǎn)未發(fā)生變化。R′反射率反映了原始光譜反射率的增減性，經(jīng)一階微分變換后，放大了原始光譜曲線的上升和下降等特征，曲線波動(dòng)劇烈。二階微分反映了原始反射率曲線的凹凸性，二階微分光譜反射率在0 周?chē)?dòng)。

圖2 土樣的原始與數(shù)學(xué)變換后的光譜反射率曲線

通過(guò)PCA 確定PLSR 建模使用的主成分個(gè)數(shù)，選擇RMSEP 最小時(shí)對(duì)應(yīng)的主成分個(gè)數(shù)，R，LGR，1/R，R′和R"的主成分個(gè)數(shù)分別為2，2，2，4 和1。表明經(jīng)簡(jiǎn)單數(shù)學(xué)變換（LGR，1/R），主成分個(gè)數(shù)未曾發(fā)生變化，經(jīng)過(guò)一階微分后，主成分個(gè)數(shù)提高為4 個(gè)，表明原始光譜中一些微弱的變化，在經(jīng)微分后得到更強(qiáng)烈的信號(hào)，導(dǎo)致主成分分析時(shí)提取更多的主成分來(lái)解釋這些變化，即一階微分變換挖掘出了更多土壤信息。二階微分反射率在0 周?chē)?dòng)，且浮動(dòng)幅度類(lèi)似，僅一個(gè)主成分即可概括絕大多數(shù)自變量信息。

Boruta 特征選擇篩選出用于RFR 建模的變量如表2 所示，R 和LGR 篩選出的特征波段除2450 nm 和2486nm 外完全一致，說(shuō)明R 進(jìn)行對(duì)數(shù)變換對(duì)識(shí)別敏感波段并未有明顯改善。進(jìn)行倒數(shù)變換后，敏感波段從19 個(gè)增加至39 個(gè)，敏感波段所處的波段范圍基本一致，均集中在1900～2500nm。經(jīng)一階微分變換后，敏感波段的數(shù)量從21 個(gè)提升至38 個(gè)，二階微分后，篩選出的敏感波段數(shù)量減少到15 個(gè)。微分變換后，敏感波段的位置發(fā)生較大改變，由未變換前短波紅外范圍擴(kuò)散至可見(jiàn)光范圍，短波紅外范圍的敏感波段減少，且分布較分散，說(shuō)明經(jīng)過(guò)微分變換凸顯了可見(jiàn)光范圍的光譜特征，挖掘出更多可見(jiàn)光波段的光譜信息。

表2 Boruta 篩選的敏感波段數(shù)量及位置

2.2 模型對(duì)比分析

結(jié)合PLSR，RFR 和不同光譜變換的波段反射率構(gòu)建的土壤C/N 預(yù)測(cè)模型評(píng)價(jià)指標(biāo)如表3 所示，不同模型和不同變換形式對(duì)模型精度有顯著影響。兩種機(jī)器學(xué)習(xí)算法中，未經(jīng)微分變換（R，LGR 和1/R）的光譜反射率建立的模型有相似精度，PLSR 模型 中R，LGR 和1/R 的R2分 別 為0.698，0.697 和0.694。表明經(jīng)對(duì)數(shù)和倒數(shù)變換的光譜反射率對(duì)預(yù)測(cè)模型精度沒(méi)有顯著提升，甚至模型精度有輕微降低，RFR 模型中也發(fā)現(xiàn)相同規(guī)律。經(jīng)微分變換的光譜反射率建立的模型精度得到提高，PLSR 中，相比于R 建立的預(yù)測(cè)模型，R′，R" 建立的預(yù)測(cè)模型R2分別提升了42.6%和33.0%，RFR 中，相比于R 建立的預(yù)測(cè)模型，R′，R"建立的預(yù)測(cè)模型R2分別提升了42.9%和42.0%，表明一階微分能挖掘潛在的土壤信息，一定程度上消除噪聲，提高預(yù)測(cè)模型精度，提高幅度相似證明PLSR 中的主成分提取和Boruta 篩選效果相似。

表3 土壤C/N 不同模型精度對(duì)比

對(duì)比PLSR 和RFR 模型，PLSR 建立的預(yù)測(cè)模型精度高于RFR 模型，PLSR 通過(guò)提取的主成分解釋波段與C/N 的線性關(guān)系，而RFR 對(duì)非線性關(guān)系有更強(qiáng)的解釋能力，PLSR 算法相比于RFR 算法更適用于建立高光譜數(shù)據(jù)與土壤C/N 之間的關(guān)系。R及其4 種變換后的光譜反射率建立的PLSR 模型比RFR 模型有不同程度地提高，R2提高23.7～27.6%，即每種變換的提高幅度差異較小。預(yù)測(cè)結(jié)果與實(shí)測(cè)值之間散點(diǎn)圖如圖3 所示，表明預(yù)測(cè)值與實(shí)測(cè)值之間有很好的相關(guān)性，預(yù)測(cè)結(jié)果是可靠的，PLSR 模型的預(yù)測(cè)值與實(shí)測(cè)值的擬合直線相比與RFR 模型的擬合直線更接近于1:1 線，表明PLSR模型預(yù)測(cè)效果比RFR 模型更好。

圖3 不同光譜變換和不同建模方法預(yù)測(cè)的土壤C/N與實(shí)測(cè)值比較

3 結(jié)論

高光譜技術(shù)依靠其數(shù)據(jù)獲取快捷，波段信息豐富等優(yōu)勢(shì)被廣泛應(yīng)用于估測(cè)土壤理化性質(zhì)，但野外環(huán)境復(fù)雜，受到土壤水分、土壤雜質(zhì)和大氣條件等影響，背景噪聲大，采集土樣實(shí)驗(yàn)室內(nèi)測(cè)定高光譜數(shù)據(jù)可以獲得較為穩(wěn)定的光譜信息，結(jié)果相對(duì)理想，相比復(fù)雜的化學(xué)測(cè)定土壤C/N 方法，又省時(shí)省力無(wú)污染。本研究通過(guò)對(duì)土樣高光譜曲線的重采樣、平均和平滑預(yù)處理，并對(duì)原始光譜反射率曲線進(jìn)行4 種數(shù)學(xué)變換，使用主成分分析、Boruta 特征篩選、PLSR 和RFR 對(duì)土壤C/N 分析研究，得到的具體結(jié)論如下：

（1）原始光譜R，對(duì)數(shù)變換LGR 和倒數(shù)變換1/R，Boruta 篩選出的敏感波段范圍相似，集中在1900～2500nm，微分變換后篩選出的敏感波段均勻分布在整個(gè)波段區(qū)間，微分變換能有效挖掘光譜信息。

（2）在PLSR 和RFR 模型中，對(duì)R 進(jìn)行對(duì)數(shù)和倒數(shù)變換并不能提高模型精度，對(duì)R 進(jìn)行微分變換，在兩種算法建立的模型精度均得到了顯著提高，R2最高提升42.9%。

（3）兩種算法均能有效預(yù)測(cè)土壤C/N，PLSR 模型精度比RFR 模型精度高23.7～27.6%，PLSR 更適用于小樣本數(shù)據(jù)集建模預(yù)測(cè)，基于PLSR 和R’構(gòu)建 的 最 優(yōu) 模 型，R2為0.995，RMSE 為0.216，MAE為0.165。