近紅外光譜結(jié)合連續(xù)投影算法檢測綜纖維素含量

熊智新 馬璞璠 梁 龍 房桂干

（1.南京林業(yè)大學(xué)林業(yè)資源高效加工利用協(xié)同創(chuàng)新中心，江蘇南京，210037；2.南京林業(yè)大學(xué)輕工與食品學(xué)院，江蘇南京，210037；3.中國林業(yè)科學(xué)研究院林產(chǎn)化學(xué)工業(yè)研究所，江蘇南京，210042）

隨著制漿造紙工業(yè)的發(fā)展，優(yōu)質(zhì)紙漿的需求日益提高，單一品種的制漿材難以維持生產(chǎn)的需求，多種制漿材混合制漿已成為必然［1-3］。制漿材組分中綜纖維素含量直接影響紙漿得率和質(zhì)量［4］，在工業(yè)生產(chǎn)中快速檢測該指標(biāo)值，對及時優(yōu)化工藝參數(shù)、提高紙漿質(zhì)量和制漿造紙工業(yè)的智能化水平有著重要意義。

近紅外光譜（Near Infrared Spectroscopy，NIRS）技術(shù)是近年來發(fā)展較快的分析技術(shù)，以其快速無損的特點結(jié)合化學(xué)計量學(xué)方法在各行業(yè)得到廣泛應(yīng)用［5］。很多學(xué)者對近紅外技術(shù)在制漿材組分檢測上做了大量的研究［2，4，6-10］，但是這些方法均使用全譜圖進(jìn)行建模分析。由于近紅外光譜主要反應(yīng)含氫基團(tuán)振動的倍頻和合頻吸收，因此光譜信息嚴(yán)重重疊，含有大量的共線性信息［11］。直接使用全譜圖分析會使模型引入大量的冗余及無關(guān)信息，對檢測過程造成不利影響。如果在建模前對全譜進(jìn)行波長選擇，篩選出對預(yù)測信息有意義的波長點或波長區(qū)間，不僅可以剔除冗余變量和無關(guān)變量，簡化模型，加快檢測速度的同時提高模型的預(yù)測能力和穩(wěn)健性，而且優(yōu)選出的波長還可以指導(dǎo)儀器設(shè)計，簡化分光系統(tǒng)，制造出高性價比的專用儀器，有利于在工業(yè)生產(chǎn)中推廣近紅外光譜技術(shù)［12］。

本課題以常用的5種制漿材總計82個樣品為研究對象，使用國產(chǎn)便攜式近紅外光譜儀采集各個樣品的光譜數(shù)據(jù)，按國標(biāo)測定其綜纖維素含量。數(shù)據(jù)集預(yù)處理后，使用連續(xù)投影算法（Successive Projections Algorithm，SPA）［13-14］選擇波長組合，建立偏最小二乘（Partial Least Square，PLS）回歸模型用于預(yù)測未知樣品的綜纖維素含量，為制漿材綜纖維素含量快速檢測提供了一種簡單、精確的方法。

1 材料與方法

1.1 樣品制備

選用中國林業(yè)科學(xué)研究院林產(chǎn)化學(xué)工業(yè)研究所提供的5種常用制漿材：包括3種闊葉木制漿材（相思木、桉木、楊木）和2種針葉木制漿材（馬尾松、杉木），共計82個樣品。原木經(jīng)去皮、切削、磨粉、過篩，選取40～60目的木粉，供后續(xù)綜纖維素含量測定和采集光譜使用。

1.2 綜纖維素含量測定

綜纖維素含量按照國家標(biāo)準(zhǔn)GB/T 2677.10—1995測定。每個樣品測定2次，取均值為其含量測定值。

1.3 樣品光譜的采集

光譜儀器采用無錫訊杰光遠(yuǎn)科技有限公司基于微機(jī)電系統(tǒng)（Micro-Electro-Mechanical System，MEMS）技術(shù)生產(chǎn)的數(shù)字微鏡器件（Digital Mirror Device，DMD）型便攜式近紅外光譜儀IAS-5000。為了便于實際生產(chǎn)中推廣應(yīng)用，該型儀器采用固定光柵以抗振，因此攜帶和安裝方便，且價格相對便宜。與實驗室分析儀器相比，該型儀器采集的光譜數(shù)據(jù)信噪比稍低，需要通過以化學(xué)計量學(xué)為基礎(chǔ)的軟件技術(shù)予以修正，以獲得更加準(zhǔn)確穩(wěn)定的模型，滿足制漿材化學(xué)組分過程分析的精度需求。

儀器基本參數(shù)如下：DMD掃描寬度12.87 nm，波長范圍900～1700 nm，分辨率10 nm，平均次數(shù)50次。采集樣品光譜時，將樣品放入測量杯中用200 g砝碼壓平以減少由于樣品分布不均勻?qū)е碌纳⑸溆绊?。每個樣品重復(fù)裝樣采集6次，取平均光譜作為該樣品的原始光譜。

1.4 光譜數(shù)據(jù)處理

為了避免異常樣本對建模過程影響，采用基于蒙特卡羅交叉驗證的異常樣本識別方法［15］對數(shù)據(jù)集中的光譜異常、組分異常及光譜和組分的共同異常樣本進(jìn)行檢驗。該方法用蒙特卡羅隨機(jī)抽樣（Monte Carlo Sampling，MCS）原理選取部分樣本作為校正集建立PLS回歸模型，以剩余樣本進(jìn)行預(yù)測，得到其預(yù)測殘差。以上步驟重復(fù)N次后，統(tǒng)計各個樣本的預(yù)測殘差均值（MEAN）和方差（STD），超出設(shè)定閾值的樣本判為異常。

剔除異常樣本后把剩下樣品按2∶1的比例劃分為校正集和預(yù)測集?？倲?shù)據(jù)集以綜纖維素含量進(jìn)行升序排列，每3個樣品為一個小組，把中間的樣本放入預(yù)測集，其余放入校正集，保證校正集樣本數(shù)據(jù)的代表性。

由于光譜儀采集的光譜信號中不僅包含有用信息，還有儀器噪聲、樣品背景以及雜散光等干擾信號，因此需要對原始光譜進(jìn)行預(yù)處理消除無用信息［11］。分別采用移動均值濾波（Move Average Filter，MAF）、多元散射校正（Multiplicative Scatter Correction，MSC）、一階導(dǎo)數(shù)（First-order Derivative，1stD）3種方法對光譜進(jìn)行預(yù)處理，根據(jù)建模驗證效果選擇合適的預(yù)處理方法。

1.5 SPA選擇特征波長

SPA是一種正向變量選擇方法，通過在向量空間中使用投影的操作來獲取具有最小共線性的變量子集。算法所選的新變量是未選擇變量在已選變量的正交子空間上的最大投影變量。波長選擇中近紅外光譜預(yù)設(shè)3個參數(shù)，即初始波長數(shù)p、波長點個數(shù)j、初始波長點k（0）。SPA的具體步驟如下：

（1）第一次迭代（N=1）開始前，任選一列波長向量xi作為初始波長點xk（0），即k（0）=i，i∈1，…，p。

（6）迭代次數(shù)N=N+1，若N＜j，返回到第二步驟循環(huán)。

最終選取的波長集合為｛k（N），N=0，…，j-1｝。算法對每次選擇的結(jié)果進(jìn)行建模預(yù)測分析，波長集合的最佳初始波長點p和最佳波長點個數(shù)j以建模驗證的最小預(yù)測均方根誤差RMSEP確定［12］。

1.6 模型評價

建模預(yù)測過程以校正均方根誤差RMSEC、校正決定系數(shù)交叉驗證均方根誤差RMSECV、交叉驗證決定系數(shù)和預(yù)測均方根誤差RMSEP、預(yù)測決定系數(shù)對模型的性能進(jìn)行綜合評價［11］。其中均方根誤差值越小，決定系數(shù)值越接近于1，模型的性能越好。

2 結(jié)果與分析

2.1 樣品原始近紅外光譜的獲取

由于光譜儀探測器在1640～1700 nm范圍內(nèi)噪聲較大，因此截取波長范圍900～1640 nm共741個波長點用于分析，得到82個木粉樣品的原始近紅外光譜如圖1所示。不同樣品的綜纖維素含量有較大差異（如表1所示），但是其近紅外譜圖總體趨勢相似，因此必須借助化學(xué)計量學(xué)的方法才能從光譜中提取定性和定量分析信息。

圖1 制漿材樣品原始近紅外光譜圖

2.2 異常樣本剔除及數(shù)據(jù)集的劃分結(jié)果

采用1.4節(jié)的異常樣本識別方法，經(jīng)過2000次MCS建模運算后，計算各樣本預(yù)測殘差的MEAN和STD繪于圖2中，這里設(shè)置閾值為3，即超出整體MEAN和STD平均值3倍的樣本判為異常［15］。從圖2中可以看出，65號樣本數(shù)據(jù)偏離主體且超出閾值限，判斷為異常樣本。

圖2 異常樣本的檢測圖

剔除異常樣本后數(shù)據(jù)集樣本個數(shù)為81個，按綜纖維素含量分布劃分54個樣本作為校正集，27個樣本作為預(yù)測集，各制漿材樣本數(shù)量及綜纖維素的含量在兩個集合中的分布如表1所示。由表1可以看出，不同制漿材在校正集和預(yù)測集中都有一定的比例，且校正集樣本的綜纖維素含量范圍大于預(yù)測集，所以校正集具有較好的代表性。

2.3 光譜的預(yù)處理結(jié)果

分別采用MAF、MSC、1stD 3種方法對數(shù)據(jù)集進(jìn)行預(yù)處理，采用留一法交叉驗證建立PLS回歸模型，以RMSEC和RMSECV以及對應(yīng)的決定系數(shù)對預(yù)處理方法建模性能進(jìn)行評價，各評價參數(shù)如表2所示。

從表2中可以看出，經(jīng)MAF預(yù)處理的光譜建模驗證效果最差，這主要是因為近紅外光譜屬于若干個不同基頻的倍頻和合頻譜帶的組合，沒有銳鋒和基線分離的峰，卻存在大量的重疊峰和肩峰（如圖1所示），平滑處理雖然可以提高信噪比，但無法消除基線漂移對提取光譜有效信息的嚴(yán)重干擾。導(dǎo)數(shù)處理則可以消除基線漂移，克服譜帶重疊，強化譜帶特征，提高了光譜分辨率和靈敏度，進(jìn)而提高模型的精度。但導(dǎo)數(shù)處理同時也放大了噪聲，尤其是本課題研究所用便攜式儀器采用光柵式分光系統(tǒng)，和實驗室常用的傅里葉變換光譜儀相比，波長采樣點相對稀疏，會使導(dǎo)數(shù)處理產(chǎn)生較大的誤差，影響模型預(yù)測結(jié)果［16］。而MSC處理，在不影響光譜本身信噪比的前提下，可以在一定程度上消除由于木粉顆粒分布不均勻及顆粒大小不同產(chǎn)生的散射影響。處理后的光譜如圖3所示。由圖3可以看出，采用MSC預(yù)處理后極大地減少了光譜基線漂移，且模型驗證效果最好。因此本課題后續(xù)研究選擇MSC作為光譜數(shù)據(jù)集的預(yù)處理方法。

表1 綜纖維素含量分布

表2 不同預(yù)處理方法建模結(jié)果

圖3 MSC預(yù)處理后的樣品光譜

2.4 波長的選擇

如前所述，近紅外光譜含有大量的重疊峰和肩峰，因此原始光譜數(shù)據(jù)的741個波長點信息重疊嚴(yán)重，且夾雜隨機(jī)噪聲多，降低了模型預(yù)測精度。利用SPA進(jìn)行波長選擇，提取出最能表征全譜圖信息的波長組合，剔除光譜中的冗余信息，既能減小建模工作量，降低模型復(fù)雜度，又可得到預(yù)測性能更好的模型。首先對MSC預(yù)處理后的校正數(shù)據(jù)集進(jìn)行波長點的迭代選擇，逐個篩選出最大投影對應(yīng)的波長點；再以每次所選波長點建立多元線性回歸模型進(jìn)行預(yù)測，最后根據(jù)預(yù)測集RMSEP的最小值確定最優(yōu)的波長點個數(shù)。不同波長點個數(shù)對應(yīng)的RMSEP值如圖4所示，最終選擇25個波長點，它們在全光譜中的位置如圖5所示。

圖4 不同波長點個數(shù)的RMSEP值

為了比較SPA選擇波長組合的效果，同時使用相關(guān)系數(shù)法、競爭性重加權(quán)自適應(yīng)選擇（Competitive Adaptive Reweighed Sampling，CARS）算法［17］進(jìn)行波長選擇。其中相關(guān)系數(shù)法設(shè)置閾值為0.42，即選擇出波長與綜纖維素含量相關(guān)系數(shù)r＞0.42的波長點184個。CARS算法迭代50次，每次隨機(jī)抽取80%的校正集樣本建模，重復(fù)500次取最優(yōu)結(jié)果，所選波長點28個，在全譜圖中分布位置也標(biāo)識在圖5中。這里需要說明的是，為了便于對比，繪圖時兩種方法所選波長點所在譜圖做了平移，因此圖中吸光度并不表示實際值。

圖5 SPA及CARS選擇的特征波長點對比

木材中的綜纖維素是纖維素和半纖維素的總稱，而C—H和O—H鍵則是這兩類碳水化合物的主要結(jié)構(gòu)。已有研究表明［18-19］，纖維素和半纖維素（以聚木糖為代表）在近紅外區(qū)的吸收峰主要有位于C—H二級倍頻區(qū)（1150～1200 nm）的1185、1212、1262 nm等，位于C—H合頻區(qū)（1300～1500 nm）的1352和1364 nm，位于C—H一級倍頻區(qū)（1600～1800 nm）附近的1588 nm，以及位于O—H一級倍頻區(qū)（1430 nm附近）的1434和1484 nm，其中C—H合頻和O—H一級倍頻譜峰信息部分重疊。為便于對比，將3種方法選擇的波長點列舉于表3。從表3可以看出，由于儀器和實驗條件不同、制漿材綜纖維素本身構(gòu)成的復(fù)雜性以及近紅外光譜固有的寬頻帶特征，本課題中3種方法選擇的波長位置和文獻(xiàn)報道存在一定偏移，但大部分挑選自以上綜纖維素特征吸收區(qū)域。其中相關(guān)系數(shù)法挑選的波長呈區(qū)間連續(xù)分布，CARS和SPA方法均呈離散分布，且SPA算法所選波長點更為均勻地分散在各特征吸收區(qū)間。文獻(xiàn)［19］也表明，1400～1600 nm范圍是不同聚合體中的O—H一級倍頻及C—H合頻吸收區(qū)域，使得這部分蘊含重疊信息較多，其中1570 nm處存在分子間鍵合羥基的弱吸收峰［11］，而這種鍵合在蔗糖、淀粉和纖維素等天然產(chǎn)物中大量存在，這也是圖5中兩種方法在1460～1600 nm區(qū)域均入選較多波長點的原因，以利于模型充分表達(dá)綜纖維素的特征吸收及校正背景吸收。

由于各官能基團(tuán)在近紅外光譜帶中重疊嚴(yán)重，以及實驗使用的便攜式儀器波長分辨率相對較低，因此很難精準(zhǔn)分析出所選波長點對應(yīng)的具體基團(tuán)信息，所以在近紅外分析技術(shù)中，一般通過建模預(yù)測精度來判斷波長選擇結(jié)果的有效性。幾種波長選擇方法結(jié)合PLS建立的最優(yōu)回歸模型性能如表4所示。

表3 不同波長選擇方法所選波長結(jié)果

表4 不同波長選擇方法PLS建模結(jié)果

比較表4各模型性能指標(biāo)可以看出，（r＞0.42）-PLS模型的RMSEP與全光譜-PLS模型相比顯著降低，說明去掉相關(guān)系數(shù)小的波長點后，建模預(yù)測效果有了一定的改善。從統(tǒng)計學(xué)理論上講，選擇相關(guān)系數(shù)大的波長，其包含的綜纖維素含量信息越多，建立的模型預(yù)測精度也相應(yīng)越高。不過由于相關(guān)系數(shù)法是基于線性統(tǒng)計方法建立的，而本課題中采用的木粉漫反射吸光度與綜纖維素以及其他組分含量并不服從線性的朗伯-比爾定律，影響了該方法波長選取的可靠性［16］。同時，所選的波長點較多（184個），且連續(xù)分布，譜帶重疊嚴(yán)重也會使模型中引入過多的冗余信息甚至噪聲信號，降低了模型的穩(wěn)健性和精度。

CARS算法優(yōu)選的波長組合雖然提高了模型性能，但是由于算法本身存在隨機(jī)因素，結(jié)果不穩(wěn)定且重現(xiàn)性低，因此，需要進(jìn)行大量的重復(fù)運算才能選擇出較好的結(jié)果。通過對該方法選出的28個波長點分布情況進(jìn)一步分析表明（如表3所示），這些波長雖然也散布在全光譜不同位置，但其選出的波長點大多分布在相鄰的位置，而SPA方法選出的25個特征波長點相對間隔較大，因此入選波長也存在較多冗余信息，影響了模型精度和穩(wěn)健性。而且CARS方法引入的波長1454 nm與木材中水分的強特征吸收峰1450 nm相近［11］，相關(guān)系數(shù)法波段區(qū)間也表明該波長和綜纖維素含量相關(guān)性不明顯。這將使基于CARS方法建立的綜纖維素分析模型（CARS-PLS）易受制漿材水分變化的干擾。

SPA-PLS模型的預(yù)測性能最好，這是因為SPA算法基于數(shù)學(xué)正交最大投影原理選擇特征波長，以消除共線性影響，剔除了含有大量無關(guān)信息的波長，選出的25個特征波長點為最小冗余組合，且能夠表征全譜圖的信息，極大地減少了建模所需的波長變量個數(shù)，同時避開了水分特征吸收區(qū)域，有利于模型在實際測量中克服制漿材水分含量波動的干擾，有利于提高模型的穩(wěn)健性和精度。因此在3種方法中，SPA選擇特征波長的建模效果最好。

3 結(jié) 論

本課題結(jié)合連續(xù)投影算法（SPA）和近紅外光譜技術(shù)對5種常用制漿材的綜纖維素含量進(jìn)行了快速檢測，通過建模實驗進(jìn)行研究及分析。

3.1 針對便攜式光譜儀采集的木粉樣品漫反射光譜數(shù)據(jù)，采用多元散射校正（MSC）進(jìn)行預(yù)處理，能較好地消除散射和平移誤差，提高分析信息提取質(zhì)量。

3.2 SPA-PLS模型的預(yù)測精度優(yōu)于對比的其他模型，預(yù)測均方根誤差RMSEP和決定系數(shù)R2p分別為0.8306和0.9801，滿足工業(yè)應(yīng)用精度需求，這為制漿材綜纖維素含量的快速檢測提供了一種更優(yōu)的方法。

3.3 SPA所選出的25個波長組合相較于CARS方法選擇出的28個特征波長組合在全譜圖上分布更分散，最大程度上消除冗余信息，其選出的波長組合也最能表征全譜圖中的綜纖維素含量信息。這對指導(dǎo)開發(fā)專用型的近紅外綜纖維素快速分析儀并在制漿造紙行業(yè)推廣具有重要意義。