基于CARS和K-S的马拉硫磷农药浓度吸收光谱预测模型研究

甄欢仪，马瑞峻，陈 瑜，孙小鹏，马创立

华南农业大学工程学院，广东 广州 510642

引 言

马拉硫磷是一种低毒高效、残效期短的有机磷杀虫剂，不仅被广泛应用于农业生产，而且其原药及其混合剂是谷物常用的储粮保护剂，但不合理地使用马拉硫磷农药会对人体和环境造成不可逆转的影响。有机磷农药的速测方法有酶联免疫法、酶抑制法和生物传感器法等。其中酶抑制法因其成本低、检测时间短和操作简单等特点，成为蔬菜农药残留快速检测最为常用的方法，但受不同的前处理方法、基质干扰方法和最低检测限高的影响，酶抑制检测的样品部分结果存在假阳性和假阴性风险[1-2]。寻找一种成本低、灵敏度高、操作简单方便和能够实现现场准确快速检测有机磷农药的方法有重要的现实意义和实际应用价值。

光谱技术是一种快速、无损检测技术，但直接使用光谱技术检测有机磷农药，存在光谱背景噪声干扰明显和吸收峰偏移不稳定、特征不明显、检测精度难以满足定量要求等问题。然而将光谱技术和化学计量学分析方法结合，对光谱数据进行去噪声预处理、提取关键变量和特征样本以及建立预测模型来压缩有用信息，可有效地解决光谱检测中存在的问题，并被不同领域的研究人员广泛应用于实际检测中[3-4]。目前采用光谱技术快速有效地定量检测水体中有机磷农药的研究鲜有报道。

本工作以马拉硫磷农药为研究对象，采用紫外/可见吸光度光谱仪获取不同浓度实验样品的光谱数据，结合竞争性自适应重加权采样(competitive adaptive reweighted sampling，CARS)变量选择算法和经典的Kennard-Stone(K-S)算法样本选择方法，建立偏最小二乘(PLS)定量分析预测模型。为验证CARS算法的变量选择上的性能，将其与常用有效的变量筛选方法蒙特卡洛无消息变量消除法(Monte Carlo-uninformative variable elimination，MC-UVE)进行系统的比较。并将K-S算法与蒙特卡洛交叉验证(Monte Carlo cross-validation，MCCV)在样本优选性能进行比较。期望通过对关键变量和更优样本的筛选得到鲁棒性好且简单可靠的定量预测模型，实现对马拉硫磷农药浓度的快速有效的定量检测。

1 实验部分

1.1 样本

马拉硫磷实验样本配制： 用1/10 000电子天平称取98%马拉硫磷标准样本0.204 1(±0.000 2) g，用少量甲醇超声溶解，然后加入甲醇准确定容至400 mL，摇匀，得到浓度为500 mg·L-1的马拉硫磷标准液。使用纯净水作为稀释剂配制不同浓度的实验样本，浓度为0.1～1.00 mg·L-1(浓度梯度为0.05 mg·L-1)，共19个样本； 1.50～10 mg·L-1(浓度梯度为0.5 mg·L-1)，共18个样本； 12～50 mg·L-1(浓度梯度为2 mg·L-1)，共20个样本； 55～100 mg·L-1(浓度梯度为5 mg·L-1)，共10个样本。合计67个浓度样本。

划分样本集前，参照文献[5]先进行异常样本剔除检查，结果表明并无异常样本，然后进行下一步样本划分。采用基于变量之间欧式(欧几里得)距离的K-S算法，按3∶1的比例在特征空间中划分出校正集(50个样本)和预测集(17个样本)。

1.2 紫外/可见光谱数据采集

使用美国Ocean Optics海阳光学公司的便携式吸光度测量仪(型号为FLAME-CHEM-UV-VIS)采集马拉硫磷的每个样本浓度的光谱数据。光谱仪的光学分辨率为1.5 nm，波长范围为200～850 nm，其间隔为0.39 nm，得到波长变量总数为1 914个。光谱数据采集软件设置积分时间为25 ms，平滑度为2，每个样本扫描20次，取平均值为最终光谱。

1.3 关键变量筛选与样本优选的方法

1.3.1 关键变量筛选方法

CARS算法基于变量回归系数绝对值|bi|，遵循“适者生存”原则，使用指数衰减函数(exponentially decreasing function，EDF)和自适应重加权采样技术(adaptive reweighted sampling，ARS)去除回归系数|bi|值相对小的波长变量，留下回归系数|bi|值权重值较大且共线性小的变量组成新的变量子集。算法的具体执行步骤可参考文献[8]。

MC-UVE算法[6]基于变量回归系数bi的稳定性值(回归系数bi均值和标准差比值的绝对值)，可得到由回归系数较稳定且含有有效信息较多的变量组成最优变量子集。该方法在许多光谱波长变量筛选研究中已被证明能有效剔除无信息变量。

1.3.2 样本优选方法

建模样本的优劣关乎整个模型的质量，影响模型预测能力。蒙特卡洛交叉验证MCCV通过多次的蒙特卡罗(MC)采样，随机将原来校正集样本划分校正子集和预测子集，建立相应的PLS模型并得到模型的预测子集浓度预测偏差(RMSEP)，最后通过设定样本的预测残差平均值和预测残差标准差的阈值选择更优样本集。K-S法基于光谱特征选取样本。

1.4 模型评价

对PLS模型的评价采用决定系数R2、均方根误差RMSE(包括校正集均方根误差RMSE和预测集均方根误差RMSEP)和预测集的RPD。一个好的模型应具有高的R2和低的RMSE。RPD是预测集的标准偏差与预测均方根误差的比值，反映了模型的分辨能力和稳健性，当RPD≥3时，模型可应用于定量分析和实际检测。本研究中的样本划分、数据预处理、PLS建模、关键波长变量选择、样本优选都基于The Unscrambler X10.4和MATLAB R2017b软件平台进行。

2 结果与讨论

2.1 光谱数据预处理

马拉硫磷不同浓度梯度下的原始吸收光谱(200～850 nm)如图1所示。

图1 马拉硫磷原始吸收光谱图Fig.1 Original absorption spectrum ofmalathion

由图1可知，由于光谱仪在其适用的波长范围两端光源能量较低，光谱基线两端出现漂移下滑现象，同时全光谱存在较大的背景噪声干扰，检测时，在200～250 nm波长范围内，光谱出现吸收峰向左偏移不稳定和低浓度与高浓度光谱曲线之间空缺的现象，另外超过750 nm光谱基线下滑现象逐渐严重，影响农药浓度在可见光部分有效的吸光度检测。研究中将200.08～750.04 nm的光谱数据分别使用基线校正、基线补偿、不同点的S-G平滑、一阶导数、二阶导数、标准化和均值中心化等数据预处理方法进行处理，进一步获取有效的分析信号以提高模型预测性能。通过对不同预处理后PLS模型结果对比分析，表明均值中心化预处理效果最好。后续的变量筛选、样本优选和模型建立研究均采用均值中心化处理后的光谱数据，如图2所示。

2.2 关键变量的选取与模型结果

由图1和图2可知，马拉硫磷在紫外波段范围内的吸收光谱的吸收峰在203.52～217.6 nm之间偏移不稳定，最大吸收波长点难以确定。而传统的标准曲线定量分析方法只能针对有明显吸收峰出现且其相应的波长点位置稳定的农药样品进行分析，使用传统的标准曲线定量分析法建立的马拉硫磷定量预测模型稳定性会降低，所以需要对马拉硫磷光谱数据采用化学计量学方法进行处理并建立定量预测分析模型。

图2 均值中心化后的马拉硫磷吸收光谱图Fig.2 Mean centering absorption spectrum of malathion

图3 竞争性自适应重加权算法选出的关键变量(a): 采样变量数的变化趋势； (b): RMSECV值的变化趋势；(c): 每个变量的回归系数路径变化趋势Fig.3 Key variables selected by CARS

(a): The changing trend of the number of sampled variables; (b): The changing trend of RMSECV value; (c): The changing trend of regression coefficits of each variables

马拉硫磷吸收光谱的每条光谱线含1 582个波长变量，波长变量数目非常多，需要在不削弱原来模型预测性能的情况下，对波长变量和校正样本集进行优选，简化模型，建立一个简单和预测能力强的马拉硫磷浓度定量预测模型。

光谱数据在执行一次CARS算法(设置MC采样为50次)后，随采样次数增加，波长变量保留数变化情况、50个不同的变量子集模型交叉验证均方根误差RMSECV值、波长变量回归系数变化路径分别如图3(a)，(b)和(c)所示。

如图3(a)所示，随着采样次数增加，由于指数衰减函数EDF的作用，波长变量的保留数目变化先快速后缓慢，在采样10次之前先大幅减少，经历了一个“粗选”过程之后再慢慢“精选”； 如图3(b)所示，RMSECV先减小，说明多重共线性或含无关信息的波长变量被剔除，模型预测能力提高，RMSECV达最小值[对应图3(c)星号垂直标示的第34次]之后逐渐增大，说明一些含有有关马拉硫磷浓度值重要信息的波长变量被剔除，导致模型的预测性能下降。因此由CARS算法在第34次采样时有最小RMSECV，获得了18个关键波长变量。

执行一次MC-UVE算法后的波长变量筛选结果如图4所示。

图4 蒙特卡洛无消息变量消除法选出的关键变量(a)： 波长变量回归系数稳定性值；(b)： 不同变量子集的PLS模型的RMSEPFig.4 Key variables selected by MC-UVE

(a)： Stability value of regression coefficient of wavelength variable； (b)： The RMSEP value corresponding to PLS model obtained by different subset of variables

1 582个波长变量的回归系数稳定性值变化情况如图4(a)所示，图中的水平虚线为稳定性阈值，该阈值是RMSEP最小的变量子集中波长变量回归系数的最小稳定性值，低于该值的波长变量被剔除； 将样正集和预测集的波长变量按稳定性值从大到小排序，每次加入10个波长变量进行PLS建模和预测，直到加入所有的波长变量，得到不同变量子集PLS模型RMSEP的结果如图4(b)所示，当波长变量数的保留数目为300时，RMSEP=1.595 0为最小值，获得最优变量子集并得到稳定性阈值，之后保留变量数目增多，RMSEP增大，说明无关信息增加，影响了预测的结果。

潜变量(latent variables，LVs)是一个非常重要的模型参数，直接影响到模型的预测性能。当潜变量选择少时，所建立的模型可能丢失较多有用的光谱信息，出现欠拟合现象，而当潜变量选择多时，模型中可能包含过多噪声，出现过拟合现象[7]。因此合理选择建模LVs是充分利用光谱信息和降低噪声的有效方法之一。通过不同潜变量PLS模型预测RMSEP合理选择各模型的潜变量，采用Full-PLS和CARS-PLS模型潜变量选择过程如图5所示。

采用全波段变量以及由上述CARS算法和MC-UVE算法各自得到的波长变量子集，并分别选择合理的LVs建立相应的马拉硫磷浓度PLS定量预测模型，综合对比分析两种变量筛选方法在关键波长变量上的选取性能，分别建立的模型结果见表1。

由图4和表1可知，运行一次MC-UVE后，得到LVs=5时，MC-UVE-PLS的minRMSEP=1.595 0。同时由图5可知，当LVs=5时，得到Full-PLS和CARS-PLS的minRMSEP分别为1.226 8和0.923 0，三个模型的LVs相同。对于MC-UVE-PLS模型，其变量数虽然从全波段1 582个变量减少到300个，但其RMSEP高于以及RPD明显低于Full-PLS模型的值，可能是因为一些回归系数较稳定但含有无关信息的变量被保留了下来，导致模型的预测性能没有提高。而CARS-PLS模型在全面提高模型的预测性能上，仅用了18个变量，占全波段变量数的1.13%，其 RMSEP就比Full-PLS模型的值降低了24.8%，且RPD值有所增加。

图5 不同潜变量的RMSEP的变化情况(a)： 不同潜变量Full-PLS的RMSEP变化情况；(b)： 不同潜变量CARS-PLS的RMSEP变化情况Fig.5 Variation of RMSEP with LVs(a)： Variation of RMSEP of Full-PLS model with LVs；(b)： Variation of RMSEP of CARS-PLS model with LVs

表1 不同波长变量集的马拉硫磷浓度PLS模型预测结果Table 1 The PLS prediction model results of malathion concentration by using the different variables set

通过对上述模型结果对比分析，表明CARS算法在马拉硫磷浓度吸收光谱关键变量的筛选性能上优于MC-UVE算法，说明CARS算法不但不会削弱模型的预测能力，还能进行关键变量的有效筛选，简化模型，提高模型稳健性。由此可知，模型的预测能力还与建模变量数和质量有关，通过提取关键变量和消除变量之间的多重共线性能够增强模型的预测能力。

2.3 样本优选与模型结果

对于建模集样本，除了要消除光谱变量之间存在较为严重的共线性之外，样本与样本之间差异性的大小也会影响模型的预测性能。合理有效地消除或者降低这种差异性，也有助于提高定量预测模型的预测性能。采用MCCV法对50个样本运行一次(设置MC采样2 500次)，以预测集(样本浓度预测)的残差均值为横坐标，残差标准差为纵坐标，得到50个样本的分布如图6所示。

将图6中预测残差平均值大于3和预测残差标准差大于0.8的6个样本剔除，被剔除的奇异样本序号分别为24，45，47，48，49和50。最后得到优选的44个校正集样本。

为了验证K-S算法样本优选方法的性能，采用MCCV作为对照方法与其进行合理有效的比较，用K-S算法也筛选44个建模样本。由MCCV法和K-S算法分别得到的44个建模样本与CARS算法筛选的18个关键变量组成新的建模集，两模型LVs的选择过程分别如图7(a,b)所示，然后分别建立相应的PLS定量预测模型，结果见表2。其中马拉硫磷浓度CARS-CCVs-PLS模型预测结果如图8所示。

图6 MCCV法50个样本分布图Fig.6 The distribution diagram of 50 samplesobtained from MCCV

由图7(a,b)可知，当LVs为5时，CARS-K-Ss-PLS和CARS-MCCVs-PLS的minRMSEP分别为0.863 4和1.026 2。

由表2和图8可知，CARS-K-Ss-PLS预测模型的RMSEC与RMSEP的值更相近，即RMSE值更稳定，且CARS-K-Ss-PLS的RPD高于CARS-MCCVs-PLS的值，说明CARS-K-Ss-PLS模型的预测能力更好，K-S法能够对校正集样本进行优选，且能够提高预测模型预测能力。MCCV法通过样本残差均值和残差标准偏差这两个统计参数来剔除样本浓度预测结果情况较差的样本，而K-S算法在特征空间中采用样本光谱为其特征变量，基于变量之间的欧式距离来均匀选取样本，充分考虑了样本光谱信息的影响。

图7 不同潜变量的RMSEP的变化情况

(a)： 不同潜变量CARS-K-Ss-PLS的RMSEP变化情况； (b)： 不同潜变量CARS-MCCVs-PLS的RMSEP变化情况

Fig.7VariationofRMSEPwithLVs

(a)： Variation of RMSEP of CARS-K-Ss-PLS model with LVs； (b)： Variation of RMSEP of CARS-MCCVs-PLS model with LVs

表2 基于关键变量和优选样本的PLS模型预测结果Table 2 PLS model prediction results based on key variables and selected effective samples

图8 CARS-K-Ss-PLS模型预测结果Fig.8 Prediction result of CARS-K-Ss-PLS model

3 结 论

(1)对马拉硫磷不同浓度梯度吸收光谱原始数据采用不同的预处理方法并建模，结果表明均值中心化预处理结果最优。

(2)采用CARS算法和MC-UVE算法分别对全波段变量进行筛选并建模，结果表明CARS-PLS模型预测结果最优，建模变量数最少，选出的18个关键波长变量比全波段的1 582个变量减少了98.87%，简化模型的同时提高了模型的预测能力。

(3)采用CARS算法筛选的18个波长变量，分别结合MCCV法和K-S算法各自优选的44个样本建模，结果表明K-S法在样本优选性能上略优，可对校正集样本进行优选，提高模型预测性能。

(4)通过对比分析各种算法建立的马拉硫磷浓度定量预测模型，结果表明采用CARS算法结合K-S法建立的CARS-K-Ss-PLS模型最优，可用于马拉硫磷浓度快速准确的定量检测。

(5)该研究为采用紫外/可见吸光度光谱法快速实时检测水体中马拉硫磷农药浓度提供了理论依据。