基于EEMD-db3的皇竹草重金属相对含量高光谱反演方法

技术领域

本发明属于植物重金属检测技术领域，具体涉及一种基于EEMD-db3的皇竹草重金属相对含量高光谱反演方法。

背景技术

皇竹草是一种营养价值高、适口性好的植物，广泛用于动物饲料的制备。由于皇竹草生物炭中特有的结构，其对土壤中的Cu、Cd、Zn等重金属离子具有良好的吸附效果，可被用来控制土壤重金属的污染。皇竹草喂给鱼、牛、羊等动物时，部分重金属元素被动物吸收，同时在动物体内积累。长期摄入会导致动物慢性中毒。过量会导致急性中毒，甚至死亡，并容易危及整个食物链的食品安全。因此有必要测定皇竹草的重金属含量是否达到饲料标准。传统的化学方法检测皇竹草中重金属含量费时费钱。高光谱技术分辨率高，可以快速获取目标连续波段的高光谱反射率信息，是一种简单、快速、无损的检测方法，广泛应用于作物、土壤成分、食物成分检测、水中漂浮物检测等定量检测中。

EEMD方法的本质是一种叠加高斯白噪声的多次经验模式分解，利用了高斯白噪声具有频率均匀分布的统计特性，通过每次加入同等幅值的不同白噪声来改变信号的极值点特性，之后对多次EMD得到的相应IMF进行总体平均来抵消加入的白噪声，从而有效抑制模态混叠的产生。db小波主要应用在离散型的小波转换，是最常使用到的小波转换，通常使用在数位信号分析、信号压缩跟噪声去除。

目前，利用光谱无损检测方法估算堆积物中重金属的含量的报道较少，利用EEMD方法在植物重金属含量检测领域的应用并不广泛。

发明内容

本发明的主要目的在于克服现有技术的缺点与不足，提出一种基于EEMD-db3的皇竹草重金属相对含量高光谱反演方法。

为了达到上述目的，本发明采用以下技术方案：

基于EEMD-db3的皇竹草重金属相对含量高光谱反演方法，包括以下步骤：

S1、对不同重金属浓度下种植的皇竹草，获取高光谱数据，测定重金属含量；

S2、采用EEMD-db3的方法对皇竹草原始光谱数据进行降噪处理，并进行特征波段选取；

S3、构建皇竹草重金属相对含量反演模型，选取最优反演模型用于皇竹草重金属相对含量高光谱反演。

进一步的，步骤S1具体为：

将鲜叶皇竹草样品通过高光谱成像仪采集光谱反射率数据，根据国际标准测定皇竹草内重金属含量；

根据国际标准测定皇竹草内重金属含量具体为：

将鲜叶皇竹草样品在70℃烘箱中杀青30min，55℃烘干48h，得到植物干样，再用粉碎机粉碎，过100目筛，置于密封塑料袋内；

采用HNO

进一步的，步骤S2包括：

采用EEMD方法分解高光谱信号得到IMF分量，采用db3降噪，然后将降噪后IMF分量与剩余分量相加重构原光谱信号数据，再对降噪后光谱数据进行一阶导、倒数对数数学变换；

db3小波母函数对降噪后皇竹草光谱数据进行分解，得到高频与低频小波系数，最后通过相关性曲线筛选显著相关波段后，采用逐次投影算法与竞争性自适应重加权算法提取特征变量，同时保留相关系数最大的波段作为特征波段。

进一步的，步骤S2中，分解得到IMF分量具体包括：

步骤一、将高斯白噪声n(t)添加到目标信号y(t)，得到一个新的信号x(t):

x(t)＝y(t)+n(t)

步骤二、按照EMD方法对信号x(t)进行分解，得到一组IMF分量与一个残余分量:

重复步骤一和步骤二i次，每次重复加入不同得白噪声n

进一步的，步骤S2中，采用db3降噪具体为：

利用高斯白噪声频谱的零均值原理，将分解出的m组IMF分量和残余分量的整体平均值作为最终结果，减少白噪声的影响：

将分解得到的IMF分量通过高通滤波器g[n]和低通滤波器h[n]进一步分解为高频滤波和低频滤波，再将低频滤波进一步分解成高频滤波和低频滤波，分解公式如下：

将低频系数和高频系数的db3小波进行重构：

进一步的，步骤S2中，将降噪后IMF分量与剩余分量相加重构原光谱信号数据，具体为：

将经db3小波降噪并重构的IMF分量与剩余分量相加重构原光谱信号数据：

进一步的，步骤S2中，对降噪后光谱数据进行一阶导、倒数对数数学变换，具体为：

对降噪后的光谱数据进行一阶导、倒数对数数学变换；

一阶导公式如下：

倒数对数公式如下：

进一步的，采用逐次投影算法提取特征变量，同时保留相关系数最大的波段作为特征波段具体为：

从一个波长开始，然后在每次迭代中加入一个新变量，直至所选变量数达到设定值N，其实现步骤表示如下：

步骤一、初始化n＝1,x

步骤二、确定未选波段变量：

步骤三、计算未选波段和初始化波段的投影映射：

Px

步骤四、确定最大投影：

步骤五、x

步骤六、n＝n+1，当n

进一步的，采用竞争性自适应重加权算法提取特征变量，同时保留相关系数最大的波段作为特征波段具体为：

竞争性自适应重加权算法将光谱各波段变量视为一个个体，筛选较强适应能力的个体作为光谱特征变量，降低原始光谱波段之前的高度共线性；由蒙特卡罗交叉验证法选择最优潜在波段变量，其中将蒙特卡罗采样次数设定为100，对采样次数进行反复迭代，筛选出特征波段，其实现步骤表示如下：

步骤一、采用蒙特卡洛采样法，每次随机选择80％的样本进入建模集，剩余的20％作为预测集建立PLS模型；记录每一次采样过程PLS模型中的回归系数的绝对值权重：

其中，|b

步骤二、利用指数衰减函数EDF强行去除回归系数绝对值权重相对较小的波长，在第i次基于MC采样建立PLS模型时，根据EDF得到保留的波长点的比例R

R

其中，μ和k是常数，按照以下两种情况计算：

情况一、在一次采样并进行相应计算时，所有的波长都参与了建模分析，因此保留的波长点的比例为1；

情况二、在最后一次采样完成并进行相应计算时，只剩下两个波长参与PLS建模，此时保留的波长点的比例为

步骤三、在每次采样时，都从上一次采样时的变量数中采用自适应加权采样ARS选择数量为Ri*n个的波长变量，进行PLS建模，计算RMSECV；

步骤四、在N次采样完成之后，CARS算法得到了N组候选的特征波长子集，以及对应的RMSECV值，选择RMSECV最小值所对应的波长变量子集为特征波长。

进一步的，步骤S3具体为：

据相关性曲线选择相关系数最大的波段作为自变量，分别构建皇竹草重金属相对含量线性、抛物线、指数和对数反演模型；

将样本随机划分为训练集和预测集，采用决定系数R

其中，y

通过比较各反演模型的决定系数R

本发明与现有技术相比，具有如下优点和有益效果：

1、本发明采用逐次投影算法与竞争性自适应重加权算法提取特征变量，同时保留相关系数最大的波段作为特征波段，通过EEMD-db3算法建立皇竹草重金属含量反演模型；验证决策系数表明，本发明方法可定量分析Cd的含量，而Zn估计模型可以更准确地定量含量；本发明的EEMD-db3预处理效果更加稳定，不输传统降噪方法。

2、本发明对重金属Zn的分析效果较好，预处理后相关系数最大值大于0.8。

附图说明

图1是本发明方法的流程图；

图2a是实施例中原始光谱曲线；

图2b是实施例中EMD分解结果；

图2c是实施例中EEMD-db3降噪后的光谱曲线。

具体实施方式

下面结合实施例及附图对本发明作进一步详细的描述，但本发明的实施方式不限于此。

实施例

本实施例在某大学生态实验场进行，在含有不同重金属浓度的实验区，以4株/m的密度种植皇竹草。试验土壤呈酸性，pH值在4.5到6.0之间。Cd浓度为0.0115～0.0650mg/kg，Cu浓度为70.25～120.65mg/kg，Zn浓度为45.40～60.38mg/kg。栽植期间，生长季第1、3个月分别施复混肥25g/株和50g/株。随机采集8株皇竹草。每株皇竹草的上部被分为5个部分：老叶、幼叶、上茎、中茎以及下茎。

如图1所示，基于EEMD-db3的皇竹草重金属相对含量高光谱反演方法，包括以下步骤：

S1、对不同重金属浓度下种植的皇竹草，获取高光谱数据，测定重金属含量；在本实施例中，具体为：

将鲜叶皇竹草样品通过HyperSIS-VNIR-QE高光谱成像仪采集光谱反射率数据，光谱波段范围为400-1000nm。

对各样本相同位置的点的反射率数据进行黑白板矫正。

测定重金属含量具体为在70℃烘箱中杀青30min，55℃烘干48h，得到植物干样，再用粉碎机粉碎，过100目筛，置于密封塑料袋内。植物重金属采用HNO

S2、采用EEMD-db3的方法对皇竹草原始光谱数据进行降噪处理，并进行特征波段选取；具体的，采用EEMD方法对皇竹草原始光谱数据进行降噪处理，再对降噪后光谱数据进行一阶导、倒数对数数学变换。即将EEMD分解得到IMF分量采用db3降噪，然后将降噪后IMF分量与剩余分量相加重构原光谱信号数据。所述db3小波母函数对降噪后皇竹草光谱数据进行分解，得到个高频与低频小波系数。最后通过相关性曲线筛选显著相关波段后，采用逐次投影算法与竞争性自适应重加权算法提取特征变量，同时保留相关系数最大的波段作为特征波段。

在本实施例中，具体为：

S21、将高斯白噪声n(t)添加到目标信号y(t)，得到一个新的信号x(t):

x(t)＝y(t)+n(t)

S22、按照EMD方法对信号x(t)进行分解，得到一组IMF分量与一个残余分量:

S23、重复步骤一和步骤二i次，每次重复加入不同得白噪声n

S24、利用高斯白噪声频谱的零均值原理，将分解出的m组IMF分量和残余分量的整体平均值作为最终结果，减少白噪声的影响：

S25、将EEDM分解得到的IMF分量通过高通滤波器g[n]和低通滤波器h[n]进一步分解为高频滤波和低频滤波，再将低频滤波进一步分解成高频滤波和低频滤波，分解公式如下：

S26、将低频系数和高频系数的db3小波进行重构：

S27、将经db3小波降噪并重构的IMF分量与剩余分量相加重构原光谱信号

S28、对降噪后的光谱数据进行一阶导、倒数对数数学变换；

一阶导公式如下：

倒数对数公式如下：

如图2a、图2b以及图2c所示，分别为皇竹草原始光谱图、EMD分解结果图以及EEMD-db3去噪后光谱曲线图，从中对比可以看出，EEMD-db3具有较好的去噪效果。

在本实施例中，采用逐次投影算法提取特征变量，同时保留相关系数最大的波段作为特征波段具体为：

从一个波长开始，然后在每次迭代中加入一个新变量，直至所选变量数达到设定值N，其实现步骤表示如下：

步骤一、初始化n＝1,x

步骤二、确定未选波段变量：

步骤三、计算未选波段和初始化波段的投影映射：

步骤四、确定最大投影：

步骤五、x

步骤六、n＝n+1，当n

在本实施例中，采用竞争性自适应重加权算法提取特征变量，同时保留相关系数最大的波段作为特征波段具体为：

竞争性自适应重加权算法将光谱各波段变量视为一个个体，筛选较强适应能力的个体作为光谱特征变量，降低原始光谱波段之前的高度共线性；由蒙特卡罗交叉验证法选择最优潜在波段变量，其中将蒙特卡罗采样次数设定为100，对采样次数进行反复迭代，筛选出特征波段，其实现步骤表示如下：

步骤一、采用蒙特卡洛采样法，每次随机选择80％的样本进入建模集，剩余的20％作为预测集建立PLS模型；记录每一次采样过程PLS模型中的回归系数的绝对值权重：

其中，|b

步骤二、利用指数衰减函数EDF强行去除回归系数绝对值权重相对较小的波长，在第i次基于MC采样建立PLS模型时，根据EDF得到保留的波长点的比例R

R

其中，μ和k是常数，按照以下两种情况计算：

情况一、在一次采样并进行相应计算时，所有的波长都参与了建模分析，因此保留的波长点的比例为1；

情况二、在最后一次采样完成并进行相应计算时，只剩下两个波长参与PLS建模，此时保留的波长点的比例为

步骤三、在每次采样时，都从上一次采样时的变量数中采用自适应加权采样ARS选择数量为Ri*n个的波长变量，进行PLS建模，计算RMSECV；

步骤四、在N次采样完成之后，CARS算法得到了N组候选的特征波长子集，以及对应的RMSECV值，选择RMSECV最小值所对应的波长变量子集为特征波长。

S3、构建皇竹草重金属相对含量反演模型，选取最优反演模型用于皇竹草重金属相对含量高光谱反演；在本实施例中，具体为：

据相关性曲线选择相关系数最大的波段作为自变量，分别构建皇竹草重金属相对含量线性、抛物线、指数和对数反演模型；

将40个样本数据随机划分为32个训练集，8个预测集；采用决定系数R

其中，y

通过比较各反演模型的决定系数R

如下表1所示，为模型类型与反演结果表。

/>

表1

还需要说明的是，在本说明书中，诸如术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。

对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其他实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。