分数阶微分Whittaker平滑器

单 鹏, 何 年, 李志刚, 吴 缀

(东北大学秦皇岛分校 控制工程学院, 河北 秦皇岛 066000)

现代光谱分析技术是一种通过分析物质的光谱来实现对物质类别及化学成分鉴定的方法,是一项结合了光谱检测、数据处理、模式识别和计算机等多种理论和方法的综合性分析技术,被广泛运用于生产生活中[1-3].因受到传感器、背景噪声、电路噪声以及随机误差等因素的影响,噪声是光谱信号探测时无法避免的现象.光滑去噪成为光谱分析中一个重要的预处理步骤,能够提高后续光谱分析算法的准确性和可靠性.光谱去噪算法得到了人们的广泛重视,产生了多种从不同角度实现光滑去噪的算法,从严格意义上可分成几类.

基于噪声与光谱信号的统计特性差异,可采用平滑滤波来消除或减少噪声,例如均值滤波、中值滤波、Kalman滤波、移动多项式平滑滤波等;基于噪声与光谱信号的频率特征差异,可采用傅里叶变换、小波变换[4]和经验模态分解[5](empirical mode decomposition,EMD)等;基于噪声与信号波形之间的差异,可采用形态滤波等.采用平滑滤波的算法需要人为设定一些参数(如移动窗口的宽度,拟合阶次),参数直接影响着算法的性能.具有代表性的算法有Savitzky-Golay(SG)多项式最小二乘法[6],SG对于一些较大的数据,计算量变大,不能快速计算出平滑序列;需要人为设定一些参数,且不能自动选择;使用宽度固定的平滑窗口,当平滑窗口宽度越大,平滑效果越好,但是光谱失真会越严重,即会导致峰畸变或者丢失.

相较于SG算法,Whittaker Smoother(WS)算法[7]计算速度快、实现简单、可以通过快速计算留一交叉验证实现自动选参.WS算法采用整数阶微分(2或3阶)对真实信号的粗糙度进行约束,并通过调节正则化参数λ来平衡粗糙度和保真度两个冲突的目标,最终找到一个适当的平滑序列.用整数阶微分表示粗糙度简单有效,但是这种表示是片面的,所包含的信息量极为有限,不能很好地描述问题.即整数阶微分表示不依赖于历史数据,缺少全局相关性,对于一些复杂连续过程整数阶往往不能很好地对其进行准确描述.

近些年来,分数阶微分研究成为热点,被广泛地应用于各个领域[8-12].在光谱分析领域,Toman等[13]提出用半阶导数对谱峰进行分辨;卢小泉等[14]采用分数阶导数处理噪声信号;Kharintsev等[15]提出一种利用分数阶导数来提取光谱参数的方法;Mocak等[16]利用分数阶微分改善信号的分辨率;李远禄等[17]提出一种基于分数阶微分的重叠峰分辨的方法;Li等[18]提出用分数阶导数结合Haar小波来分辨重叠的洛伦兹峰;徐继刚等[19]利用分数阶导数处理油气红外光谱数据用来提高多元校正模型分析精度;张东等[20]和田安红等[21]将分数阶微分应用在盐渍土高光谱数据预处理中,表明分数阶微分对干旱荒漠区土壤光谱数据进行潜在信息挖掘是有效的;本文通过分数阶微分来增强对真实光谱信号的粗糙度表示,进一步提高平滑效果.本文将RL和GL分数阶微分定义引入WS算法中,解决了在WS算法实际应用中存在的相关问题.实验采用偏最小二乘法对Marzipan在多种不同仪器下的光谱进行定量分析,预测糖分和水分.实验结果表明分数阶WS算法比整数阶更加灵活有效,光谱定量分析精度更高.

1 算法原理

1.1 Whittaker Smoother算法原理

WS算法:对于一列长度为m,等间距的噪声序列为y,通过WS得到平滑序列为z.对于这两个序列,有以下要求:①对原始数据的保真度.通常用函数s=∑i(yi-zi)2表示.②序列z的粗糙度.这里用z的一阶微分来表示，Δzi=zi-zi-1,对其进行平方求和,表示z的粗糙度:R=∑i(Δzi)2.然后,平衡两个相互冲突的目标,这两个目标的平衡组合记为Q=S+λR,λ是一个用户选择的数.惩罚最小二乘的思想是找出使Q最小的一个z的序列.记为

(1)

为了避免大量繁琐的代数运算,这里引入矩阵和向量,目标函数表示为

Q=|y-z|2+λ|Ddz|2.

(2)

其中:y,z为列向量;Dd为(m-1)×m的1阶微分矩阵(用户也可以选其他整数阶阶次).对式(2)求导,令其等于零,可得序列z的表达式:

(3)

其中E为单位阵.

其中，Sλ满足SλI=I,I为单位向量.

(4)

其中:yi,zi分别为序列y,z的第i个点;Sii(λ)为平滑矩阵第i行第i列对应的元素.即留一交叉验证误差Scv代数和矩阵形式分别表示为[22]

(5)

(6)

其中，g为平滑矩阵Sλ对角线上元素构成的对角矩阵.

1.2 分数阶微分定义

目前常用的主要有三种分数阶微分定义:Riemann-Liouvile,Grumwald-Letnikov,Caputo.本文选用前两种定义,同样为避免大量繁杂的代数运算,选择矩阵形式表示分数阶微分.

1.2.1 Riemann-Liouvile定义

Kilicman等[23]给出Riemann-Liouvile微分运算矩阵:

(7)

其中:Γ(α)为Gamma函数;

ξ1=1,ξp=pα+1-2(p-1)α+1+(p-2)α+1,

(p=2,3,…,m-i+1).

1.2.2 Grumwald-Letnikov定义

Grumwald-Letnikov分数阶微分定义为

(8)

式中:α为阶数;h为微分步长;t与a分别为微分的上、下限;Γ(x)为Gamma函数.设一元函数f(x)定义域x∈[a,t]，令h=1,n=[t-a]，则由定义可推导出函数f(x)分数阶微分的差值表达式:

(9)

1.3 Fractional Differential Whittaker Smoother(FDWS)

1.3.1 算法原理

序列z的粗糙度R用分数阶微分来代替,根据上文给出的分数阶微分矩阵可以写出各点的微分r=Dαz，粗糙度为各点微分的平方和.目标函数Q=|y-z|2+λ|Dαz|2,序列z的表达式为

记RL分数阶微分实现为FDWS_RL,GL分数阶微分实现为FDWS_GL.

1.3.2 交叉验证

在一定的区间内,通过交叉验证误差(选择相应的参数:正则化系数λ、微分阶次α)可以选出一条较为理想的光滑信号,而且在一定程度上还可以进行人工调参,用户可根据需求选择.

2 实 验

2.1 实验条件

实验在Windows 10操作系统下实现;主要配置:Core I5 CPU主频2.50 GHz,8 GB RAM;编程在Matlab,2019a版下实现.

2.2 实验数据

1) NMR数据:该数据集是斯坦福大学开发的小波工具包(Wavelab toolbox)里所采用的经典数据集(见http://statweb.stanford.edu/～wavelab/),Eilers[7]曾对其进行光滑处理,验证WS算法有效性.本文验证对尖锐峰信号的提取是否有效.

2) Marzipan数据集:该数据集由7种光谱仪测量的32个Marzipan样品的近红外(NIR)或红外(IR)光谱组成,用于水分和糖分浓度的预测(见http://www.models.life.ku.dk/Marzipan).本文选用5种光谱仪器的光谱数据进行光滑处理,验证分数阶WS算法对后续光谱分析精度的影响.

2.3 FDWS_RL,FDWS_GL对比实验

对光谱数据进行光滑处理的目的在于去除噪声的同时能够尽可能多地保留真实信息.整数阶WS算法在处理一些含有尖窄峰的信号时(如NMR光谱),峰值信号丢失过多(见图2放大图).与之相反,分数阶在保证同样去噪效果下,可以通过调节阶次,来实现对真实尖峰信息的保留.通过图2观察得到,在一定的阶次范围内,峰值信号的保留与微分矩阵的阶次成正比,阶次越高,保留的有用信息越多.此外,RL,GL两种计算方式实验结果差别不大.最后固定正则化系数(λ=10),分别计算各阶次交叉验证均方根误差,结果如表1所示.观察发现,均方根误差与阶次大致呈反比,实验结果再次表明,在一定范围内,随着阶次越高,保留的真实信息越多.

图1核磁共振谱的RL平滑

图2 核磁共振谱的不同阶次的平滑

表1 交叉验证均方根误差

数据预处理主要是为了后续的光谱定量分析服务,为了进一步验证分数阶微分对后续光谱定量分析精度的影响,本文对Marzipan数据集进行实验.根据KS[24]算法将样本按3∶1比例划分,24个样本作为训练集,8个样本作为测试集,采用偏最小二乘回归方法(PLS)[25]分别建立原始光谱、原始光谱+FDWS_RL(含整数阶)处理、原始光谱+FDWS_GL处理后的PLS模型.PLS采用5折交叉验证,从1～15中选出最佳潜变量个数.测试集样本的均方根预测误差结果(见表2,表3)表明:5种仪器的光谱在光滑预处理后预测精度都有所提高,其中Bomem,IA260仪器的光谱,经过FDWS_RL光滑预处理后对糖分的预测精度显著提高,两种光谱均方根预测误差分别由1.475,3.579降到0.926,0.488;Nit仪器下光谱对水分的均方根预测误差由0.423降到0.278.FDWS_GL实验结果与之相似.可以发现预测均方根误差最小值所对应的阶次均为分数阶,反映出分数阶WS算法比整数阶更具有应用价值.从理论上讲,分数阶WS阶次连续可调,同时涵盖了整数阶,因此光滑效果永远不会低于整数阶WS算法,最差保持原有整数阶WS算法的光滑效果.实际实验结果同样表明,用分数阶微分作光滑处理时,能够有效地提取有用信号,相较于整数阶WS算法可以进一步提高后期光谱分析算法的准确性.

表2 PLS预测水分、糖分均方根误差(RL)

表3 PLS预测水分、糖分均方根误差(GL)

3 结 论

1) 本文算法采用分数阶微分表示光谱的粗糙度,使算法更加灵活,表示更加丰富.

2) 分数阶WS算法比整数阶WS算法能保留更多的真实信息.在核磁共振谱的实验中得到证明.

3) 分数阶WS算法可以进一步提高后续光谱定量分析的精度,同时表明对原始光谱数据进行光滑处理可以增强有效信号的相对强度,提高光谱分析算法的准确性.在Marzipan红外光谱实验中得到证明.