近年来,由于人口急剧增长、工业迅猛发展,大量未经处理的工业固体废弃物(主要来自采掘业、有色金属冶炼、化学原料及化学制品等)在土壤表面堆放和倾倒,有害废水向土壤中渗透,农药的不科学使用以及大气中的有害气体及飘尘随雨水降落到土壤等,导致土壤重金属铜(Cu)、铅(Pb)、汞(Hg)、镉(Cd)等污染问题日趋严重[1-2]。重金属在土壤中不能被降解,其超限后会污染农作物并严重影响农产品品质,且可以通过食物链进入人体,破坏人体正常的生理功能,进而对人体造成危害[3-5]。因此,如何快速准确鉴别作物受重金属污染的元素类别以及监测其所受污染程度已成为目前的研究热点。常规的生化检测方法主要是通过采集作物样本来检测其中重金属元素含量、浓度及一些生化特征参数等,但这种方法存在费时费力、工作量大等缺点,且无法满足无损、快速、便捷等污染监测需求。高光谱遥感具有光谱信息丰富且光谱连续、波段多且范围窄、从可见光到近红外有数十甚至数百个波段、理化特性与参量反演能力强等特点,因此已成为现今环境、农作物等重金属污染监测的热点研究内容之一[6-7]。
重金属在作物叶片内积累会导致叶细胞畸形,使细胞质膜的选择透性、组成和结构均发生改变,而细胞结构改变、成分分解会使光在叶片内传输时出现不连续的折射,从而引起整个叶片到冠层的光谱变化。研究发现,重金属胁迫下作物的电磁波反射光谱变异能反映受重金属污染后植物的色素、细胞结构和含水量等生化参数的特征信息,并能响应出重金属污染元素类别[8]。因此,研究作物光谱变化特征和污染信息提取的光谱技术至关重要。基于AR模型的功率谱估计是现代功率谱估计中最常用的方法,而利用AR模型时必须计算出AR模型的参数,其中Burg算法被广泛应用。近年来,Burg算法被应用在很多领域,如罗治军等[9]将Burg算法运用到经验小波变换(Empirical wavelet transform,EWT)中,使得EWT算法能在功率谱上找到明确的模态边界;王路晗等[10]用Burg算法对风切边雷达信号进行处理;李佳等[11]利用Burg算法获得车流量的雷达信号信息,对双波束车流量的速度进行检测;钱莹晶等[12]采用Burg算法较为准确地估计出脉冲编码雷达中目标物体的多普勒频偏;Iman等[13]将Burg算法用于测量患者心率等;但Burg算法在高光谱信号处理与分析方面鲜有报道。本文尝试将Burg算法用于玉米叶片的高光谱数据处理,以甄别玉米受重金属污染的元素类别和所受污染程度。
本文将以不同浓度Cu和Pb污染胁迫下的盆栽玉米及其植株叶片光谱作为研究对象,采用Burg算法绘制玉米叶片光谱的功率谱曲线。考虑到该算法适应于短序列,实验将首先选取玉米的紫谷、绿峰和红边3个特征区间,并对特征区间的高光谱数据进行一阶光谱微分(Spectral derivative,SD)和奇异值分解(Singular value decomposition,SVD)处理,以消除噪声等。然后结合Burg算法对预处理后的光谱信号进行功率谱估计,通过曲线特征对受重金属Cu和Pb胁迫生长的玉米叶片进行污染元素种类辨别和污染程度诊断,同时基于不同年份同等条件的Cu、Pb污染玉米叶片光谱数据检验该模型的稳定性与可行性。
1 材料与方法 1.1 实验设计 1.1.1 玉米植株培养实验选用有底漏的花盆培植“密糯8号”玉米种子,重金属Cu和Pb以不同浓度CuSO4∙5H2O和Pb(NO3)2溶液的形式翻土加入到花盆中,待玉米出苗后向盆中喷洒适量且等量的NH4NO3、KH2PO4和KNO3营养液。实验设计Cu、Pb污染浓度为0、50、200 mg·kg-1和400 mg·kg-1,记为CK(0)、Cu(50)、Pb(50)、Cu(200)、Pb(200)、Cu(400)、Pb(400),每个浓度设置3组平行试验,玉米培育期间保持各盆栽生长条件一致。
1.1.2 玉米叶片光谱测量使用SVC HR-1024I型地物光谱仪采集玉米叶片反射光谱。采集时选择50 W卤素灯作为光源,探头视场角设置为4°,并垂直于玉米叶片表面5 cm处,采集输出的光谱反射率并用平面板对光谱进行标准化。测量时每株采摘3片代表性叶片[玉米顶部刚长出不久的新叶片(New,N),底部完全成熟的叶片(Old,O)和位于中部的叶片(Middle,M)]分别进行光谱测量,测量过程中将玉米叶片放置于黑色塑料袋包裹的硬纸片上以排除其他因素的干扰。每盆获得3组数据,去除异常光谱值后取平均值作为研究数据。
1.1.3 Cu和Pb含量测定收集采集光谱数据后的玉米叶片,对其进行冲洗(去除灰尘杂质)、干燥(烘干水分),用高纯硝酸、高氯酸消化(排除其他离子的干扰并使Cu、Pb以离子的形式存在)等预处理,Cu和Pb含量的测定依据《硅酸岩石化学分析方法》(GB/T 14506.30—2010)。叶片中Cu含量的测定采用电感耦合等离子体发射光谱仪进行定量分析,叶片中Pb含量的测定采用WFX-120原子吸收分光光度计[14]。对不同浓度Cu、Pb处理下玉米叶片中的老(O)、中(M)、新(N)3种叶片分别测量3次后取平均值作为该浓度下叶片中Cu、Pb的含量。实验所测得的玉米叶片中Cu和Pb的含量如表 1所示。
光谱微分(SD)变换是一种常用的光谱特征分析方法,光谱微分变换可在一定程度上限制低频噪声对目标光谱的影响,并能有效削弱大气效应等的干扰[15],进而增强光谱曲线的细微变化。光谱一阶微分的计算公式为:
(1) |
式中:λi-1、λi和λi+1为玉米叶片光谱数据相邻波段的波长;ρ′(λi)是波长λi的一阶微分值;ρ(λi+1)和ρ(λi-1)为波长λi+1和λi-1对应的玉米叶片光谱反射率。
1.3 奇异值分解奇异值分解(SVD)是线性代数中一种矩阵分解技术,在信号处理方面有重要应用。SVD分解是谱分析理论在任意矩阵上的推广[16]。
假设原矩阵为Am×n(m < n),经过奇异值分解后为:
(2) |
式中:Um×m和Vn×n为单位正交矩阵,其对应的列向量和行向量为单位向量且互相正交,都是矩阵A对应的特征向量;Σm×n为对角矩阵,Σm×n=diag(σ1,σ2,σ3,…,σm),其元素σi为按降序排列矩阵A的奇异值,它的对角元素即为Am×n的奇异值,Σm×n中除对角线外,其余元素均为0。相关研究表明,矩阵Σm×n中的奇异值即为原始数据集Am×n中的全部重要特征,其他特征大多为多余的或是噪声。
就光谱信号的SVD处理而言,设实验所测得的玉米叶片光谱数据矩阵为矩阵A,则矩阵A中存在m× m的正交矩阵U=[u1,u2,…,um]和n×n的正交矩阵V= [v1,v2,…,vn],则:
(3) |
式中:σi为玉米光谱数据A的奇异值;ui为相应于奇异值σi的左奇异向量;vi为相应于奇异值σi的右奇异向量,且满足:
(4) |
(5) |
式中:i=1,2,…,m;λi是方阵AAT的特征值;λi′是方阵ATA的特征值。根据式(4)和式(5)可求得矩阵U和V。
(6) |
(7) |
(8) |
根据式(6)、式(7)、式(8)得到玉米光谱数据A奇异值σi为:
(9) |
AR模型是现代谱估计中最常用的一种功率谱估计模型,而AR模型在功率谱估计时,必须计算出模型的参数。AR模型参数求解以线性预测理论为基础,发展出许多算法,主要有Yule⁃Walker法、协方差法、Burg法和改进的协方差法,Burg算法以其良好的分辨率和稳定性得到广泛应用[11]。
Burg算法为一种递推算法,不需要计算序列的自相关,是直接通过原始序列求出模型参数,进而求得功率谱。其基本思想是根据线性预测器的前向和后向预测功率的总均方误差之和最小的准则直接从观测数据来估计反射系数,然后利用Lenvinson - Durbin算法的递推公式求出AR模型参数[17]。其具体算法如下:
(1)定义前向和后向误差分别为ef(n)和eb(n),预测误差功率为E。初始化前向、后向误差及预测误差功率:
(10) |
式中:x(n)为玉米叶片光谱信号序列,n=0,1,…,N-1。
(11) |
(2)计算反射系数:
(12) |
式中:m=1,2,…,p,p为模型的阶数;efm-1 (n)表示m-1阶前向预测误差;ebm - 1 (n)为m-1阶后向预测误差。
(3)根据Lenvinson递推算法计算模型参数及预测误差功率:
(13) |
式中:k=1,2,…,m-1;am(k)表示第m阶模型的第k个系数。
(14) |
(15) |
(4)递推高一阶前、后向预测误差:
(16) |
(17) |
把m更新为m+1,重复步骤(2)至步骤(4),直到预测误差功率Em满足要求。
(5)得到参数后,通过公式(18)计算玉米叶片光谱信号的功率谱密度为:
(18) |
式中:P为功率谱密度;ak为AR模型参数;σ2为序列的方差。
2 结果与讨论 2.1 SD-SVD-Burg算法的应用与讨论实验所获得的Cu和Pb胁迫下玉米叶片光谱如图 1所示。从曲线看出,玉米叶片的反射光谱曲线整体呈现先上升后下降的趋势,在波长700~2 500 nm间健康叶片的反射波峰强度高于受重金属元素胁迫的叶片光谱,但受Cu、Pb胁迫叶片的光谱曲线差异较小,仅从光谱曲线无法辨别玉米受重金属元素的污染类别与污染程度。根据植物光谱的光谱特性,将Cu和Pb胁迫下玉米叶片光谱划分为紫谷、绿峰和红边3个特征区间作为研究对象,各特征区间对应的范围如表 2所示。波长380~450 nm出现一个吸收谷即紫谷,主要由栅栏组织中的光合色素强烈吸收所引起,由于相对地反射绿光,因此健康植物的视觉效果是绿色的。波长500~600 nm出现一个小反射峰即绿峰,是描述色素状态和健康状况的重要指示波段。波长670~780 nm形成陡坡即红边,是植物反射光谱中最具特征的指示波段。重金属会对植物的色素、细胞结构和含水量等生化参数产生影响,因此选择特征波段来反映健康植物和受重金属污染的植物是可行的。
实验分别对Cu与Pb胁迫梯度下各特征区间的光谱经式(1)处理,得到一阶微分光谱,然后分别对各特征区间的微分光谱进行SVD处理,获得各特征区间的光谱奇异值。因Burg算法适应于短序列信号,经光谱一阶微分以及SVD处理后数据在去除噪声的同时有效压缩了数据。然后利用Burg算法绘制功率谱曲线,首先输入序列x(n),即各特征区间经预处理后的玉米叶片光谱以矩阵形式作为一个信号序列,接着设置模型阶数,阶数是Burg算法的重要参数,使用不同的阶数对算法的结果影响也不同。因此,实验分别选择阶数p为3、4、5进行仿真,通过Lenvinson递推算法求得模型参数,利用式(18)求得功率谱密度。最终健康玉米叶片及Cu、Pb胁迫下紫谷、绿峰、红边的信号处理结果如图 2~图 5所示。
分析图 2和图 3可知,在紫谷区间,健康玉米叶片的功率谱曲线在阶数p=3时,曲线较平缓,无波峰波谷,而受Cu、Pb污染的玉米叶片均出现1个波峰,由此可定性分析玉米叶片是否已受Cu或Pb污染。从图 3看出,在p=5时,受Cu污染的曲线表现为单波峰,而受Pb污染的曲线表现为双波峰,因此通过此特征可快速定性分析玉米所受重金属污染的元素类别(Cu或Pb)。随着Cu和Pb胁迫浓度的增加,受Cu或Pb污染的曲线波峰均出现左移现象,因此可通过曲线波峰的位置来判断玉米叶片受重金属离子的胁迫程度。在p=4时,受Cu污染的曲线波峰较陡且随着胁迫浓度的增加,波峰逐渐有变缓的趋势,而受Pb污染的曲线波峰较缓,随着胁迫浓度的增加波峰逐渐变缓,因此在p=4时也可通过波峰坡度来判断玉米叶片所受重金属污染的类别及污染程度。
分析图 2和图 4可知,在绿峰区间,健康玉米叶片在p=3时,曲线表现为单波峰,而在Cu、Pb胁迫下曲线均出现3个波峰,因此可通过波峰个数来判断玉米叶片是否受到重金属元素污染。从图 4可知,当p=5时,受Cu污染的曲线出现4个波峰,受Pb污染的曲线为3个波峰,由此可辨别玉米所受重金属污染的元素类别(Cu或Pb)。当p=4或5时,观察到Cu或Pb胁迫浓度越大,波峰越陡,因此可根据波峰的坡度来粗略诊断玉米受重金属离子的污染程度。
分析图 2和图 5可知,在红边区间,健康玉米叶片在p=5时,仅出现一个波峰,而叶片受Cu、Pb胁迫下的曲线有2~3个波峰,可由此来诊断玉米叶片是否受到重金属元素的污染。从图 5看出,受Cu与受Pb胁迫下的功率谱曲线未表现出明显差异,无法通过曲线特征来甄别玉米所受重金属元素的污染类别;随着重金属污染程度的增加,曲线的波峰有变缓的趋势,因此可根据波峰的坡度来判断玉米所受重金属离子的污染程度。
由图 2~图 5可知,在紫谷和绿峰区间通过功率谱密度曲线图,均能较好地定性分析玉米是否受Cu、Pb污染并能判断出受哪种元素污染,利用曲线的波峰变化能诊断玉米所受重金属离子的污染程度;在红边区间,可辨别出玉米是否受到重金属元素污染以及所受污染程度,但未能诊断出污染元素的类别。
为进一步从定量的角度去检测紫谷、绿峰和红边区间的功率谱曲线对重金属Cu与Pb胁迫程度的监测效果,分别统计了在阶数p=5时,不同胁迫浓度下功率谱曲线的平均功率与玉米叶片中Cu、Pb含量的相关系数,结果如表 3所示。由表可知,Cu胁迫下紫谷、绿峰及红边3个区间的功率谱曲线平均功率与玉米叶片中Cu含量的相关系数均在0.9以上,表现出较强的相关性,说明对重金属Cu胁迫程度的监测达到很好的效果;Pb胁迫下在绿峰区间功率谱密度曲线平均功率与玉米叶片中Pb含量的相关系数达到0.995 8,效果理想。
为验证SD-SVD-Burg算法在监测重金属污染信息方面的稳定性和普适性,选取不同年份同等实验条件的Cu、Pb污染玉米叶片光谱数据作为验证组对该模型进行可行性检验。根据表 3,选取与Cu和Pb相关性最高的绿峰区间作为研究对象,验证组统计的Cu与Pb胁迫下绿峰功率谱曲线平均功率与玉米叶片中Cu与Pb含量的相关系数如表 4所示。由表 4可知,在绿峰区间,功率谱曲线的平均功率与玉米叶片中Cu含量相关系数为0.921 3,与Pb含量的相关系数达0.991 5,均表现出较强的相关性,因此该算法在区分作物重金属污染类别和监测污染程度方面具有一定的有效性和可行性。
(1)在紫谷与绿峰区间,健康玉米叶片与受Cu、Pb污染的叶片波峰个数均不相同,因此可根据曲线的波峰个数定性分析玉米是否受到污染,并能判断出玉米所受重金属污染的元素类别(Cu或Pb)。利用曲线的波峰位置及波峰坡度变换可定性诊断玉米所受重金属元素的污染程度。
(2)通过统计不同胁迫浓度下功率谱曲线的平均功率与玉米叶片中Cu、Pb含量的相关系数,得出绿峰曲线的平均功率可较好地定量分析Cu、Pb重金属污染程度。利用不同年份的玉米叶片光谱数据对该算法进行稳定性验证,表明SD-SVD-Burg算法可以作为监测玉米是否受Cu、Pb污染的一种参考方法。
[1] |
Battsengel E, Murayama T, Fukushi K, et al. Ecological and human health risk assessment of heavy metal pollution in the soil of the ger district in Ulaanbaatar, Mongolia[J]. International Journal of Environmental Research and Public Health, 2020, 17(13): 4668. DOI:10.3390/ijerph17134668 |
[2] |
任加国, 王彬, 师华定, 等. 沱江上源支流土壤重金属污染空间相关性及变异解析[J]. 农业环境科学学报, 2020, 39(3): 530-541. REN Jia-guo, WANG Bin, SHI Hua-ding, et al. Spatial correlation and variation analysis of soil heavy metals contamination in upper source tribu butary of Tuojiang River, China[J]. Journal of Agro-Environmental Science, 2020, 39(3): 530-541. |
[3] |
陈文轩, 李茜, 王珍, 等. 中国农田土壤重金属空间分布特征及污染评价[J]. 环境科学, 2020, 41(6): 2822-2833. CHEN Wen-xuan, LI Qian, WANG Zhen, et al. Spatial distribution characteristics and pollution evaluation of heavy metals in arable land soil of China[J]. Environmental Sciences, 2020, 41(6): 2822-2833. |
[4] |
任传义, 程军勇, 张延平, 等. 竹笋地土壤重金属污染潜在生态风险及食用笋健康风险评价[J]. 农业环境科学学报, 2017, 36(5): 855-862. REN Chuan-yi, CHENG Jun-yong, ZHANG Yan-ping, et al. Assessment of heavy metals potential ecological hazards of soil and health risk of bamboo shoots[J]. Journal of Agro-Environmental Science, 2017, 36(5): 855-862. |
[5] |
张云霞, 周浪, 肖乃川, 等. 鬼针草(Bidens pilosa L.)对镉污染农田的修复潜力[J]. 生态学报, 2020, 40(16): 5805-5813. ZHANG Yun-Xia, ZHOU Lang, XIAO Nai-chuan, et al. Remediation potential of B. pilosa L. in cadmium-contaminated farmland[J]. Acta Ecologica Sinica, 2020, 40(16): 5805-5813. |
[6] |
杨可明, 张伟, 付萍杰, 等. 铜污染土壤光谱特征及弱差信息的LHPSD甄测模型[J]. 光谱学与光谱分析, 2019, 39(7): 2228-2236. YANG Ke-ming, ZHANG Wei, FU Ping-jie, et al. The LH-PSD analysis model of Cu contaminated soil spectral characteristics and weak characteristic information[J]. Spectroscopy and Spectral Analysis, 2019, 39(7): 2228-2236. |
[7] |
陈功伟, 赵思颖, 倪才英. 高光谱监测技术在重金属污染土壤上的应用[J]. 中国科学院大学学报, 2019, 36(4): 560-566. CHEN Gong-wei, ZHAO Si-ying, NI Cai-ying. Hyperspectral monitoring of soil contaminated by heavy metals[J]. Journal of the University of Chinese Academy of Sciences, 2019, 36(4): 560-566. |
[8] |
史超, 黄超, 李书, 等. 矾山铜钼矿区重金属与覆盖植被胁迫相关性光谱特征[J]. 地质科技通报, 2020, 39(3): 202-210. SHI Chao, HUANG Chao, LI Shu, et al. Spectral characteristics and correlation of heavy metal and vegetation cover stress in Fanshan copper-molybdenum[J]. Geological Science and Technology Bulletin, 2020, 39(3): 202-210. |
[9] |
罗治军, 刘涛, 张玉玲, 等. 一种改进型AR-EWT算法及其在振动模态分解中的应用研究[J]. 机械设计, 2019, 36(增刊1): 1-6. LUO Zhi-jun, LIU Tao, ZHANG Yu-ling, et al. A revised AR-EWT and its application on the vibration mode decomposition[J]. Journal of Machine Design, 2019, 36(Suppl1): 1-6. |
[10] |
王路晗, 张兴敢. 基于Burg改进算法的风切变风速估计[J]. 南京大学学报(自然科学), 2016, 52(5): 925-931. WANG Lu-han, ZHANG Xing-gan. Speed estimation of shear wind based on the improved Burg algorithm[J]. Journal of Nanjing University(Natural Science), 2016, 52(5): 925-931. |
[11] |
李佳, 曹林, 王东峰, 等. 基于小波降噪的雷达时频信号互相关测速算法[J]. 现代电子技术, 2020, 43(9): 14-20. LI Jia, CAO Lin, WANG Dong-feng, et al. Wavelet denoising based cross-correlation velocity detection algorithm of radar time-frequency signal[J]. Modern Electronic Technology, 2020, 43(9): 14-20. |
[12] |
钱莹晶, 周群. Burg频谱估计算法的硬件加速方法研究[J]. 电子测量与仪器学报, 2015, 29(9): 1382-1390. QIAN Ying-jing, ZHOU Qun. Research on hardware acceleration method for Burg spectrum estimation algorithm[J]. Journal of Electronic Measurement and Instrumentation, 2015, 29(9): 1382-1390. |
[13] |
Iman R, Yeon-Mo Y, Ki M. Blood volume pulse extraction for Noncontact heart rate measurement by digital camera using singular value decomposition and Burg algorithm[J]. Energies, 2018, 11(5): 1076. DOI:10.3390/en11051076 |
[14] |
丘雪红, 曹莉, 韩日畴. 冬虫夏草与产地植物及土壤的重金属含量测定[J]. 环境昆虫学报, 2020, 42(2): 267-273. QIU Xue-hong, CAO Li, HAN Ri-chou. Determination of heavy metals in Ophiocordyceps sinensis and its habitat plants and habitat soil[J]. Journal of Environmental Entomology, 2020, 42(2): 267-273. |
[15] |
徐凯健, 田庆久, 徐念旭, 等. 基于时序NDVI与光谱微分变换的森林优势树种识别[J]. 光谱学与光谱分析, 2019, 39(12): 3794-3800. XU Kai-jian, TIAN Qing-jiu, XU Nian-xu, et al. Classifying forest dominant trees species based on high dimensional Time-Series NDVI data and differential transform methods[J]. Spectroscopy and Spectral Analysis, 2019, 39(12): 3794-3800. |
[16] |
郭明军, 李伟光, 杨期江, 等. 基于SVD原理的PCA特征频率提取算法及其应用[J]. 华南理工大学学报(自然科学版), 2020, 48(1): 1-9. GUO Ming-jun, LI Wei-guang, YANG Qi-jiang, et al. PCA feature frequency extraction algorithm based on SVD principle and its application[J]. Journal of South China University of Technology(Natural Science edition), 2020, 48(1): 1-9. |
[17] |
Gal A K, Maria F, Paolo B, et al. Burg matrix divergence based MultiMetric learning[J]. Frontiers in Artificial Intelligence and Applications, 2016, 285. |