基于Savitzky-Golay平滑插值的烟草近红外光谱模型转移

王 迪，冯伟华，郭军伟，王 锐，刘惠民，宗国浩，刘绍锋，王永胜，赵 乐

中国烟草总公司郑州烟草研究院，郑州高新技术产业开发区枫杨街2号 450001

近红外光谱分析技术是一种结合仪器科学、化学计量学等学科对样品定性或定量的二次分析技术，以无损、快速、无污染等优点而受到越来越多研究人员的青睐，已经在食品、药物、石油、农业及烟草等领域得到广泛应用［1-3］。为了提高近红外预测模型的适用性，往往通过模型转移方法将主机模型应用到子机仪器上以避免重复建模［4］，常用的模型转移算法如截距/斜率校正（S/B）［5］、光谱空间变换（SST）［6］、直接标准化（DS）［7］等均需要主机光谱与子机光谱的波数一致才能够计算。在日常分析过程中发现，即使是相同品牌的傅里叶近红外光谱仪，在更新换代之后，由于激光器的变化，导致其光谱成像波数不一致，造成模型转移算法失效，而插值方法作为有效解决这一类似问题的方法之一，常被用在模型转移前处理过程中，使子机光谱与主机光谱变换一致后再进行模型转移。

插值［8］是离散函数逼近的重要方法，根据函数在已知有限点处的取值，利用插值方法能够估算出函数在其他未知点处的近似值。通过研究Zero、Slinear、Quadratic、Cubic、Nearest［9-11］5种插值方法，分别对光谱数据进行插值及Savitzky-Golay［12］（SG）平滑处理。而SG作为近红外光谱经典平滑预处理方法被广泛采用，本项目组在前期的研究过程中曾使用该方法对光谱预处理后实现了建模应用。在此基础上，本研究中提出了一种新的基于SG平滑的插值方法（SG-Inter）对光谱进行变换，从而达到将主机、子机两种仪器光谱补齐的目的；通过利用Zero、Slinear、Quadratic、Cubic、Nearest以及SG-Inter方法对光谱进行插值处理得到新的子机光谱数据后，再经过光谱空间变换（SST）将光谱转移后预测烟草总植物碱、总糖、总氮、还原糖、氯、钾、淀粉、新植二烯和绿原酸共9种具有代表性的烟草化学指标的质量分数，根据预测结果统计、分析、评价不同插值方法对模型转移效果的影响。

1 材料与方法

1.1 材料、试剂和仪器

收集2020年全国各省级中烟工业有限责任公司的初烤烟叶样品，共800个；采用Kennard-Stone（KS）方法［13］筛选出200个代表性样品，分别采集其近红外光谱。

主机MPA（1代）和子机TANGO（2代）光谱仪（德国Bruker公司）；ZM200型粉碎机（德国Retsch公司）；BSA124S型电子天平（感量0.000 1 g，德国Satorious公司）；AA3连续流动分析仪（德国BRAN+RUBBE公司）。

1.2 方法

1.2.1 样品处理与化学指标检测

参考行业标准［14］干燥处理样品，直至可用手指捻碎。将样品通过粉碎机粉碎研磨，并过0.250 mm（60目）分样筛，混匀后装入密封袋。分别采用行业标准，无相关标准方法的采用文献方法，测定样品中的总植物碱［15］、总糖和还原糖［16］、总氮［17］、氯［18］、钾［19］、淀粉［20］、新植二烯［21］以及绿原酸［22］的质量分数。

1.2.2 光谱采集与预处理

将仪器的激光能量范围设置为4 000～10 000 cm-1，分辨率设置为8 cm-1，扫描次数设置为64次，分别用主机MPA（1代）和子机TANGO（2代）在相同的条件下采集近红外光谱，采用Savitzky-Golay平滑（窗口17，2次多项式，1阶导数）进行光谱预处理。

1.2.3 主机模型

选择项目组在前期研究中所构建的模型作为烟草近红外分析主机模型，用于本研究中模型转移效果评价。

1.2.4 插值方法的研究与应用

采用Zero、Slinear、Quadratic、Cubic、Nearest 5种插值方法分别对原始光谱以及SG平滑后的光谱数据进行插值处理，将子机数据与主机数据补齐，随后利用光谱空间变换（SST）方法对插值后的光谱进行处理，并预测模型转移后烟叶中总植物碱、总糖、总氮、还原糖、氯、钾、淀粉、新植二烯和绿原酸的质量分数。

1.2.5 基于Savitzky-Golay的平滑插值方法

基于Savitzky-Golay平滑插值（SG-Inter）的模型转移，即采用平滑处理与插值同步计算来进行模型转移。本研究中，SG-Inter（窗口大小为17，1阶导数，2次多项式）随滑动窗口在子机光谱上的移动，求得窗口内子机光谱数据的2阶多项式并1阶求导，把子机光谱SG平滑点的最邻近主机光谱波长点值代入到函数中求得插值。把插值变换融合到数据预处理过程中，保留原始光谱的变化趋势，将波数为1 456的子机光谱平滑插值到与主机光谱相同的1 555个波长点，子机与主机光谱波数调整一致后进行模型转移并预测总植物碱、总糖、总氮、还原糖、氯、钾、淀粉、新植二烯和绿原酸的质量分数。具体计算原理如下：

Savitzky-Golay平滑是一种卷积滑动窗口的加权平均算法，设滤波窗口的宽度w=2i+1，i是半窗宽度，i=1，2，3，…，n；x代表数据点在窗口内的相对位置，其取值为［-i，…，0，…，i］，数据点所在位置对应的函数值为P（x）。根据窗口内的数据点，构造一个n阶多项式进行拟合得到f（x）表达式，见公式（1）［12］：

式中：cn0，cn1，…，cnn分别代表拟合函数f（x）中拟合数据点的系数，经过最小二乘拟合，得到残差E的表达式，见公式（2）［12］：

目标是使残差E趋于最小，将公式（2）中各项系数导数εz分别设置为0，z=（0，1，2，3，…，n），得到公式（3）［12］：

将公式（3）化简得到公式（4）：

当滑动窗口大小与平滑阶数固定时，将待拟合窗口［P（-i），…，P（0），…，P（i）］内数据带入公式（4），可求得多项式系数列表［cn0，cn1，…，cnn］T。如图1所示，设平滑窗口w=5，每次求解窗口内第w/2=3个位置的平滑多项式（实心圆点为原始信号点，空心圆点为待插值信号点，实心方点为待插值信号点的最邻近原始信号点），随着卷积核窗口以步长为1向前平移，Savitzky-Golay平滑方法将重新拟合卷积核窗口平移后的平滑多项式，并求得窗口中心点所在位置的平滑后数值。本研究中将已知光谱波长点的对应数值xi与最邻近待预测波长点对应数值xj进行等比例缩放得到比例系数m=xi/xj，根据已知波长点值的位置yi，求得待预测波长点的相对位置yj=yi·m，将yj代入到平滑多项式f（x）中求解，得到待预测波长点的数据值。经过卷积核窗口的滑动，最终得到平滑在线插值后的光谱数据列表［f（y0），f（y1），…，f（yn）］，再经过一系列数据预处理操作后，将转移后光谱数据列表应用于主机模型进行化学成分预测分析。

图1 平滑插值（SG-Inter）原理图（w=5）Fig.1 Smoothing interpolation theory（Window=5）

1.2.6 光谱转移及化学成分预测

以MPA（1代）扫描的标准样品光谱为主机光谱，MPA（2代）扫描的标准样品光谱为子机光谱，将子机光谱平滑求导后插值处理，使子机光谱与主机光谱的横坐标光谱波数调整一致，并按照光谱空间变换方法（SST）［6］将调整后的子机光谱转移至主机光谱，使用主机模型对转移后的子机光谱进行化学指标分析预测，待得到指标预测结果后再统计分析。

1.2.7 数据处理方法

使用计算机编程语言Python实现本研究中的各种模型计算。

2 结果与分析

2.1 插值方法对光谱的直接影响

2.1.1 原始光谱直接插值后的谱图分析

使用Zero、Slinear、Quadratic、Cubic、Nearest 5种插值方法对子机原始光谱进行插值处理，并与原始光谱对比得到图2。可知，原始光谱能够与这几种插值方法处理后的光谱数据基本重合。说明插值方法能够解决子机光谱与主机光谱波数不一致的问题。

图2 光谱直接插值谱图对比Fig.2 Comparison of spectra after direct interpolation

模型转移过程中通常需要对光谱进行平滑求导预处理以消除仪器噪音的影响，将光谱平滑处理对比后得到图3。从图3中能明显看出不同方法插值后的光谱与原始光谱平滑处理后差异较大，通过对几处差异明显的出峰点局部放大后，可以清晰、直观地看出各插值方法处理后的光谱与原始光谱的差异。由于原始光谱直接插值时只考虑了插值区间内的数据处理，没有将数据点左右两侧数据值变化趋势考虑在内，但光谱SG平滑预处理时，需要对每个平滑窗口的表达式进行1阶求导，即光谱直接插值的结果与原始光谱在谱图出峰点位置的凸函数表达不一致，导致平滑后的数据存在明显差异。

图3 插值光谱平滑求导（1阶）预处理后的谱图对比Fig.3 Comparison of spectra after interpolation and 1st derivative smoothing

2.1.2 光谱平滑后再插值处理的谱图分析

利用Zero、Slinear、Quadratic、Cubic、Nearest 5种插值方法对经过SG平滑后的子机光谱进行插值，同时利用本研究中提出的基于Savitzky-Golay平滑的插值方法（SG-Inter）对光谱进行处理后与原始光谱平滑后对比得到图4。可知，平滑后光谱数据插值与原始光谱平滑后的数据能够基本重合。说明平滑后再插值可以解决子机光谱与主机光谱波数不一致而导致的模型转移问题。进一步将出峰点光谱局部放大后，发现经过Cubic、Quadratic及SG-Inter插值方法处理后的数据能够较好地与原始光谱相重合，而其他几类方法如Zero、Nearest等均存在较为明显的波动。

图4 平滑求导（1阶）后再插值的光谱对比Fig.4 Comparison of spectra after 1st derivative smoothing and interpolation

由于光谱的平滑过程已经对原始数据进行降噪处理，结合谱图变化，从光谱形态变化的角度分析后不难发现，插值方法能够在此基础上充分获取数据点前后的变化趋势，更易于寻找插值函数多项式使子机光谱与主机光谱的波数调整一致。因此，先平滑再插值的效果要优于先插值再平滑的效果，并有效支撑模型转移和化学指标定量预测。

2.2 插值方法对模型转移预测结果的影响

使用不同方法插值后的光谱进行模型转移及指标预测，研究插值方法对指标预测的定量影响。为了进一步对子机光谱的预测结果进行统计分析，选择各指标的均方根误差（Root mean square error，RMSE）、模型决定系数（R2）、相对分析误差［23］（Residual predictive deviation，RPD）作为评价指标，分别对模型转移效果进行分析比较。

2.2.1 基于原始光谱直接插值的模型转移预测结果分析

利用Zero、Slinear、Quadratic、Cubic、Nearest 5种插值方法对子机原始光谱直接插值，通过对插值后的光谱进行空间变换方法（SST）［6］实现光谱转移，将各指标预测结果与主机模型预测结果统计分析后得到表1。可知，基于原始光谱直接插值的模型转移及指标预测结果的RMSE 远大于主机模型交叉验证结果。说明模型的预测精度距离主机模型还有一定差距。同样，基于光谱直接插值的转移模型的R2普遍较小，说明模型与预测指标之间的关联关系还不够紧密。此外，通常认为，若RPD＜1.4，表明所建模型预测结果不可靠；若1.4≤RPD≤2.0，表明模型的预测结果可以接受；若RPD＞2.0，则表明模型的预测准确性很高，能够用于模型分析［23］。从表1 可知，整体上，基于Cubic、Quadratic 等非线性插值方法处理后的各指标预测结果比其他插值后的预测效果好，模型的RPD 均大于1.4，在可接受范围内，但相对于主机模型的预测性能还有较大的差距。

表1 原始光谱直接插值下9种烟草化学指标的RMSE、R2、RPD值Tab.1 RMSE，R2 and RPD values of nine tobacco chemical indexes after interpolation

就总植物碱、总糖、总氮、还原糖、氯、钾、淀粉、新植二烯和绿原酸的质量分数而言，绿原酸的模型预测能力最差，且预测误差最大。对于模型转移方法，基于子机原始光谱直接插值的不同方法整体预测效果均比较差，转移模型决定系数均较小，预测误差相对较大，尤其是基于Slinear、Quadratic、Cubic和Nearest插值方法处理的光谱，模型转移预测能力均处于同一水平；经过Zero 插值处理的光谱模型转移预测效果最差。

2.2.2 基于光谱平滑数据插值的模型转移预测结果

将子机原始光谱SG平滑预处理后，利用Zero、Slinear、Cubic、Nearest和Quadratic方法以及本研究中提出的基于Savitzky-Golay的平滑插值（SG-Inter）方法对数据进行变换，各化学指标定量预测的统计分析结果见表2。

表2 光谱平滑后再插值下的RMSE、R2、RPD值Tab.2 RMSE，R2 and RPD values after smoothing and interpolation

可知，对于总糖和氯来说，SG-Inter插值处理的光谱预测误差明显优于其他插值方法；对于总植物碱、总氮、钾、还原糖、淀粉和新植二烯来说，无论是基于平滑后再Cubic、Quadratic等插值方法还是平滑插值（SG-Inter）的模型转移预测误差RMSE均较小且差别不大，但基于Zero插值方法的上述9种化学指标模型转移预测结果相对不稳定且误差较大；对于多酚类物质绿原酸，所建主机模型决定系数不是很高，对不同插值方法处理后的光谱进行模型转移，化学指标预测结果误差较大，且模型不够稳健，后续仍需要对主机模型进行优化，对子机模型转移进行更深入的研究。综合来看，尽管与主机模型相比，这9种化学指标的子机模型预测RMSE均大于主机模型交叉验证的RMSE，且外部验证误差稍大于主机模型内部交叉验证误差还可以被理解，但仍然有提升的空间；而对于总植物碱、总糖、总氮、还原糖和新植二烯来说，SG-Inter插值和Nearest、Cubic、Quadratic等插值方法处理后的光谱整体上模型预测能力与主机模型相当且差距不大，且均较为良好；但对于氯、钾、绿原酸和淀粉来说，与主机模型预测结果的R2相比仍有一定差距。

综上，针对本研究中所实测的200个代表性样品，基于所提出的SG-Inter插值方法处理子机光谱数据，在经过光谱空间变换（SST）实现子机光谱转移后，对9种指标进行主机模型预测，结果表明：除绿原酸外，其余8种指标的模型转移预测结果的RMSE均相对最小，转移模型的R2与RPD在多种指标中的评价性能均良好、略优或与其他插值方法效果相当。表明SG-Inter方法可以顺利支撑模型转移及预测工作的开展，而且该方法将光谱平滑预处理过程与插值计算结合在一起，在保证模型预测能力和精度的前提下，能够大大提升模型转移的工作效率。而对于绿原酸等其他指标的模型转移工作需要进一步从近红外光谱形成机理着手分析数据，减少SG-Inter方法处理过程中造成的光谱精度损失，从而提升指标预测的准确度。

3 结论

①使用插值的方法能够有效解决子机光谱与主机光谱波数及波长点不一致而导致模型转移中主机模型失效的问题。②子机光谱直接插值后与原始光谱能相对重合，但平滑处理后的光谱差异显著，且模型转移预测结果误差较大；子机光谱SG平滑处理后再插值的光谱与原始光谱能相对重合，且模型转移预测效果良好。③基于Savitzky-Golay平滑插值（SG-Inter）的模型转移对多种指标质量分数的预测效果略优于Cubic、Quadratic、Slinear等5种插值方法，指标预测误差相对较小，模型预测能力较好，且计算流程清晰，原理简单，能够将平滑处理与插值分析集为一体并同步高效运算，可减少预处理环节的分析误差，提高模型转移的工作效率。