柠檬酸钠(C₆H₅Na₃O₇, 比利时ACROS公司), 氯金酸(HAuCl₄· 4H₂O, 阿达玛斯试剂公司), 和2-硝基苯甲酸(C₇H₅NO₄, 北京伊诺凯科技有限公司), 芘(C₁₆H₁₀, 国药集团化学试剂公司)。 实验中所用的己烷、 氯仿和乙酸乙酯均为分析纯, 羧甲基纤维素(CMC, 上海梯希爱化工贸易有限公司)和硅藻土。 食用油是从当地超市购买的, 所有实验均使用超纯水。

用于制备金纳米颗粒(Au NPs)的玻璃器皿首先在王水(HNO₃/HCl, 1∶ 3, V/V)中浸泡20 min, 然后用超纯水清洗。 Au NPs 通过柠檬酸钠还原法制备。 将100 mL HAuCl₄ (0.583 2 mmol· L^-1) 水溶液在搅拌条件下加热回流。 然后, 将4.0 mL柠檬酸钠(1%)快速加入回流的溶液中, 继续回流30 min, 冷却至室温。 制备完成的球形金纳米粒子, 直径范围为50至60 nm。

将一毫升食用油与固定浓度的芘溶液点在硅藻土板上, 距离薄层色谱板底部1.5 cm。 在室温条件下, 将其垂直放入流动相中进行样品分离, 然后使用紫外灯(254 nm)与碘色法追踪分析物的位置, 计算保留因子。 然后将3 μ L浓缩30倍金纳米粒子沉积在分析物斑点上。 利用Origin2017作图。 为了现场快速检测, 我们使用785 nm的激发激光的便携式拉曼光谱仪(BWS465-iRman; B& W-Tek, USA), 激光功率为30 mW, 积分时间2 s, 将拉曼光谱数据保存到csv文件中。

由于采集到的原始光谱图上存在很多噪声, 因此, 对光谱数据进行进一步的分析之前, 需要对光谱数据进行基线矫正和降噪等处理。 通过对比Savitzky-Golay平滑滤波、 小波变换和自适应迭代加权惩罚最小二乘等方法, 其实验结果为自适应迭代加权惩罚最小二乘预处理的方法结果最好, 如图1所示。 该方法是在惩罚最小二乘的基础上, 通过在迭代过程中自适应地调节拟合基线与原始信号之间残差平方和的权重, 快速灵活地找到不规则变化的基线并予以扣除。 随机选取9条原始数据使用自适应迭代加权惩罚最小二乘进行预处理, 结果如图2所示。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 预处理算法结果比较图Fig.1 Comparison chart of preprocessing algorithm results

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 光谱预处理图Fig.2 Spectral preprocessing diagram

拉曼光谱进行物质的定量分析是一个热门研究问题, 以往采用的方法无法完善地拟合化学物质之间的非线性关系。 为此本文方法充分利用PCA提取拉曼光谱非线性特征和BP神经网络预测能力强的优点, 提出一种ML-PCA-BP的拉曼光谱预测模型, 具体算法流程如下:

(1)对拉曼光谱的特征峰进行分峰拟合, 采用高斯函数作为基础算法, 拟合获取其峰的关键信息如高度、 半高宽、 面积, 根据拟合程度, 修改高斯拟合函数中参数初始值或者重新选择特征峰对应的起点与终点;

(2)对以拉曼位移作为特征变量的拉曼强度值与高斯拟合获得的参数组成的数据进行归一化处理;

(3)使用主成分分析的方法处理归一化后的数据, 进行降维提取关键参数;

(4)将训练样本的关键参数输入到BP神经网络进行训练, 采用梯度下降算法, 并初始化权重值与偏置值, 为了提高泛化能力采用L2正则化优化BP神经网络; 并且可以根据预测值的误差修改主成分个数、 网络中的参数, 获得最优网络结构, 最终进行浓度预测。

如图2显示了不同浓度的拉曼光谱, 从图中可以明显地看出拉曼位移范围[571.74~598.31]与[1 214.50~1 249.43]有明显的波形, 主要是由于芳环的伸缩模式和CH面内弯曲振动而造成的, 因此属于芘的特征峰。

为提高浓度的预测精度, 从图形中获取特征峰的一些关键参数。 光谱图中每个谱峰都具有明确的物理意义且携带不同的物理信息, 在[571.74~598.31]与[1 214.50~1 249.43]范围内是芘的特征峰。 在且考虑量子力学修正可知拉曼散射强度与样品单位体积内的散射分子数(样品浓度)成正比。 式(1)中采用高斯函数对定位的两个谱峰进行分峰拟合, 由式(1), 式(2)可获得各不同频移位置处谱峰的面积、 半高宽和峰高等参数, 得到特征量参数如表1。

表1

Table 1

表1(Table 1)

表1 光谱图高斯分峰拟合特征参数 Table 1 Characteristic parameters of Gaussian peak division fitting of spectrogram 序号 峰1 峰2 A1 FWHM1 S1 A2 FWHM2 S2 1 124.325 7 12.656 1 4 910.73 5 0.000 0 0.000 0 0.000 2 119.442 5 12.863 6 4 773.723 0.000 0 0.000 0 0.000 3 493.242 5 15.040 4 7 892.856 0.000 0 0.000 0 0.000 4 345.855 3 8.692 1 3 200.509 240.675 3 17.686 3 4530.114 5 208.257 7 10.170 0 2 254.869 171.501 2 19.081 6 3482.765 6 724.547 9 8.805 6 6 792.444 378.714 6 17.222 6 6 941.474 7 1 023.286 0 8.308 3 9 051.223 685.171 1 18.817 4 13 721.440 8 1 387.343 0 8.906 9 13 155.530 986.704 3 14.659 0 15 393.210 9 1 626.462 0 8.531 7 14 774.130 615.577 6 15.424 8 10 105.450 10 2 034.161 0 9.104 2 19 716.720 869.742 3 14.338 7 13 272.450 … … … … … … … 60 14 045.193 8 8.769 6 131 133.650 10 978.590 0 12.894 8 15 0692.634

表1 光谱图高斯分峰拟合特征参数 Table 1 Characteristic parameters of Gaussian peak division fitting of spectrogram

y=y0+Ae-(z-xc)22(FWHM2ln4)2(1)

s=∫-∞+∞y0+Ae-(z-xc)22(FWHM2ln4)2dz=AFWHM2ln42π(2)

其中, y为拉曼散射强度; z为拉曼频移; A为拉曼谱峰高度; x_c为谱峰中心偏移量; FWHM为谱峰半高宽。 得到全部特征量参数如表1。

主成分分析(PCA)主要对数据的协方差矩阵进行特征分解, 从而得到相应的特征值与特征向量, 在实现数据维数缩减的同时保留数据的主体信息。 研究中, 由于特征向量较多, 为了加快运算速率, 所以利用PCA对SERS图范围为[571.74~598.31]与[1 214.50~1 249.43]的特征峰数据与拟合特征峰获得的高度、 半高宽、 面积等参数进行主成分分析。 提取出关键的特征信息, 降低维度, 从图3中可以发现, 即保留3个主成分以后其贡献率达到99.9%。 主成分选取的个数越多, 包含的信息越全面。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 主成分方差分析图Fig.3 Principal component analysis of variance

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 神经网络结构图Fig.4 Neural network structure diagram

BP神经网络模型主要需要确定以下参数: 神经网络的层数, 输入层、 输出层、 隐藏层的节点个数, 各层神经元的激活函数, 训练方法。

(1)参数初始化: 权重值与偏置值使用正态分布生成初始值;

(2)输入层: 设置5个神经元, 为PCA算法获取的5个特征参数;

(3)隐藏层: 设置10层, 分别为300, 270, 240, 210, 180, 150, 120, 90, 60, 30;

(4)输出层: 一个神经元代表浓度;

(5)激活函数: 激活函数采用relu, 输入层到隐藏层与隐藏层内部均采用relu函数, 隐藏层到输出层采用全连接, 无激活函数;

(6)迭代优化器: adam函数, 动量项的梯度下降法;

(7)损失函数: 均方方差损失函数(MSE);

(8)正则化: L2正则化, 主要进行权重值的修改。

神经网络输入层表示为

x=[x1x2x3xi…xt] i=1, 2, …, t(3)

式(3)中, i=1, 2, …, t, t为输入层神经元的个数, x为主成分分析获得的特征向量。

随之隐藏层第一层的的表达式为

h1=[h11h21h31…hj1…ht1](4)

设 wij(1)为输入层第j个神经元第到隐藏层第一层第i个神经元之间的连接权重; bi(1)为隐藏层第1层神经元的偏置, 那么由此可以推出

hi1=f∑j=1qwij(1)xj+b(1)(5)

以及输出层的公式为

y=∑j=1λwij(k)hjk-1+bk(6)

式(6)中, k为输出层的层数, wij(k)为从k-1层第j个神经元与k层第i个神经元之间的连接权重 , hjk-1代表第k-1层的第j个神经元, b^(k)为第k层第i个神经元的偏置, λ 代表k-1层的神经元个数。

m个训练样本的损失函数为

E=12m∑i=1m(di-yi)2(7)

式(7)中, d为真实值, y为预测值, i=1, 2, …m样本个数。

紧接着进行L2正则化的损失函数

EL2=E+λ2m∑ww2(8)

式(8)中, λ 是一个超参数, 范围[0, 1], w则是我们训练的深井网络中每一层的权重矩阵。

根据损失函数反向修改权重的表达式

wij(n)=wik(n)-α∂EL2∂wij(n)(9)

进一步可以推导出

wij(n)=wij(n)-α∂E∂wij(n)-αλmwij(n)(10)

式(10)中, α 为学习率。

根据损失函数反向修改偏置量的表达式

bi(n)=bi(n)-α∂Ebi(n)(11)

(1)采用了薄层色谱和表面增强拉曼散射(SERS)光谱相结合的方法已成功地用于从食用油样品中识别多环芳烃, 拉曼光谱数据保存到csv文件中;

(2)使用自适应迭代加权惩罚最小二乘算法进行光谱预处理;

(3)将光谱数据进行分峰拟合获取高度、 面积、 半高宽等关键参数;

(4)数据集分成测试集与训练集, 数据集可以描述为一个X阶矩阵(n× d), 包含n个长度为d的光谱, 其中d是在[571.74~598.31]与[1 214.50~1 249.43]之间测量到的光谱点与两个峰对应的高度、 面积、 半高宽等参数, 选取10条数据作为测试集;

(5)将高斯拟合的参数与两个峰的特征变量进行归一化, 再进行主成分分析, 依据结果误差依次进行调整主成分个数;

(6)最后将主成分分析获得的特征参数输入基于L2正则化BP神经网络进行训练, 获得最终的预测浓度。

最终基于ML-PCA-BP模型预测出来的芘浓度的预测值与真实值的均方根误差为0.31, 说明真实值与预测值两者的数值十分接近, 预测的浓度相对精准。

为了评价模型的性能, 分别使用相对决定系数R²和均方根误差(RMSE)对模型的预测精度进行了分析。

RMSE=1n∑i=1n(di-yi)2(12)

R2=1-∑i=1n(di-yi)2∑i=1n(di-y-)2(13)

式(13)中, y-为真实浓度的平均值。

从图5可以看出, ML-PCA-BP模型预测值与真实值的函数表达式为y=0.97x+0.68; 偏最小二乘的表达式为y=0.58x+6.64; LR(linear regression)的函数表达式为y=0.45x+9.28。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 测试集芘预测浓度与真实浓度之间的关系Fig.5 Relationship between predicted pyrene concentration and actual pyrene concentration in test set

根据图5和表2中, ML-PCA-BP, PLSR, LR三种算法的对比, 其中ML-PCA-BP的决定系数R²最接近1, 拟合的线性回归方程最优, 预测值与真实值的斜率接近为1, 从图6中和表2中可以发现ML-PCA-BP的预测值与真实值距离最为相近, 其次是PLSR, LR的预测值与真实值相对较远, ML-PCA-BP均方根误差最小为0.31, PLSR算法的均方根误差为1.85, LR的均方根误差为2.28。 从以上的图表中则可以看出ML-PCA-BP算法的预测最为精准。

表2

Table 2

表2(Table 2)

表2 模型对比表 Table 2 Model comparison 模型 R2 RMSE ML-PCA-BP 0.99 0.31 PLSR 0.58 1.85 LR 0.26 2.28

表2 模型对比表 Table 2 Model comparison

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 真实值与实测值对比Fig.6 Comparison between real value and measured value

ML-PCA-BP模型相比于LR算法能够较好提取光谱的非线性特征, 相比于PLSR算法具有较好的泛化能力与较强的预测能力, 因此ML-PCA-BP模型可以更加准确快速地检测出多环芳烃的浓度。