基于CARS和K-S的馬拉硫磷農(nóng)藥濃度吸收光譜預(yù)測模型研究

甄歡儀，馬瑞峻，陳 瑜，孫小鵬，馬創(chuàng)立

華南農(nóng)業(yè)大學(xué)工程學(xué)院，廣東 廣州 510642

引 言

馬拉硫磷是一種低毒高效、殘效期短的有機磷殺蟲劑，不僅被廣泛應(yīng)用于農(nóng)業(yè)生產(chǎn)，而且其原藥及其混合劑是谷物常用的儲糧保護劑，但不合理地使用馬拉硫磷農(nóng)藥會對人體和環(huán)境造成不可逆轉(zhuǎn)的影響。有機磷農(nóng)藥的速測方法有酶聯(lián)免疫法、酶抑制法和生物傳感器法等。其中酶抑制法因其成本低、檢測時間短和操作簡單等特點，成為蔬菜農(nóng)藥殘留快速檢測最為常用的方法，但受不同的前處理方法、基質(zhì)干擾方法和最低檢測限高的影響，酶抑制檢測的樣品部分結(jié)果存在假陽性和假陰性風險[1-2]。尋找一種成本低、靈敏度高、操作簡單方便和能夠?qū)崿F(xiàn)現(xiàn)場準確快速檢測有機磷農(nóng)藥的方法有重要的現(xiàn)實意義和實際應(yīng)用價值。

光譜技術(shù)是一種快速、無損檢測技術(shù)，但直接使用光譜技術(shù)檢測有機磷農(nóng)藥，存在光譜背景噪聲干擾明顯和吸收峰偏移不穩(wěn)定、特征不明顯、檢測精度難以滿足定量要求等問題。然而將光譜技術(shù)和化學(xué)計量學(xué)分析方法結(jié)合，對光譜數(shù)據(jù)進行去噪聲預(yù)處理、提取關(guān)鍵變量和特征樣本以及建立預(yù)測模型來壓縮有用信息，可有效地解決光譜檢測中存在的問題，并被不同領(lǐng)域的研究人員廣泛應(yīng)用于實際檢測中[3-4]。目前采用光譜技術(shù)快速有效地定量檢測水體中有機磷農(nóng)藥的研究鮮有報道。

本工作以馬拉硫磷農(nóng)藥為研究對象，采用紫外/可見吸光度光譜儀獲取不同濃度實驗樣品的光譜數(shù)據(jù)，結(jié)合競爭性自適應(yīng)重加權(quán)采樣(competitive adaptive reweighted sampling，CARS)變量選擇算法和經(jīng)典的Kennard-Stone(K-S)算法樣本選擇方法，建立偏最小二乘(PLS)定量分析預(yù)測模型。為驗證CARS算法的變量選擇上的性能，將其與常用有效的變量篩選方法蒙特卡洛無消息變量消除法(Monte Carlo-uninformative variable elimination，MC-UVE)進行系統(tǒng)的比較。并將K-S算法與蒙特卡洛交叉驗證(Monte Carlo cross-validation，MCCV)在樣本優(yōu)選性能進行比較。期望通過對關(guān)鍵變量和更優(yōu)樣本的篩選得到魯棒性好且簡單可靠的定量預(yù)測模型，實現(xiàn)對馬拉硫磷農(nóng)藥濃度的快速有效的定量檢測。

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]先進行異常樣本剔除檢查，結(jié)果表明并無異常樣本，然后進行下一步樣本劃分。采用基于變量之間歐式(歐幾里得)距離的K-S算法，按3∶1的比例在特征空間中劃分出校正集(50個樣本)和預(yù)測集(17個樣本)。

1.2 紫外/可見光譜數(shù)據(jù)采集

使用美國Ocean Optics海陽光學(xué)公司的便攜式吸光度測量儀(型號為FLAME-CHEM-UV-VIS)采集馬拉硫磷的每個樣本濃度的光譜數(shù)據(jù)。光譜儀的光學(xué)分辨率為1.5 nm，波長范圍為200～850 nm，其間隔為0.39 nm，得到波長變量總數(shù)為1 914個。光譜數(shù)據(jù)采集軟件設(shè)置積分時間為25 ms，平滑度為2，每個樣本掃描20次，取平均值為最終光譜。

1.3 關(guān)鍵變量篩選與樣本優(yōu)選的方法

1.3.1 關(guān)鍵變量篩選方法

CARS算法基于變量回歸系數(shù)絕對值|bi|，遵循“適者生存”原則，使用指數(shù)衰減函數(shù)(exponentially decreasing function，EDF)和自適應(yīng)重加權(quán)采樣技術(shù)(adaptive reweighted sampling，ARS)去除回歸系數(shù)|bi|值相對小的波長變量，留下回歸系數(shù)|bi|值權(quán)重值較大且共線性小的變量組成新的變量子集。算法的具體執(zhí)行步驟可參考文獻[8]。

MC-UVE算法[6]基于變量回歸系數(shù)bi的穩(wěn)定性值(回歸系數(shù)bi均值和標準差比值的絕對值)，可得到由回歸系數(shù)較穩(wěn)定且含有有效信息較多的變量組成最優(yōu)變量子集。該方法在許多光譜波長變量篩選研究中已被證明能有效剔除無信息變量。

1.3.2 樣本優(yōu)選方法

建模樣本的優(yōu)劣關(guān)乎整個模型的質(zhì)量，影響模型預(yù)測能力。蒙特卡洛交叉驗證MCCV通過多次的蒙特卡羅(MC)采樣，隨機將原來校正集樣本劃分校正子集和預(yù)測子集，建立相應(yīng)的PLS模型并得到模型的預(yù)測子集濃度預(yù)測偏差(RMSEP)，最后通過設(shè)定樣本的預(yù)測殘差平均值和預(yù)測殘差標準差的閾值選擇更優(yōu)樣本集。K-S法基于光譜特征選取樣本。

1.4 模型評價

對PLS模型的評價采用決定系數(shù)R2、均方根誤差RMSE(包括校正集均方根誤差RMSE和預(yù)測集均方根誤差RMSEP)和預(yù)測集的RPD。一個好的模型應(yīng)具有高的R2和低的RMSE。RPD是預(yù)測集的標準偏差與預(yù)測均方根誤差的比值，反映了模型的分辨能力和穩(wěn)健性，當RPD≥3時，模型可應(yīng)用于定量分析和實際檢測。本研究中的樣本劃分、數(shù)據(jù)預(yù)處理、PLS建模、關(guān)鍵波長變量選擇、樣本優(yōu)選都基于The Unscrambler X10.4和MATLAB R2017b軟件平臺進行。

2 結(jié)果與討論

2.1 光譜數(shù)據(jù)預(yù)處理

馬拉硫磷不同濃度梯度下的原始吸收光譜(200～850 nm)如圖1所示。

圖1 馬拉硫磷原始吸收光譜圖Fig.1 Original absorption spectrum ofmalathion

由圖1可知，由于光譜儀在其適用的波長范圍兩端光源能量較低，光譜基線兩端出現(xiàn)漂移下滑現(xiàn)象，同時全光譜存在較大的背景噪聲干擾，檢測時，在200～250 nm波長范圍內(nèi)，光譜出現(xiàn)吸收峰向左偏移不穩(wěn)定和低濃度與高濃度光譜曲線之間空缺的現(xiàn)象，另外超過750 nm光譜基線下滑現(xiàn)象逐漸嚴重，影響農(nóng)藥濃度在可見光部分有效的吸光度檢測。研究中將200.08～750.04 nm的光譜數(shù)據(jù)分別使用基線校正、基線補償、不同點的S-G平滑、一階導(dǎo)數(shù)、二階導(dǎo)數(shù)、標準化和均值中心化等數(shù)據(jù)預(yù)處理方法進行處理，進一步獲取有效的分析信號以提高模型預(yù)測性能。通過對不同預(yù)處理后PLS模型結(jié)果對比分析，表明均值中心化預(yù)處理效果最好。后續(xù)的變量篩選、樣本優(yōu)選和模型建立研究均采用均值中心化處理后的光譜數(shù)據(jù)，如圖2所示。

2.2 關(guān)鍵變量的選取與模型結(jié)果

由圖1和圖2可知，馬拉硫磷在紫外波段范圍內(nèi)的吸收光譜的吸收峰在203.52～217.6 nm之間偏移不穩(wěn)定，最大吸收波長點難以確定。而傳統(tǒng)的標準曲線定量分析方法只能針對有明顯吸收峰出現(xiàn)且其相應(yīng)的波長點位置穩(wěn)定的農(nóng)藥樣品進行分析，使用傳統(tǒng)的標準曲線定量分析法建立的馬拉硫磷定量預(yù)測模型穩(wěn)定性會降低，所以需要對馬拉硫磷光譜數(shù)據(jù)采用化學(xué)計量學(xué)方法進行處理并建立定量預(yù)測分析模型。

圖2 均值中心化后的馬拉硫磷吸收光譜圖Fig.2 Mean centering absorption spectrum of malathion

圖3 競爭性自適應(yīng)重加權(quán)算法選出的關(guān)鍵變量(a): 采樣變量數(shù)的變化趨勢； (b): RMSECV值的變化趨勢；(c): 每個變量的回歸系數(shù)路徑變化趨勢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個波長變量，波長變量數(shù)目非常多，需要在不削弱原來模型預(yù)測性能的情況下，對波長變量和校正樣本集進行優(yōu)選，簡化模型，建立一個簡單和預(yù)測能力強的馬拉硫磷濃度定量預(yù)測模型。

光譜數(shù)據(jù)在執(zhí)行一次CARS算法(設(shè)置MC采樣為50次)后，隨采樣次數(shù)增加，波長變量保留數(shù)變化情況、50個不同的變量子集模型交叉驗證均方根誤差RMSECV值、波長變量回歸系數(shù)變化路徑分別如圖3(a)，(b)和(c)所示。

如圖3(a)所示，隨著采樣次數(shù)增加，由于指數(shù)衰減函數(shù)EDF的作用，波長變量的保留數(shù)目變化先快速后緩慢，在采樣10次之前先大幅減少，經(jīng)歷了一個“粗選”過程之后再慢慢“精選”； 如圖3(b)所示，RMSECV先減小，說明多重共線性或含無關(guān)信息的波長變量被剔除，模型預(yù)測能力提高，RMSECV達最小值[對應(yīng)圖3(c)星號垂直標示的第34次]之后逐漸增大，說明一些含有有關(guān)馬拉硫磷濃度值重要信息的波長變量被剔除，導(dǎo)致模型的預(yù)測性能下降。因此由CARS算法在第34次采樣時有最小RMSECV，獲得了18個關(guān)鍵波長變量。

執(zhí)行一次MC-UVE算法后的波長變量篩選結(jié)果如圖4所示。

圖4 蒙特卡洛無消息變量消除法選出的關(guān)鍵變量(a)： 波長變量回歸系數(shù)穩(wěn)定性值；(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個波長變量的回歸系數(shù)穩(wěn)定性值變化情況如圖4(a)所示，圖中的水平虛線為穩(wěn)定性閾值，該閾值是RMSEP最小的變量子集中波長變量回歸系數(shù)的最小穩(wěn)定性值，低于該值的波長變量被剔除； 將樣正集和預(yù)測集的波長變量按穩(wěn)定性值從大到小排序，每次加入10個波長變量進行PLS建模和預(yù)測，直到加入所有的波長變量，得到不同變量子集PLS模型RMSEP的結(jié)果如圖4(b)所示，當波長變量數(shù)的保留數(shù)目為300時，RMSEP=1.595 0為最小值，獲得最優(yōu)變量子集并得到穩(wěn)定性閾值，之后保留變量數(shù)目增多，RMSEP增大，說明無關(guān)信息增加，影響了預(yù)測的結(jié)果。

潛變量(latent variables，LVs)是一個非常重要的模型參數(shù)，直接影響到模型的預(yù)測性能。當潛變量選擇少時，所建立的模型可能丟失較多有用的光譜信息，出現(xiàn)欠擬合現(xiàn)象，而當潛變量選擇多時，模型中可能包含過多噪聲，出現(xiàn)過擬合現(xiàn)象[7]。因此合理選擇建模LVs是充分利用光譜信息和降低噪聲的有效方法之一。通過不同潛變量PLS模型預(yù)測RMSEP合理選擇各模型的潛變量，采用Full-PLS和CARS-PLS模型潛變量選擇過程如圖5所示。

采用全波段變量以及由上述CARS算法和MC-UVE算法各自得到的波長變量子集，并分別選擇合理的LVs建立相應(yīng)的馬拉硫磷濃度PLS定量預(yù)測模型，綜合對比分析兩種變量篩選方法在關(guān)鍵波長變量上的選取性能，分別建立的模型結(jié)果見表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模型，其變量數(shù)雖然從全波段1 582個變量減少到300個，但其RMSEP高于以及RPD明顯低于Full-PLS模型的值，可能是因為一些回歸系數(shù)較穩(wěn)定但含有無關(guān)信息的變量被保留了下來，導(dǎo)致模型的預(yù)測性能沒有提高。而CARS-PLS模型在全面提高模型的預(yù)測性能上，僅用了18個變量，占全波段變量數(shù)的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模型預(yù)測結(jié)果Table 1 The PLS prediction model results of malathion concentration by using the different variables set

通過對上述模型結(jié)果對比分析，表明CARS算法在馬拉硫磷濃度吸收光譜關(guān)鍵變量的篩選性能上優(yōu)于MC-UVE算法，說明CARS算法不但不會削弱模型的預(yù)測能力，還能進行關(guān)鍵變量的有效篩選，簡化模型，提高模型穩(wěn)健性。由此可知，模型的預(yù)測能力還與建模變量數(shù)和質(zhì)量有關(guān)，通過提取關(guān)鍵變量和消除變量之間的多重共線性能夠增強模型的預(yù)測能力。

2.3 樣本優(yōu)選與模型結(jié)果

對于建模集樣本，除了要消除光譜變量之間存在較為嚴重的共線性之外，樣本與樣本之間差異性的大小也會影響模型的預(yù)測性能。合理有效地消除或者降低這種差異性，也有助于提高定量預(yù)測模型的預(yù)測性能。采用MCCV法對50個樣本運行一次(設(shè)置MC采樣2 500次)，以預(yù)測集(樣本濃度預(yù)測)的殘差均值為橫坐標，殘差標準差為縱坐標，得到50個樣本的分布如圖6所示。

將圖6中預(yù)測殘差平均值大于3和預(yù)測殘差標準差大于0.8的6個樣本剔除，被剔除的奇異樣本序號分別為24，45，47，48，49和50。最后得到優(yōu)選的44個校正集樣本。

為了驗證K-S算法樣本優(yōu)選方法的性能，采用MCCV作為對照方法與其進行合理有效的比較，用K-S算法也篩選44個建模樣本。由MCCV法和K-S算法分別得到的44個建模樣本與CARS算法篩選的18個關(guān)鍵變量組成新的建模集，兩模型LVs的選擇過程分別如圖7(a,b)所示，然后分別建立相應(yīng)的PLS定量預(yù)測模型，結(jié)果見表2。其中馬拉硫磷濃度CARS-CCVs-PLS模型預(yù)測結(jié)果如圖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預(yù)測模型的RMSEC與RMSEP的值更相近，即RMSE值更穩(wěn)定，且CARS-K-Ss-PLS的RPD高于CARS-MCCVs-PLS的值，說明CARS-K-Ss-PLS模型的預(yù)測能力更好，K-S法能夠?qū)πＵ瘶颖具M行優(yōu)選，且能夠提高預(yù)測模型預(yù)測能力。MCCV法通過樣本殘差均值和殘差標準偏差這兩個統(tǒng)計參數(shù)來剔除樣本濃度預(yù)測結(jié)果情況較差的樣本，而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 基于關(guān)鍵變量和優(yōu)選樣本的PLS模型預(yù)測結(jié)果Table 2 PLS model prediction results based on key variables and selected effective samples

圖8 CARS-K-Ss-PLS模型預(yù)測結(jié)果Fig.8 Prediction result of CARS-K-Ss-PLS model

3 結(jié) 論

(1)對馬拉硫磷不同濃度梯度吸收光譜原始數(shù)據(jù)采用不同的預(yù)處理方法并建模，結(jié)果表明均值中心化預(yù)處理結(jié)果最優(yōu)。

(2)采用CARS算法和MC-UVE算法分別對全波段變量進行篩選并建模，結(jié)果表明CARS-PLS模型預(yù)測結(jié)果最優(yōu)，建模變量數(shù)最少，選出的18個關(guān)鍵波長變量比全波段的1 582個變量減少了98.87%，簡化模型的同時提高了模型的預(yù)測能力。

(3)采用CARS算法篩選的18個波長變量，分別結(jié)合MCCV法和K-S算法各自優(yōu)選的44個樣本建模，結(jié)果表明K-S法在樣本優(yōu)選性能上略優(yōu)，可對校正集樣本進行優(yōu)選，提高模型預(yù)測性能。

(4)通過對比分析各種算法建立的馬拉硫磷濃度定量預(yù)測模型，結(jié)果表明采用CARS算法結(jié)合K-S法建立的CARS-K-Ss-PLS模型最優(yōu)，可用于馬拉硫磷濃度快速準確的定量檢測。

(5)該研究為采用紫外/可見吸光度光譜法快速實時檢測水體中馬拉硫磷農(nóng)藥濃度提供了理論依據(jù)。