李冠穩(wěn)，高小紅，肖能文，肖云飛

(1.青海師范大學 地理科學學院，青海 西寧 810008；2.青海省自然地理與環(huán)境過程重點實驗室，青海 西寧 810008；3.中國環(huán)境科學研究院，北京 100012)

1 引言

土壤有機質(zhì)是土壤肥力和土壤質(zhì)量的重要指標，是農(nóng)業(yè)土壤最重要的參數(shù)之一[1]?？焖?、精準地掌握有機質(zhì)含量的空間變化，是精準農(nóng)業(yè)實施和農(nóng)業(yè)可持續(xù)發(fā)展的重要內(nèi)容[2]。使用傳統(tǒng)化學方法測量土壤有機質(zhì)含量，分析過程周期長、成本高、一次只能檢測一個項目，且對環(huán)境有一定污染，很難大規(guī)模推廣使用[3]?？梢?近紅外(Visible and near infrared，Vis-NIR)光譜分析技術能夠快速、大范圍地重復獲取同一區(qū)域的土壤信息，逐漸成為土壤屬性信息快速與長期監(jiān)測的重要手段之一；且Vis-NIR光譜分辨率高、波段信息豐富，這使得可見-近紅外光譜分析技術在土壤有機質(zhì)的預測分析中表現(xiàn)出巨大的研究潛力[4-7]。但在實際應用中，NIR光譜區(qū)域是由含氫基團的倍頻和合頻吸收峰組成，光譜信息重疊嚴重，篩選出土壤有機質(zhì)的光譜響應波段是簡化模型和提高模型預測能力的關鍵。

特征波長選擇是可見-近紅外光譜研究的一個重要步驟，己經(jīng)引起了越來越多學者的關注[8-9]。李艷坤等[10]基于集群策略和UVE技術，并進一步結合小波變換，得到了更為簡約的模型，提高了PLS模型的預測穩(wěn)定性能。劉珂等[11]通過一致性策略和連續(xù)投影算法結合從全譜波長中選出的一系列波長子集，然后分別基于這些波長子集建立模型，取得了較為滿意的預測效果。林志丹等[12]應用SPA和GA進行波長優(yōu)化，并建立土壤有機質(zhì)Vis-NIR估算模型，結果顯示，對原始光譜進行特征波長優(yōu)選能夠顯著提高模型的精度。競爭性自適應重加權算法(Competitive adaptive reweighted sampling，CARS)是由梁逸曾團隊開發(fā)的一種特征波長變量選擇算法，以偏最小二乘模型中回歸系數(shù)絕對值大小確定最優(yōu)變量子集[13]，而穩(wěn)定競爭性自適應重加權算法(Stbility CARS，sCARS)以變量的穩(wěn)定性為衡量指標，并延續(xù)CARS算法的變量篩選流程，被證明是一種較優(yōu)的特征變量選擇方法[14]。如劉國富等[15]基于sCARS策略挑選NIR光譜區(qū)域特征變量，變量選擇的穩(wěn)定性和準確性都得到了增強，提高了模型精度，預測均方根誤差和相關系數(shù)分別為0.054 3和0.990 8。丁泊洋等[16]采用sCARS算法挑選特征變量建立多元校正模型，預測相關系數(shù)RP為0.978 1，具有較好的預測能力。然而張曉羽等[14]、劉國富等[15]和胡靜等[16]均是利用sCARS方法篩選特征變量，并建立線性的偏最小二乘回歸(Partial Least Squares Regression，PLSR)模型，與非線性的隨機森林(Random forest，RF)建模方法結合的并不多見。與PLSR模型相比，RF模型魯棒性更好，對異常值和噪聲的敏感度更低。

因此，本研究基于青海省湟水流域401個表層土壤的Vis-NIR光譜，應用sCARS方法進行特征波長變量篩選，建立較為簡潔、穩(wěn)定性更好的PLSR和RF模型，并與CARS、IRIV、SPA和GA方法的PLSR和RF模型結果進行比較，探索sCARS算法結合RF模型快速估測土壤有機質(zhì)含量的可行性，為土地質(zhì)量評價和高空間分辨率數(shù)字化土壤制圖提供數(shù)據(jù)支持。

2 材料與方法

2.1 土壤光譜數(shù)據(jù)采集與預處理

我們于2015、2016年10—11月期間，采集青海省湟水流域表層土壤(0～20 cm)共428個土壤樣品，土壤類型主要為栗鈣土、黑鈣土、灰鈣土、山地草甸土、高山草甸土以及灰褐土；并于室內(nèi)自然風干，研磨，過100目篩。有機質(zhì)含量采用重鉻酸鉀-外加熱法測定。使用美國ASD FieldSpec 4光譜儀采集土壤Vis-NIR光譜數(shù)據(jù)。于暗室內(nèi)將過篩的土壤樣品倒入涂黑的盛樣器皿中，減少了外界雜散光的影響，提高光譜質(zhì)量。盛樣器皿直徑為10 cm、高度為1.5 cm。光源為光譜儀配套的75 W鹵素燈，天頂角為30°，距樣品表面45 cm，光線幾乎是平行入射到樣品上，減少了由于土壤顆粒分布不均勻所造成的陰影影響。儀器光纖探頭視場角為25°，垂直向下距樣品表面10 cm處，探頭接收土壤光譜的區(qū)域直徑為5 cm，小于盛樣器皿的直徑，這樣既能避免外界雜散光的影響，又能使光纖探頭接收到的信號均為土壤樣品的反射光譜信息。儀器預熱30 min之后進行白板定標，每個土壤樣品采集4個方向(間隔90°)共20條光譜曲線，為減少測量時土壤樣品光譜各向異性的影響，取20條光譜曲線的算術平均值作為該土壤樣品的實際反射光譜數(shù)據(jù)[17]。土壤樣品最終光譜曲線如圖1(a)所示。剔除原始光譜中噪聲較大的波段(350～400 nm和2 401～2 500 nm)，并聯(lián)合使用多元散射校正(Multiplicative scatter correction，MSC)、中值濾波(Median filter，MF)和一階微分(1st derivative)對原始光譜進行預處理。圖1(b)為經(jīng)MSC-MF-1st Der預處理后的光譜曲線，從圖中可以看出，原始光譜經(jīng)預處理后，不同有機質(zhì)含量光譜曲線等級特征不再明顯，有效地消除了基線漂移及其他背景的干擾，光譜曲線的細節(jié)特征更加突出。

圖1 土壤樣品原始光譜(a)及預處理光譜(b)反射率曲線Fig.1 Raw(a)and pretreatment spectral(b)reflectance curve of soil samples

2.2 方 法

2.2.1 穩(wěn)定競爭性自適應重加權采樣算法(sCARS)

矩陣XN×P為所測樣本光譜數(shù)據(jù)，N為樣本數(shù)量，P為變量數(shù)。sCARS算法具體步驟為：

(1)計算每個波長變量的穩(wěn)定性值cj，cj定義如公式(1)：

(1)

(2)使用強制波長選擇和自適應性重加權采樣方法(ARS)篩選出變量穩(wěn)定性值較大的組成一個變量子集，篩選出的變量數(shù)占全波段的比率由指數(shù)衰減函數(shù)(Exponential decay function，EDF)計算。

(3)重復步驟(1)～(2)形成循環(huán)，最終得到K個變量子集，建立PLSR模型，然后采用十折交叉驗證對這些變量子集進行評估，RMSECV值最小時對應的變量子集作為最后的特征變量子集，K為sCARS算法的循環(huán)次數(shù)。

2.2.2 隨機森林(Random forest，RF)

RF模型是一種分層非參數(shù)方法，融合了隨機特征選擇和Bagging算法兩大機器學習技術，與傳統(tǒng)的分類器算法相比，不但能較好地容忍異常值和噪聲，而且能同時處理連續(xù)型和離散型數(shù)據(jù)[18]。RF模型建模步驟如下：

(1)利用bootstrap重抽樣技術從原始訓練集N中有放回地重復隨機抽取k個樣本生成新的訓練樣本集合；

(3)每棵樹最大限度地生長，使每個節(jié)點的不純度達到最小，不做任何修剪；

(4)生成多棵樹以形成隨機森林，利用隨機森林分類器對新的數(shù)據(jù)進行判別與分類，分類的結果是由樹分類器的投票數(shù)決定的。

2.3 建模樣本集劃分

異常樣本的存在會對模型的性能產(chǎn)生嚴重的干擾，因此在光譜建模分析之前有必要對異常樣本進行識別與剔除[19]。采用主成分分析結合馬氏距離法剔除異常樣本，共剔除異常樣本27個，最終用于分析的土壤樣本共401個。

將異常樣本剔除后的401個土樣按有機質(zhì)含量從高到低排序，按2∶1的比例劃分校正集和驗證集樣本。表1為校正集和驗證集土壤有機質(zhì)含量統(tǒng)計表。校正集中土壤有機質(zhì)含量范圍為4.86～148.74 g·kg-1，平均值為32.47 g·kg-1;驗證集有機質(zhì)含量范圍為8.26～133.56 g·kg-1，平均值為32.16 g·kg-1。濃度梯度法所劃分的校正集樣本組分含量涵蓋了預測集樣本組分含量，避免了過多的“特殊”樣本劃分為建模集，這樣建立的模型能夠更好地預測未知樣本。

表1 校正集和驗證集土壤有機質(zhì)含量統(tǒng)計表Tab.1 Soil organic matter content statistics of calibration sets and validation sets g·kg-1

2.4 模型精度評價

采用Chang等[20]給出的評判等級，當RPD小于1.4時，表明模型不具備估算能力；當RPD大于等于1.4小于2時，表明模型可對樣本進行粗略估算，且可以通過改進模型方法提高模型的預測能力；當RPD大于等于2時，表明模型可以較好地對樣本進行估算。

3 結果與討論

3.1 特征變量選擇

3.1.1 sCARS算法特征變量選擇

sCARS算法以變量穩(wěn)定性作為變量選擇衡量指標，增強了變量選擇的穩(wěn)定性，并延續(xù)CARS算法變量篩選流程。圖2為采用sCARS算法挑選特征變量過程圖，從圖2(a)中可以看出，隨著sCARS算法迭代次數(shù)的增加，所保留的波長數(shù)量逐漸減少，且減少速度由快到慢，表明sCARS算法挑選特征波長變量過程中具有“粗選”和“精選”兩個階段，且“粗選”和“精選”兩個階段存在轉(zhuǎn)折點。圖2(b)為十折交叉驗證RMSECV值變化趨勢圖，可以得知，隨著運行次數(shù)的增加，RMSECV值呈先由大到小再由小到大的變化趨勢。當運行次數(shù)為27次時，RMSECV值最小，表明在1～27次變量篩選運行過程中，剔除了與土壤有機質(zhì)含量相關性較小的波長，對建模結果影響不大；而27次之后RMSECV值開始上升，可能是由于刪除了與土壤有機質(zhì)含量相關的變量導致RMSECV值增大，模型效果變差。結合圖2(c)回歸系數(shù)路徑變化圖可以發(fā)現(xiàn)，當運行次數(shù)為27次時，RMSECV值最小，即選擇的特征波長子集最佳，共選擇51個特征變量，僅占總變量數(shù)的2.55%。圖3為sCARS算法挑選的51個特征變量在一條光譜曲線上的分布情況。

圖2 sCARS算法變量篩選流程Fig.2 Variable selection process by sCARS method

圖3 sCARS方法挑選的特征變量分布圖Fig.3 Distribution map of characteristic variables selected by sCARS method

3.1.2 CARS、IRIV、SPA、GA算法特征變量選擇

CARS算法利用指數(shù)衰減函數(shù)和自適應重加權技術優(yōu)選出偏最小二乘模型中回歸系數(shù)絕對值大的變量點，去除權重值較小的點，再基于十折交叉驗證，選出均方根誤差最小的變量子集，確定為最優(yōu)變量組合。本研究基于CARS算法共選擇59個特征變量，占全部變量數(shù)的2.95%。CARS算法的優(yōu)點是速度快，最終選出的特征變量的化學意義也比較容易解釋，但其選擇的特征變量不穩(wěn)定。

IRIV算法是由中南大學梁逸曾教授課題組提出的一種基于模型集群分析策略的波長選擇算法[21]，將信息變量分為強信息變量、弱信息變量、干擾變量和無信息變量。IRIV由隨機子集生成、子集模型建立、模型參數(shù)分析三個環(huán)節(jié)構成，相對于一般的波長選擇算法，IRIV算法具有在波長選擇時呈現(xiàn)出軟收縮的特點，因此一般能更為穩(wěn)妥地保留有效波長，但其缺點是計算量較大，因此應用受到限制[22]。本研究基于IRIV算法保留的強信息和弱信息變量數(shù)為63個，占全部變量的3.15%。

SPA算法是一種新興的波長選擇算法[23]，其原理為基于連續(xù)投影策略選擇與某一點波長線性相關最小的波長構成一個波長子集，重復上述操作，直至全部波長點選擇完畢；然后基于這些波長子集建立模型，根據(jù)模型精度進而挑選出最優(yōu)的波長子集。本研究采用SPA算法共選擇出5個最優(yōu)特征變量，占全部變量的0.25%，分別為1 361，1 758，1 909，2 049，2 213 nm。SPA算法可以盡可能地消除波長變量間共線性的影響，提高特征變量的選擇能力，但其缺點是在挑選特征變量過程中傾向于選擇共線性較小的變量點而不是有效變量點，因此該算法選擇特征變量也不穩(wěn)定。

GA算法是一種通過模擬自然進化過程搜索最優(yōu)解的方法[24]。借鑒生物的自然選擇和遺傳機理，遺傳算法主要通過編碼、種群初始化、適應度函數(shù)、遺傳操作和終止條件等步驟優(yōu)化選擇。GA算法具有全局最優(yōu)、易實現(xiàn)等特點，成為目前最為常用的一種波長選擇算法。但同時由于隨機選擇初始種群，選擇、交叉和變異都具有很強的隨機性，因此不能保證每個波長選擇結果的一致性，故本研究擬采用多次(10次)運行GA算法，選取特征變量篩選結果中出現(xiàn)頻率較高的波長，最終作為特征波長用于構建模型，按該方法從原始光譜中共選取186個特征波長變量，占全部變量的9.3%。

圖4為CARS、IRIV、SPA、GA算法挑選的特征變量在一條光譜曲線上的分布。從圖3和圖4中可以看出，5種變量篩選方法挑選的特征波長變量主要分布在1 900～2 400 nm的近紅外光譜區(qū)域，其中sCARS、CARS、IRIV、GA法篩選的特征變量在可見-近紅外光譜區(qū)域均有分布，而SPA算法挑選的特征變量較分散地分布于近紅外光譜區(qū)域內(nèi)，可見光區(qū)域均未被選擇。

圖4 CARS(a)、IRIV(b)、SPA(c)和GA(d)算法篩選特征變量分布圖。Fig.4 Distribution map of characteristic variables selected by CARS(a),IRIV(b),SPA(c)and GA(d)method.

3.3 PLSR建模

表2 不同變量篩選方法PLSR建模精度Tab.2 Accuracies of PLSR modeling with different variable selection methods

圖5 sCARS-PLSR模型預測值和實測值散點圖Fig.5 Scatter diagram of predicted and measured values for the sCARS-PLSR model

3.4 RF建模

圖6為sCARS-RF模型校準集和驗證集樣本實測值和預測值的散點圖。從圖中可以看出，sCARS-RF模型校正集和驗證集數(shù)據(jù)點均較為均勻地分布在1∶1直線的兩側(cè)，達到了較高的預測水平，這與上述分析一致。

表3 不同變量篩選方法RF建模精度Tab.3 Accuracies of RF modeling with different variable selection methods

圖6 sCARS-RF模型預測值和實測值散點圖Fig.6 Scatter diagram of predicted and measured values for the sCARS-RF model

4 討論

PLRS模型中，sCARS算法模型精度高于CARS、IRIV、GA、SPA和全波段；RF模型中，基于5種變量選擇算法模型精度與全波段模型精度相差不大，但其構建模型的變量數(shù)卻顯著減少，大大提高了建模效率。對原始光譜進行特征變量篩選，在保證模型精度的同時大大降低了模型的復雜度?；贑ARS、GA和SPA算法挑選的特征變量建模，雖能簡化模型，但變量選擇的穩(wěn)定性較差，挑選的特征變量不總是能反映屬性信息。IRIV算法雖能較穩(wěn)妥地保留有效波長，但其缺點是計算量較大，因此應用受到限制。sCARS算法以變量的穩(wěn)定性作為衡量指標，變量選擇分“粗選”和“精選”兩個階段，既提高了變量選擇效率，又增加了變量選擇的穩(wěn)定性和準確性。但需注意的是，RF模型的精度并沒有像PLSR模型通過應用sCARS算法挑選特征變量而大大增加，且sCARS-PLSR模型精度仍然不如全譜RF模型，這可能是由于RF模型在Vis-NIR光譜數(shù)據(jù)分析中考慮到大量非線性關系，在PLSR模型與變量選擇方法的任何組合中都沒有觀察到這個特征，這一結果也支持了上述的討論，對土壤有機質(zhì)含量的Vis-NIR光譜分析應該采用非線性校準方法以獲得最佳預測效果。sCARS算法挑選的特征變量包含了土壤有機質(zhì)含量最有效的信息，可以代替RF模型的全部原始光譜。

5 結論

以青海省湟水流域401個土壤樣本的有機質(zhì)含量為研究對象，應用sCARS、CARS、IRIV、SPA和GA算法從全波段光譜數(shù)據(jù)中篩選特征變量，分別建立基于特征波段和全波段的PLSR和RF預測模型，取得了較好的預測效果。主要研究結論如下：

(2)RF模型的預測效果優(yōu)于PLSR模型。與采用全波段建模相比，使用特征變量建立PLSR模型，模型精度均有提高；采用特征變量構建RF模型對模型預測精度提高幫助不明顯，但其構建模型的變量數(shù)卻顯著減少，大大提高了建模效率。對全波段進行特征變量篩選，在保證模型精度的同時大大降低了模型的復雜度。

(3)sCARS算法以變量穩(wěn)定性作為變量選擇衡量指標，有效克服了CARS、IRIV、SPA和GA算法的不足，既增強了變量選擇的穩(wěn)定性和準確性，又提高了變量選擇效率，與RF模型結合可實現(xiàn)土壤有機質(zhì)含量快速、無損、精準估測。