基于背景扣除提取峰面积的拉曼光谱定量分析方法及系统

本发明属于光谱检测及光谱分析领域，具体涉及一种基于背景扣除提取峰面积的拉曼光谱定量分析方法及系统。

背景技术：

拉曼光谱检测技术是光谱检测技术的一种，可用于物质的浓度监测，因其具有光源相对简单、对物质识别能力好、样品不需要复杂处理等特点而获得广泛应用。现有便携式拉曼光谱设备可直接用于无任何处理的样品，使得其广泛应用在探矿、毒品检测、食品检测等多种领域。

便携式拉曼光谱监测系统一般需设置相应的数据分析算法才可检测出物质的信息。当前的拉曼光谱定量分析算法一般包括去背景和定量分析两步：第一步采用特定的自动化算法，将拉曼光谱中的荧光、杂散光等成分祛除，获得比较纯粹的拉曼光谱信号；第二步对拉曼信号进行定量分析，具体为进行谱峰识别，提取峰面积，然后与光谱数据库进行对比，找出相似物质，并进行多元线性拟合，此外，也可以用机器学习算法。

但是，现有背景扣除的主要判断标准是拉曼光谱的谱峰形状，相比于荧光及其他成分，拉曼峰光谱成分较窄，因而一般的背景扣除工作是根据光谱形状分辨，扣除光谱中较宽、变化平缓的部分，保留尖锐、变化剧烈的部分，如多项式拟合、自适应迭代重加权惩罚最小二乘法等方法。但是，根据光谱形状不同、拟合背景方法不同，最终得到的背景曲线也会不同，这些不同之处会影响定量分析的准确性，导致检测结果产生较大的偏差，进而导致拉曼光谱定量分析系统的精度下降。因此，现有拉曼光谱定量分析算法第二步的结果严重依赖第一步的结果，若第一步产生较大的误差，第二步会产生有更大的偏差，使得检测结果产生较大的偏差。

技术实现要素：

本发明的目的是解决现有拉曼光谱分析方法存在结果差异较大，精度较低的问题，提供一种基于背景扣除提取峰面积的拉曼光谱定量分析方法及系统，该方法和系统可获取更高的拉曼光谱定量分析精度。

为实现上述目的，本发明通过以下技术方案来实现；

一种基于背景扣除提取峰面积的拉曼光谱定量分析方法，包括以下步骤：

步骤一、将样品分为两组，分别为校准样品和待测样品，所述校准样品和待测样品均包括目标物质和参考物质；

步骤二、采集经过校准样品和待测样品后的拉曼光信号，并将拉曼光信号处理为数字信号，所述数字信息包括波长信息和光谱强度信息；

步骤三、对校准样品的目标物质信号进行去奇点噪声处理；

3.1)读取校准样品中目标物质的光谱强度信息y1(i,j)，y1(i,j)为数组，i,j为数组元素序号，i范围从1到n，n为单个光谱强度信息的长度；j范围从1到m1，m1为校准样品采集的光谱强度信息的个数；

3.2)设定奇点阈值k；

3.3)对步骤3.1)的光谱强度信息去奇点噪声处理，得到处理后的光谱强度信息y2(i,j)；

y2(1,j)＝y1(1,j)j＝1、2、……m1；

y2(n,j)＝y1(n,j)j＝1、2、……m1；

i＝2、3、.....、n-1，j＝1、2、……m1

步骤四、设置峰面积函数；

4.1)查询校准样品的目标物质的拉曼峰信息，确定拉曼峰的峰尖a和两侧峰谷b、c的位置；

4.2)根据波长信息x(i)、峰尖a和两侧峰谷b、c的位置，获取峰尖的数组元素序号ai和峰谷的数组元素序号bi、ci；其中，ai为|x(i)-a|最小时，i的值；bi为|x(i)-b|最小时，i的值；ci为|x(i)-c|最小时，i的值；

4.3)根据光谱强度信息y2(i,j)、峰尖的数组元素序号ai和峰谷的数组元素序号bi、ci，获取目标物质的峰尖信号ay、峰谷信号by和峰谷信号cy；其中，ay＝y2(ai，j)，by＝y2(bi，j)，cy＝y2(ci，j)；j＝1……m1；

4.4)根据波长信息x(i)、峰尖的数组元素序号ai和峰谷的数组元素序号bi、ci，获取每条光谱在峰尖a的拉曼位移信息ax和在峰谷b、c的拉曼位移信息bx、cx；其中，ax＝x(ai)，bx＝x(bi)，cx＝x(ci)；

4.5)获取目标物质的峰面积函数s(k)；

s(k)＝((by-ay)×k+(by-cy)×(1-k))×(cx-ax)/2

其中，k为校准样品中目标物质峰面积动态参数；

步骤五、重复步骤三和步骤四的数据处理过程，对校准样品的参考物质信号进行相同的处理，获取参考物质的峰面积函数s1(k1)；其中，k1是校准样品中参考物质拉曼峰峰面积的动态参数；

步骤六、获取扫描动态参数k、k1；

6.1)获取校准样品的标签浓度x(j)和目标物质的相对光谱强度y(j)；

y(j)＝s(k)/(s1(k1))；j＝1……m1；

6.2)对校准样品的标签浓度x(j)和目标物质的相对光谱强度y(j)进行线性拟合，获得拟合的相关系数r²；

6.3)不断改变k和k1，记录r²最大时的k和k1值，以及线性拟合所得的斜率c1和截距c2；

步骤七、对待测样品的拉曼光信号进行处理；

重复步骤三至步骤五的数据处理过程，对待测样品的拉曼光信号进行处理，得到待测样品的峰面积函数s’(k’)和s1’(k1’)；

其中，s’(k’)为待测样品中目标物质的峰面积函数；k’为待测样品中目标物质的拉曼峰峰面积的动态参数；s1’(k1’)为待测样品中参考物质的峰面积函数；k1’为待测样品中参考物质的拉曼峰峰面积的动态参数；

步骤八、计算待测样品的相对光谱强度y’(j)；

y’(j)＝s’(k’)/s1’(k1’)j＝1……m2；

其中，m2为待测样品采集的光谱强度信息的个数；

步骤九、将步骤六得到的k和k1值代入待测样品的相对光谱强度y’(j)，得到待测样品的相对光谱强度y’(j)＝s’(k)/s1’(k1)；

步骤十、根据步骤六得到的斜率c1和截距c2，计算待测样品浓度x’(j)＝(y’(j)-c2)/c1，并输出待测样品浓度x’(j)。

进一步地，步骤3.2)中，k为正整数，且不大于n/2。

同时，本发明还提供一种基于背景扣除提取峰面积的拉曼光谱定量分析系统，包括激光光源、拉曼探头、光栅光谱仪和分析模块；激光光源产生激光，激光进入拉曼探头；拉曼探头将激光作用于样品，并收集经过样品的拉曼光，随后将拉曼光传输给光栅光谱仪；光栅光谱仪检测拉曼探头传入的拉曼光，将拉曼光信号处理为数字信号；分析模块对光栅光谱仪输出的信号进行处理，获取样品的数据；所述分析模块被处理器执行时实现上述方法的骤三至步骤十。

进一步地，所述激光光源为单波长窄线宽激光器。

与现有技术相比，本发明方法具有以下技术效果：

1.本发明方法和系统在拉曼峰峰面积提取中，使用了可优化参数，然后根据参考样本的结果，拟合获得优化参数的值，相当于根据参考样本，对每个拉曼峰的峰面积计算方法进行了优化，该方式比统一的拉曼峰峰面积计算方法更精细，获得的结果也就具有更高进度。

2.本发明可以在不改动普通拉曼光谱定量分析系统硬件的基础上，可大幅度提高定量分析精度。

附图说明

图1为本发明基于背景扣除提取峰面积的拉曼光谱定量分析系统示意图；

图2为本发明实施例中设置峰面积函数的示意图；

图3为本发明实施例中得到的光谱图；

图4为本发明实施例中目标物质葡萄糖的拉曼峰的峰尖和左右峰谷的位置示意图；

图5为本发明实施例中参考物质水的拉曼峰峰尖a和左右峰谷b、c的位置示意图图；

图6为本发明方法得到的拟合结果示意图；

图7为现有方法得到的拟合结果示意图。

具体实施方式

以下结合附图和具体实施例对本发明的内容作进一步详细描述。

本发明提供一种基于背景扣除提取峰面积的拉曼光谱定量分析方法及系统，从而提高了拉曼光谱定量分析的精度。该系统和方法中，首先获取拉曼光谱的峰面积，根据峰的左右边界点设置动态参数；然后扫描动态参数，根据标定物质浓度进行线性拟合，得到最佳的动态参数；然后根据最佳动态参数，对未标定结果进行分析，获得定量分析结果。经过实验验证，在最佳动态参数下，本系统和方法得到的拟合精度，优于传统拉曼光谱定量分析方法。

如图1所示，本发明基于背景扣除提取峰面积的拉曼光谱定量分析系统主要由四个部分组成，分别是激光光源、拉曼探头、光栅光谱仪和分析模块。激光光源产生拉曼泵浦光，传输到拉曼探头。拉曼探头接收激光光源产生的泵浦光，照射到样品上，并将样品产生的拉曼光收集起来，传输到光谱仪。光栅光谱仪接收拉曼探头传输而来的拉曼光，将其转化为数字化的光谱信息，并传输给实验室计算机。分析模块接收光栅光谱仪传输过来的光谱信息，经过定量分析算法的处理，给出定量分析结果，并输出。

本发明基于背景扣除提取峰面积的拉曼光谱定量分析方法主要包括去奇点噪声、设置峰面积函数、扫描动态参数、输出结果等步骤，其具体过程如下：

步骤一、将样品分为两组，分别为校准样品和待测样品，校准样品和待测样品均包括目标物质和参考物质，目标物质为样品中待分析的成分，参考物为样品中选定的成分保持不变的物质；

步骤二、采集经过校准样品和待测样品后的拉曼光信号，并将拉曼光信号处理为数字信号，该数字信息包括波长信息和光谱强度信息；

步骤三、对步骤二获取的校准样品中目标物质的拉曼光信号进行去奇点噪声处理；

3.1)读取校准样品中目标物质的光谱强度信息y1(i,j)，y1(i,j)为数组，i,j为数组元素序号，i范围从1到n，n为单个光谱强度信息的长度；j范围从1到m1，m1为校准样品采集的光谱强度信息的个数；

3.2)设定奇点阈值k，k为正整数，且不大于n/2；

3.3)对步骤3.1)的光谱强度信息去奇点噪声处理，得到处理后的光谱强度信息y2(i,j)；

y2(1,j)＝y1(1,j)j＝1、2、……m1；

y2(n,j)＝y1(n,j)j＝1、2、……m1；

i＝2、3、.....、n-1，j＝1、2、……m1

步骤四、设置峰面积函数；

4.1)查询校准样品中目标物质的拉曼峰信息，确定拉曼峰的峰尖a和两侧峰谷b、c的位置；

4.2)根据波长信息x(i)、峰尖a和两侧峰谷b、c的位置，获取峰尖的数组元素序号ai和峰谷的数组元素序号bi、ci；其中，ai为|x(i)-a|最小时，i的值；bi为|x(i)-b|最小时，i的值；ci为|x(i)-c|最小时，i的值；

4.3)根据光谱强度信息y2(i,j)、峰尖的数组元素序号ai和峰谷的数组元素序号bi、ci，获取校准样品的峰尖信号ay、峰谷信号by和峰谷信号cy；其中，ay＝y2(ai，j)，by＝y2(bi，j)，cy＝y2(ci，j)；j＝1……m1；

4.4)根据波长信息x(i)、峰尖的数组元素序号ai和峰谷的数组元素序号bi、ci，获取每条光谱在峰尖a的拉曼位移信息ax，两侧峰谷b、c的拉曼位移信息bx、cx；其中，ax＝x(ai)，bx＝x(bi)，cx＝x(ci)；

4.4)获取校准样品中目标物质峰面积函数s(k)；

s(k)＝((by-ay)×k+(by-cy)×(1-k))×(cx-ax)/2；

其中，k为校准样品中目标物质峰面积动态参数；

步骤五，参照步骤三和步骤四，对校准样品中的参考物质信号进行处理，获取校准样品中参考物质的峰面积表示s1(k1)；其中，k1是校准样品中参考物质拉曼峰峰面积的动态参数；

步骤六、获取扫描动态参数k、k1；

扫描动态参数是根据校准数据获取最优的峰面积参数；

6.1)获取校准样品的标签浓度x(j)，获取目标物质的相对光谱强度y(j)；

y(j)＝s(k)/(s1(k1))j＝1……m1；

6.2)对校准样品的标签浓度x(j)和目标物质的相对光谱强度y(j)进行线性拟合，获得拟合的相关系数r²；

6.3)不断改变k和k1，记录r²最大时的k和k1值，及线性拟合所得的斜率c1和截距c2；

步骤七、对待测样品的拉曼光信号进行处理；

参照步骤三至步骤五，对待测样品的拉曼光信号进行处理，得到待测样品的峰面积函数s’(k’)和s1’(k1’)；

其中，s’(k’)为待测样品中目标物质的峰面积函数；k’为待测样品中目标物质的拉曼峰峰面积的动态参数；s1’(k1’)为待测样品中参考物质的峰面积函数；k1’为待测样品中参考物质的拉曼峰峰面积的动态参数；

步骤八、计算待测样品的相对光谱强度y’(j)；

y’(j)＝s’(k’)/s1’(k1’)j＝1……m2；

其中，m2为待测样品采集的光谱强度信息的个数；

步骤九、将步骤六得到的k和k1值代入待测样品的相对光谱强度y’(j)，得到待测样品的相对光谱强度y’(j)＝s’(k)/s1’(k1)；

步骤十、根据步骤六得到的斜率c1和截距c2，计算待测样品浓度x’(j)＝(y’(j)-c2)/c1，并输出待测样品浓度x’(j)。

下面以具体的实施例对本发明方法和系统进行描述和效果验证。

本发明使用的拉曼光谱采集系统为光栅式拉曼光谱采集系统，主要由激光光源、拉曼探头、光栅光谱仪、计算机组成。其中，激光光源要求为单波长窄线宽激光器，光栅光谱仪要求具有比较优秀的光谱线宽分辨率。激光器产生激光，传入拉曼探头。样品置于拉曼探头中，与激光相互作用，产生拉曼光。拉曼探头将产生的拉曼光收集起来，传输给光谱仪。光谱仪将传入的拉曼光进行分析，获得数字化的拉曼光谱信息(包含拉曼位移信息和每个拉曼位移对应的拉曼光谱强度)。这些信息传入实验室计算机进行处理，经过特定的算法，最终获得定量分析结果输出。

本实例中，使用的激光器是785nm的窄线宽连续输出半导体激光器，输出功率为200mw。使用的光谱仪中光栅为全息衍射光栅，探测器为日本滨松的ccd。光谱仪的光谱探测范围为光谱分辨能力为3cm^-1，探测波长为795.16-919.5nm，对应拉曼位移为162.8-1863.6cm^-1。

本实例中，对制成的样品分组，将其中一部分作为校准样本，另一部分作为待测样本。使用的校准样品和待测样品为葡萄糖明胶牛血清混合物。先将葡萄糖溶于牛血清，配置成不同浓度的溶液，然后和明胶混合并加热，最终形成葡萄糖明胶牛血清混合物，作为样品。实验室计算机进行处理时，需要输入校准样品的待测物浓度，进行特定的算法处理，才能给出待测样品的结果。

步骤一、将样品分为两组，分别为校准样品和待测样品；

步骤二、获取所有样品的拉曼光谱后，将其按照校准样品和待测样品区分，并将校准样品运用算法进行以下步骤的处理，获得目标物质(葡萄糖)和标定物质(水)的拉曼峰动态参数；

步骤三、对步骤二获取的校准样品中目标物质的拉曼光信号进行去奇点噪声处理；

3.1)读取校准组的光谱强度数据y(i)，i＝1—1024；

3.2)设定奇点阈值k＝3；

3.3)对步骤3.1)的光谱强度信息去奇点噪声处理，得到处理后的光谱强度信息y2(i,j)；

步骤四、设置峰面积函数；

设置峰面积函数为本方法的核心部分，其原理如下：

在图2中，拉曼峰的峰尖为b，两侧峰谷为a和c，按照多项式拟合算法，此拉曼光谱的背景，是a、c之间的一条柔和曲线；按照自适应迭代重加权惩罚最小二乘法，此拉曼光谱的背景，是a、c之间的一条直线。因此算法的选择，以及具体计算时的细节，都会导致背景曲线的变化；而真实的拉曼光谱背景曲线，应该是a、c之间的一条曲线，但是这条曲线的形状，是未知的，因此，本方法设置一个参数，用来表示由于曲线形状导致的背景变化；

a、b、c之间的面积，可以用三角形近似表示，三角形的面积，为底×高/2；底＝cx-ax，而高，因为a、c之间的背景曲线形状，会有所不同，不过大致应该在(by-ay)和(by-cy)之间。因此，设置比例参数k，将高表示为：(by-ay)*k+(by-cy)*(1-k)；

最终面积表示为：s(k)＝((by-ay)*k+(by-cy)*(1-k))*(cx-ax)/2

步骤五、参照步骤三和步骤四，对校准样品中的参考物质信号进行处理，获取校准样品中参考物质的峰面积表示s1(k1)；其中，k1是校准样品中参考物质拉曼峰峰面积的动态参数；

步骤六、获取扫描动态参数k、k1；

扫描动态参数是根据校准数据获取最优的峰面积参数。

本实例中，用葡萄糖拉曼峰1128cm^-1代表目标物质，用水峰1650cm^-1作为参考信号，采集获得的样品拉曼光谱如图3所示：

在1128cm^-1位置和1650cm^-1位置，分别获得其对应的a、b、c，如图4和图5所示：

以校准葡萄糖浓度为x，以s1128(k)/s1650(k1)为y，进行线性拟合；

扫描k和k1，获得最优的线性拟合结果，记录k和k1和线性拟合的斜率c1和截距c2；

步骤七、对待测样品的拉曼光信号进行处理；

参照步骤三至步骤五，对待测样品的拉曼光信号进行处理，得到待测样品的峰面积函数s’(k’)和s1’(k1’)；其中，s’(k’)为待测样品中目标物质的峰面积函数；k’为待测样品中目标物质的拉曼峰峰面积的动态参数；s1’(k1’)为待测样品中参考物质的峰面积函数；k1’为待测样品中参考物质的拉曼峰峰面积的动态参数；

步骤八、计算待测样品的相对光谱强度y’(j)；

y’(j)＝s’(k’)/s1’(k1’)；j＝1……m2；

其中，m2为待测样品采集的光谱强度信息的个数；

步骤九、根据k和k1，计算对应的s1128(k)和s1650(k1)，并计算对应的y；

步骤十、根据步骤六得到的斜率c1和截距c2，计算待测样品浓度x’(j)＝(y’(j)-c2)/c1，并输出待测样品浓度x’(j)。

作为对照组，使用自适应迭代重加权惩罚最小二乘法去背景，然后获取1128cm^-1峰面积，并用1650cm^-1峰面积作归一化。然后使用偏最小二乘法进行拟合，结果如图7所示。

从图6和图7的结果中可以看出，采用本发明的拟合结果r²为0.9997，使用传统方法的r²为0.97，本方法比传统方法优秀。同时，从图中可以看出，本方法所获得的拉曼光谱强度，更加贴近拟合曲线。