近年来，由于人口急剧增长、工业迅猛发展，大量未经处理的工业固体废弃物（主要来自采掘业、有色金属冶炼、化学原料及化学制品等）在土壤表面堆放和倾倒，有害废水向土壤中渗透，农药的不科学使用以及大气中的有害气体及飘尘随雨水降落到土壤等，导致土壤重金属铜（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以不同浓度CuSO₄∙5H₂O和Pb（NO₃）₂溶液的形式翻土加入到花盆中，待玉米出苗后向盆中喷洒适量且等量的NH₄NO₃、KH₂PO₄和KNO₃营养液。实验设计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所示。

1.2 一阶光谱微分

光谱微分（SD）变换是一种常用的光谱特征分析方法，光谱微分变换可在一定程度上限制低频噪声对目标光谱的影响，并能有效削弱大气效应等的干扰^[15]，进而增强光谱曲线的细微变化。光谱一阶微分的计算公式为：

(1)

式中：λ_i-1、λ_i和λ_i+1为玉米叶片光谱数据相邻波段的波长；ρ′（λ_i）是波长λ_i的一阶微分值；ρ（λ_i+1）和ρ（λ_i-1）为波长λ_i+1和λ_i-1对应的玉米叶片光谱反射率。

1.3 奇异值分解

奇异值分解（SVD）是线性代数中一种矩阵分解技术，在信号处理方面有重要应用。SVD分解是谱分析理论在任意矩阵上的推广^[16]。

假设原矩阵为A_m×n（m < n），经过奇异值分解后为：

(2)

式中：U_m×m和V_n×n为单位正交矩阵，其对应的列向量和行向量为单位向量且互相正交，都是矩阵A对应的特征向量；Σ_m×n为对角矩阵，Σ_m×n=diag（σ₁，σ₂，σ₃，…，σ_m），其元素σ_i为按降序排列矩阵A的奇异值，它的对角元素即为A_m×n的奇异值，Σ_m×n中除对角线外，其余元素均为0。相关研究表明，矩阵Σ_m×n中的奇异值即为原始数据集A_m×n中的全部重要特征，其他特征大多为多余的或是噪声。

就光谱信号的SVD处理而言，设实验所测得的玉米叶片光谱数据矩阵为矩阵A，则矩阵A中存在m× m的正交矩阵U=[u₁，u₂，…，u_m]和n×n的正交矩阵V= [v₁，v₂，…，v_n]，则：

(3)

式中：σ_i为玉米光谱数据A的奇异值；u_i为相应于奇异值σ_i的左奇异向量；v_i为相应于奇异值σ_i的右奇异向量，且满足：

(4)

(5)

式中：i=1，2，…，m；λ_i是方阵AA^T的特征值；λ_i′是方阵A^TA的特征值。根据式（4）和式（5）可求得矩阵U和V。

(6)

(7)

(8)

根据式（6）、式（7）、式（8）得到玉米光谱数据A奇异值σ_i为：

(9)

1.4 基于AR模型的Burg算法

AR模型是现代谱估计中最常用的一种功率谱估计模型，而AR模型在功率谱估计时，必须计算出模型的参数。AR模型参数求解以线性预测理论为基础，发展出许多算法，主要有Yule⁃Walker法、协方差法、Burg法和改进的协方差法，Burg算法以其良好的分辨率和稳定性得到广泛应用^[11]。

Burg算法为一种递推算法，不需要计算序列的自相关，是直接通过原始序列求出模型参数，进而求得功率谱。其基本思想是根据线性预测器的前向和后向预测功率的总均方误差之和最小的准则直接从观测数据来估计反射系数，然后利用Lenvinson - Durbin算法的递推公式求出AR模型参数^[17]。其具体算法如下：

（1）定义前向和后向误差分别为e^f（n）和e^b（n），预测误差功率为E。初始化前向、后向误差及预测误差功率：

(10)

式中：x（n）为玉米叶片光谱信号序列，n=0，1，…，N-1。

(11)

（2）计算反射系数：

(12)

式中：m=1，2，…，p，p为模型的阶数；e^f_m-1 (n)表示m-1阶前向预测误差；e^b_{m - 1} (n)为m-1阶后向预测误差。

（3）根据Lenvinson递推算法计算模型参数及预测误差功率：

(13)

式中：k=1，2，…，m-1；a_m（k）表示第m阶模型的第k个系数。

(14)

(15)

（4）递推高一阶前、后向预测误差：

(16)

(17)

把m更新为m+1，重复步骤（2）至步骤（4），直到预测误差功率E_m满足要求。

（5）得到参数后，通过公式（18）计算玉米叶片光谱信号的功率谱密度为：

(18)

式中：P为功率谱密度；a_k为AR模型参数；σ²为序列的方差。

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，效果理想。

2.2 算法稳定性验证与讨论

为验证SD-SVD-Burg算法在监测重金属污染信息方面的稳定性和普适性，选取不同年份同等实验条件的Cu、Pb污染玉米叶片光谱数据作为验证组对该模型进行可行性检验。根据表 3，选取与Cu和Pb相关性最高的绿峰区间作为研究对象，验证组统计的Cu与Pb胁迫下绿峰功率谱曲线平均功率与玉米叶片中Cu与Pb含量的相关系数如表 4所示。由表 4可知，在绿峰区间，功率谱曲线的平均功率与玉米叶片中Cu含量相关系数为0.921 3，与Pb含量的相关系数达0.991 5，均表现出较强的相关性，因此该算法在区分作物重金属污染类别和监测污染程度方面具有一定的有效性和可行性。

3 结论

（1）在紫谷与绿峰区间，健康玉米叶片与受Cu、Pb污染的叶片波峰个数均不相同，因此可根据曲线的波峰个数定性分析玉米是否受到污染，并能判断出玉米所受重金属污染的元素类别（Cu或Pb）。利用曲线的波峰位置及波峰坡度变换可定性诊断玉米所受重金属元素的污染程度。

（2）通过统计不同胁迫浓度下功率谱曲线的平均功率与玉米叶片中Cu、Pb含量的相关系数，得出绿峰曲线的平均功率可较好地定量分析Cu、Pb重金属污染程度。利用不同年份的玉米叶片光谱数据对该算法进行稳定性验证，表明SD-SVD-Burg算法可以作为监测玉米是否受Cu、Pb污染的一种参考方法。