CN110162740B - 一种用于光谱分辨增强的逆矩阵迭代反卷积方法 - Google Patents

Abstract

本发明公开了一种用于光谱分辨增强的逆矩阵迭代反卷积方法，包括如下步骤：步骤一、序列卷积和卷积方阵；步骤二、卷积方阵累积乘与卷积核函数峰展宽；步骤三、峰分辨增强。本发明通过原始函数与反卷积矩阵乘，可实现峰宽收窄目的，更进一步的，提出了反卷积单位矩阵构造方法，从而实现预期精度的反卷积计算效果，计算过程快速、可控，结果稳定、准确，适用范围广，可用于拉曼、红外等分子光谱，也可用于质谱、核磁共振、XRD、XRF等具有对称峰型的其他光谱的分辨增强，实施例以拉曼光谱为例，给出了实施步骤和实现效果。

Description

一种用于光谱分辨增强的逆矩阵迭代反卷积方法

技术领域

本发明涉及光谱分辨增强领域，特别涉及一种用于光谱分辨增强的逆矩阵迭代反卷积方法。

背景技术

光谱，例如拉曼、红外、XRD、XRF和NMR，都对分辨率有较高的要求。仪器价格随其分辨性能显著上升，因此，一直以来都期望能在信号处理上找到一种稳定、准确且广泛的分辨增强算法。

分辨率下降的本质，是由于仪器或者其他因素与真实信号卷积的结果，在仪器方面例如狭缝、光源或激发源的単色性、分光器和检测器性能等，其他因素包括温度、压力、布朗运动等随机展宽。表达在函数上，即观测结果是真实信号与多种卷积核函数的复合卷积的结果。理论上认为，只要找到合适的卷积核函数，常规反卷积方法，可提高分辨率。

目前，最稳定有效的反卷积方法是快速傅里叶(FFT)反卷积方法，以此基础，还发展出了各类方法，此类方法的缺点是需要定义准确的核函数，但找到具备普遍性的核函数难度大，直接影响了方法的广泛性和易用性。现实数据影响因素多，核函数复杂，实际上要提出普遍适用的统一的核函数非常困难，要解决各种仪器在不同场合的分辨增强问题，同时满足稳定、准确和广泛的要求也是非常困难的。

针对核函数寻找困难，发展起来的盲反卷积方法可以有效解决核函数寻找问题，可用于如图形图像分辨增强这类结果非严格一致的反卷积，但该方法同时又带来结果重复性问题，而难以满足需要严格一致、结果可重复的光谱要求。

由于常规的卷积核函数，如高斯函数、劳伦兹函数等都具有自身性质经过卷积后仍然保持，只是峰宽扩展的属性，而沃伊特函数是不同比例的高斯和劳伦兹函数的卷积，这些函数几乎已经包含了拉曼、红外、NMR等具有对称峰型的光谱的卷积核函数，针对光谱的核函数复杂主要是因为高斯函数、洛伦兹函数占比和展宽程度不同导致的。

因此，本发明提出了一种新的卷积和反卷积实现方法，用于光谱的分辨增强。通过构造一个包含卷积核函数微小单元的卷积方阵，代替卷积核函数直接定义，该卷积方阵及其逆分别对应了卷积计算和反卷积计算。通过构造的卷积方阵组合乘和迭代乘，可以方便地近似出与实际峰包含的各种卷积核函数，稳定、准确和广泛地应用于光谱分辨增强。

发明内容

为解决上述现有技术存在的问题，本发明的目的在于提供一种用于光谱分辨增强的逆矩阵迭代反卷积方法，用于光谱的分辨增强。通过构造一个包含卷积核函数微小单元的卷积方阵，代替卷积核函数直接定义，该卷积方阵及其逆分别对应了卷积计算和反卷积计算。通过构造的卷积方阵组合乘和迭代乘，可以方便地近似出与实际峰包含的各种卷积核函数，稳定、准确和广泛地应用于光谱分辨增强。

为达到上述目的，本发明的技术方案为：

一种用于光谱分辨增强的逆矩阵迭代反卷积方法，包括如下步骤：

步骤一、序列卷积和卷积方阵

对于序列f(n)和g(n)的卷积结果为：

如果f(n)是一个含有m个值的光谱序列，g(n)序列截断为2m-1 个元素，式(1)表达为：

改写成矩阵形式：

式(3)中，保留其中可计算部分，即：

式(4)中F的元素重新编号，有：

考察式(4a)，如果卷积核g是有限峰宽的对称函数，其中g(m) 是峰值，而峰值前后b个元素值趋于0；同时，f序列的最前和最后的b个元素值也趋于零；那么，式(4)可完成计算，计算结果与式 (2)相等。即得到结论，本发明提出的式(4)是一种通过构造对称卷积核方阵完成卷积计算的新方法。用矩阵乘的形式表达为：

F＝f·G (5)

进一步的，

f＝F·inv(G) (6)

式(6)的含义是，已知序列F与卷积方阵的逆相乘，得到去卷积结果f。

步骤二、卷积方阵累积乘与卷积核函数峰展宽

以高斯函数为例，对于

记：

G⁽¹⁾(y)＝Convolve^{[g(x)，g(x)，x，y]} (8)

有：

即，高斯函数自身卷积后仍为高斯函数，峰宽增加至卷积次数0 的

倍。

根据该结论，矩阵G连乘，记为：

只要单位G峰宽足够小，只需调整n，就可以逼近至需要精度的卷积和反卷积要求。而不同的卷积核函数，只需要替换G中的序列 g(n)即可。

步骤三、峰分辨增强

峰宽是影响分辨的主要因素，通过缩小峰宽(半高峰宽，FWHM) 达到分辨增强，使得相互重叠覆盖的峰得到辨识。式(9)表明正向卷积导致峰变宽，而逆向可实现峰收窄，根据式(6)和式(10)，可构造出相应的反卷积核函数单位矩阵，继而通过迭代计算得到需要精度的反卷积核函数矩阵。

进一步的，该方法的具体步骤为：

步骤一、生成单位卷积和反卷积矩阵

1)输入需要处理的光谱序列值峰f(n)，前后保留或增加若干个不含出峰的平坦数据点；

2)根据f(n)中元素数量m，以及f中峰的性质，确定反卷积所用的分布函数；

3)根据计算精度要求，确定分布函数峰宽(半峰宽)，可以在 0.1到1之间选择；

4)生成g(n)序列，其中元素数量为2m-1，峰值位置为第m个元素；将g(n)序列值放在第1行，依次向后平移，生成第2行、……，直至第m行，平移中缺失元素用0或者‘NaN’替代，得到尺寸为 m×(3m-2)的矩阵M；

5)从M中，切出第m列至2m-1列，得到一个m×m的方阵，该方阵即为单位卷积矩阵G；

6)对G求逆，即得到单位反卷积矩阵IG；

生成的单位卷积和反卷积矩阵是斜对角对称矩阵，针对不同元素数量光谱序列，不需要每次重复生成，首先生成一个元素数量大的单位阵，后续需要用到时，从大矩阵中截取相应尺寸的方阵；

步骤二、应用反卷积矩阵增强光谱分辨率

1)选择需要分辨增强的光谱峰；对于拉曼或红外光谱，由于不同波数上峰展宽情况不一致，为了保证分辨增强稳定、准确，从完整图谱上选择需要分辨增强的部分，同时尽可能保留一部分头尾基线；

2)选择反卷积核函数类型；根据峰型特点和宽度，判断展宽主要影响因素后选择；如果无法判断或者缺乏先验知识，也可以随机或者枚举试算，考察分辨率增强结果；

3)依据分峰精度的要求，定义单位卷积矩阵的序列峰宽，生成单位反卷积矩阵；

4)输入的原始峰F，按照f＝F·inv(G⁽ⁿ⁾)，或者f＝F·IG⁽ⁿ⁾，其中IG⁽ⁿ⁾与G⁽ⁿ⁾类似，即：

计算去卷积后的分辨增强峰f；

5)不断增加n进行迭代，峰宽降至重叠峰可满足识别要求时可停止计算，或者继续计算，信号的信噪比达到用户可容忍极限；

需要注意的，分辨率与信号的信噪比关系密切，分辨率提升实际是信号频率从低频向高频的转换，由于分辨率提高，原来隐含在信号中的低频噪声，也会同时表现为显著的高频噪声，导致信噪比恶化；为降低信噪比限制，应该选择信噪比好的信号进行分辨增强；本发明的去卷积分辨增强是渐进过程，因此过程中可以视情况，穿插平滑降噪，改善输出结果的信噪比。

进一步的，所述步骤一的步骤1)中，前后各大于10个数据点，以便于在最终结果中舍弃。

进一步的，所述步骤一的步骤2)中，所述分布函数具体为高斯分布、洛伦兹分布。

进一步的，所述步骤一的步骤6)中，具体为：

已生成元素数量为m×m的大矩阵G(m×m)，如果只需要p×p个元素(p＜m)，从G(m×m)中截取G(p×p)即可；G矩阵及其子集可以写成：

相对于现有技术，本发明的有益效果为：

本发明通过原始函数与反卷积矩阵乘，可实现峰宽收窄目的，更进一步的，提出了反卷积单位矩阵构造方法，从而实现预期精度的反卷积计算效果，计算过程快速、可控，结果稳定、准确，适用范围广，可用于拉曼、红外等分子光谱，也可用于质谱、核磁共振、XRD、XRF等具有对称峰型的其他光谱的分辨增强，实施例以拉曼光谱为例，给出了实施步骤和实现效果。

附图说明

图1为甲苯的拉曼光谱；

图2为900-1100cm^-1波数段甲苯拉曼光谱；

图3为应用洛伦兹核函数反卷积矩阵迭代光谱分辩增强效果图；

图4为应用高斯核函数反卷积矩阵迭代光谱分辩增强效果图；

图5为CCl₄的原始峰和分辨增强后的对比图；

图6为长焦距大型拉曼光谱仪的测量结果

具体实施方式

下面结合附图和具体实施方式对本发明技术方案做进一步详细描述：

实施步骤

1、生成单位卷积和反卷积矩阵

1)输入需要处理的光谱序列值峰f(n)，前后保留或增加若干个不含出峰的平坦数据，通常要求前后各大于10个数据点，以便于在最终结果中舍弃。

2)根据f(n)中元素数量m，以及f中峰的性质，确定反卷积所用的分布函数，例如高斯分布、洛伦兹分布等；

3)根据计算精度要求，确定分布函数峰宽(半峰宽)，可以在 0.1到1之间选择；

4)生成g(n)序列，其中元素数量为2m-1，峰值位置为第m个元素。将g(n)序列值放在第1行，依次向后平移，生成第2行、......，直至第m行，平移中缺失元素用0或者‘NaN’替代，得到尺寸为 m×(3m-2)的矩阵M。

5)从M中，切出第m列至2m-1列，得到一个m×m的方阵，该方阵即为单位卷积矩阵G。

6)对G求逆，即得到单位反卷积矩阵IG。

需要指出的，生成的单位卷积和反卷积矩阵是斜对角对称矩阵，也就是说，针对不同元素数量光谱序列，不需要每次重复生成，可以首先生成一个元素数量较大的单位阵，后续需要用到时，从大矩阵中截取相应尺寸的方阵即可。

例如，已生成元素数量为m×m的大矩阵G(m×m)，如果只需要p×p 个元素(p＜m)，从G(m×m)中截取G(p×p)即可。G矩阵及其子集可以写成：

2、应用反卷积矩阵增强光谱分辨率

1)选择需要分辨增强的光谱峰。对于拉曼或红外光谱，由于不同波数上峰展宽情况不一致，为了保证分辨增强稳定、准确，从完整图谱上选择需要分辨增强的部分，同时尽可能保留一部分头尾基线。

2)选择反卷积核函数类型。根据峰型特点和宽度，判断展宽主要影响因素后选择。如果无法判断或者缺乏先验知识，也可以随机或者枚举试算，考察分辨率增强结果。

3)依据分峰精度的要求，定义单位卷积矩阵的序列峰宽，生成单位反卷积矩阵。

4)输入的原始峰F，按照f＝F·inv(G⁽ⁿ⁾)，或者f＝F·IG⁽ⁿ⁾，其中IG⁽ⁿ⁾ 与G⁽ⁿ⁾ 类似，即：

计算去卷积后的分辨增强峰f。/>

5)不断增加n进行迭代，峰宽降至重叠峰可满足识别要求时可停止计算，或者继续计算，信号的信噪比达到用户可容忍极限。

需要注意的，分辨率与信号的信噪比关系密切，分辨率提升实际是信号频率从低频向高频的转换，由于分辨率提高，原来隐含在信号中的低频噪声，也会同时表现为显著的高频噪声，导致信噪比恶化。为降低信噪比限制，应该选择信噪比好的信号进行分辨增强；本发明的去卷积分辨增强是渐进过程，因此过程中可以视情况，穿插平滑降噪，改善输出结果的信噪比。

实施例

以甲苯的拉曼光谱作为实施例。

图1是甲苯的拉曼光谱。

1.选择900-1100cm^-1波数段甲苯拉曼光谱F，为反卷积分辨增强示例，如图2。

2.此段光谱的像素点数为121，即f(n)中元素数量m为121。可用洛伦兹函数反卷积将Voigt峰中的洛伦兹函数贡献消除；

3.根据计算精度要求，确定洛伦兹峰的单位峰宽0.1；

4.按照洛伦兹函数生成g(n)序列，其中元素数量为2m-1，峰值中心位为第m个元素。g(n)序列值放在第1行，依次向后平移，生成第2行、……，直至第m行，平移中缺失元素用0或者‘NaN’替代，得到尺寸为m×(3m-2)的矩阵M。

5.从M中，切出第m列至2m-1列，得到一个m×m的方阵，该方阵即为单位卷积矩阵G。

6.对G求逆，即得到单位反卷积矩阵IG。

7.输入的原始峰F，按照f＝F·inv(G⁽ⁿ⁾)，或者f＝F·IG⁽ⁿ⁾，其中IG⁽ⁿ⁾ 与G⁽ⁿ⁾ 类似，即：

计算去卷积后的分辨增强峰f。不断增加n进行迭代。

图3分别是原始峰、迭代110次、200次和400次的效果。999 波数处的原始峰和110次洛伦兹反卷积处理后的半峰宽(FWHM)分别是6.09和4.79，1027波数处的原始峰和处理后半峰宽(FWHM) 分别是6.51和5.13。

图3中，迭代110次时峰型未显示失真。以洛伦兹峰继续反卷积，峰型出现失真，表明在110迭代时，洛伦兹峰成分已经消除完毕。

8.继续选择反卷积核函数为高斯函数进行反卷积光谱分辩增强。高斯反卷积矩阵生成步骤同前3-7步，其中相关参数改为，g(n)改为高斯函数，单位峰宽0.5，迭代65000次，分辨增强结果如图4。图中两个主峰的FWHM分别为2.63和3.07，与原始峰相比缩减超过一倍，更重要的是原来重叠在999波数的主峰内的弱峰得到了清晰辨识。

本发明提出的分辨增强，峰的性质可以良好保留，可以显著增强重叠峰的识别能力。

为了证实方法的准确性和有效性，还对采用分辨率为6cm^-1普通拉曼光谱仪测量CCl₄光谱中重叠的多重峰做了辨识。采用的反卷积操作参数为：选择洛伦兹峰宽0.1，迭代100次，然后，选择高斯峰宽0.5，迭代35000次。

图5是CCl₄的原始峰和分辨增强后的对比，图6是长焦距大型拉曼光谱仪的测量结果。

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何不经过创造性劳动想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求书所限定的保护范围为准。

Claims (4)

1.一种用于光谱分辨增强的逆矩阵迭代反卷积方法，包括如下步骤：

步骤一、序列卷积和卷积方阵

对于序列f(n)和g(n)的卷积结果为：

如果f(n)是一个含有m个值的光谱序列，g(n)序列截断为2m-1个元素，式(1)表达为：

改写成矩阵形式：

式(3)中，保留其中可计算部分，即：

式(4)中F的元素重新编号，有：

考察式(4a)，如果卷积核g是有限峰宽的对称函数，其中g(m)是峰值，而峰值前后b个元素值趋于0；同时，f序列的最前和最后的b个元素值也趋于零；那么，式(4)可完成计算，计算结果与式(2)相等；即得到结论式(4)是一种通过构造对称卷积核方阵完成卷积计算的方法；用矩阵乘的形式表达为：

F＝f·G (5)

进一步的，

f＝F·inv(G) (6)

式(6)的含义是，已知序列F与卷积方阵的逆相乘，得到去卷积结果f；

步骤二、卷积方阵累积乘与卷积核函数峰展宽

以高斯函数为例，对于

记：

G⁽¹⁾(y)＝Convolve^{[g(x)，g(x)，x，y]} (8)

有：

即，高斯函数自身卷积后仍为高斯函数，峰宽增加至卷积次数n的

倍；

根据该结论，矩阵G连乘，记为：

只要单位G峰宽足够小，只需调整n，就可以逼近至需要精度的卷积和反卷积要求；而不同的卷积核函数，只需要替换G中的序列g(n)即可；

步骤三、峰分辨增强

峰宽是影响分辨的主要因素，通过缩小峰宽，达到分辨增强，使得相互重叠覆盖的峰得到辨识；式(9)表明正向卷积导致峰变宽，而逆向可实现峰收窄，根据式(6)和式(10)，可构造出相应的反卷积核函数单位矩阵，继而通过迭代计算得到需要精度的反卷积核函数矩阵；

该方法的具体步骤为：

步骤一、生成单位卷积和反卷积矩阵

1)输入需要处理的光谱序列值峰f(n)，前后保留或增加若干个不含出峰的平坦数据点；

2)根据f(n)中元素数量m，以及f中峰的性质，确定反卷积所用的分布函数；

3)根据计算精度要求，确定分布函数峰宽或半峰宽，在0.1到1之间选择；

4)生成g(n)序列，其中元素数量为2m-1，峰值位置为第m个元素；将g(n)序列值放在第1行，依次向后平移，生成第2行、……，直至第m行，平移中缺失元素用0或者‘NaN’替代，得到尺寸为m×(3m-2)的矩阵M；

5)从M中，切出第m列至2m-1列，得到一个m×m的方阵，该方阵即为单位卷积矩阵G；

6)对G求逆，即得到单位反卷积矩阵IG；

生成的单位卷积和反卷积矩阵是斜对角对称矩阵，针对不同元素数量光谱序列，不需要每次重复生成，首先生成一个元素数量大的单位阵，后续需要用到时，从大矩阵中截取相应尺寸的方阵；

步骤二、应用反卷积矩阵增强光谱分辨率

1)选择需要分辨增强的光谱峰；对于拉曼或红外光谱，由于不同波数上峰展宽情况不一致，为了保证分辨增强稳定、准确，从完整图谱上选择需要分辨增强的部分，同时尽可能保留一部分头尾基线；

2)选择反卷积核函数类型；根据峰型特点和宽度，判断展宽主要影响因素后选择；如果无法判断或者缺乏先验知识，则选用随机或者枚举试算，考察分辨率增强结果；

3)依据分峰精度的要求，定义单位卷积矩阵的序列峰宽，生成单位反卷积矩阵；

4)输入的原始峰F，按照f＝F·inv(G⁽ⁿ⁾)，或者f＝F·IG⁽ⁿ⁾，其中IG⁽ⁿ⁾与G⁽ⁿ⁾类似，即：

计算去卷积后的分辨增强峰f；

5)不断增加n进行迭代，峰宽降至重叠峰可满足识别要求时可停止计算，或者继续计算，信号的信噪比达到用户可容忍极限；

需要注意的，分辨率与信号的信噪比关系密切，分辨率提升实际是信号频率从低频向高频的转换，由于分辨率提高，原来隐含在信号中的低频噪声，也会同时表现为显著的高频噪声，导致信噪比恶化；为降低信噪比限制，应该选择信噪比好的信号进行分辨增强；去卷积分辨增强是渐进过程，因此过程中视情况，穿插平滑降噪，改善输出结果的信噪比。

2.根据权利要求1所述的方法，其特征在于，所述步骤一的步骤1)中，前后各大于10个数据点，以便于在最终结果中舍弃。

3.根据权利要求1所述的方法，其特征在于，所述步骤一的步骤2)中，所述分布函数具体为高斯分布、洛伦兹分布。

4.根据权利要求1所述的方法，其特征在于，所述步骤一的步骤6)中，具体为：

已生成元素数量为m×m的大矩阵G(m×m)，如果只需要p×p个元素，其中p<m，从G(m×m)中截取G(p×p)即可；G矩阵及其子集写成：

/>

Images