楊可明 李 燕 程 鳳 高 鵬 張 超

(中國礦業(yè)大學(xué)(北京)煤炭資源與安全開采國家重點實驗室， 北京 100083)

0 引言

近年來，礦山開采、污水排放、礦物肥料使用等造成的土壤重金屬污染問題日益嚴(yán)重，給人類生產(chǎn)生活帶來嚴(yán)重影響[1]。銅(Cu)是土壤重金屬污染的主要元素之一，Cu在土壤中的累積對植被生長發(fā)育產(chǎn)生嚴(yán)重影響[2-3]。同時，通過食物鏈進(jìn)入人體，將破壞人體神經(jīng)系統(tǒng)和免疫系統(tǒng)[4]。因此，植被重金屬污染的實時、有效監(jiān)測顯得尤為重要。研究表明，受重金屬污染的植被光譜曲線會發(fā)生畸變[5-7]，因此很多學(xué)者利用高光譜遙感技術(shù)獲取重金屬脅迫下的植被光譜信息，旨在通過監(jiān)測植被光譜特征變化來診斷重金屬污染程度。

電磁波在傳播過程中不可避免受到大氣等影響而產(chǎn)生噪聲，噪聲對光譜數(shù)據(jù)中污染信息的解讀造成極大干擾[8]，尋找一種有效提取作物污染弱信息的方法成為高光譜遙感重金屬污染監(jiān)測領(lǐng)域的研究熱點。近年來越來越多的特征提取方法被引入到高光譜遙感領(lǐng)域，主要有小波變換[9-10](Wavelet transform，WT)、分形特征提取法[11-12](Fractal dimension，F(xiàn)D)、經(jīng)驗?zāi)B(tài)分解[13](Empirical mode decomposition，EMD)等。上述方法在一定程度上均可較好地提取出光譜變異信息，但仍存有一定局限性。WT難以選擇合適的小波基函數(shù)，F(xiàn)D受到信號長度的影響，EMD分解過程中存在模態(tài)混疊現(xiàn)象。變分模態(tài)分解(Variational mode decomposition，VMD)是近年來提出的一種新的自適應(yīng)性信號分解方法，該方法解決了EMD方法的問題，并具有很好的噪聲魯棒性[14]，目前被廣泛應(yīng)用于信號分析及故障特征提取等方面[15-16]，而在高光譜遙感領(lǐng)域鮮見報道。

本文引入VMD理論，并結(jié)合多尺度熵(Multiscale entropy，MSE)構(gòu)建一種基于變分模態(tài)分解的MSE光譜弱信息提取及定量描述模型(VMD-MSE)。利用VMD對受Cu2+污染的玉米光譜數(shù)據(jù)進(jìn)行分解獲得模態(tài)分量，以MSE方法對其進(jìn)行定量描述獲得VMD-MSE模型值(VM)，探討VM值與玉米葉片中Cu2+含量的相關(guān)關(guān)系，并基于VM值構(gòu)建玉米葉片中Cu2+含量預(yù)測模型，以實現(xiàn)玉米重金屬污染程度的診斷，為農(nóng)作物重金屬污染監(jiān)測研究提供參考。

1 理論與方法

1.1 變分模態(tài)分解(VMD)

VMD分解過程實質(zhì)為變分問題的求解過程，主要分為變分約束問題的建立和求解兩部分，最終將信號分解為數(shù)個有限帶寬的固有模態(tài)分量[17-18]。該過程包含自適應(yīng)濾波組，表現(xiàn)了較好的噪聲魯棒性。

(1)變分約束問題的建立

假設(shè)每個“模態(tài)”是具有中心頻率的有限帶寬，變分問題描述為尋求K個模態(tài)函數(shù)uk(t)(k=1，2，…，K)，使得每個模態(tài)的估計帶寬之和最小，約束條件為各模態(tài)之和等于輸入信號f。具體步驟如下：①每個模態(tài)通過Hilbert變換計算與之相關(guān)的解析信號。②對于每個模態(tài)，通過加入指數(shù)項調(diào)整各自估計的中心頻率，把模態(tài)的頻譜變換到基帶上。③通過對解調(diào)信號進(jìn)行H1高斯平滑對帶寬進(jìn)行估計。④得到一個變分約束問題

(1)

其中

(u1，u2，…，uK)=u

(ω1，ω2,…，ωK)=ω

式中δ(t)——單位脈沖函數(shù)

u——分解得到的K個模態(tài)分量

ω——各模態(tài)分量的中心頻率

f——輸入的原始信號

最后求解該問題。

(2)變分約束問題的求解

①引入二次懲罰因子α和拉格朗日算子λ(t)，將約束性變分問題變?yōu)榉羌s束性變分問題，構(gòu)成擴展的拉格朗日表達(dá)式為

L({uk},{ωk},λ)=

(2)

式中α——二次懲罰因子，在高斯噪聲存在的情況下可以保證信號的重構(gòu)準(zhǔn)確度

λ(t)——拉格朗日算子，用來保持約束條件的嚴(yán)格性

②利用ADMM算法迭代搜索求取上述擴展的拉格朗日函數(shù)的鞍點，即式(1)約束變分模型的最優(yōu)解，其中解得模態(tài)分量uk及中心頻率ωk分別為

(3)

(4)

1.2 多尺度熵(MSE)

MSE是基于樣本熵(Sample entropy, SpEn)的一種時間序列復(fù)雜性的度量方法，用于從不同尺度度量時間序列的復(fù)雜程度[19]。假設(shè)原始數(shù)據(jù)為X={x1，x2，…，xN}，則MSE的具體計算步驟如下[20]：

(1)給定嵌入維數(shù)m、相似容限r(nóng)及尺度因子τ=1，2，…，τmax，建立新的粗粒序列

(5)

對于每個τ，原序列被分為N/τ個長度為τ的粗粒序列。

(2)計算樣本熵

①由步驟(1)所得的粗?；蛄薪M成m維向量

Y(i)=(yi(τ),yi+1(τ),…，yi+m-1(τ))
(1≤i≤N-m)

(6)

②計算Y(i)、Y(j)間的距離

(7)

(8)

(9)

④將維數(shù)增至m+1，重復(fù)①～③，得到Bm+1(r)。

⑤此序列的樣本熵為

(10)

當(dāng)N取有限值時，式(10)表示為

(11)

利用式(6)～(11)計算每個尺度序列的SampnEn，即可得到

MSE(X)=SampnEn(y(τ),m,r)

(12)

由以上計算過程可知，MSE由τ、m和r3個參數(shù)決定。經(jīng)過多次實驗分析，選取τmax=5，m=2，r=0.15δ(其中δ為X的標(biāo)準(zhǔn)差)時，光譜奇異信息規(guī)律性最顯著，能夠較好地反映光譜信號中的弱信息特征。

1.3 VMD-MSE模型

基于不同濃度Cu2+脅迫梯度下的玉米葉片光譜形態(tài)仍極為相似的特點，難以利用其進(jìn)行有效的污染診斷，提出了一種VMD-MSE光譜特征提取及定量描述模型：將各脅迫梯度下的玉米葉片光譜經(jīng)VMD分解，選取能夠表征光譜奇異特征的模態(tài)函數(shù)uk，計算其多尺度熵，實現(xiàn)光譜變異特征的定量描述。算法流程如圖1所示。

圖1 VMD-MSE模型流程Fig.1 Flow chart of VMD-MSE model

2 實驗數(shù)據(jù)及模型應(yīng)用

2.1 實驗數(shù)據(jù)獲取

(1)植株培養(yǎng)：選用有底漏的盆缽進(jìn)行“中糯1號”玉米種子培育。設(shè)置脅迫濃度為0(ck)、100、200、300、400、700、800、900 μg/g的CuSO4·5H2O溶液，將其翻土加入到玉米實驗盆缽中，每梯度設(shè)置3組平行實驗，共24盆。2016年5月6日對玉米種子進(jìn)行催芽處理，出苗后向盆栽中添加NH4NO3、KH2PO4和KNO3營養(yǎng)液。玉米培育期間定期通風(fēng)與澆水，保持適宜的培育溫度與濕度。

(2)光譜數(shù)據(jù)采集：2016年7月17日對玉米葉片反射光譜進(jìn)行測量。在50 W鹵素?zé)艄庠凑丈湎拢瑢⒐庾V范圍為350～2 500 nm的SVC HR-1024I型地物光譜儀探頭視場角設(shè)置為4°并垂直于葉片進(jìn)行光譜采集，采集的光譜使用平面白板進(jìn)行標(biāo)準(zhǔn)化。選取每盆玉米植株的老、中、新3種代表性葉片進(jìn)行光譜測試，每盆獲得3組數(shù)據(jù)。各脅迫梯度光譜由3組平行實驗的9條光譜數(shù)據(jù)求均值所得。各脅迫梯度下的玉米葉片光譜曲線如圖2所示。

圖2 不同濃度Cu2+脅迫下玉米葉片光譜數(shù)據(jù)Fig.2 Spectral data of corn leaves stressed by different Cu2+ concentrations

(3)Cu2+含量測定：2016年9月16日對采集過光譜數(shù)據(jù)的葉片進(jìn)行沖洗、干燥、粉碎等樣品預(yù)處理，再經(jīng)高純硝酸和高氯酸消化處理后用WFX-120型原子吸收分光光度計對葉片中Cu2+含量進(jìn)行測定，每梯度測量3次后取平均值作為該梯度葉片中的Cu2+含量，測量結(jié)果如表1所示。隨著脅迫梯度的增大，葉片中Cu2+含量的總體趨勢表現(xiàn)為先增大后減小，原因可能是隨著土壤中Cu2+增加，玉米對其吸收逐漸增大；而當(dāng)土壤中Cu2+濃度過高時，對植被根系產(chǎn)生毒害作用，造成根部受損，對Cu2+吸收逐漸降低。本文選取脅迫梯度為0、100、200、300、400 μg/g的5組Cu2+含量數(shù)據(jù)用于建立回歸模型，脅迫梯度為700、800、900 μg/g的3組數(shù)據(jù)用于模型驗證。

2.2 VMD-MSE光譜弱信息探測

以0、100 μg/g脅迫濃度下的玉米葉片光譜數(shù)據(jù)為例，分析VMD-MSE模型的特征提取與描述

表1 不同脅迫濃度下的玉米葉片中Cu2+含量Tab.1 Cu2+ contents in corn leaves stressed by different copper concentrations μg/g

過程：

(1)對葉片光譜進(jìn)行分解模態(tài)數(shù)K為2～7的VMD處理，獲得第K模態(tài)函數(shù)uk，如圖3所示。由植被反射光譜特性曲線可知，在960 nm和1 900 nm附近，植被光譜受水分和大氣的影響，出現(xiàn)較大干擾信息。由圖3可知，光譜信號噪聲隨著分解次數(shù)的增大而增大。K=3時，960 nm和1900 nm處噪聲較弱，且全波段奇異信息顯著；從K=4開始，1 000 nm和1 900 nm處噪聲開始出現(xiàn)，且奇異信息逐漸弱化。因此選擇K=3既實現(xiàn)模態(tài)分離，又有效抑制噪聲，可有效提取玉米銅污染脅迫信息。

圖3 VMD分解結(jié)果Fig.3 Results of decomposition of VMD

(2)計算由VMD分解得到的第3模態(tài)函數(shù)u3的多尺度熵，即為VMD-MSE的模型值(記為VM)，不同尺度下的模型值記為VMi(i=1，2，…，τmax)。

3 結(jié)果與分析

3.1 VMD-MSE污染監(jiān)測

對0、100、200、300、400 μg/g脅迫濃度下的玉米葉片光譜應(yīng)用VMD-MSE模型進(jìn)行污染信息提取與測度。各脅迫濃度下的玉米葉片光譜的VM值如圖4所示。不同脅迫濃度下，VM曲線變化趨勢有所差異，但總體而言，曲線變化趨勢相似，即隨著τ的增大而增大。在所有尺度因子上，VM值與作物受污染程度呈負(fù)相關(guān)關(guān)系，即受污染越嚴(yán)重，VM值越小。顯然，各脅迫濃度下的光譜細(xì)節(jié)特征通過不同尺度VM(即VMi)得到體現(xiàn)，能更好地揭示光譜曲線畸變特征。

圖4 不同脅迫濃度下玉米葉片光譜的VM值Fig.4 VM values of corn leaf spectrum under different stress concentrations

圖5 不同脅迫梯度下各尺度VM值與玉米葉片中Cu2+含量相關(guān)關(guān)系Fig.5 Relationships between VM on five scales and Cu2+ contents in corn leaves stressed by different Cu2+ concentrations

不同脅迫濃度下，VMi與葉片中Cu2+含量的相關(guān)變化如圖5所示，兩者呈現(xiàn)相反的變化趨勢。隨著脅迫濃度的增加，葉片中Cu2+含量先增加后減小，相應(yīng)地VMi表現(xiàn)為先減小后增大的趨勢，兩者相關(guān)系數(shù)絕對值均大于0.90，其中VM1與Cu2+相關(guān)性最強，相關(guān)系數(shù)為-0.968。據(jù)此可認(rèn)為，各尺度VM均可在一定程度上診斷玉米受銅脅迫程度，其中VM1效果最優(yōu)。

3.2 VM值與玉米銅含量的回歸分析

為了進(jìn)一步確定VMD-MSE模型值與玉米銅脅迫的定量關(guān)系，以不同脅迫濃度下玉米葉片中Cu2+含量為因變量y，以VMi為自變量x進(jìn)行回歸分析，VMi值與銅含量的線性回歸方程、決定系數(shù)(R2)和均方根誤差(RMSE)如表2所示。VM1～VM5的回歸決定系數(shù)均大于0.88，且顯著性水平P<0.05。整體上，從VM1到VM4，R2逐漸減小，RMSE逐漸增大，這可能是由于隨時間序列復(fù)雜度的增加，樣本熵估計準(zhǔn)確率下降所致。

為了驗證上述模型的可靠性與準(zhǔn)確度，選用700、800、900 μg/g 3個脅迫濃度下的樣本數(shù)據(jù)對模型進(jìn)行檢驗，圖6描述了玉米葉片中Cu2+含量實測值與模型預(yù)測值的相關(guān)變化。各預(yù)測模型的Cu2+含量預(yù)測值與實測數(shù)據(jù)具有很好一致性。VM1和VM2回歸模型的預(yù)測結(jié)果較優(yōu)，相關(guān)系數(shù)大于0.70，其中VM1模型預(yù)測結(jié)果最優(yōu)，相關(guān)系數(shù)為0.992；VM3～VM5回歸模型預(yù)測效果較差，相關(guān)系數(shù)小于0.650，可能是由于隨著污染程度的增加，光譜信號的奇異信息隨機性增加，且隨度量尺度的增大，熵估計準(zhǔn)確率下降，使VM3～VM5模型穩(wěn)定性及適用性變差。因此，通過VMD-MSE模型值VM1可以建立玉米Cu污染程度的預(yù)測模型，從而實現(xiàn)監(jiān)測玉米Cu污染狀況的目的。

表2 基于不同尺度VM值的玉米葉片中Cu2+含量回歸模型Tab.2 Regression models for VM on different scales with Cu2+ content in corn leaves

圖6 玉米葉片中Cu2+含量實測值與預(yù)測值的相關(guān)關(guān)系Fig.6 Relationships between measured and predicted values of Cu2+ contents in corn leaves

4 結(jié)論

(1)將VMD運用到高光譜信息提取中，從相似光譜信號中挖掘有利于診斷重金屬Cu污染的畸變信息。結(jié)果表明，對原始光譜數(shù)據(jù)進(jìn)行3次VMD分解后，可有效抑制噪聲，實現(xiàn)模態(tài)分離，從而提取出隱藏在光譜噪聲中的污染弱信息。

(2)將MSE方法應(yīng)用于VMD分解結(jié)果，揭示了污染信息不同尺度的內(nèi)在差異。

(3)基于VMD和MSE構(gòu)建的VMD-MSE模型可有效診斷作物銅污染程度，模型值VM與玉米葉片中Cu2+含量呈現(xiàn)顯著負(fù)相關(guān)，兩者回歸精度達(dá)到0.88以上，其中VM1模型預(yù)測效果最優(yōu)。本研究構(gòu)建的預(yù)測模型能在一定程度上預(yù)測玉米受污染程度，但具有一定局限性，即預(yù)測精度隨污染程度增大而降低，這種情況可通過增加樣本數(shù)量對模型進(jìn)行優(yōu)化而得到改善。