基于DE算法和NURBS曲線的ICP-AES光譜基線校正方法研究

廉小親，陳彥銘，王宇喬，劉 鈺

1.北京工商大學(xué)人工智能學(xué)院，北京 100048 2.北京工商大學(xué)中國輕工業(yè)工業(yè)互聯(lián)網(wǎng)與大數(shù)據(jù)重點(diǎn)實(shí)驗(yàn)室，北京 100048

引 言

電感耦合等離子體原子發(fā)射光譜法(inductively coupled plasma atomic emission spectrometry，ICP-AES)是一種以電感耦合等離子體為激發(fā)光源的原子發(fā)射光譜分析技術(shù)，具有多元素同時(shí)檢測、分析速度快以及準(zhǔn)確度高等特點(diǎn)，被廣泛應(yīng)用于稀土、冶金、化工、無機(jī)材料和醫(yī)學(xué)等眾多領(lǐng)域[1-2]。然而在ICP-AES測量過程中，受儀器暗電流和光源雜散光的影響，測量所得的光譜通常存在一定程度的基線漂移現(xiàn)象，造成元素含量定量分析誤差，因此在測量過程中往往需要對光譜進(jìn)行基線校正處理[3-4]。

目前較為常見的基線校正方法包括迭代多項(xiàng)式擬合法[5-6]、導(dǎo)數(shù)法[7-8]、移動(dòng)窗口平滑法[9-10]、小波變換法[11]、分段多項(xiàng)式擬合法和分段樣條擬合法[12-14]等。迭代多項(xiàng)式擬合法的實(shí)現(xiàn)步驟較為簡單，但是對于不同的光譜信號，多項(xiàng)式擬合階數(shù)往往需要論證確定，難以建立一種通用性較強(qiáng)的光譜基線擬合模型。導(dǎo)數(shù)法和移動(dòng)窗口平滑法實(shí)現(xiàn)原理同樣相對簡單，但二者均會(huì)在時(shí)域內(nèi)降低光譜的分辨率，影響光譜的峰形，不利于后期的定量分析。小波變換也能夠有效去除光譜基線，但需要針對不同的光譜信號建立相應(yīng)的小波基函數(shù)，因此算法的運(yùn)算量相對較高，應(yīng)用效率也相對較低。分段多項(xiàng)式擬合法在一定程度上克服了迭代多項(xiàng)式擬合法中階數(shù)難以確定的問題，但是相比于分段樣條擬合法，分段多項(xiàng)式擬合法有時(shí)無法保證全波段范圍下的擬合誤差，臨界處擬合基線的平滑性也較差。樣條擬合法根據(jù)內(nèi)接點(diǎn)的分布方式可分為均勻樣條擬合和非均勻樣條擬合，均勻樣條擬合由于內(nèi)接點(diǎn)的分布缺乏目的性，當(dāng)光譜基線漂移較為劇烈時(shí)，最終的擬合效果和校正效果較為一般；而非均勻樣條擬合則彌補(bǔ)了均勻樣條擬合的不足之處，能夠擬合出更為理想的光譜基線，但非均勻樣條擬合過程中通常需要針對不同的光譜基線選定合適的內(nèi)接點(diǎn)序列和控制點(diǎn)序列，因此會(huì)引入大量的人力成本，且難以保證每次選取的內(nèi)接點(diǎn)序列和控制點(diǎn)序列較為適合。

針對以上問題，本文提出一種基于差分進(jìn)化算法(differential evolution，DE)和非均勻B樣條(non uniform rational B-spline，NURBS)曲線模型的ICP-AES光譜基線校正方法，通過對比光譜信號基線校正前和校正后的特征分布情況，基于NURBS曲線模型建立光譜基線校正評價(jià)函數(shù)，并利用差分進(jìn)化算法進(jìn)行求解，最終擬合出相應(yīng)的光譜基線；該方法不僅能夠充分發(fā)揮非均勻B樣條曲線在光譜基線擬合問題上的優(yōu)勢，也能夠利用差分進(jìn)化算法全局尋優(yōu)的特點(diǎn)自行計(jì)算出合適的內(nèi)接點(diǎn)序列和控制點(diǎn)序列，有效地克服了非均勻B樣條曲線擬合過程中內(nèi)接點(diǎn)序列和控制點(diǎn)序列選取困難的問題，同時(shí)也能夠獲得更好的基線校正效果。

1 實(shí)驗(yàn)部分

1.1 ICP-AES光譜基線校正原理

本文提出的基于差分進(jìn)化算法和非均勻B樣條曲線的ICP-AES光譜基線校正方法流程如圖1所示。原始光譜數(shù)據(jù)首先通過高斯濾波消除部分高頻噪聲信號，其次篩選光譜信號中的極小值點(diǎn)，以極小值點(diǎn)集合的相對標(biāo)準(zhǔn)偏差為度量準(zhǔn)則、以NURBS曲線模型的控制點(diǎn)序列和內(nèi)接點(diǎn)序列為尋優(yōu)參數(shù)，建立光譜基線校正評價(jià)函數(shù)，然后利用差分進(jìn)化算法求解該評價(jià)函數(shù)的最優(yōu)解，即使得極小值點(diǎn)集合相對標(biāo)準(zhǔn)偏差最小時(shí)的控制點(diǎn)序列和內(nèi)接點(diǎn)序列，最后，通過NURBS曲線模型擬合該控制點(diǎn)序列和內(nèi)接點(diǎn)序列對應(yīng)的曲線作為光譜基線，以此實(shí)現(xiàn)光譜基線校正。

圖1 基線校正流程

1.2 光譜預(yù)處理

在ICP-AES測量過程中，受雜散光影響，光譜中包含大量噪聲信號，嚴(yán)重影響基線擬合過程中波谷點(diǎn)的選擇，因此在進(jìn)行基線擬合前首先需進(jìn)行預(yù)處理以消除噪聲。本文采用高斯濾波對ICP-AES光譜進(jìn)行預(yù)處理。高斯濾波是數(shù)字信號處理中常用的濾波方式，其本質(zhì)為原始信號與高斯卷積核的離散卷積運(yùn)算，如式(1)所示。

(1)

式(1)中，f可視為原始信號，g可視為高斯卷積核。由于ICP-AES光譜信號基本服從高斯分布，因此通過高斯卷積核進(jìn)行卷積能夠盡可能地保留譜線信息，避免降低光譜分辨率；另一方面，本文對實(shí)測的ICP-AES光譜數(shù)據(jù)進(jìn)行分析，抽樣統(tǒng)計(jì)不同波段中的噪聲信號分布特征，擬合噪聲數(shù)據(jù)的概率密度函數(shù)以及概率分布，如圖2和圖3所示，根據(jù)統(tǒng)計(jì)結(jié)果可知，ICP-AES光譜中噪聲信號概率密度函數(shù)基本服從高斯分布，因此認(rèn)為通過高斯濾波能夠較為有效地去除噪聲信號。

圖2 光譜噪聲信號概率密度

圖3 光譜噪聲信號概率分布

1.3 建立光譜基線模型

本文采用非均勻B樣條曲線模型進(jìn)行基線擬合。非均勻B樣條曲線常用于自由曲線和曲面的構(gòu)造，具有良好的低階光滑性，因此同樣適用于擬合光譜基線。記數(shù)據(jù)點(diǎn)序列集合為X，X={(x1|u,x1|v),(x2|u,x2|v),(x3|u,x3|v), …,(xn|u,xn|v)}，采用節(jié)點(diǎn)符號t將數(shù)據(jù)點(diǎn)序列X進(jìn)行劃分，并記節(jié)點(diǎn)序列集合為T，T=[t-k+1,t-k+2, …,t-1,t0,t1, …,tn,tn+1,tn+2, …,tn+k]，則節(jié)點(diǎn)序列滿足式(2)關(guān)系

t-k+1≤t-k+2≤…≤t-1≤t0

(2)

式(2)中，t1到tn稱為內(nèi)節(jié)點(diǎn)，其余稱之為外節(jié)點(diǎn)，n為內(nèi)節(jié)點(diǎn)數(shù)，k為B樣條曲線的階數(shù)；若內(nèi)節(jié)點(diǎn)均勻分布，則最終構(gòu)成的曲線稱之為均勻B樣條曲線，反之則為非均勻B樣條曲線。

對于k階的B樣條曲線，其表達(dá)式如式(3)所示。

(3)

(4)

在計(jì)算過程中，控制點(diǎn)序列和曲線坐標(biāo)點(diǎn)序列通常均以坐標(biāo)的形式進(jìn)行運(yùn)算，因此式(3)也可寫為

(5)

同理，可記控制點(diǎn)序列為C={(c0|u,c0|v),(c1|u,c1|v), …,(cn|u,cn|v)}。

1.4 建立ICP-AES光譜基線校正評價(jià)函數(shù)

根據(jù)1.3節(jié)NURBS曲線的定義可知，控制點(diǎn)序列C和內(nèi)接點(diǎn)序列T直接決定了NURBS曲線的分布，因此若采用NURBS曲線模型擬合光譜基線，則需要尋找一組合適的控制點(diǎn)序列C和內(nèi)接點(diǎn)序列T，若控制點(diǎn)序列C和內(nèi)接點(diǎn)序列T的選取不恰當(dāng)，將會(huì)直接影響最終基線校正的結(jié)果。

依靠專家經(jīng)驗(yàn)和光譜分布情況進(jìn)行控制點(diǎn)序列C和內(nèi)接點(diǎn)序列T的選取一方面效率較低，另一方面其選擇結(jié)果也難以保證最終光譜基線校正的效果，因此本文建立了ICP-AES光譜基線校正評價(jià)函數(shù)，該評價(jià)函數(shù)基于NURBS曲線建立，并以光譜基線校正前后極小值的分布特征為評價(jià)指標(biāo)，通過求解評價(jià)函數(shù)的全局最優(yōu)解即可得到合適的控制點(diǎn)序列C和內(nèi)接點(diǎn)序列T，有效改善了NURBS曲線的應(yīng)用效率，同時(shí)也提高了控制點(diǎn)序列C和內(nèi)接點(diǎn)序列T選取的準(zhǔn)確性。

基線校正前和基線校正后，光譜及其極小值點(diǎn)分布差異如圖4和圖5所示。

圖4 基線校正前后光譜信號

圖5 基線校正前后光譜極小值

根據(jù)基線校正前和校正后的光譜信號可知，光譜數(shù)據(jù)點(diǎn)中的極小值序列的分布發(fā)生了明顯變化，顯然，基線校正前，光譜數(shù)據(jù)點(diǎn)中極小值序列的離散程度更高，基線校正后，光譜數(shù)據(jù)點(diǎn)中的極小值序列的離散程度更低，分布更加接近極小值序列的算術(shù)平均值。

根據(jù)標(biāo)準(zhǔn)偏差的數(shù)學(xué)定義可知，標(biāo)準(zhǔn)偏差常用于度量數(shù)據(jù)分布的離散程度，用以衡量數(shù)據(jù)值偏離算術(shù)平均值的程度，因此可以利用光譜數(shù)據(jù)點(diǎn)中極小值序列的標(biāo)準(zhǔn)偏差值作為光譜基線校正程度的度量準(zhǔn)測，以此來建立光譜基線校正評價(jià)函數(shù)，達(dá)到基線校正的目的。

設(shè)濾波處理后且未校正基線的光譜強(qiáng)度值服從模型I′(x)，基線校正后的光譜強(qiáng)度值服從模型I(x)，本文采用非均勻B樣條曲線擬合光譜基線，因此設(shè)基線模型為B(C,T)；由于光譜強(qiáng)度值的疊加可視為線性疊加，因此上述模型滿足如式(6)所示的關(guān)系：

I(x)=I′(x)-B(C,T)

(6)

設(shè)光譜中極小值序列集合為P，P={(p0,I′(p0)),(p1,I′(p1)), …,(pn,I′(pn))}，并記PX={p0,p1, …,pn}，PY={I′(p0),I′(p1), …,I′(pn)}，經(jīng)基線校正后，集合PX對應(yīng)的校正強(qiáng)度值集合為IY，IY={I(p0),I(p1), …,I(pn)}。

以校正強(qiáng)度值集合IY的標(biāo)準(zhǔn)偏差為度量準(zhǔn)測，建立評價(jià)函數(shù)，如式(7)所示。

F(C,T)=STD({I(p0)，I(p1)，…，I(pn)})=

STD({I′(p0)-yv|yu=p0，I′(p1)-yv|yu=p1，…，

I′(pn)-yv|yu=pn})

(7)

式(7)中，yu和yv定義與式(5)中定義相同，均表示NURBS曲線橫、縱坐標(biāo)，STD表示計(jì)算集合的標(biāo)準(zhǔn)偏差。因此，光譜基線校正問題可以轉(zhuǎn)換為，確定一組合適的C、T取值，使得F(C,T)取值達(dá)到最小，即，針對F(C,T)進(jìn)行全局尋優(yōu)的問題。而在全局尋優(yōu)方面，差分進(jìn)化算法則能夠快速準(zhǔn)確的達(dá)到上述目標(biāo)。

1.5 求解評價(jià)函數(shù)全局最優(yōu)解

根據(jù)1.4節(jié)內(nèi)容可知，本文提出的ICP-AES光譜基線校正方法可以轉(zhuǎn)換為求解評價(jià)函數(shù)的全局最優(yōu)解。遺傳算法在全局尋優(yōu)問題上往往能夠獲得較為理想的效果，本文采用基于遺傳算法改進(jìn)的差分進(jìn)化算法求解式(7)的全局最優(yōu)解。

差分進(jìn)化算法主要包括種群初始化、變異、交叉和選擇四個(gè)步驟。種群初始化后，循環(huán)執(zhí)行變異、交叉和選擇操作，直到達(dá)到最大進(jìn)化代數(shù)；在每一輪進(jìn)化中對種群所有個(gè)體向量執(zhí)行變異操作，得到變異向量；然后將變異向量與個(gè)體向量進(jìn)行交叉操作，得到交叉向量；最后對比交叉向量與個(gè)體向量的適應(yīng)度，選取適應(yīng)度較小的向量作為下一代種群的個(gè)體向量。差分進(jìn)化算法流程如圖6所示。

圖6 差分進(jìn)化算法流程圖

(1)種群初始化

記種群個(gè)體向量為X=[x1,x2, …,xD]，其中D為個(gè)體向量維度；設(shè)定種群個(gè)體向量中各參數(shù)搜索范圍，并生成最大個(gè)體向量和最小個(gè)體向量，如式(8)所示。

(8)

通過初始化生成NP個(gè)向量，表示種群的所有個(gè)體，則第G代種群中第i個(gè)個(gè)體向量可表示為

(9)

個(gè)體向量中第j個(gè)參數(shù)的生成規(guī)則如式(10)所示。

j=1，2，…，D

(10)

(2)個(gè)體適應(yīng)度評價(jià)

將當(dāng)前進(jìn)化代數(shù)內(nèi)所有的個(gè)體向量Xi, G代入適應(yīng)度函數(shù)F(·)，輸出的結(jié)果F(Xi, G)為個(gè)體的適應(yīng)度值，適應(yīng)度值越小，則相應(yīng)的個(gè)體越優(yōu)。

(3)變異

每一輪進(jìn)化中，對當(dāng)前進(jìn)化代數(shù)G下的所有的個(gè)體向量Xi, G執(zhí)行變異操作得到變異向量Vi, G。本文采用DE/rand/1/bin變異算子，其表達(dá)式如式(11)所示。

Vi,G=Xir1，G+SF(Xir2，G-Xir3，G)

(11)

式(11)中，ir1，ir2和ir3為當(dāng)前進(jìn)化代數(shù)內(nèi)種群中的三個(gè)隨機(jī)個(gè)體向量索引，且必須滿足ir1≠ir2≠ir3，SF為[0, 1]內(nèi)的常數(shù)，表示縮放因子。

(4)交叉

將所有個(gè)體向量Xi, G與對應(yīng)的變異向量Vi, G的信息成分交叉互換，得到交叉向量Ui, G。交叉規(guī)則如式(12)所示。

(12)

式(12)中，CR為(0, 1)內(nèi)的常數(shù)，表示交叉概率，jrand為[1,D]內(nèi)的隨機(jī)整數(shù)，表示隨機(jī)交叉參數(shù)索引。

(5)選擇

將所有個(gè)體向量Xi, G和交叉向量Ui, G代入適應(yīng)度函數(shù)進(jìn)行比較，選擇適應(yīng)度較小的向量作為下一代進(jìn)化的新個(gè)體向量，如式(13)所示。

(13)

2 結(jié)果與討論

2.1 測試方案

2.1.1 測試數(shù)據(jù)來源

本文搭建一臺(tái)CCD型掃描式ICP-AES光譜儀，并利用該光譜儀掃描一組樣品溶液對應(yīng)的ICP-AES光譜信號，該溶液所含元素及各元素對應(yīng)的特征波長如表1所示。

表1 樣品溶液所含元素及特征波長

選取其中的294.547～297.400 nm的光譜作為原始樣本并進(jìn)行歸一化處理，如圖7所示；對原始樣本進(jìn)行高斯濾波處理，消除高頻噪聲，并作為最終的基線校正測試樣本，如圖8所示。

圖7 原始樣本光譜信號

圖8 測試樣本光譜信號

2.1.2 基線校正模型特征參數(shù)取值范圍

由1.4節(jié)和1.5節(jié)內(nèi)容可知，利用差分進(jìn)化算法求解本文所建立的光譜基線校正評價(jià)函數(shù)時(shí)，首先要確定出評價(jià)函數(shù)中特征參數(shù)的取值范圍，即控制點(diǎn)序列C和內(nèi)接點(diǎn)序列T的取值范圍。

由于控制點(diǎn)序列C和內(nèi)接點(diǎn)序列T的取值范圍決定了評價(jià)函數(shù)全局最優(yōu)解的分布區(qū)間，進(jìn)而決定了NURBS曲線的趨勢，因此對最終基線擬合結(jié)果的影響較大。本文為保證基線擬合結(jié)果的可靠性，充分利用差分進(jìn)化算法全局尋優(yōu)的能力，以及本文建立的評價(jià)函數(shù)的優(yōu)點(diǎn)，設(shè)定內(nèi)接點(diǎn)序列T和控制點(diǎn)序列C的取值范圍為有效取值區(qū)間內(nèi)的最大范圍。對于NURBS曲線，內(nèi)接點(diǎn)序列T的有效取值區(qū)間為[0, 1]，故本文繼續(xù)采用該取值區(qū)間作為差分進(jìn)化算法的搜索范圍；而對于控制點(diǎn)序列C，本文則以樣本光譜信號的最小值和最大值分別作為取值區(qū)間的臨界值。

2.2 測試結(jié)果

2.2.1 差分進(jìn)化算法計(jì)算結(jié)果

設(shè)定差分進(jìn)化算法初始條件，取NP=100，G=60，進(jìn)行迭代運(yùn)算，最終計(jì)算完成后，控制點(diǎn)序列C分布如圖9所示，內(nèi)接點(diǎn)序列T分布結(jié)果如圖10所示。

圖9 控制點(diǎn)序列C分布

圖10 內(nèi)接點(diǎn)序列T分布

2.2.2 NURBS曲線擬合結(jié)果

根據(jù)計(jì)算出的控制點(diǎn)序列C和內(nèi)接點(diǎn)序列T擬合對應(yīng)的NURBS曲線作為測試樣本的光譜基線，并對樣本進(jìn)行基線校正，得到校正后的光譜，如圖11所示。

圖11 測試樣本光譜、擬合基線和基線校正后光譜

2.2.3 不同噪聲情況下基線擬合結(jié)果對比

在本文提出的基于差分進(jìn)化算法和NURBS曲線模型的ICP-AES光譜基線校正方法中，基線擬合效果和基線校正效果一定程度上受噪聲分布影響；當(dāng)噪聲未充分消除時(shí)，測試樣本光譜信號中極小值點(diǎn)序列的密度將會(huì)增大，導(dǎo)致最終在擬合基線的過程中存在一定程度的過擬合現(xiàn)象，如圖12所示。

圖12 噪聲未充分消除時(shí)基線擬合情況及基線校正情況

而通過高斯濾波或其他濾波方式消除噪聲時(shí)，若濾波效果過于激進(jìn)，則一定程度上可以避免過擬合現(xiàn)象，但同時(shí)則可能會(huì)導(dǎo)致測試樣本光譜信號時(shí)域分辨率下降，波峰波谷出現(xiàn)失真，最終影響基線擬合效果和定性定量分析效果。

因此在實(shí)驗(yàn)過程中，若存在噪聲干擾較為復(fù)雜的情況，可以通過專家判定方式，視情況選定極小值序列P的取值，并建立評價(jià)函數(shù)進(jìn)行求解，從而獲得較好的基線校正和基線擬合結(jié)果。

3 結(jié) 論

提出了一種基于差分進(jìn)化算法和NURBS曲線模型的ICP-AES光譜基線校正方法，該方法依據(jù)基線校正前和基線校正后光譜信號的特征變化以及NURBS曲線模型建立評價(jià)函數(shù)，并通過差分進(jìn)化算法進(jìn)行求解，計(jì)算出適當(dāng)?shù)膬?nèi)接點(diǎn)序列T和控制點(diǎn)序列C，進(jìn)而擬合出合適的光譜基線，達(dá)到基線校正的目的。

本文提出的方法不僅能夠充分利用NURBS曲線的自身優(yōu)勢擬合出不同類型的光譜基線，同時(shí)也能夠利用差分進(jìn)化算法全局尋優(yōu)的特點(diǎn)計(jì)算出合適的內(nèi)接點(diǎn)序列T和控制點(diǎn)序列C，極大的提高了NURBS曲線的應(yīng)用效率，也改善了NURBS曲線應(yīng)用于光譜基線校正領(lǐng)域時(shí)內(nèi)接點(diǎn)序列T和控制點(diǎn)序列C難以確定的局限性。但本文所提出的方法仍具有一定的不足之處，例如當(dāng)光譜信號中噪聲較為復(fù)雜時(shí)，可能會(huì)出現(xiàn)基線過擬合現(xiàn)象，此時(shí)仍然會(huì)引入部分人為工作量，因此在后續(xù)的優(yōu)化中，將針對以上問題進(jìn)行重點(diǎn)研究，進(jìn)一步提高NURBS曲線在光譜基線校正過程中的應(yīng)用效率。