CN109839362A

CN109839362A - 基于渐进式去噪技术的红外光谱定量分析方法

Info

Publication number: CN109839362A
Application number: CN201910211897.5A
Authority: CN
Inventors: 彭思龙; 段开红; 谢琼; 田瑞华; 许鹏程; 张秀娟; 雅梅; 刘彦宏; 李晓娟
Original assignee: Inner Mongolia Boda Quality Inspection Technology Co Ltd; Beijing Zhongke Zhicao Technology Co Ltd
Current assignee: Inner Mongolia Boda Quality Inspection Technology Co Ltd; Beijing Zhongke Zhicao Technology Co Ltd
Priority date: 2019-03-20
Filing date: 2019-03-20
Publication date: 2019-06-04
Anticipated expiration: 2039-03-20
Also published as: CN109839362B

Abstract

本发明属于分析化学技术领域，旨在解决现有的稀疏表示算法只进行一次去噪，无法完全去除噪声的问题。本发明提供了一种基于渐进式去噪技术的红外光谱定量分析方法，包括：构建目标函数；将增广拉格朗日函数引入目标函数；初始化稀疏表示系数、拉格朗日乘子和惩罚系数和迭代数；更新；利用检验函数进行检验，判断是否满足检验条件；如果满足检验条件，则停止迭代；输出稀疏表示系数；应用硬阈值函数对稀疏表示系数α进行去噪处理；应用保留的系数重建红外光谱数据，应用设定指标评价定量分析效果，如果定量分析效果满足设定标准，则停止迭代；输出定量分析方法模型。本发明能够以渐进式的方式来逐次消除噪声，提高定量分析模型的准确性和鲁棒性。

Description

基于渐进式去噪技术的红外光谱定量分析方法

技术领域

本发明属于分析化学技术领域，具体提供一种基于渐进式去噪技术的红外光谱定量分析方法。

背景技术

FTIR光谱分析法(即傅立叶变换红外光谱分析方法)是一种基于样本光谱，对化学品进行定性和定量分析的快速无损检测方法。光谱是由叠加在背景噪声上的一系列峰组成，有效的信息主要集中在峰值中，FTIR光谱中不可避免地存在噪声，噪声会掩盖光谱中的有用信息，从而降低定量分析模型的准确性和鲁棒性。因此，在线FTIR定量分析之前需要进行光谱去噪。

现有的去噪算法有：平均法，SavitzkyGolay平滑法。平均法需要采集多组光谱数据进行平均，耗时长，不满足在线分析的实时性；SavitzkyGolay平滑法随着窗口宽度的增大，光谱的峰值会变得平滑，导致一部分有效信息会丢失。稀疏表示法是目前常用的一种方法，通过稀疏解重构谱图从而达到去噪的目的。但是由于噪声的随机性，我们不能实现噪声的稀疏表示，现有的稀疏表示算法只进行一次去噪，并不能完全去除噪声。

因此，本领域需要一种基于渐进式去噪技术的红外光谱定量分析方法来解决上述问题。

发明内容

为了解决现有技术中的上述问题，即为了解决现有的稀疏表示算法只进行一次去噪，并不能完全去除噪声的问题，本发明提供了一种基于渐进式去噪技术的红外光谱定量分析方法，所述红外光谱定量分析方法包括：

S1：构建目标函数；

S2：将增广拉格朗日函数引入所述目标函数；

S3：初始化所述增广拉格朗日函数的稀疏表示系数α、拉格朗日乘子μ和惩罚系数λ以及迭代数k；

S4：更新α^(k+1)、λ^(k+1)和μ^(k+1)；

S5：利用检验函数进行检验，并判断是否满足检验条件；

S6：如果满足所述检验条件，则停止迭代；如果不满足所述检验条件，则使迭代数k＝k+1，并返回S4；

S7：输出稀疏表示系数α；

S8：应用硬阈值函数对稀疏表示系数α进行去噪处理；

S9:应用保留的系数重建红外光谱数据，应用设定指标评价定量分析效果，如果定量分析效果满足设定标准，则停止迭代；如果定量分析效果不满足所述设定标准，则返回S3；

S10:输出定量分析方法模型。

在上述红外光谱定量分析方法的优选技术方案中，步骤S1具体包括：

应用L1范数作为稀疏约束，将所述目标函数构建为：

s.t.x-Φα＝0

其中，x为测量光谱，Φ为字典，||α||₁为每个系数的绝对值之和。

在上述红外光谱定量分析方法的优选技术方案中，步骤S2具体包括：

将下列增广拉格朗日函数引入所述目标函数：

对α和μ求偏导数，具体为：

令得到：

在上述红外光谱定量分析方法的优选技术方案中，步骤S3具体包括：设定α＝0、μ＝0、λ＝0、k＝0。

在上述红外光谱定量分析方法的优选技术方案中，步骤S4具体包括：

将α^(k+1)更新为：

将λ^(k+1)更新为：

将μ^(k+1)更新为：

μ^(k+1)＝μ^(k)+(x-Φα^(k))

其中，x为测量光谱，Φ为字典。

在上述红外光谱定量分析方法的优选技术方案中，所述检验函数为：

其中，x为测量光谱，Φ为字典。

在上述红外光谱定量分析方法的优选技术方案中，所述检验条件为：

其中，∈为噪声级。

在上述红外光谱定量分析方法的优选技术方案中，步骤S8具体包括：

应用硬阈值函数将稀疏表示系数α中的参数的绝对值从小到大依次排列；

去除每次预定百分比p％内较小的系数，硬阈值T表示为：

T＝p％

在上述红外光谱定量分析方法的优选技术方案中，步骤S9具体包括：

应用保留的系数重建红外光谱数据；

将重建的红外光谱数据划分为训练集和测试集；

应用偏最小二乘法交叉验证，用均方根误差RMSECV指标评价定量分析效果；

如果RMSECV小于设定的阈值，则停止迭代；否则返回S3。

本领域技术人员能够理解的是，在本发明的优选技术方案中，通过构建目标函数并将增广拉格朗日函数引入目标函数，再通过参数初始化以及更新α^(k+1)、λ^(k+1)和μ^(k+1)，然后迭代直到满足检验条件，输出稀疏表示系数α，再应用硬阈值函数去除稀疏表示系数α中较小的系数并应用保留下来的系数重建红外光谱数据，然后对定量分析效果进行分析，迭代直到满足设定标准然后输出定量分析方法模型，通过这样的设置，即通过重构和迭代，也就是以渐进式的方式来逐次消除噪声，从而提高定量分析模型的准确性和鲁棒性。

附图说明

图1是不同γ值对于的GGD；

图2是软阈值函数曲线和硬阈值函数曲线；

图3是本发明的PSD-PLS算法训练过程的流程图；

图4是玉米数据集的原始光谱；

图5是验证参数RMSEP比对表；

图6是验证参数R2比对表；

具体实施方式

下面参照附图来描述本发明的优选实施方式。本领域技术人员应当理解的是，这些实施方式仅仅用于解释本发明的技术原理，并非旨在限制本发明的保护范围。

基于背景技术指出的现有的稀疏表示算法只进行一次去噪，并不能完全去除噪声的问题，本发明提供了一种基于渐进式去噪技术的红外光谱定量分析方法，旨在采用渐进式去噪，提高光谱数据的去噪效果，提高定量分析模型的准确性和鲁棒性。

具体地，每个测量光谱x是由纯谱s和噪声n组成的，因此可表示为：

x＝s+n (1)

对于稀疏表示法，如果字典用Φ表示，则纯谱s可表示为：

s＝Φα (2)

其中，α是稀疏表示系数。本发明中，分别用||α||₀，||α||₁和||α||₂分别表示α的非零元素的个数，系数的绝对值和和欧几里得距离。由于构造的字典是冗余的，因此α不是唯一的，选择用L1范数作为稀疏约束，该方法表示为：

s.t.x-Φα＝0

在本发明中，字典的构造可以采用广义高斯分布(GGD)，需要说明的是，广义高斯分布(GGD)在图像处理和信号处理中有着广泛的应用，其形状参数不同，分布的形状也不同，高斯分布和拉普拉斯分布是GGD的两种特殊情况。鉴于丰富的形状分布，使用GGD重建光谱信号。字典Φ由代表光谱的GGD线型的列构成。GGD的概率密度函数可以表示为：

其中，

Γ(.)是γ函数，公式(4)中，υ，σ，γ，β分别代表GGD的均值，标准差，形状系数和尺度系数。衰减率由形状参数γ决定，γ值越小，则衰减速度越快。如图1所示：

当γ＝2，β＝2，广义高斯分布可表示如下：

即高斯分布。

当γ＝1，广义高斯分布可表示如下：

即拉普拉斯分布。

在字典中，基本函数为高斯分布、拉普拉斯分布和其他GGD线性图。随着峰位置参数，线型图宽度参数β和形状参数γ变化，可以构造出不同的字典。

在本发明中，提供一种改进的ALM快速求解算法，公式(3)描述的优化问题的约束条件和目标函数，可以写成一个标准形式，引入拉格朗日函数，并表示为：

其中的μ和λ分别代表拉格朗日乘子和惩罚系数。为了得到最优解，分别用拉格朗日函数对α和μ求偏导数，具体为：

公式(10)中，sign()表示符号函数，令公式(10)等于0，可以得到：

拉格朗日乘子μ，可表示为：

μ^(k+1)＝μ^(k)+(x-Φα^(k)) (13)

其中，k为迭代数。除了更新μ，惩罚系数λ也需要更新，在该算法中，λ是自适应地更新，以达到快速收敛。该算法是收敛的。且α^(k)的收敛速度约为1/2。参数λ的自适应更新表示为：

在上述中，该算法可以描述如下：

S1：构建目标函数；

S2：将增广拉格朗日函数引入目标函数；

S3：初始化增广拉格朗日函数的稀疏表示系数α、拉格朗日乘子μ和惩罚系数λ以及迭代数k；

S4：更新α^(k+1)、λ^(k+1)和μ^(k+1)；

S5：利用检验函数进行检验，并判断是否满足检验条件；

S7：输出稀疏表示系数α；

具体地，上述的步骤S4具体包括：

将α^(k+1)更新为：

将λ^(k+1)更新为：

将μ^(k+1)更新为：

μ^(k+1)＝μ^(k)+(x-Φα^(k))

其中，检验函数为：

检验条件为：

其中，∈为噪声级。

也就是说，如果则停止迭代；否则令k＝k+1，返回并重复S4。

本发明可以应用PSD算法，由于噪声的随机性，因此不能实现噪声的稀疏表示，所以稀疏表示参数α需要从上述的ALM-NPR算法中得出，并用于降噪。硬阈值函数和软阈值函数在小波降噪处理中广泛使用，相应的曲线在图2中表示。

在软阈值函数中，经处理参数和原始参数间存在着恒定偏差，导致了原始光谱和重构光谱间也存在偏差，所以，提出了硬阈值函数，将参数的绝对值从小到大依次排列，然后去掉每次预定百分比p％内较小的系数，硬阈值T表示如下：

T＝p％ (15)

即在上述的步骤S7之后，本发明的红外光谱定量分析方法还包括：S8：应用硬阈值函数对稀疏表示系数α进行去噪处理。需要说明的是，在稀疏表示系数α中的所有系数中，除了几个较大的系数，其余的较小的系数都非常小，并且这些较小的系数的值远小于这几个较大的系数，因此，本算法中就是去掉这些明显很小的系数，从而利用这几个较大的系数来重构光谱。

本发明基于ALM-NPR算法和硬阈值函数，可实现渐进式光谱去噪。PSD算法的主要目的是舍弃在ALM-NPR算法渐进式去噪过程中产生的较小系数，保留的系数用来重构光谱。通过重构和迭代，能够逐次消除噪声。

结合PLS算法，可以对PSD算法进行评价。在上述的步骤S8之后，本发明的红外光谱定量分析方法还包括：S9：应用保留的系数重建红外光谱数据，应用设定指标评价定量分析效果，如果定量分析效果满足设定标准，则停止迭代；如果定量分析效果不满足所述设定标准，则返回S3；S10:输出定量分析方法模型。在一种可能的情形中，将实验数据集分为训练集和验证集，在训练过程中，每次舍弃p％的参数，其余的参数用来重构光谱，并用PLS算法进行评价。训练集的交叉验证均方根误差(RMSECV)可用来评价PSD-PLS算法，当RMSECV小于设定的阈值时，训练停止，可得到PLS模型和PSD算法的参数。本领域技术人员可以在实际应用中灵活地设置上述设定的阈值，只要通过该设定的阈值能够判定是否需要停止训练即可。然后，基于训练得到的PSD-PLS模型，可得到在线FTIR光谱的定量分析模型，详细过程可以参照图3。当然，除了采用上述的交叉验证均方根误差(RMSECV)来评价PSD-PLS算法，还可以采用其他的方式来评价PSD-PLS算法，例如：决定系数R²。

为了测试本发明的PSD-PLS算法，可以采用一组真实的数据进行测试，并将计算结果跟SavitzkyGolay平滑(SGS)和小波等去噪方法得到的结果进行对比。在预处理过程中，将校正后的光谱用于定量分析，用RMSEP(预测均方根误差，用于评价模型对外部样本的预测能力)和R²(决定系数，用于评价模型的精确度和回归效果)来评价PSD-PLS算法，算法用MATLAB 2014a版本编写，在3.20GHz的个人计算机上运行。

实验数据选用玉米数据集，该数据由FTIR光谱仪MP5测量得到，波长范围为1100-2498nm，间隔为2nm(700通道)，每个样品成分包括水分、油、蛋白质和淀粉值。样本光谱图如图4所示。样本数量一共有80个，将前80％划分为训练集，剩下的20％样本作为验证集，超参数用10折交叉验证方法得出。在图4中，可以看出，原始光谱是光滑的，噪声较小。为了展示PSD-PLS法在噪声FTIR光谱定量分析中的效果，在80个样本光谱中增加了高斯分布噪声，高斯噪声的平均值和标准差分别是0和光谱强度的0.1％，选择T＝1％作为硬阈值函数的设定百分比。

玉米的训练模型是基于水分、油、蛋白质和淀粉对光谱的响应。为了评价预处理方法，引入RMSEP和R²来计算验证集的准确性和鲁棒性。交叉验证完成后，得到SGS算法的最优参数(水分、油、蛋白质和淀粉的窗宽)。预处理完成后，校正光谱用于定量分析，可与未经预处理的验证结果(如图5和图6所示中NO一行)做对比，比较了原始光谱和叠加高斯噪声的光谱，结果在如图5所示的RMSEP表和如图6所示的R²表中展示，表中标*的代表叠加高斯噪声的光谱。

从图5和图6可以得出，对于玉米数据集的原始数据，PSD^-PLS跟SGS、小波等去噪法相比，计算得到的RMSEP和R²略有改进。当原始数据叠加了高斯噪声后，用PSD-PLS法计算得到的RMSEP和R²，明显优于其他方法。特别是玉米油成分这组数据，用PLS回归模型计算得到的R²是0.5998，因此不能用PLS模型进行预测。但是，用PSD-PLS算法，R²提高到0.8355。结果表明FTIR光谱不能一次去噪成功，渐进式光谱去噪在定量分析中是很有必要的。

通过上述实验结果表明，本发明的PSD-PLS方法成功地解决了在线FTIR光谱定量分析问题。基于PSD-PLS算法，FTIR光谱法可广泛地应用到在线定量分析中。

至此，已经结合附图所示的优选实施方式描述了本发明的技术方案，但是，本领域技术人员容易理解的是，本发明的保护范围显然不局限于这些具体实施方式。在不偏离本发明的原理的前提下，本领域技术人员可以对相关技术特征作出等同的更改或替换，这些更改或替换之后的技术方案都将落入本发明的保护范围之内。

Claims

1.一种基于渐进式去噪技术的红外光谱定量分析方法，其特征在于，所述红外光谱定量分析方法包括：

S1：构建目标函数；

S2：将增广拉格朗日函数引入所述目标函数；

S4：更新α^(k+1)、λ^(k+1)和μ^(k+1)；

S5：利用检验函数进行检验，并判断是否满足检验条件；

S7：输出稀疏表示系数α；

S8：应用硬阈值函数对稀疏表示系数α进行去噪处理；

S10:输出定量分析方法模型。

2.根据权利要求1所述的红外光谱定量分析方法，其特征在于，步骤S1具体包括：

应用L1范数作为稀疏约束，将所述目标函数构建为：

s.t.x-Φα＝0

3.根据权利要求1所述的红外光谱定量分析方法，其特征在于，步骤S2具体包括：

将下列增广拉格朗日函数引入所述目标函数：

对α和μ求偏导数，具体为：

令得到：

4.根据权利要求1所述的红外光谱定量分析方法，其特征在于，步骤S3具体包括：

设定α＝0、μ＝0、λ＝0、k＝0。

5.根据权利要求1所述的红外光谱定量分析方法，其特征在于，步骤S4具体包括：

将α^(k+1)更新为：

将λ^(k+1)更新为：

将μ^(k+1)更新为：

μ^(k+1)＝μ^(k)+(x-Φα^(k))

其中，x为测量光谱，Φ为字典。

6.根据权利要求1所述的红外光谱定量分析方法，其特征在于，

所述检验函数为：

其中，x为测量光谱，Φ为字典。

7.根据权利要求6所述的红外光谱定量分析方法，其特征在于，

所述检验条件为：

其中，∈为噪声级。

8.根据权利要求1所述的红外光谱定量分析方法，其特征在于，步骤S8具体包括：

去除每次预定百分比p％内较小的系数，硬阈值T表示为：

T＝p％。

9.根据权利要求1所述的红外光谱定量分析方法，其特征在于，步骤S9具体包括：

应用保留的系数重建红外光谱数据；

将重建的红外光谱数据划分为训练集和测试集；

应用偏最小二乘法进行交叉验证，用均方根误差RMSECV指标评价定量分析效果；

如果RMSECV小于设定的阈值，则停止迭代；否则返回S3。