不同pH值土壤中铅含量的太赫兹光谱反演建模研究

李 超，李 斌，张丽琼，叶大鹏，郑书河

1. 福建农林大学机电工程学院，福建 福州 350012 2. 北京农业信息技术研究中心，北京 100097 3. 农业部农业遥感机理与定量遥感重点实验室，北京 100097

引 言

重金属具有生物毒性、 不可降解性及在生物体中的累积性，通过废水、 大气沉降持续释放至周围环境的表层土壤中，不断累积的重金属在土壤中通过有机体转化为有机复合物，最终对自然生态系统和人类健康造成威胁[1-2]。铅作为环境中典型的重金属污染元素之一，在进入土壤后，容易通过溶解、 沉淀、 络合和吸附等方式以5种(可交换态、 碳酸盐结合态、 铁锰氧化物结合态、 有机结合态和残渣态)不同的化学结合态富集于土壤中[3]。因此土壤中铅含量的检测对于防控土壤重金属污染具有重要意义。

目前利用光谱反演土壤重金属含量的研究主要集中在重金属种类、 重金属分布区域等方面。相关研究人员以河滩平原[4]、 矿区[5-6]、 农田[7]为例，利用近红外光谱和高光谱探究了土壤光谱与Cd，Cu，As，Ni和Pb等重金属含量之间的关系，通过光谱反演建立了含量预测模型。尽管上述土壤重金属含量光谱反演建模都取得了不错的预测效果，然而对于利用太赫兹对不同pH土壤中重金属含量反演的研究却少有报道。相关研究[8-9]表明，土壤光谱信息与重金属化学结合态的组分具有一定相关性，而pH值是影响土壤重金属化学结合态的重要因素之一。因此，有必要探索并建立重金属在不同pH土壤中的最佳反演模型，以实现土壤重金属含量的准确预测。

李斌[10]等开展了土壤中铅含量的初步研究，发现土壤样品中铅含量与对应的太赫兹吸收谱之间存在一定的对应关系，因而得出了利用太赫兹光谱技术对土壤重金属含量进行测定具有可行性的结论。以此为理论依据，本文尝试利用太赫兹光谱技术建立不同pH土壤中铅含量的反演预测模型。首先采集土壤样品的太赫兹光谱，并对光谱做一定预处理，然后对预处理后的光谱数据，利用SPA选取光谱的特征频率，最后分别采用PLS，SVM和BPNN建立铅含量的反演预测模型，并通过比较模型效果，确定铅在不同pH土壤中的最佳预测模型。

1 实验部分

1.1 样品制备

制备样品所需的纯净土壤采集于北京市农林科学院院内实验田，采集0～20 cm的表层土壤，约2.5 kg。取部分采集到的土壤样品送至北京中科联宏检测技术公司进行土壤的主要成分检测，根据国家标准《土壤环境质量农用地土壤污染风险管控标准》[11]检测土壤样品中的重金属元素是否已经超标。经检测后，采集的土壤样品中铅含量背景值为29.476 mg·kg-1，符合国家标准。

为模拟铅离子自然状态下进入不同pH的土壤，利用NaOH溶液和HCl溶液分别配制pH为5.5，7.0，8.5的土壤，选用0.1 mol·L-1的Pb(NO3)2标准溶液作为铅源加入土壤中。每种pH条件下按铅含量为50～1 000 mg·kg-1，梯度为50 mg·kg-1配制含铅土壤样品20组，每组2个，共计120个。另配制5个不含铅的土壤样品用作比较和分析铅进入土壤前后的太赫兹光谱曲线变化。含铅土壤经过干燥、 研磨后称取200 mg含铅土壤样品，使用手动液压型压片机(Specac GS15011, 英国)进行压片，压力为3.5 T，压制时间3 min。由于土壤样品本身特性，在压片过程中容易出现松散现象。为解决这一问题，在土壤中加入少量聚乙烯粉末实现样品制备工艺的优化。最终成型的压片，厚度约为1.0～1.1 mm，直径约为13 mm。

1.2 光谱数据采集与预处理

实验采用德国Menlo Systems公司的太赫兹时域光谱系统TERA K15。为避免空气中的水蒸气对太赫兹波的强烈吸收所产生的干扰，将太赫兹波的产生与探测装置置于透明亚克力密闭箱内，并连续充入干燥氮气，保证测试环境相对湿度小于5%，环境温度在23 ℃左右。将压片固定在样品架上, 然后置于THz-TDS系统中进行扫描以采集0～80 ps的时域光谱数据，每一个压片取不同部位3次测量的平均值作为该样品的光谱数据。最后将时域光谱数据输入到PC端的Teralyzer软件程序处理，即可从时域光谱数据中得到相应的吸收光谱数据。

样品的光谱数据在首端和高频段的信噪比较低，因此本文只取0.075～2 THz波段的数据用于后续分析。应用MSC、 基线校正和Savitzky-Golay平滑对光谱数据进行预处理，减少噪声、 样品粒度和光程变化等因素对光谱产生的影响。其中，采用基线倾斜的方法对吸收光谱曲线进行基线校正，采用三次多项式5点平滑的方法对吸收光谱曲线进行Savitzky-Golay平滑。

1.3 模型建立与评价参数

PLS是一种结合了多元线性回归分析、 典型相关分析和主成分分析三种算法优点的化学计量分析法，广泛应用于光谱分析中[12]。在本研究的PLS建模中，根据样品数量和维度，初始设定15个主因子数，采用留一交互验证法对校正集进行交互验证，取交互验证均方根误差(RMSECV)最小值对应的主因子个数作为建模的最佳主因子数。

SVM是一种有监督的学习方法，用于数据分析和模式识别[13]。本研究以径向基函数(RBF)为核函数，采用粒子群算法计算c和g，初始设置种群的粒子为20，学习因子c1和c2分别为1.5和1.7；设定惯性权重λ初始值为0.9，终止迭代次数为200。粒子群算法在粒子的适应度值趋于稳定或达到最大迭代次数后，寻优结束，得出最佳的SVM模型参数c和g。

BPNN是一种基于误差反向传播算法的非线性多层前馈神经网络，它包括输入层、 隐含层和输出层[14]。本研究初始设置输入层节点为样本的特征维度，而输出层节点、 目标误差、 最大迭代次数和学习率等参数分别为10，100，200和0.01，隐含层节点数根据经验公式(1)选取，根据网络训练效果确定最优的隐含层节点数。当网络训练达到目标误差或最大迭代次数时，网络训练停止。

(1)

式(1)中，m为隐含层节点数，n为输入层节点数，l为输出层节点数，α为1～10之间的常数。

在Matlab 2016a(MathWorks, USA)软件中完成对原始光谱的预处理、 特征频率的筛选、 PLS，SVM和BPNN模型的建立。在建立模型时，每种pH条件下的40个样本，按照3∶1的比例划分成校正集和预测集。本文选用校正集相关系数(Rc)，校正集均方根误差(RMSEC)，预测集相关系数(Rp)，预测集均方根误差(RMSEP)和剩余预测偏差(RPD)作为铅含量预测效果的评价指标。

2 结果与讨论

2.1 光谱特征及预处理

重金属铅进入土壤前后的吸收光谱曲线如图1(a)所示。可以观察到，铅含量为100 mg·kg-1的土壤样品与纯净的土壤样品相比，吸收系数增大且在1.75 THz处尤为明显，这说明铅在进入土壤后能对太赫兹光谱产生一定影响。样品的原始吸收曲线经过MSC、 基线校正和Savitzky-Golay平滑等预处理后的光谱曲线如图1(b)所示。从图中可以看出，经过预处理后的光谱曲线，无效噪声得到抑制，样品之间的分辨信息得到增强。且不同样品之间整体表现为吸收系数随铅含量的增加而逐渐增加。

图1 样品吸收曲线(a)：铅进入土壤前后的吸收曲线；(b)：预处理后的吸收曲线Fig.1 Absorption curves of samples(a)：Absorption curves of lead before and after entering soil；(b)：Absorption curves after pretreatment

2.2 全光谱的PLS建模分析

全光谱的PLS模型最佳主因子个数和预测结果如表1所示。pH 8.5的样品模型预测结果最好，Rc，RMSEC，Rp和RMSEP分别为0.991 3，27.51 mg·kg-1，0.983 9和33.35 mg·kg-1，RPD为6.85。pH 7.0的样品相关系数和RPD结果较好，但均方根误差偏大。而pH 5.5的样品模型预测结果较差，Rc，RMSEC，Rp和RMSEP分别为0.688 9，133.53 mg·kg-1，0.604 3和164.88 mg·kg-1，RPD仅为1.25。结果表明pH 7.0和pH 5.5的样品模型预测结果需要进一步提高。

表1 全光谱的PLS建模及预测结果Table 1 The PLS modeling and prediction results of full spectrum

2.3 SPA特征频率选择

实验中发现，0.075～2.0 THz的光谱数据量较大且数据之间存在共线性和大量冗余的问题，容易造成模型不稳定，预测效果较差。因此需要对原始数据进行降维以提取特征变量，减少冗余以及共线性数据的影响[11]。连续投影算法(SPA)是一种使矢量空间共线性最小化的前向变量选择算法，在光谱分析中有广泛的应用。本文对样品的光谱数据进行特征频率选择，选取性能好的特征频率，将进一步提高模型分析的准确性。

如图2所示，以SPA算法对样品光谱数据进行筛选，得到了特征频率的分布图，根据对应编号最终筛选出特征频率。pH 8.5的样品筛选出11个特征频率(0.2，0.419，0.6，0.788，0.888，1，1.206，1.306，1.406，1.469和1.544 THz)。pH 7.0和pH 5.5的样品，分别筛选了10个特征频率(0.144，0.581，0.775，1.038，1.094，1.206，1.325，1.425，1.531，1.569 THz)和13个特征频率(0.081，0.219，1.356，1.431，1.513，1.588，1.656，1.731，1.8，1.844，1.875，1.925和1.95 THz)。

图2 土壤样本的特征频率分布图(a)：pH 8.5；(b)：pH 7.0；(c)：pH 5.5Fig.2 Distribution of variables selected plot by SPA(a)：pH 8.5；(b)：pH 7.0；(c)：pH 5.5

2.4 基于特征频率的预测模型

基于SPA选择的特征频率建立了SPA-PLS，SPA-SVM和SPA-BPNN的预测模型，其中SPA-PLS的最佳主因子个数，SPA-SVM的最优c和g，SPA-BPNN的最优隐含层节点数以及三种模型的预测结果如表2所示。可以看出，在经过SPA选择特征频率后的建模效果普遍比全光谱效果好，这表明SPA算法选择的特征频率不仅含有原始光谱的有效信息，能代替原始光谱进行建模，而且减少了数据量，达到简化模型并提高模型的精度以及稳健性的目的。其中，pH 8.5的样品模型预测结果中，效果最好的为SPA-PLS模型，Rc，Rp，RMSEC，RMSEP和RPD分别为0.997 7，0.994 6，14.52 mg·kg-1，22.70 mg·kg-1和9.63。pH 7.0的样品模型预测结果较全光谱有了很大提升，其中效果最好的为SPA-SVM模型，Rc，Rp，RMSEC，RMSEP和RPD分别为0.996 2，0.975 7，20.25 mg·kg-1，33.04 mg·kg-1和4.56；pH 5.5样品SPA-PLS模型中的预测效果相比于全光谱虽有所提升，但RMSEP仍然高于样品的梯度值50 mg·kg-1，说明pH 5.5条件下的样品数据不适合类似于PLS的线性预测模型。但在非线性的预测模型SPA-SVM和SPA-BPNN中，得到了不错的预测效果，其中预测效果最好的是SPA-BPNN模型，Rc，Rp，RMSEC，RMSEP和RPD分别为0.968 7，0.974 4，48.83 mg·kg-1，55.03 mg·kg-1和4.44。基于同一种模型作比较，还可以得出3种pH条件下的样品的预测效果高低排序依次为pH 8.5>pH 7.0>pH 5.5，可能是由于铅在进入碱性土壤后，更容易与土壤产生络合反应，形成影响太赫兹光谱曲线的化学结合态，而具体是哪种形态铅影响太赫兹光谱的吸收，目前尚且不清楚，有待后续研究。其中三种pH条件下样品的最佳模型预测图如图3所示。

表2 基于特征频率的建模及预测结果Table 2 The modeling and prediction results based on sensitive frequencies

图3 三种pH条件下样品的最佳模型预测图(a)：pH 8.5；(b)：pH 7.0；(c)：pH 5.5Fig.3 The best prediction results(a)：pH 8.5；(b)：pH 7.0；(c)：pH 5.5

3 结 论

应用太赫兹光谱技术对不同pH土壤中铅含量的最佳反演预测模型进行了探索性研究。制备了3种不同pH的含铅土壤样品，采集样品的太赫兹光谱并做预处理。对预处理后的光谱数据进行了特征频率选择和建模分析，并通过模型评价参数对不同模型进行比较。实验结果表明：(1)在经过SPA选择特征频率后的建模效果普遍比全光谱的效果好；(2)其中pH 8.5的样品最佳预测模型为SPA-PLS，Rc，Rp，RMSEC，RMSEP和RPD分别为0.997 7，0.994 6，14.52 mg·kg-1，22.70 mg·kg-1和9.63；(3)pH 7.0的样品最佳预测模型为SPA-SVM，Rc，Rp，RMSEC，RMSEP和RPD分别为0.996 2，0.975 7，20.25 mg·kg-1，33.04 mg·kg-1和4.56；(4)pH 5.5的样品最佳预测模型为SPA-BPNN，Rc，Rp，RMSEC，RMSEP和RPD分别为0.968 7，0.974 4，48.83 mg·kg-1，55.03 mg·kg-1和4.44。(5)建模效果高低排序依次为pH 8.5>pH 7.0>pH 5.5。本研究为不同pH土壤中铅含量的反演预测提供了一种新思路，但土壤取样具有区域局限性，重金属类型单一，且样本数量也有待进一步增加。在将来的研究工作中将进一步完善实验的不足之处，探索更加精准的土壤铅及其他重金属含量的反演预测模型。