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

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

1 引言

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

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

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

2 材料與方法

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

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

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

2.2 方 法

2.2.1 穩(wěn)定競(jìng)爭(zhēng)性自適應(yīng)重加權(quán)采樣算法(sCARS)

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

(1)計(jì)算每個(gè)波長(zhǎng)變量的穩(wěn)定性值cj，cj定義如公式(1)：

(1)

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

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

2.2.2 隨機(jī)森林(Random forest，RF)

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

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

(3)每棵樹最大限度地生長(zhǎng)，使每個(gè)節(jié)點(diǎn)的不純度達(dá)到最小，不做任何修剪；

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

2.3 建模樣本集劃分

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

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

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

2.4 模型精度評(píng)價(jià)

采用Chang等[20]給出的評(píng)判等級(jí)，當(dāng)RPD小于1.4時(shí)，表明模型不具備估算能力；當(dāng)RPD大于等于1.4小于2時(shí)，表明模型可對(duì)樣本進(jìn)行粗略估算，且可以通過改進(jìn)模型方法提高模型的預(yù)測(cè)能力；當(dāng)RPD大于等于2時(shí)，表明模型可以較好地對(duì)樣本進(jìn)行估算。

3 結(jié)果與討論

3.1 特征變量選擇

3.1.1 sCARS算法特征變量選擇

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

圖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īng)重加權(quán)技術(shù)優(yōu)選出偏最小二乘模型中回歸系數(shù)絕對(duì)值大的變量點(diǎn)，去除權(quán)重值較小的點(diǎn)，再基于十折交叉驗(yàn)證，選出均方根誤差最小的變量子集，確定為最優(yōu)變量組合。本研究基于CARS算法共選擇59個(gè)特征變量，占全部變量數(shù)的2.95%。CARS算法的優(yōu)點(diǎn)是速度快，最終選出的特征變量的化學(xué)意義也比較容易解釋，但其選擇的特征變量不穩(wěn)定。

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

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

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

圖4為CARS、IRIV、SPA、GA算法挑選的特征變量在一條光譜曲線上的分布。從圖3和圖4中可以看出，5種變量篩選方法挑選的特征波長(zhǎng)變量主要分布在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模型預(yù)測(cè)值和實(shí)測(cè)值散點(diǎn)圖Fig.5 Scatter diagram of predicted and measured values for the sCARS-PLSR model

3.4 RF建模

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

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

圖6 sCARS-RF模型預(yù)測(cè)值和實(shí)測(cè)值散點(diǎn)圖Fig.6 Scatter diagram of predicted and measured values for the sCARS-RF model

4 討論

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

5 結(jié)論

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

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

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