范賢光, 王海濤, 王 昕, 許英杰, 王秀芬, 闕 靖
廈門大學(xué)物理與機電工程學(xué)院, 福建 廈門 361005
基于非均勻B樣條的拉曼光譜基線校正算法
范賢光, 王海濤, 王 昕*, 許英杰, 王秀芬, 闕 靖
廈門大學(xué)物理與機電工程學(xué)院, 福建 廈門 361005
基線校正是一種常用的消除光譜熒光干擾的方法, 是拉曼光譜數(shù)據(jù)處理的必要步驟之一。 傳統(tǒng)的多項式擬合基線校正算法, 簡單且易于實現(xiàn), 但是擬合階次難以確定, 靈活性較差。 使用非均勻B樣條代替多項式進行擬合, 在保留原有算法優(yōu)點的基礎(chǔ)上, 利用原始拉曼譜圖的峰位置信息自適應(yīng)地確定非均勻B樣條的節(jié)點向量, 然后以固定階次擬合光譜基線。 B樣條自身具有分段光滑的特性, 而計算樣條節(jié)點的節(jié)點向量自適應(yīng)選取算法中的峰位置信息通過使用兩次具有不同母函數(shù)的連續(xù)小波變換(continuous wavelet transform, CWT)來獲取, 既加強了原始光譜數(shù)據(jù)與B樣條算法本身的聯(lián)系, 也克服了傳統(tǒng)多項式擬合的不足。 為了驗證本文算法的有效性, 選取了甲基對硫磷和某品牌菜籽油兩種被測物進行實驗, 并使用該算法進行了基線校正, 并與兩種其他的基線校正算法與進行了對比。 實驗結(jié)果表明, 該方法利用固定的擬合階次就能達到較好的校正效果, 所需要的參數(shù)較少, 校正結(jié)果不會出現(xiàn)過擬合或欠擬合的現(xiàn)象, 是一種有效的拉曼光譜基線校正算法。
拉曼光譜; 基線校正; 非均勻B樣條; 節(jié)點向量; 峰位置
拉曼光譜(Raman spectroscopy)作為一種鑒定物質(zhì)結(jié)構(gòu)的分析測試手段, 廣泛應(yīng)用于材料、 化工、 石油、 高分子、 生物、 環(huán)保、 地質(zhì)等領(lǐng)域[1]。 由于熒光干擾等因素的存在, 光譜往往會發(fā)生較大漂移[2], 因此, 需要通過一定的方法對拉曼光譜進行基線校正。
目前, 主要的基線校正的方法有多項式擬合[3]、 小波變換[4]、 求導(dǎo)[5]等方法。 多項式擬合是最常用的基線校正算法, 通常使用最小二乘法, 通過循環(huán)迭代獲得拉曼光譜的基線。 此方法易于實現(xiàn), 但精度不高, 當(dāng)譜峰較多時, 擬合階次難以確定, 容易出現(xiàn)過擬合, 并且由于沒有考慮局部性, 當(dāng)某些數(shù)據(jù)點漂移時, 會影響整體擬合效果。
小波變換法利用小波變換對原始光譜進行分解, 對于分解到小波域的信號, 采取低頻置零, 高頻閾值過濾的操作, 截取到的有用信號再進行重構(gòu)。 此方法的前提條件是, 基線相對于原信號是變化緩慢的低頻信號, 基線和光譜特征峰的頻率有明顯界限。 此外, 小波函數(shù)的選取和變換尺度都較難確定, 且對校正效果的好壞有直接影響。 求導(dǎo)的方法遵循小波變換法“信號分解-信號處理-信號重構(gòu)”這樣的處理原則, 但是校正基線后, 會改變原始光譜的形狀。
針對以上問題, 本文提出一種基于非均勻B樣條的基線校正方法, 利用小波變換尋峰算法自適應(yīng)地確定非均勻B樣條的節(jié)點向量, 并使用固定階次的非均勻B樣條, 通過循環(huán)迭代擬合基線, 從而達到較好的基線校正效果。
1.1 非均勻B樣條擬合
B樣條方法以其低階光滑的特性, 被廣泛的應(yīng)用于自由曲線曲面的造型中[6-7]。 B樣條的曲線方程為
(1)
其中,di(i=0, 1, …,n)表示第i個控制系數(shù), B樣條曲線的控制多邊形就是由這些控制點依次連接而成。Ni, k(u)(i=0, 1, …,n)表示第i個k次(k+1階)B樣條基函數(shù), 定義為
(2)
其中,ui稱為節(jié)點, 它是節(jié)點向量U={u0,u1, …,un}這個單調(diào)非減集合的一個元素。 定義域被節(jié)點細分為一個個的節(jié)點區(qū)間。 如果這些節(jié)點區(qū)間是均勻的, 則稱使用這些節(jié)點的B樣條曲線為均勻B樣條曲線, 否則, 稱之為非均勻B樣條曲線。 計算一條k次B樣條曲線, 由方程(1)可知, 共需n+1個控制系數(shù)di(i=0, 1, …,n)以及n+1個k次B樣條基函數(shù)。 控制系數(shù)一般通過最小二乘法獲得[8], B樣條基函數(shù)則利用節(jié)點向量通過式(2)獲得。 因此, 節(jié)點向量對于計算整個B樣條曲線有著非常重要的作用, 節(jié)點的分布直接影響B(tài)樣條算法的擬合效果。
1.2 節(jié)點向量自適應(yīng)選取方法
拉曼光譜基線變化較大的位置一般發(fā)生在拉曼光譜的譜峰位置附近, 因此本文初步選取峰位置點作為B樣條算法的節(jié)點。 文獻[9]提出了一種利用連續(xù)小波變換實現(xiàn)拉曼譜峰識別的方法, 本文利用該方法獲得拉曼光譜的譜峰位置, 進而獲得非均勻B樣條節(jié)點向量。 小波變換是一種利用小波把信號從時域轉(zhuǎn)換到小波域的手段, 素有“數(shù)學(xué)顯微鏡”的美譽[10-11]。 連續(xù)小波變換的數(shù)學(xué)表達式可以定義為
(3)
其中,s(t)是原信號,a為尺度因子,b為時間平移因子,Ψa, b(t)是小波母函數(shù),C(a,b)為小波變換系數(shù), 表示不同尺度不同時間平移下小波函數(shù)的線性組合。 小波系數(shù)表示了選定小波母函數(shù)與處于特定時移b和特定尺度a下與原信號的相似度。 譜峰識別的算法中, 選取墨西哥帽小波作為小波母函數(shù)
(4)
使用墨西哥帽小波作為小波母函數(shù), 對拉曼信號進行連續(xù)小波變換后, 在每個尺度a上的小波系數(shù), 靠近峰中心位置附近時, 會出現(xiàn)模極大值。 并且這個模極大值隨尺度的變化而變化, 當(dāng)小波變換尺度與原譜峰峰形最匹配時, 模極大值達到最大。 將小波系數(shù)中的模極大值連接起來, 則形成如圖1所示的脊線。 該圖的橫坐標為拉曼位移, 縱坐標為連續(xù)小波變換尺度。 那些又長并且由較大的模極大值連接成的脊線就是較大的譜峰的位置。 由于噪聲的存在, 最初檢測到的譜峰會混雜有假峰, 因此, 需要對假峰進行剔除處理。 文獻[12]介紹了譜峰位置檢測的大量準則, 如信噪比、 脊線、 模極大值和峰寬等。 本文中的峰位置檢測按照基于信噪比和脊線的方法, 其中信噪比用來剔除假峰, 脊線長度閾值用來提取較大的峰位置。
獲得所有譜峰位置后, 其峰位置點P=[p0,p1, …,pn]還需進行進一步的處理, 才能作為非均勻B樣條的節(jié)點。 這是由于為保證B樣條擬合的效果, 其節(jié)點的分布不宜過密, 也不宜過稀。 過密的節(jié)點分布, 一方面增加擬合的計算量, 另一方面容易導(dǎo)致過擬合; 而過稀疏的節(jié)點分布, 容易導(dǎo)致欠擬合。 為保證節(jié)點的合理分布, 設(shè)定兩個參數(shù)d1和d2, 對峰位置P=[p0,p1, …,pn]進行處理, 從而確定B樣條節(jié)點向量。
Fig.1 Identified ridge lines based on CWT coefficient image
節(jié)點向量的獲取方法可分為以下幾個步驟:
(1) 使用墨西哥帽小波作為母函數(shù)對原拉曼信號進行連續(xù)小波變換;
(2) 連接不同尺度下的模極大值序列, 構(gòu)建脊線, 具體參考文獻[9];
(3) 根據(jù)參考文獻[13]的方法剔除假峰, 得到最終的峰位置向量P=[p0,p1, …,pn];
(5) 處理后的峰值點向量P賦值給U, 作為非均勻B樣條的節(jié)點向量。
1.3 非均勻B樣條基線校正的實現(xiàn)
本文利用改進的尋峰算法確定非均勻B樣條的節(jié)點向量, 然后利用非均勻B樣條循環(huán)逼近曲線基線, 原始光譜扣除基線后, 即得到基線校正后的曲線。 非均勻B樣條基線校正算法的原理可由圖2表示, 整個算法的基本步驟如下:
(1) 采樣獲得m個點構(gòu)成的原始光譜數(shù)據(jù)S0, 令y0=S0;
(2) 使用1.2節(jié)節(jié)點確定算法計算非均勻B樣條節(jié)點向量U;
(3) 對y0使用4階(3次)非均勻B樣條擬合, 得到y(tǒng)n作為初次擬合的基線;
(4)yn與y0逐點比較, 取較小的點賦值給yn;
(5)yn與y0的相對誤差若小于迭代閾值ζ, 則判定yn為光譜曲線的基線, 否則, 令y0=yn, 返回步驟3繼續(xù)執(zhí)行, 直到滿足條件為止;
(6) 校正后的光譜數(shù)據(jù)Scorrect=S0-yn。
Fig.2 Process of baseline correction by Non-Uniform B-spline fitting
2.1 材料及儀器
本文所選取的被測物質(zhì)是濃度為10 ppm的甲基對硫磷和濃度為10 ppm某品牌菜籽油。 所使用的實驗儀器為Ocean Optics公司的QE65Pro。
2.2 結(jié)果分析
利用本文算法對甲基對硫磷(parathion-methyl, PM)和菜籽油(colza oil, CO)進行處理, 兩種物質(zhì)的原始拉曼譜圖和擬合出的基線如圖3所示。 其中, 峰值檢測算法參數(shù)為: 小波變換尺度scales=1∶70, 剔除脊線的參數(shù)SNR=5和ridgeLength=10, 得到的峰值點分別為
PPM=[481, 614, 730, 828, 913, 1 006, 1 076, 1 196, 1 275, 1 357, 1 433, 1 505, 1 646, 2 370]
PCO=[390, 840, 972, 1 093, 1 255, 1 292, 1 436, 1 655, 1 736, 1 869, 1 889, 1 941, 1 975, 2 025, 2 087, 2 146, 2 183, 2 235, 2 342, 2 403]
選取的循環(huán)迭代閾值ζ=0.5, 非均勻B樣條擬合階數(shù)為4。 按照上述給定參數(shù), 根據(jù)本文算法擬合出的甲基對硫磷和菜籽油的拉曼譜圖的基線如圖3所示。 由圖3(a)可以看出, 本文所測的甲基對硫磷的拉曼譜圖在1 000~1 500 cm-1附近峰的分布較密, 拉曼位移從1 700~2 500 cm-1沒有明顯的拉曼峰。 根據(jù)峰值檢測算法計算出的甲基對硫磷的峰值點PPM, 峰值較密附近存在多達7個峰值點, 1 700~2 500 cm-1僅有1個峰值點, 因此PPM不能直接作為B樣條節(jié)點向量, 需要對PPM進行進一步處理。
根據(jù)節(jié)點向量自適應(yīng)選取方法得到參數(shù)d1=100,d2=250, 則處理后的非均勻B樣條節(jié)點為
UPM=[481, 614, 730, 913, 1 196, 1 357, 1 505, 1 646, 1 846, 2 180, 2 370]
UCO=[390, 615, 840, 972, 1 093, 1 255, 1 436, 1 655, 1 869, 1 975, 2 087, 2 161, 2 235, 2 342]
處理后的B樣條節(jié)點向量UPM, 節(jié)點分布合理, 得到的校正基線也較為平滑, 沒有在峰位置過密的1 000~1 500 cm-1區(qū)域出現(xiàn)過擬合, 也沒有在峰值點較少的1 700~2 500 cm-1上出現(xiàn)欠擬合。 校正后的甲基對硫磷的拉曼光譜圖如圖4(a)所示。 本文算法在進行基線校正時, 并沒有對原始拉曼數(shù)據(jù)進行平滑去噪等預(yù)處理操作, 圖3(b)為實驗所測的菜籽油的拉曼譜圖, 箭頭標示處存在高頻噪聲, 對比圖4(b)菜籽油基線校正后的拉曼譜圖, 可以發(fā)現(xiàn)菜籽油的譜圖的高頻噪聲處基線擬合依然平滑, 校正后的譜圖也沒有發(fā)生過擬合的現(xiàn)象, 證明了本算法的有效性。
Fig.3 Raman spectroscopy and its baseline by Non-Uniform B-spline fitting
選取傳統(tǒng)多項式擬合算法與改進的懲罰最小二乘算法[13]與本文算法進行比較。 圖5給出了使用5階多項式擬合的基線校正算法所得出的菜籽油拉曼譜圖及其校正基線。 所
Fig.4 Raman spectroscopy after baseline correction by Non-Uniform B-spline fitting
Fig.5 Raman spectroscopy and its baseline by polynomial fitting of colza oil
測得的菜籽油拉曼譜圖漂移較大, 且存在高頻噪聲, 利用5階多項式擬合菜籽油基線時, 800 cm-1處開始的基線出現(xiàn)明顯的過擬合現(xiàn)象, 1 800 cm-1處又出現(xiàn)欠擬合現(xiàn)象。 圖6為利用改進的懲罰最小二乘法獲得的甲基對硫磷的校正基線及其校正后的拉曼譜圖。 文獻[13]將該方法做成了R語言開源包baselineWavelet。 該方法先利用兩個小波變換獲取峰的信息, 然后在峰位置處使用直線連接峰首尾作為基線, 非峰部分使用懲罰最小二乘法擬合基線, 拼接后得到最終基線。 此方法計算出的基線與原始譜圖結(jié)合緊密, 不易出現(xiàn)過擬合的現(xiàn)象。 但是, 由于基線是兩種方法拼接而成, 該方法的基線整體不夠平滑, 如圖6(a)所示。 比較圖4(a)和圖6(b), 本文算法與改進的懲罰最小二乘算法校正效果差別不大, 但是, 此方法相較于本文算法, 默認參數(shù)眾多, 兩次小波變換計算量也較大。
Fig.6 Raman spectroscopy of parathion-methyl
(a): Raman spectroscopy and its baseline; (b): Raman spectroscopy after baseline correction by baseline wavelet
綜上所述, 采用基于非均勻B樣條的基線校正方法, 充分利用了B樣條的低階光滑特性, 克服了傳統(tǒng)多項式基線校正算法擬合階數(shù)難確定的缺陷, 在基線偏移嚴重和存在高頻噪聲時, 仍能夠擬合出光滑的基線, 并沒有出現(xiàn)欠擬合和過擬合的現(xiàn)象, 對于存在基線漂移的拉曼光譜信號, 達到了較好的基線校正效果。
提出了一種基于非均勻B樣條的拉曼光譜基線校正方法, 先利用基于小波的譜峰識別算法獲得拉曼譜圖的峰位置, 然后對這些峰位置點進行處理, 從而確定出非均勻B樣條的節(jié)點向量, 最后利用B樣條以固定的階次循環(huán)迭代擬合拉曼譜圖的基線。 與傳統(tǒng)的多項式擬合算法相比, 本文的自適應(yīng)節(jié)點確定算法充分利用了不同拉曼譜圖的信息, 讓擬合出的基線與原始拉曼譜圖結(jié)合更緊密; 同時, 本文的方法能夠以固定的階次擬合光譜基線, 并且能夠有效地避免過擬合和欠擬合。 因此, 本文所提出的方法, 能夠較好地實現(xiàn)拉曼光譜的基線校正, 并且調(diào)整參數(shù)少、 易于實現(xiàn), 可以作為一種有效的拉曼光譜基線校正算法。
[1] CHU Xiao-li(褚小立). Molecular Spectroscopy Analytical Technology Combined with Chemometrics and Its Applications(化學(xué)計量學(xué)方法與分子光譜分析技術(shù)). Beijing: Chemical Industry Press(北京: 化學(xué)工業(yè)出版社), 2011. 311.
[2] McCreer R L. Raman Spectroscopy for Chemical Analysis. Wiley-Interscience, 2000.
[3] FENG Xin-wei, ZHU Zhong-liang, SHEN Meng-jie, et al. Computer and Applied Chemsity, 2009, 26(6): 759.
[4] Jagtiani A V, Sawant R, Carletta J. Measurement Science and Technology, 2008, 19(6): 1.
[5] O’Grady A, Dennis AC, Denvir D. Analytical Chemistry, 2001, 73: 2058.
[6] Piegl L, Tiller W. The NURBS Book(非均勻有理B樣條). Translated by ZHAO Gang, MU Guo-wang, WANG La-zhu(趙 罡, 穆國旺, 王拉柱, 譯). Beijing: Tsinghua University Press(北京: 清華大學(xué)出版社), 2010. 60.
[7] GUO Jian-feng, ZHU Chang-qing(郭建峰, 朱長青). Engineering of Surveying and Mapping(測繪工程), 2000, 9(4): 25.
[8] SHI Fa-zhong(施法中). CAGD & NURBS(計算機輔助幾何設(shè)計與非均勻有理B樣條). Beijing: Higher Education Press(北京: 高等教育出版社), 2013. 303.
[9] Du P, Kibbe W A, Lin S M. Bioinformatics, 2006, 22: 2059.
[10] Daubechies I. Ten Lectures on Wavelets(小波十講). Translated by LI Jian-ping(李建平, 譯). Beijing: National Defense Industry Press(北京: 國防工業(yè)出版社), 2011.
[11] SUN Yan-kui(孫延奎). Wavelet Analysis and Applications. Beijing: China Machine Press, 2005.
[12] Yang C, He Z Y, Yu W C. BMC Bioinformatics, 2009. 10.
[13] Zhang Z M, Chen S, Liang Y Z, et al. Journal of Raman Spectroscopy, 2010, 41: 659.
*Corresponding author
Baseline Correction Algorithm for Raman Spectroscopy Based on Non-Uniform B-Spline
FAN Xian-guang, WANG Hai-tao, WANG Xin*, XU Ying-jie, WANG Xiu-fen, QUE Jing
School of Physics and Mechanical & Electrical Engineering, Xiamen University, Xiamen 361005, China
As one of the necessary steps for data processing of Raman spectroscopy, baseline correction is commonly used to eliminate the interference of fluorescence spectra. The traditional baseline correction algorithm based on polynomial fitting is simple and easy to implement, but its flexibility is poor due to the uncertain fitting order. In this paper, instead of using polynomial fitting, non-uniform B-spline is proposed to overcome the shortcomings of the traditional method. Based on the advantages of the traditional algorithm, the node vector of non-uniform B-spline is fixed adaptively using the peak position of the original Raman spectrum, and then the baseline is fitted with the fixed order. In order to verify this algorithm, the Raman spectra of parathion-methyl and colza oil are detected and their baselines are corrected using this algorithm, the result is made comparison with two other baseline correction algorithms. The experimental results show that the effect of baseline correction is improved by using this algorithm with a fixed fitting order and less parameters, and there is no over or under fitting phenomenon. Therefore, non-uniform B-spline is proved to be an effective baseline correction algorithm of Raman spectroscopy.
Raman spectroscopy; Baseline correction; Non-uniform B-spline; Node vector; Peak position
Dec. 23, 2014; accepted Apr. 18, 2015)
2014-12-23,
2015-04-18
國家重大科學(xué)儀器設(shè)備開發(fā)專項(2011YQ03012417)資助
范賢光, 1980年生, 廈門大學(xué)物理與機電工程學(xué)院機電工程系副教授 e-mail: fanxg@xmu.edu.cn *通訊聯(lián)系人 e-mail: xinwang@xmu.edu.cn
O657.3
A
10.3964/j.issn.1000-0593(2016)03-0724-05