儀器噪聲在近紅外光譜偏最小二乘模型內(nèi)的傳播效應(yīng)

張若秋，杜一平

(華東理工大學(xué) 上海市功能性材料化學(xué)重點實驗室 化學(xué)與分子工程學(xué)院，上海 200237)

近紅外(Near-infrared,NIR)光譜分析是近年來最引人關(guān)注的分析技術(shù)之一，在農(nóng)業(yè)、制藥、化工等眾多領(lǐng)域得到廣泛應(yīng)用[1-2]。目前NIR技術(shù)的發(fā)展重點是應(yīng)用，簡便、快速、無損、準(zhǔn)確等優(yōu)點正是其受重視并被廣泛應(yīng)用于國民經(jīng)濟發(fā)展各個領(lǐng)域的本質(zhì)原因。但在使用NIR技術(shù)的基層單位，技術(shù)隊伍的研究能力往往不如研發(fā)單位，更不及高校和專門的科研機構(gòu)。而NIR技術(shù)使用和推廣的難點是近紅外光譜分析模型的建立需要一些專門的技術(shù)，尤其是依賴于化學(xué)計量學(xué)。因此，模型的建立和維護是近紅外光譜技術(shù)發(fā)展的一個突出難點和瓶頸。

偏最小二乘(Partial least squares,PLS)是近紅外光譜多元校正領(lǐng)域使用最為廣泛的算法之一[3-4]。在光譜多元校正的實際應(yīng)用中，有很多因素會影響該算法的建模過程和模型的預(yù)測能力，其中包括光譜儀器噪聲。然而，目前絕大多數(shù)研究都針對如何使用濾波器或平滑算法減少光譜噪聲[5-6]，鮮有針對儀器噪聲如何影響偏最小二乘建模過程及其預(yù)測能力的研究。本文重點闡述和討論儀器噪聲怎樣通過第一個隱變量的計算被引入模型中，并隨著后續(xù)隱變量不斷在模型中傳遞和累積從而對整個建模過程產(chǎn)生影響，可為今后進(jìn)一步研究在建模中減少或抑制噪聲在模型內(nèi)的傳播和放大提供理論依據(jù)。

1 儀器噪聲對偏最小二乘建模的影響

樣本的測量光譜用X表示，組分含量用y表示，多元校正關(guān)系模型為y=Xβ+e?？紤]測量光譜存在誤差，將X表示為：

X=X0+E

(1)

其中，X0為不含測量誤差的真光譜，E表示儀器噪聲，假設(shè)矩陣中的每個元素均滿足均值為0且方差相同的高斯分布。類似地，含量真值用y0表示(本文只考慮光譜存在誤差，含量無誤差，即y0=y)。

1.1 偏最小二乘算法實現(xiàn)

PLS的主流實現(xiàn)方式有非線性迭代偏最小二乘(NIPALS)[7]和SIMPLS兩種，為方便理論推導(dǎo)本文使用SIMPLS，下文給出SIMPLS的偽代碼[8]。SIMPLS算法流程如下：

Fora=1,…,A;

a=1:對S0進(jìn)行奇異值分解;

a>1:對(I-P(PTP)-1PT)S0進(jìn)行奇異值分解;

計算得到權(quán)重向量r=第一個左奇異向量;

計算得到得分向量t=X0r;

將r、t和p分別存入R、T和P;

End

無噪光譜矩陣X0含有p個變量和n個樣本，濃度向量y0含有n個樣本。R、T和P分別代表權(quán)重矩陣、得分矩陣和載荷矩陣，而βpls則為偏最小二乘模型的回歸系數(shù)向量。為方便，本文使用LSV1(S0)表示S0的第一個左奇異向量。

1.2 儀器噪聲傳遞

當(dāng)PLS模型的隱變量數(shù)量(nLVs)為k個時，在不同隱變量下所形成的空間和由Krylov序列形成的空間滿足如下等式[5]：

{Rk}={S0,XTXS0,(XTX)2S0,…,(XTX)k-1S0}

(2)

{Pk}={XTXS0,(XTX)2S0,…,(XTX)kS0}

(3)

以上兩個等式表明協(xié)方差矩陣S0(當(dāng)預(yù)測指標(biāo)只有1個時是協(xié)方差向量)和XTX完全確定了βpls。前者確定模型的初始方向，而后者則影響每個隱變量的實際計算。當(dāng)光譜數(shù)據(jù)中包含了儀器噪聲時，S0和pk的計算可改寫為如下形式：

(4)

(5)

因此，通過上述公式發(fā)現(xiàn)偏最小二乘模型將從第一個隱變量開始發(fā)生偏離真模型的現(xiàn)象，噪聲同時被引入到模型中，并在隨后隱變量的計算過程中被不斷傳遞和累積：

(6)

Rk=[LSV1(S1),LSV1(S2),…,LSV1(Sk)]

(7)

βpls=RkRkTS0

(8)

上述噪聲傳遞過程同時也會對偏最小二乘模型的預(yù)測能力產(chǎn)生影響。預(yù)測殘差可表示為：

(9)

(10)

(11)

(12)

圖1 噪聲傳遞和偏最小二乘預(yù)測誤差關(guān)系示意圖Fig.1 The scheme of the relationship between noise propagation and prediction error

噪聲傳遞和偏最小二乘的預(yù)測誤差關(guān)系如圖1所示，不同噪聲水平下的最優(yōu)隱變量數(shù)也將發(fā)生改變，噪聲水平越高，最優(yōu)隱變量的數(shù)值越小。這說明噪聲傳遞不僅影響偏最小二乘模型的預(yù)測能力，還會給模型的選擇帶來影響。

2 數(shù)據(jù)集與計算方法

2.1 數(shù)據(jù)集模擬

本文用模擬數(shù)據(jù)對上述討論進(jìn)行驗證。模擬數(shù)據(jù)集包含了2 000個建模樣本和2 000個獨立預(yù)測樣本，其變量數(shù)為1 000。數(shù)據(jù)集中包含1個待測組分和40個干擾組分，共41個組分。每個模擬組分的純光譜信號均滿足一個隨機均值和方差的高斯分布。建模集中待測組分和干擾組分的濃度值服從0到1的均勻分布，而在預(yù)測集中服從0.05到0.95的均勻分布。模擬光譜由濃度矩陣和純光譜信號矩陣相乘得到。

對于每個固定的光譜信號xi(i=1,2,…,n)均產(chǎn)生一個服從標(biāo)準(zhǔn)高斯分布的噪聲向量ei，噪聲水平的大小定義為ei的最大值和xi的比值。更高的噪聲水平意味著更低的信噪比(S/N)，反之亦然。3個不同水平的噪聲(0.3%、0.75%和1%)用于給無噪光譜數(shù)據(jù)添加人為噪聲。

2.2 指標(biāo)與計算軟件評價

評價所建立的NIR模型通常考慮模型誤差和相關(guān)系數(shù)兩個指標(biāo)，本文只考慮模型誤差指標(biāo)。建模時一般將數(shù)據(jù)集劃分為校正集(C)、交互檢驗集(CV)和預(yù)測集(P)，因此模型誤差也包括校正誤差、交互檢驗誤差和預(yù)測誤差。模型誤差用均方根誤差(Root mean squared error,RMSE)表示：

(13)

已有不少研究表明蒙特-卡洛采樣(Monte-carlo sampling,MCS)能夠降低模型過擬合的風(fēng)險[9-10]，因此本文采用基于MCS的方法，即采樣誤差分布分析的交互檢驗(Cross-validation based on sampling error profile analysis,SEPA-CV)[11]確定偏最小二乘的最優(yōu)隱變量大小。MCS的采樣數(shù)和建模樣本比例分別設(shè)定為2 000和0.8。所有計算均使用MATLAB完成(Version 2010a,The MathWorks,USA)。

為準(zhǔn)確量化偏最小二乘模型的噪聲水平，采用Durbin-Watson(DW)[12-16]統(tǒng)計量評價噪聲對模型的影響程度。DW統(tǒng)計量由下式計算得到：

(14)

xi和xi-1是p維向量中的兩個連續(xù)元素。如果整個向量中連續(xù)元素間的相關(guān)性很弱，例如該向量中每個元素均為隨機變量，則DW統(tǒng)計量的大小將收斂于2。一個擁有至少100個元素的隨機向量的DW統(tǒng)計量的95%置信區(qū)間在1.7和2.3之間。

3 結(jié)果與討論

3.1 噪聲傳遞對偏最小二乘模型回歸向量的影響

首先使用DW統(tǒng)計量對3種不同噪聲水平的光譜所建立的偏最小二乘模型的權(quán)重向量r、載荷向量p和回歸系數(shù)向量βpls進(jìn)行噪聲程度的估計，其結(jié)果見圖2。

由圖2可知，偏最小二乘模型中這3種向量的DW統(tǒng)計量均隨著隱變量數(shù)的增加而增加，這反映了噪聲在偏最小二乘模型中的傳遞和累積現(xiàn)象，與“1.2”理論推導(dǎo)的結(jié)果相符合。對于無噪光譜所建立的PLS模型，在1～41個隱變量范圍內(nèi)，r、p和βpls的DW統(tǒng)計量遠(yuǎn)小于有噪光譜建模得到的各中間向量的DW值。r和p的DW統(tǒng)計量在第42個隱變量處突然增至2左右，說明此時這兩個向量已無法包含任何有用的光譜信息，因為模擬數(shù)據(jù)集的真實組分?jǐn)?shù)是41。然而βpls的DW在第42個隱變量處卻增加不明顯，這是因為βpls是由全部權(quán)重向量r的線性組合，其DW的大小應(yīng)該由這些權(quán)重向量的噪聲程度共同決定。對于包含了儀器噪聲的光譜所建立的PLS模型來說，在前25個隱變量下r、p和βpls的DW與無噪光譜PLS模型相比其對應(yīng)的向量DW相差不大，但在第25個隱變量之后它們的差異顯著增加，這是因為噪聲的傳遞和累積的顯著程度會隨著隱變量數(shù)的增加而增加。

圖2 模擬數(shù)據(jù)集在隱變量1～42下建立的偏最小二乘模型的r(A)、p(B)和βpls(C)的DW統(tǒng)計量Fig.2 DW values of r(A),p(B) and βpls(C) from simulated dataset with different noise levels at nLVs from 1 to 42

3.2 噪聲傳遞對模型預(yù)測能力的影響

圖3 模擬數(shù)據(jù)集中不同噪聲水平光譜所建立 模型的預(yù)測誤差示意圖Fig.3 The prediction error of PLS models built by simulated dataset with different noise levels

采用上述模擬數(shù)據(jù)驗證噪聲傳遞對PLS模型預(yù)測能力的影響，結(jié)果見圖3和表1。Mp為無噪光譜所建立的偏最小二乘模型的RSS，隨著隱變量的增加，模型復(fù)雜度不斷增加，Mp則不斷降低，直至隱變量數(shù)量和模型的真實組分?jǐn)?shù)(41)達(dá)到一致，此時模型的復(fù)雜度最優(yōu)。使用無噪光譜建立的PLS模型的預(yù)測誤差并不包含Np項，因為Np是由儀器噪聲產(chǎn)生。由圖3可知，Np的大小受兩方面影響，一是其隨著PLS隱變量數(shù)量的增加而增加，二是其隨著噪聲水平的增加而增加，這與“1.2”的理論分析一致。由圖4和表1可知，Np的存在和增大會使模型的預(yù)測能力下降。

表1 模擬數(shù)據(jù)集的交互檢驗和RMSEP結(jié)果Table 1 Results of cross-validation and RMSEP to simulated dataset

圖4和表1描述了模擬數(shù)據(jù)集交互檢驗和RMSEP的結(jié)果。對于無噪光譜建立的模型，最優(yōu)隱變量數(shù)量和期望的情況一致(均為41)，且在該隱變量下的RMSECV和RMSEP均很接近于0(1.190 2×10-6和2.774 5×10-6)，此時該模型是這套模擬數(shù)據(jù)集的真實模型。隨著光譜的噪聲水平增加，模型的RMSECV和RMSEP相應(yīng)增加，但更重要的是模型的最優(yōu)隱變量也不斷發(fā)生變化(例如：對于噪聲程度為0.3%和1%的光譜來說，最優(yōu)隱變量分別為33和31)。此外，不同噪聲水平下所建立模型的RMSECV和RMSEP的差異隨著隱變量的增加而增加，這證明了“1.2”的理論推導(dǎo)，即噪聲會在PLS模型中不斷傳遞和累積。由圖4還能發(fā)現(xiàn)，無噪光譜建立的偏最小二乘模型的預(yù)測誤差在隱變量超過41后會突然顯著增加，這意味著模型的過擬合。

圖4 模擬數(shù)據(jù)集包含不同水平噪聲的光譜所建立模型的RMSECV(A)和RMSEP(B)Fig.4 RMSECV(A) and RMSEP(B) values at nLVs from 1 to 42 for simulated dataset

圖5 第一個隱變量下PLS和PoLiSh模型的回歸系數(shù) 及二者差值Fig.5 Regression coefficients of PLS and PoLiSh and their differences

3.3 去噪算法對噪聲傳遞的驗證

采用PoLiSh算法對人為添加1%噪聲的模擬數(shù)據(jù)進(jìn)行噪聲傳遞的驗證。PoLiSh算法的原理為，在使用NIPALS計算每一個隱變量下權(quán)重向量時，對其使用Savizky-Golay平滑[17]以消除噪聲的積累。本文中Savizky-Golay平滑窗口大小為21，擬合階數(shù)為2，使用PoLiSh所建立模型的最優(yōu)隱變量數(shù)為35，相應(yīng)RMSEP值為0.076 0。與表1中引入了1%噪聲的光譜所建立的PLS模型相比，RMSEP降低且最優(yōu)隱變量的數(shù)目增加，說明PoLiSh一定程度上減少了噪聲的影響。圖5反映了第一個隱變量的部分PLS模型和PoLiSh模型，可以明顯看出使用包含儀器噪聲的光譜建立的PLS模型從第一個隱變量開始包含噪聲，而PoLiSh算法建立的模型可在一定程度上減弱噪聲的影響。

4 結(jié) 論

本文闡述了儀器噪聲如何從第一個隱變量的計算開始被引入偏最小二乘模型中并通過隨后的隱變量計算被不斷傳遞和累積，通過對偏最小二乘的計算過程相關(guān)公式的推導(dǎo)和對一套模擬光譜數(shù)據(jù)集的詳細(xì)研究，總結(jié)并論證了與無噪的光譜相比，噪聲的引入將使得偏最小二乘模型的預(yù)測能力變差且模型最優(yōu)隱變量數(shù)減少，且光譜中包含的噪聲水平越高，則其對模型預(yù)測能力和最優(yōu)隱變量的改變越大的結(jié)論。