潘高威，史晉芳，邱 榮，王慧麗，萬(wàn) 情，張志威，王 凱

（1.西南科技大學(xué) 極端條件物質(zhì)特性聯(lián)合實(shí)驗(yàn)室，四川 綿陽(yáng) 621010；2.西南科技大學(xué) 制造科學(xué)與工程學(xué)院，四川 綿陽(yáng)621010；3.陸軍勤務(wù)學(xué)院 教研保障中心，重慶 401331）

引言

激光誘導(dǎo)擊穿光譜（LIBS）技術(shù)是上世紀(jì)六十年代發(fā)展起來(lái)的一種原子發(fā)射光譜檢測(cè)技術(shù)。LIBS技術(shù)同時(shí)具有無(wú)需專(zhuān)門(mén)制樣、對(duì)樣品破壞小、靈敏度高、檢測(cè)限低等特點(diǎn)，在工業(yè)現(xiàn)場(chǎng)監(jiān)測(cè)、環(huán)境污染檢測(cè)、深空深海探測(cè)、文物鑒定等領(lǐng)域的研究熱度不斷上升[1]。中階梯光柵光譜儀相較于傳統(tǒng)光譜儀，具有分辨率高、衍射效率高的特點(diǎn)[2]，是LIBS檢測(cè)技術(shù)的重要組成部分。但是，在利用中階梯光柵光譜儀作為激光誘導(dǎo)擊穿光譜分析的數(shù)據(jù)采集儀器時(shí)，為了增強(qiáng)弱峰信號(hào)，原始光譜數(shù)據(jù)中會(huì)存在周期性拱形基線(xiàn)，這種具有明顯區(qū)域分布的拱形基線(xiàn)主要來(lái)源于激光等離子體的連續(xù)白光背景輻射和中階梯光譜儀的周期性衍射效率差異，這不僅增加了光譜數(shù)據(jù)預(yù)處理的難度，也減少了定量分析的精度和準(zhǔn)確性。因此為了減少這些特殊的基線(xiàn)對(duì)光譜數(shù)據(jù)的影響，基線(xiàn)校正成了LIBS數(shù)據(jù)處理中的一個(gè)重要環(huán)節(jié)。

常用的基線(xiàn)校正方法主要有多項(xiàng)式擬合[3-4]、小波變換[5-7]、懲罰最小二乘法[8-11]等。在多項(xiàng)式擬合中，利用最小二乘法對(duì)人為主觀(guān)選擇的擬合點(diǎn)進(jìn)行多項(xiàng)式擬合，易導(dǎo)致過(guò)擬合或者欠擬合。小波變換通過(guò)將光譜信號(hào)分解成高頻和低頻2種信號(hào)。在低頻信號(hào)中刪去基線(xiàn)，針對(duì)不同種類(lèi)的基線(xiàn)時(shí)，需要重新對(duì)小波函數(shù)和分解層數(shù)進(jìn)行選擇，該方法也會(huì)產(chǎn)生無(wú)意義的負(fù)值，在處理復(fù)雜基線(xiàn)時(shí)效果并不理想。2003年，Eilser提出了一種基于whittaker平滑器原理的非對(duì)稱(chēng)權(quán)重懲罰最小二乘法基線(xiàn)擬合算法（AsLS）[12]，并將其成功運(yùn)用于飛行時(shí)間質(zhì)譜分析的基線(xiàn)校正。趙恒等人將“局部對(duì)稱(chēng)加權(quán)”與“非對(duì)稱(chēng)加權(quán)”結(jié)合，提出了一種基于局部對(duì)稱(chēng)權(quán)重加權(quán)懲罰最小二乘法（LSRPLS）[13]，解決了非對(duì)稱(chēng)懲罰最小二乘法在噪聲情況下擬合拉曼光譜基線(xiàn)時(shí)產(chǎn)生的局部抬升問(wèn)題。但是由于LISB數(shù)據(jù)量大且光譜數(shù)據(jù)復(fù)雜，在將其他光譜基線(xiàn)校正方法運(yùn)用在LIBS數(shù)據(jù)上時(shí)，會(huì)同時(shí)存在局部的欠擬合和過(guò)擬合現(xiàn)象。因此本文針對(duì)該周期性拱形基線(xiàn)的校正問(wèn)題，提出了一種改進(jìn)的懲罰最小二乘方法，在保留光譜數(shù)據(jù)有效信息的基礎(chǔ)上，進(jìn)行了LIBS光譜數(shù)據(jù)的基線(xiàn)校正。

1 改進(jìn)算法原理

1.1 算法原理

LIBS光譜數(shù)據(jù)的基線(xiàn)主要由2部分組成，一部分是由等離子在激發(fā)后所帶來(lái)的連續(xù)背景導(dǎo)致的光譜數(shù)據(jù)強(qiáng)度的整體抬升；另一部分就是中階梯光柵光譜儀特有的周期性拱形突起所帶來(lái)的局部變化。因此，要求基線(xiàn)擬合算法能在擬合連續(xù)背景的同時(shí)，不會(huì)忽視掉這一系列拱形突起。

LIBS光譜數(shù)據(jù)的基線(xiàn)主要來(lái)源于激光等離子體的連續(xù)白光背景輻射，由于中階梯光譜儀的多級(jí)衍射和波長(zhǎng)間的衍射效率差異，導(dǎo)致連續(xù)白光背景輻射基線(xiàn)被多級(jí)衍射調(diào)制從而形成周期性拱形基線(xiàn)。這種周期性拱形基線(xiàn)的存在，會(huì)增加光譜數(shù)據(jù)預(yù)處理的難度，影響LIBS定量分析的精度。因此，進(jìn)行基線(xiàn)擬合時(shí)必須考慮該因素。

假設(shè)有一個(gè)數(shù)據(jù)長(zhǎng)度為 N的經(jīng)min-max標(biāo)準(zhǔn)化后的光譜數(shù)據(jù)y，設(shè)其擬合基線(xiàn)為z，利用 z的二階差分來(lái)定義擬合基線(xiàn)的粗糙度：

式中：R表示粗糙度； Δ2zi表示二階差分，Δ2zi=zi-2zi-1+zi-2，i∈[3,N] 。利用光譜強(qiáng)度 y 與擬合基線(xiàn) z的差值來(lái)定義擬合線(xiàn)對(duì)于原始數(shù)據(jù)的保真度S （即擬合線(xiàn)對(duì)于原始數(shù)據(jù)的偏離程度）：

式中：ωi為權(quán)重系數(shù)。ωi會(huì)隨著迭代的進(jìn)行而發(fā)生改變，迭代的本質(zhì)是將譜峰視為突變信號(hào)而基線(xiàn)視作連續(xù)信號(hào)，突變信號(hào)會(huì)在一次次迭代中逐漸減少，其對(duì)應(yīng)的權(quán)重系數(shù)也會(huì)減小，即權(quán)重的“非對(duì)稱(chēng)性”。ωi的取值方式為

式中：p 為一個(gè)極小值，通常使 0.001≤p≤0.1。當(dāng)yi高 于基線(xiàn) zi時(shí)，則說(shuō)明這部分?jǐn)?shù)據(jù)為譜峰，將其權(quán)值設(shè)置為p；反之則認(rèn)為這部分?jǐn)?shù)據(jù)屬于基線(xiàn)部分，權(quán)值設(shè)置為1 -p，保留這部分?jǐn)?shù)據(jù)。但是由于LIBS光譜數(shù)據(jù)量大，使用原始方法對(duì) ωi進(jìn)行賦值會(huì)導(dǎo)致在無(wú)譜峰位置處產(chǎn)生欠擬合，因此引入新的權(quán)值系數(shù)賦值方式：

然而擬合基線(xiàn) z越光滑，則對(duì)理想基線(xiàn)的保真度就越低。為了同時(shí)滿(mǎn)足這2個(gè)方面的要求，令R作為懲罰項(xiàng)，引入代價(jià)函數(shù)Q =S+λR。λ為懲罰項(xiàng)系數(shù)，取值在1 02～109之間。代價(jià)函數(shù)的矩陣形式為

式中：W 是主對(duì)角線(xiàn)為ω 的稀疏矩陣； D為二階差分矩陣 Dz=Δ2zi。因此求基線(xiàn)問(wèn)題就轉(zhuǎn)變成了尋找一個(gè)能夠使Q達(dá)到最小值的z。為了使代價(jià)函數(shù)Q 最小化，令其對(duì) z求偏導(dǎo)后等于零，得到：

此外，在用非對(duì)稱(chēng)懲罰最小二乘法擬合LIBS光譜數(shù)據(jù)時(shí)，基線(xiàn)會(huì)在譜峰處有所抬升，且譜峰峰值越大，此處的基線(xiàn)受影響越大。因此引入懲罰系數(shù)數(shù)組α，在峰值處施加一個(gè)較大的懲罰系數(shù)，對(duì)基線(xiàn)則施加一個(gè)較小的懲罰系數(shù)，使擬合基線(xiàn)在譜峰位置處更加平滑，從而減小譜峰峰值對(duì)基線(xiàn)的影響。懲罰系數(shù)數(shù)組 α的取值方式如（7）式：

式中：dyi=yi+1-yi。α的數(shù)值根據(jù)標(biāo)準(zhǔn)化后光譜數(shù)據(jù)的一階導(dǎo)數(shù)確定，在譜峰數(shù)據(jù)變化率大的位置施加更大的懲罰系數(shù)。

在引入了α和新的ω的取值方式后，（6）式變?yōu)椋?）式：

1.2 算法實(shí)現(xiàn)

本文算法流程如圖1所示，步驟如下：

圖1 改進(jìn)后的算法流程圖Fig.1 Flow chart of improved algorithm

1）輸入預(yù)先經(jīng)過(guò)min-max標(biāo)準(zhǔn)化后光譜強(qiáng)度部分?jǐn)?shù)據(jù)y；

2）初始化權(quán)重?cái)?shù)組 ω 和懲罰數(shù)組α；

3）計(jì)算標(biāo)準(zhǔn)差 σ ，根據(jù)（8）式，計(jì)算擬合基線(xiàn)z；

4）根據(jù)光譜數(shù)據(jù) y與擬合基線(xiàn) z的差值，對(duì)權(quán)重?cái)?shù)組 ω和懲罰系數(shù)數(shù)組 α中的元素進(jìn)行重新賦值；重復(fù)步驟3和步驟4，直至 ω 中的元素不再發(fā)生明顯變化；

5）迭代終止，將光譜數(shù)據(jù) y與擬合基線(xiàn)z相減得到去基線(xiàn)后的光譜數(shù)據(jù)。

2 仿真與實(shí)驗(yàn)

2.1 仿真

在實(shí)際的光譜檢測(cè)中無(wú)法準(zhǔn)確地得知其理論基線(xiàn)，因此需要通過(guò)仿真數(shù)據(jù)來(lái)模擬基線(xiàn)校正的過(guò)程，將改進(jìn)后的基線(xiàn)校正方法與改進(jìn)的非對(duì)稱(chēng)懲罰最小二乘法（IAsLS）[14]、懲罰最小二乘法（PLS）、非對(duì)稱(chēng)加權(quán)懲罰最小二乘法（asPLS）[8]和自適應(yīng)迭代重加權(quán)懲罰最小二乘法（airPLS）[9]4種基線(xiàn)校正方法在仿真光譜數(shù)據(jù)上進(jìn)行比較。LIBS光譜數(shù)據(jù)包含眾多重疊峰與強(qiáng)弱不一的峰，為了盡可能還原真實(shí)光譜數(shù)據(jù)，模擬光譜數(shù)據(jù)由120個(gè)位置、半峰寬各不相同的洛倫茲峰、理論基線(xiàn)和信噪比為40的噪聲組成。洛倫茲峰計(jì)算公式如（9）式：

式中：N為譜峰個(gè)數(shù)；Ai為光譜譜峰面積；hi為譜峰半高寬；ri為譜峰中心位置；x為光譜數(shù)據(jù)波長(zhǎng)范圍 [200,711.5] （nm）。

圖2為不同基線(xiàn)擬合方法對(duì)模擬光譜進(jìn)行擬合的結(jié)果，結(jié)合圖2（b）和圖2（c）可以看到，改進(jìn)后的算法在短波段和長(zhǎng)波段對(duì)模擬光譜均有良好的擬合效果，并在峰值處的基線(xiàn)抬升問(wèn)題也得到了緩解。為驗(yàn)證擬合基線(xiàn)的效果，計(jì)算去基線(xiàn)后譜峰峰值與無(wú)基線(xiàn)時(shí)譜峰峰值之間的均方根誤差Rpeak[15]以及擬合基線(xiàn)和理論基線(xiàn)之間的均方根誤差Rbase[16]，計(jì)算結(jié)果見(jiàn)表1。從表1可知，在處理模擬光譜時(shí)，本文算法所估計(jì)的基線(xiàn)與理論基線(xiàn)更為接近，且能在一定程度上保留原始譜峰信息。

表1 基線(xiàn)校正前后譜峰峰值的均方根誤差Rpeak以及擬合基線(xiàn)和理論基線(xiàn)之間的均方根誤差RbaseTable 1 Rpeak of peak intensity before and after baseline correction and Rbase between fitting baseline and theoretical baseline

圖2 模擬光譜基線(xiàn)擬合結(jié)果Fig.2 Baseline fitting results of simulated spectrum

2.2 實(shí)驗(yàn)

2.2.1 實(shí)驗(yàn)裝置與樣品配置

實(shí)驗(yàn)中采用的是共線(xiàn)雙脈沖激發(fā)方式，其示意圖如圖3所示，利用DG645控制激光出光的先后順序，第1束激光擊打在樣品表面產(chǎn)生等離子體，第2束激光則是對(duì)等離子體進(jìn)行再加熱，利用光譜儀采集等離子體在冷卻時(shí)發(fā)出的光譜。使用的光譜儀為L(zhǎng)TB公司的ARYELLE200型中階梯光譜儀，曝光時(shí)間設(shè)置為4 500 ms，延遲時(shí)間設(shè)置為3 μs，。制備樣品時(shí)，將標(biāo)準(zhǔn)土壤(GBW07386，GSS-32)或二氧化硅粉末(S5631-500)、PE微粉(HDPE，18180)和分析純化合物硝酸鈾酰，在瑪瑙研缽中研磨均勻并烘干，土壤和PE微粉的比例為3:1，共配置6個(gè)濃度的標(biāo)準(zhǔn)樣品。2束激光分別為波長(zhǎng)1 064 nm、能量90 mJ與波長(zhǎng)355 nm、能量50 mJ的激光。

圖3 LIBS實(shí)驗(yàn)示意圖Fig.3 Schematic diagram of LIBS experiment

2.2.2 繪制定標(biāo)曲線(xiàn)

為驗(yàn)證本文算法實(shí)際應(yīng)用的可行性，將其與4種不同基線(xiàn)校正方法進(jìn)行比較。圖4是對(duì)整體光譜數(shù)據(jù)去基線(xiàn)后，選取UII 454.36 nm特征譜線(xiàn)作為分析對(duì)象繪制的定標(biāo)曲線(xiàn)，可以看出，5種去基線(xiàn)方法下的定標(biāo)曲線(xiàn)模型均有良好的性能。再結(jié)合表2可知，本文算法得到的定標(biāo)曲線(xiàn)模型相關(guān)系數(shù)R2為0.997 2，且均方根誤差只有347.128 6，小于其他4種去基線(xiàn)方法得到的模型，即改進(jìn)后的基線(xiàn)校正方法處理后的樣本數(shù)據(jù)偏離定標(biāo)曲線(xiàn)的程度更小，說(shuō)明在實(shí)際使用定標(biāo)曲線(xiàn)進(jìn)行定量分析時(shí)，本文方法能使去基線(xiàn)后的定標(biāo)曲線(xiàn)模型依然有足夠的可靠性，更適合處理LIBS光譜數(shù)據(jù)的基線(xiàn)。

表2 定標(biāo)曲線(xiàn)模型的相關(guān)系數(shù)和均方根誤差Table 2 Correlation coefficients and RMSE of calibration curve model

圖4 采用不同基線(xiàn)校正方法后的定標(biāo)曲線(xiàn)Fig.4 Calibration curves with different baseline correction methods

3 結(jié)論

本文針對(duì)激光誘導(dǎo)擊穿光譜檢測(cè)技術(shù)中中階梯光柵光譜儀所引起的基線(xiàn)局部拱形突起問(wèn)題，提出一種改進(jìn)的懲罰最小二乘基線(xiàn)擬合方法。該方法通過(guò)在譜峰區(qū)局部增大懲罰函數(shù)來(lái)減小基線(xiàn)擬合方法在譜峰區(qū)域中的局部突起現(xiàn)象，同時(shí)還優(yōu)化了權(quán)值數(shù)組的判別標(biāo)準(zhǔn)，避免了基線(xiàn)擬合方法在基線(xiàn)區(qū)域出現(xiàn)“欠擬合”現(xiàn)象，并在與其他四種去基線(xiàn)方法的仿真實(shí)驗(yàn)中得到了驗(yàn)證。通過(guò)對(duì)實(shí)際數(shù)據(jù)進(jìn)行基線(xiàn)校正和繪制定標(biāo)曲線(xiàn)，也證明了本文算法在實(shí)際應(yīng)用中能提升定標(biāo)曲線(xiàn)的回歸效果，有助于后續(xù)對(duì)元素進(jìn)行定量分析。