王雪映, 劉適搏, 朱繼偉, 馬婷婷

1. 中國科學(xué)院長春光學(xué)精密機械與物理研究所, 吉林 長春 130033

2. 中國科學(xué)院大學(xué), 北京 100049

3. 長春長光格瑞光電技術(shù)有限公司, 吉林 長春 130102

4. 吉林省水文水資源局, 吉林 長春 130022

引 言

化學(xué)需氧量(chemical oxygen demand, COD)是水樣在一定條件下, 氧化1 L水樣中還原性物質(zhì)所消耗氧化劑的量, 以mg·L-1表示, 它反映了水中受還原性物質(zhì)的污染程度[1], 是水質(zhì)評價的重要指標(biāo), 其數(shù)值越大, 代表水質(zhì)越差, 受污染越嚴(yán)重[2]。 目前國際社會普遍公認(rèn)的COD檢測的經(jīng)典標(biāo)準(zhǔn)方法為重鉻酸鹽法, 此外還有高錳酸鉀法、 分光光度法、 快速消解法等檢測方法。 這些方法雖檢測精度較高, 但存在操作繁瑣、 使用有毒試劑、 易造成二次污染等缺點。

高光譜技術(shù)是水質(zhì)檢測技術(shù)發(fā)展的新趨勢, 具有檢測快捷、 無需使用有毒試劑、 無二次污染等優(yōu)點。 基于高光譜數(shù)據(jù)反演COD濃度早有研究, 林劍遠(yuǎn)等對嘉興市主城區(qū)河網(wǎng)的COD等水質(zhì)參數(shù)進(jìn)行了反演, 基于航空和水表高光譜遙感數(shù)據(jù), 對水質(zhì)采樣化驗數(shù)據(jù)和水表反射率進(jìn)行相關(guān)性分析, 計算出最佳波段組合并進(jìn)行回歸分析和擬合, 研究表明COD的最優(yōu)擬合模型為二次多項式模型, 其決定系數(shù)R2和均方差RMSE分別為0.74和2.79[3]。 埃及的Elsayed等基于高光譜數(shù)據(jù)對卡龍湖中的多個水質(zhì)參數(shù)進(jìn)行了評估, 使用手持光譜儀采集水體光譜反射率, 利用已發(fā)布的光譜反射率指數(shù)、 新兩波段光譜反射率指數(shù)和新三波段光譜反射率指數(shù)以及人工神經(jīng)網(wǎng)絡(luò)建立模型, 發(fā)現(xiàn)新三波段光譜指數(shù)與COD具有中等相關(guān)性,R2在0.52～0.64之間[4]。 董月群等以廣東省鶴山市沙坪河為研究區(qū)域, 通過無人機高光譜成像遙感數(shù)據(jù)對城市河網(wǎng)的COD等水質(zhì)參數(shù)進(jìn)行反演研究, 通過分析水表反射率0階導(dǎo)數(shù)、 一階導(dǎo)數(shù)、 二階導(dǎo)數(shù)與水質(zhì)參數(shù)濃度的相關(guān)性選取相關(guān)性最高的波段組合, 并用統(tǒng)計回歸的方法建立模型, 研究結(jié)果顯示一階導(dǎo)數(shù)對應(yīng)的線性擬合模型具有較好的反演精度,R2介于0.85～0.96之間[5]。

上述對地表水COD濃度的反演, 多是采用無人機載高光譜成像儀于室外高空拍攝或在室外利用手持光譜儀獲得光譜數(shù)據(jù), 目前少有在室內(nèi)環(huán)境下采集地表水高光譜數(shù)據(jù)并進(jìn)行COD反演的研究。 為了填補這一研究領(lǐng)域的空白, 拓展高光譜水質(zhì)分析技術(shù)的使用范圍, 本研究以地表水樣本理化分析的實測值結(jié)果為參考值, 在實驗室內(nèi)拍攝并提取出其高光譜數(shù)據(jù), 以標(biāo)準(zhǔn)漫反射板為參照計算相應(yīng)的光譜反射率, 通過Pearson相關(guān)性分析提取出原始光譜和導(dǎo)數(shù)光譜數(shù)據(jù)中與COD濃度相關(guān)性較高的特征譜數(shù)據(jù), 利用粒子群算法優(yōu)化的最小二乘支持向量機(particle swarm optimization-least squares support vector machines, PSO-LSSVM)算法, 基于原始光譜、 一階導(dǎo)數(shù)全譜、 二階導(dǎo)數(shù)全譜、 一階導(dǎo)數(shù)特征譜、 二階導(dǎo)數(shù)特征譜5種光譜數(shù)據(jù)建立不同的COD反演模型, 通過決定系數(shù)(coefficient of determination,R2)、 均方差(root mean square error, RMSE)和相對偏差(relative percent deviation, RPD)分析比較這幾種模型的預(yù)測精度和可靠性, 為開發(fā)實驗室環(huán)境下高光譜水質(zhì)分析設(shè)備進(jìn)行技術(shù)儲備。

1 實驗部分

1.1 樣品

地表水樣本來源于吉林省內(nèi)的新立城水庫、 太平池水庫、 松花江等流域的23個站點, 水樣COD濃度的實測值以相應(yīng)的理化分析結(jié)果為參考。 剔除具有極端實測值的樣本后, 數(shù)據(jù)集樣本總量為129個, COD濃度范圍為1～60 mg·L-1, 其中有56個樣本濃度在10～20 mg·L-1范圍內(nèi)。 采用隨機抽樣的方法抽取32個樣本做驗證集, 剩余的97個樣本點做訓(xùn)練集, 使得訓(xùn)練集和驗證集的樣本個數(shù)比約為3∶1, 具體劃分情況如表1所示。

表1 地表水樣本COD濃度統(tǒng)計分析

1.2 數(shù)據(jù)采集

檢測裝置如圖1所示, 實驗儀器主要有高精密電驅(qū)位移臺及其控制箱、 高光譜成像光譜儀、 計算機等; 實驗材料及用品主要有地表水樣品、 標(biāo)準(zhǔn)漫反射板等。 高光譜成像光譜儀采用的是長春長光格瑞光電技術(shù)有限公司自主研發(fā)的GRS-400PG推掃型成像光譜儀, 分辨率約2.4 nm, 光譜范圍為400～1 000 nm; 實驗光源選用四個25 W鎢絲燈泡置于成像光譜儀前端鏡頭四周, 以更好地模擬室外太陽光的均勻性; 盛樣皿口徑為12 cm, 容積為400 mL; 選用的標(biāo)準(zhǔn)漫反射板的反射率為50%, 尺寸為13 cm×13 cm。

圖1 檢測裝置實物圖

數(shù)據(jù)采集時, 每個水樣的取樣量為250 mL, 將高光譜成像光譜儀固定于位移臺垂直上方約50 cm處, 蓋上遮光布, 驅(qū)動位移臺勻速移動, 同步拍攝每個水樣及標(biāo)準(zhǔn)漫反射板的高光譜圖像, 重復(fù)采集10次以減少測試誤差。

1.3 數(shù)據(jù)預(yù)處理

1.3.1 反射率計算

根據(jù)實驗室內(nèi)拍攝的高光譜數(shù)據(jù)計算水體光譜反射率公式為

(1)

式(1)中,R(λ)為不同波長的水體光譜反射率, 下標(biāo)λ表示波長。DN水為探測器輸出的代表水體輻射光譜強度的DN值數(shù)據(jù),DNboard(λ)為同步獲得的50%標(biāo)準(zhǔn)漫反射板的DN值數(shù)據(jù),Rboard(λ)為50%標(biāo)準(zhǔn)漫反射板的反射率[6]。

1.3.2 導(dǎo)數(shù)法處理

為了降低背景噪聲干擾, 提升光譜數(shù)據(jù)對COD濃度的敏感度, 采用導(dǎo)數(shù)法進(jìn)行處理。 參考文獻(xiàn)[7-8]處理導(dǎo)數(shù)光譜的方法, 計算導(dǎo)數(shù)光譜的公式為

(2)

(3)

式(2)和式(3)中,R,R1st,R2nd分別為原始光譜反射率值、 一階導(dǎo)數(shù)光譜反射率值、 二階導(dǎo)數(shù)光譜反射率值;λ為波長。

1.3.3 特征譜提取

高光譜數(shù)據(jù)具有波段多、 數(shù)據(jù)量大的特點, 提取特征譜段可以提高計算速度, 減少冗余噪聲影響。 相關(guān)性分析是判斷兩個或多個變量之間統(tǒng)計學(xué)關(guān)聯(lián)的方法, 相關(guān)系數(shù)是用以反映變量之間相關(guān)關(guān)系密切程度的統(tǒng)計指標(biāo)。 采用Pearson相關(guān)性分析方法, 對特征譜段進(jìn)行分析并提取, 其相關(guān)系數(shù)計算公式為

(4)

Pearson相關(guān)系數(shù)的范圍為[-1, 1], 其值越接近于1或-1, 表明相關(guān)性越強; 相關(guān)系數(shù)越接近于0, 表明相關(guān)性越弱。

1.4 反演模型

最小二乘支持向量機(least squares support vector machines, LSSVM)是在支持向量機(support vector machine, SVM)中引入最小二乘線性系統(tǒng)作損失函數(shù)的機器學(xué)習(xí)方法, 用二次損失函數(shù)取代SVM中的不敏感損失函數(shù), 通過構(gòu)造損失函數(shù)將SVM的二次尋優(yōu)變?yōu)榍蠼饩€性方程。 該方法簡化了運算復(fù)雜性, 提升了運算速度, 通過結(jié)構(gòu)風(fēng)險最小化原理可提高模型泛化能力, 可以有效解決傳統(tǒng)學(xué)習(xí)方法處理有限樣本數(shù)據(jù)、 非線性、 高維數(shù)等問題的困難, 越來越廣泛地被應(yīng)用于非線性復(fù)雜系統(tǒng)的建模中。 徑向基核函數(shù)(radial basis function, RBF)應(yīng)用廣泛、 適用性強、 具有較寬收斂域, 對非線性系統(tǒng)具有很好的控制性能, 幾乎任何情況(維數(shù)過高過低、 樣本數(shù)據(jù)過多過少等)下, 均能表現(xiàn)出較理想的效果[9], 因此本研究選用RBF核函數(shù)作為模型的核函數(shù)。

在LSSVM模型中, 正則化參數(shù)和核函數(shù)參數(shù)的選擇對構(gòu)造回歸函數(shù)是至關(guān)重要的, 直接影響著模型的學(xué)習(xí)和泛化能力。 優(yōu)化模型參數(shù)的常用方法有遺傳算法(genetic algorithm, GA)[10]、 粒子群算法(particle swarm optimization, PSO)[11-13]、 蝙蝠算法(bat algorithm, BA)[14]和灰狼算法(grey wolf optimizer, GWO)[15]等。 本研究選用了PSO算法來優(yōu)化尋參, 它是一種模擬鳥群覓食的群體智能算法, 將只有位置和速度兩個屬性的粒子模擬為鳥, 每個粒子在空間中的極值作為潛在的最優(yōu)解; 所有粒子根據(jù)評價函數(shù)決定的適應(yīng)值不斷調(diào)整各自的位置和速度, 以尋找最優(yōu)解。

設(shè)置粒子種群數(shù)量為20, 最高迭代次數(shù)為50, 正則化參數(shù)γ和核函數(shù)參數(shù)σ的取值范圍均設(shè)置為[0.01, 100], 采用K折交叉驗證對數(shù)據(jù)集進(jìn)行訓(xùn)練調(diào)整,K取值為5, 用交叉驗證的均方差RMSE作為評價函數(shù)的適應(yīng)值。 經(jīng)過PSO算法計算, 可得到正則化參數(shù)γ和核函數(shù)參數(shù)σ的最優(yōu)參數(shù)值, 將其代入LSSVM算法中, 可以使LSSVM算法有較強的適應(yīng)能力和較好的預(yù)測精度。

1.5 模型評價

采用決定系數(shù)R2、 均方差RMSE和相對偏差RPD三個指標(biāo)用于模型對比和評價

(5)

(6)

(7)

2 結(jié)果與討論

2.1 水體反射率

圖2(a)為樣本水體的原始光譜曲線, 水體在波長400～620 nm處反射率較高, 在波長750～900 nm處反射率迅速下降。 使用導(dǎo)數(shù)法對原始光譜進(jìn)行一階導(dǎo)數(shù)、 二階導(dǎo)數(shù)處理, 所得的地表水光譜反射率數(shù)據(jù)分布如圖2(b)、 (c)所示。 對比分析可得: 導(dǎo)數(shù)處理有效強化了光譜曲線的特征, 呈現(xiàn)出較多的波峰、 波谷; 但仍無法直觀判斷出各波段水體光譜反射率與地表水樣COD濃度間的明顯線性關(guān)系。

圖2 水體光譜反射率

2.2 相關(guān)性分析與特征譜提取

以訓(xùn)練集樣本數(shù)據(jù)為相關(guān)性分析對象, 對水樣的光譜反射率數(shù)據(jù)和COD濃度實測值進(jìn)行相關(guān)性分析, 所得Pearson相關(guān)系數(shù)如圖3所示。 未經(jīng)處理的原始光譜反射率數(shù)據(jù)的相關(guān)性較小, 整體變化趨勢為相關(guān)性隨波長的增大而增大, 在波長872 nm處相關(guān)性最大, 相關(guān)系數(shù)為-0.309 8。 一階導(dǎo)數(shù)、 二階導(dǎo)數(shù)光譜反射率數(shù)據(jù)存在相關(guān)性較大的波段, 無相關(guān)系數(shù)隨波長變化的趨勢。 其中, 一階導(dǎo)數(shù)光譜在波長436 nm處相關(guān)性最大, 相關(guān)系數(shù)為-0.496 8; 二階導(dǎo)數(shù)光譜的相關(guān)系數(shù)波動很大, 和一階導(dǎo)光譜相比有更多明顯的峰值, 在波長855 nm處相關(guān)性最大, 相關(guān)系數(shù)為0.483 5。 以上結(jié)果說明, 原始光譜對COD濃度變化的敏感度較低, 在本研究中導(dǎo)數(shù)光譜與地表水的COD濃度之間的相關(guān)性相對更高, 利用導(dǎo)數(shù)光譜數(shù)據(jù)建?？赡芫哂懈叩念A(yù)測精度。

圖3 訓(xùn)練集相關(guān)系數(shù)

基于上述相關(guān)性分析結(jié)果, 分別選出相關(guān)系數(shù)絕對值大于0.35的波段作為特征譜段, 對應(yīng)的特征譜相關(guān)系數(shù)分布情況如圖4所示。 一階導(dǎo)數(shù)光譜相關(guān)性較好的波長共有38個, 主要分布在在648～811 nm處; 二階導(dǎo)光譜相關(guān)性較好的波長共有24個, 主要分布在在725～811 nm處。

圖4 特征譜相關(guān)系數(shù)

2.3 建模反演結(jié)果及對比

根據(jù)Pearson相關(guān)性分析結(jié)果提取特征譜數(shù)據(jù), 利用PSO-LSSVM算法基于訓(xùn)練集數(shù)據(jù)建立COD濃度反演模型, 然后將驗證集數(shù)據(jù)代入所建模型中加以驗證以確定最優(yōu)反演模型。 由于PSO-LSSVM模型的運行結(jié)果具有不確定性, 因此每組模型運行10次后取平均值來作為最終結(jié)果。

在采用全譜數(shù)據(jù)建模時, 對訓(xùn)練集原始光譜數(shù)據(jù)、 一階導(dǎo)數(shù)光譜數(shù)據(jù)、 二階導(dǎo)數(shù)光譜數(shù)據(jù)三類數(shù)據(jù), 取其波長400～900 nm內(nèi)全譜數(shù)據(jù)數(shù)據(jù), 利用PSO-LSSVM算法進(jìn)行建模, 訓(xùn)練集和驗證集的預(yù)測結(jié)果分別如圖5和圖6所示。

圖5 利用全譜數(shù)據(jù)建模的訓(xùn)練集預(yù)測結(jié)果

圖6 利用全譜數(shù)據(jù)建模的驗證集預(yù)測結(jié)果

在采用特征譜數(shù)據(jù)建模時, 因在對原始光譜數(shù)據(jù)進(jìn)行相關(guān)性分析時未能找到特征波長, 故只對訓(xùn)練集一階導(dǎo)數(shù)特征光譜數(shù)據(jù)、 二階導(dǎo)數(shù)特征光譜數(shù)據(jù)進(jìn)行建模與分析。 根據(jù)相關(guān)性分析結(jié)果, 以相關(guān)系數(shù)絕對值大于0.35為篩選準(zhǔn)則確定特征光譜, 利用PSO-LSSVM算法進(jìn)行建模, 訓(xùn)練集和驗證集的預(yù)測結(jié)果分別如圖7和圖8所示。

圖7 利用特征譜數(shù)據(jù)建模的訓(xùn)練集預(yù)測結(jié)果

圖8 利用特征譜數(shù)據(jù)建模的驗證集預(yù)測結(jié)果

表2給出了利用全譜數(shù)據(jù)和特征譜數(shù)據(jù)建模時訓(xùn)練集和驗證集預(yù)測結(jié)果的R2、 RMSE和RPD。 訓(xùn)練集的交叉驗證結(jié)果中, 利用一階全譜、 一階特征譜、 二階特征譜數(shù)據(jù)建模的R2均大于0.98, RMSE均小于1.5, RPD均大于9, 遠(yuǎn)遠(yuǎn)優(yōu)于利用原始譜和二階全譜數(shù)據(jù)建模的結(jié)果, 說明前三種模型具有較強的適應(yīng)能力和較好的預(yù)測精度。 驗證集的預(yù)測結(jié)果中, 利用原始光譜建模的R2和RPD較低, RMSE較高, 表明其預(yù)測效果較差; 利用導(dǎo)數(shù)光譜數(shù)據(jù)建模的R2和RPD較高, RMSE較低, 表明其具有較好的預(yù)測效果; 其中利用一階導(dǎo)數(shù)特征譜建模時預(yù)測結(jié)果最好,R2和RPD最大, RMSE最小。 同時, 訓(xùn)練集的交叉驗證結(jié)果和驗證集的預(yù)測結(jié)果都出現(xiàn)了經(jīng)同階導(dǎo)數(shù)預(yù)處理后利用特征譜數(shù)據(jù)建模結(jié)果優(yōu)于利用全譜數(shù)據(jù)建模結(jié)果的情況。

表2 地表水COD濃度的預(yù)測結(jié)果評價

實驗結(jié)果表明, 利用原始光譜數(shù)據(jù)建模的模型預(yù)測能力較差, 而利用導(dǎo)數(shù)光譜數(shù)據(jù)進(jìn)行建?？梢蕴岣吣Ｐ偷念A(yù)測精度; 另外, 對于導(dǎo)數(shù)光譜數(shù)據(jù), 提取特征譜數(shù)據(jù)所建立的模型比利用全譜數(shù)據(jù)建立的模型有更好地預(yù)測效果。

3 結(jié) 論

以吉林省的地表水樣品為研究對象, 分析樣品的高光譜數(shù)據(jù)和COD濃度實測值之間的相關(guān)性, 利用PSO-LSSVM算法建立地表水COD濃度的反演模型, 通過對相關(guān)反演方法的研究及對預(yù)測結(jié)果的分析, 得出如下結(jié)論:

(1)地表水COD濃度與原始光譜反射率的相關(guān)性較弱, 無明顯特征波長。 經(jīng)過導(dǎo)數(shù)法預(yù)處理后, 地表水COD濃度與光譜反射率的相關(guān)性明顯增強。

(2)利用特征譜數(shù)據(jù)建模的結(jié)果要優(yōu)于利用全譜數(shù)據(jù)建模的結(jié)果, 其中利用一階導(dǎo)數(shù)特征譜建模反演地表水COD濃度的效果最好, 精度最高, 訓(xùn)練集的R2、 RMSE和RPD分別為0.993 4、 1.063 8和12.198 5, 驗證集的R2、 RMSE和RPD分別為0.856 7、 3.822 9和2.641 4。

(3)在室內(nèi)基于高光譜數(shù)據(jù)反演地表水COD濃度的應(yīng)用具有可行性, 驗證了所提出的基于PSO-LSSVM算法、 導(dǎo)數(shù)處理法、 Pearson相關(guān)性分析法建立的地表水COD濃度反演模型, 能達(dá)到較好反演精度。 并且, 未來隨地表水樣本數(shù)據(jù)量的不斷積累, 在反演模型的適用濃度、 反演精度、 模型的穩(wěn)健性等方面都將得到提升, 日后也將進(jìn)一步探究基于高光譜數(shù)據(jù)和PSO-LSSVM算法估測其他水質(zhì)參數(shù)的可能性。