Improved iterative smoothing for correcting baseline in atomic emission spectra fitting algorithm

原子发射光谱是一种在元素分析领域广泛应用的光谱技术。它基于原子在受激条件下从激发态返回基态时发射特定频率的光谱现象，通过测量和分析这些发射谱线来确定样品中的元素种类及其相对含量。每种元素会在特定频率上拥有发射谱线。原子发射光谱技术同时具备多项优势，如分析速度快、检测限低、高灵敏度以及低成本，因此在冶金^[1]、矿物学^[2]、地质学^[3]等众多领域广泛应用。然而，用光谱仪器获得的光谱信号存在各种干扰因素，其中包括样品中基体元素对有效光谱信号的干扰以及杂散光和暗电流干扰^[4]。光谱基线的增大不仅会降低信背比，还可能影响元素后续定量分析的准确性，这会严重削弱原子发射光谱的分析能力。因此，进行基线校正是提取真实的原子发射光谱信号和实现背景消除的关键步骤。

常用的基线校正方法有迭代多项式拟合^[5-6]、小波变换^[7-9]、惩罚最小二乘^[10-12]、非对称最小二乘^[13-15]、连续背景极小点筛选分段拟合^[16]等方法。非对称最小二乘算法方法虽然简单，但容易出现欠拟合现象。多项式拟合法^[4]是通过最小二乘法确定多项式系数，来拟合光谱背景基线，但存在拟合基线高估的情况。小波变换^[8]是进行光谱信号拆解，有助于分离低频和高频信号成分，但计算量和难度较大。加权惩罚最小二乘^[10]是通过给谱峰设置小的权重，将峰值区域降低到基线，但谱峰会为基线估计带来干扰。由于原子发射光谱的数据量庞大且复杂，应用其它光谱背景校正方法时，会存在基线抬升、过拟合、欠拟合、稳定性差、不平滑等现象。

针对上述方法存在的问题，本文中提出了一种平滑迭代模型的原子发射光谱基线校正方法，通过窗口划分原始光谱区间，选取每个区间中一定数量的较小点，将较小点进行线性插值补全点数，通过平滑迭代滤波对插值后的数据处理，对比迭代前后各点斜率变化率，最终拟合出相应的光谱基线。该方法不仅能削弱因窗口宽度问题选取较小点出现异常点的问题，也能够利用平滑滤波实现异常点的优化，能获得很好的基线校正效果。

基于最小二乘原理的Stavitzky-Golay(SG)滤波^[17]是一种在移动窗口内进行加权多项式平滑。在原始光谱数据序列x_i处左右各取M个数据，利用多项式拟合所获得的点预估出x_i的数据。公式如下：

式中：α_k为滤波的多项式系数；ε为拟合曲线与原始光谱的残差；i为序列号；p₀为原始光谱数据；k为多项式的阶次；p(x_i)为平滑后的预估值。设置ε取最小值作为边界条件，求解可得多项式最佳系数。ε越小，则原数据拟合程度越好。平滑窗口在光谱中移动，每次移动步长为1。考虑到本文中是现实原子发射光谱的基线拟合，由于SG滤波的特点首尾端点处窗口长度不足，在用SG滤波的方法计算首尾两端处要具有光滑性，所以在首尾处要设定窗口一端长度不足时，将另一端的数据对称来补全窗口缺失部分的数据。

用SG滤波平滑迭代拟合光谱基线时，通过窗口取的较小点中会存在远高于基线的异常值点，异常点相对于基线是高频信号，则可视其为突变信号，当进行平滑时，异常点位置就会逐渐下降，但同时紧邻的光谱会逐渐升高。

SG滤波具有谱线平滑的特点，使用其拟合原子发射光谱基线的同时，又能对基线实现一定的平滑作用。故此，本文作者提出了改进的平滑迭代并用斜率变化率判定的拟合基线校正算法。

采用区间选取较小点，能够一定程度上预估出光谱基线的大致位置，但对于谱峰重叠较多的光谱来讲，区间选取较小点时，极有可能选取到重叠峰上的较小点，故对较小点要进行一定的处理，确保偏离基线的点能够被处理。通过SG滤波平滑迭代可以将偏离基线的异常点校正到基线位置，如图 1所示。图中点表示初始窗口选取的较小点，包括远离基线的异常点，通过对较小点形成的初始基线进行平滑，随着滤波次数增加异常点逐渐逼近基线。

Figure 1.  Spectrum smoothing iteration process for small points

由于异常点相对基线来讲是突变信号，在多次迭代后会逐渐逼近基线，可以通过对比迭代前后对应点斜率变化率来判断确定基线受异常点的影响情况，定义：

式中：K_j(x_i)为第j次滤波后点x_i处的斜率；j代表迭代次数；ΔK为滤波前后斜率的变化率。

因此，本文中探讨了一种选取区间内较小值点的确定基线点、对异常点进行迭代平滑，并利用迭代前后各点处斜率变化率来判断拟合校正基线的方法，其流程图如图 2所示。通过采用SG滤波平滑算法，连续扩大SG滤波窗口，进行多轮迭代，逐渐趋向于光谱的基线。在原始光谱经过基线校正处理后，可获得校正后的光谱。光谱基线校正信号处理的步骤如下所述：(a)定义窗口。定义一个长度为N的窗口，窗口在原始光谱R中滑动，选取窗口中b个较小点，窗口滑动时不重叠，所有较小点进行线性插值，以生成平滑数据B；(b)进行SG滤波。对基线数据集合B进行SG滤波，用F表示滤波后的数据集合；(c)计算斜率变化率。使用1阶差分法计算数据集合F中各数据点的斜率变化率；(d)概率统计。对斜率变化率进行概率统计，以确定何时终止滤波。用差分法计算滤波前后各点的斜率和滤波前后各点的斜率变化率，设定一个阈值δ，如果超过T%的数据点的斜率变化率不超过阈值δ，则滤波结束国；(e)迭代。如果不满足概率统计条件，选择数据集合F和原始数据集合B中相应位置的较小值，生成新的数据集合R′，并增加SG滤波的窗口长度；(f)得到最优基线估计结果。重复步骤(b)~步骤(e)，在斜率达到步骤(d)中的条件时，迭代结束。

Figure 2.  SG iterative filtering and fitting method flowchart

在上述基线校正算法中引入N，b，δ，T4个调节参数。通过大量实验获取各类型合金样品的原子发射光谱，互补金属氧化物半导体(complementary metal oxide semiconductor, CMOS)的像素点个数为2048时，光谱谱峰信号点的个数在15~25个左右，窗口的取值N应大于谱峰信号点的个数，范围为30~40。b的取值不宜过大或过小，过大会在异常点会增多，增加后续平滑迭代次数，计算量增大；过小会在丢失部分有用的基线信息，后续拟合基线会存在偏差，故b的取值应为窗口长度的1/3~1/2之间。δ取决于迭代滤波逐渐增大的窗口宽度，窗口每次迭代增加的宽度越大，δ相对应也要增大，δ的取值范围为：0.0001~0.0004。D的取值需和δ配合。

这一方法的核心思想是通过SG平滑迭代滤波来逐步降低噪声和异常点的影响，以获得准确的光谱基线估计。同时，通过线性插值和斜率比较来处理光谱数据中的峰重叠和异常点，以确保最终光谱基线的估计结果准确可靠，这种方法适用于处理复杂的光谱数据。

模拟光谱包含有15个高斯特征峰和信噪比为40 dB的噪声基线构成，包含2000个数据点。高斯特征峰计算公式^[18]如：

式中：A为高斯函数峰值强度；μ为峰值位置的横坐标；λ为波长；σ是标准差表征控制峰的宽度。

综合考虑简便程度以及算法效果，较小点选取的窗口大小为32，每个窗口中选择较小的12个点，SG滤波初始窗口长度为20，每次迭代窗口加15，设定一个斜率变化率的阈值δ=0.2%，如果超过90%的数据点的斜率变化率不超过阈值0.2%，迭代终止输出最佳基线。

仿真光谱以及采用的非对称最小二乘^[19](asymmetric least squares，ALS)、model-free方法^[20]、minima方法^[21]和本文中的方法得到的拟合基线如图 3所示。

Figure 3.  Simulated spectrum baseline fitting result

由图 3可以看出，在拟合无峰处基线时，图中方法都对基线有较好的拟合效果，ALS方法出现了一定的基线高估情况；在拟合谱峰处基线时，model-free方法会在重叠峰处存在凸起，拟合准确度不高；在拟合小峰处基线时，minima方法偶尔会取到小峰处的点，导致基线稳定性降低。在3种不同类型的光谱中，本文中算法均获得了拟合准确的基线，具有较好的准确性、稳定性及平滑性。

为了更具体地比较不同方法对模拟光谱信号拟合的效果，采用不同方法预估出的基线，然后求取基线相对标准偏差(relative standard deviation, RSD)^[22]。为了更详细地评估不同方法在拟合模拟光谱方面的性能，首先使用各方法拟合出的模拟光谱的基线，然后计算基线的相对标准偏差。较小的相对标准偏差表明经校正后的基线更加平稳，拟合效果更为优越，各方法拟合基线信号的相对标准偏差如表 1所示。由表中可以看出，与其它3种方法相比，本文中算法在模拟光谱数据中得到的RSD最小为8.25%，相对于minima拟合降低了25.14%，其在拟合基线时波动最小，拟合效果最好。

根据表 1所示数据，本研究方法相对于其它3种基线拟合方法，在模拟光谱中获得了最小RSD。该结果表明，本文中所提出的算法有效地实现了光谱基线的平滑，降低了波动，同时获得了最佳的拟合效果。因此，该方法可有效拟合原子发射光谱的基线状况，从而提高后续定量分析的准确性。

实验中采用自搭的光电激发光谱平台(激发光源电压为24 V，光源放电频率为400 Hz)以及自搭的激光诱导击穿光谱(laser-induced breakdown spectroscopy, LIBS)平台(激光器为Nd ∶ YAG脉冲激光器，波长为1064 nm，能量为50 mJ，频率为1 Hz)选取圆柱型低碳钢、镁合金、锌合金的标准样品。预处理利用车床将样品表面氧化层去除，避免影响后续定量分析结果。完成样品测试，获得原始原子发射光谱数据R。

对3种方法与本文中的算法进行了比较。选择不同特征的原子发射光谱进行拟合基线对比验证算法，分别为光电激发的低碳钢标样光谱、光电激发的镁合金标样光谱、激光激发锌合金标样光谱，结果如图 4所示。

Figure 4.  Application of various algorithms in fitting baselines for practical atomic emission spectroscopy

图 4a中minima方法对基线的拟合效果较好，但在谱峰附近拟合了周围的小峰；图 4b中信号区域基线曲率突变或者出现阶跃时，model-free方法基线拟合过低，拟合后的信号仍存在基线部分；图 4c中ALS方法得到的基线，在谱峰出存在凸起，基线拟合过高。当光谱区域基线曲率变化不大时，每种方法的的拟合结果很相近，当基线附近出现跃迁或小峰时，本文中算法较其它算法都更稳定、准确。综上所述，本文中算法在谱峰密度大的光谱中拟合基线，相比与其它方法不会存在基线高估的情况，同时也能保持较好的平滑性；在光谱信号区域曲率突变时，本文中算法较其它算法不会出现欠拟合的情况，能更准确、稳定得到拟合基线；在光谱较为稳定的区域，本文中算法表现出更小的波动性。

分别使用上述4种方法对3种原始光谱进行基线校正，校正后的光谱如图 5所示。通过比较本文中算法处理后的光谱底部较为平稳，保留各谱峰的强度较高；minima方法的整体光谱波动最大。

Figure 5.  Effect of various algorithms in removing baselines from practical atomic emission spectra

为了验证不同算法对元素定量分析的影响，将MgAlZn系镁合金标样中的锌元素作为定量分析对象，先用汞灯对光谱进行波长定标，再分别使用ALS、minima、model-free和本文中算法对镁合金标样光谱进行基线校正。图 6为对整体光谱数据进行基线校正后获取的有效光谱。选取255.794 nm特征谱线作为分析对象绘制的定量曲线，可以看出，4种去基线方法下的定量曲线模型都有良好的定量分析能力。

Figure 6.  Calibration curves with different baseline correction

表 2为统计的各算法的锌元素定标曲线模型的相关系数R²和预测均方根误差(root mean square error prediction, RMSEP)，可用来评价定量分析的准确度。相关系数越高证明定量分析模型越精准，本文中算法获得最大的相关性为99.85%，大于其它3种去基线方法，即本文中的算法基线校正方法处理后的样本数据的定标曲线更精准。预测均方根误差越小，说明定量分析的误差越小，由表 2可清晰看出，本文中算法的RMSEP仅为0.5912，相比minima方法预测均方根误差降低了80.96%。

针对原子发射光谱存在的背景干扰强的问题，提出了一种改进的平滑迭代的光谱信号基线拟合方法。该方法通过窗口选取区间中较小点，进行插值计算出初始基线，通过SG迭代逐渐降低初始基线中异常值点对基线的影响，对比平滑迭代前后各点斜率变化率，判断基线中异常点的影响，直到迭代到异常点对基线的影响基本忽略不计，获得到最终的基线。本文中的算法不仅避免了拟合方法在谱峰区间出现的过拟合现象，相对于其它3种方法相比，还拥有较高的稳定性和平滑性。模拟光谱基线中得到基线的RSD最小为8.25%，相对标准偏差仅为minima拟合的24.7%，为model-free拟合的70.6%，为ALS拟合的85.8%。在实际光谱中，用4种方法拟合MgAlZn系镁合金标样光谱基线，将锌元素建立定量分析模型，本文中的算法获得最大的相关性为99.85%，大于其它3种去基线方法，且预测均方根误差最小为0.5912，比minima降低了80.96%。因此，本文中提出的算法能较好拟合原子发射光谱的基线，并能为光谱信号的定量分析提供更准确的信息，可以作为一种光谱基线拟合的有效手段。