HHT时频分析土壤光谱的重金属铜离子污染信息提取模型
杨可明1, 汪国平1,2, 付萍杰1, 张伟1, 王晓峰1
1. 中国矿业大学(北京)地球科学与测绘工程学院, 北京 100083
2. 北京航天宏图信息技术股份有限公司, 北京 100195

作者简介: 杨可明, 1969年生, 中国矿业大学(北京)地球科学与测绘工程学院教授 e-mail: ykm69@163.com

摘要

土壤中不同浓度Cu2+含量映射到土壤光谱上的信息量十分微弱, 并且这些高光谱数据中也存在着难以避免的噪声, 因而本研究的关键是如何在土壤光谱复杂的噪声环境中提取微弱Cu2+信息。 经验模态分解算法(EMD)能够有效去除高光谱数据中的噪声, 且EMD是Hilbert变换对“非线性非稳定”信号时频分析的前提, 当引入Huang变换后, 可利用Hilbert-Huang变换(HHT)模型时频分析高光谱数据以实现降噪处理与信息提取。 通过时频的HHT分析不同浓度Cu2+污染下的土壤光谱, 完成从原始光谱经EMD分解出各本征模态函数(IMF)分量的包络线、 调制信号和频谱等曲线中挖掘土壤光谱的Cu2+污染信息。 研究结果表明, 相同浓度Cu2+污染时的土壤光谱HHT时频分析结果相同, 不同浓度时则不同, 所以也可依据IMF分量反演土壤Cu2+含量。 因此, 高光谱数据的HHT时频分析能为土壤光谱的信息挖掘、 光谱诊断和Cu2+含量反演等提供一种新的方法和思路。

关键词: 重金属铜污染; 土壤光谱; 经验模态分解; 希尔伯特-黄变换; 弱信息探测; 时频分析
中图分类号:TP7 文献标志码:A
A Model on Extracting the Pollution Information of Heavy Metal Copper Ion Based on the Soil Spectra Analyzed by HHT in Time-Frequency
YANG Ke-ming1, WANG Guo-ping1,2, FU Ping-jie1, ZHANG Wei1, WANG Xiao-feng1
1. College of Geoscience and Surveying Engineering, China University of Mining & Technology (Beijing), Beijing 100083, China;
2. Beijing PIESAT Information Technology Co., Ltd., Beijing 100195, China;
Abstract

The amount of information that the different Cu2+ contents in soil are mapped to the soil spectra is very weak, and the noises in the hyperspectral data are also very difficult to avoid. Thus the sky point of the research would be how to extract the weak Cu2+ information from the complex noise environment of soil spectra. The empirical mode decomposition (EMD) algorithm can effectively remove the noise in hyperspectral data. Further more, the EMD is the premise of Hilbert transform that is a kind of time-frequency analysis on nonlinear and unstable signal. Therefore the time-frequency analysis algorithm on Hilbert-Huang transform (HHT) could be used for de-noise processing and information extracting of the soil spectra after the Huang transform is introduced. In this paper, through the HHT was applied in time-frequency to analyze the soil spectra polluted by different Cu2+ concentrations, and the information mining of Cu2+ pollution in soil spectra was achieved based on the curves such as envelope, modulation signal and frequency spectrum of each intrinsic mode function (IMF) component that could be decomposed out from the original soil spectra by EMD. The study results showed that, the analyzing results of soil spectra were identical based on the time-frequency HHT when the soil was polluted by same Cu2+ concentrations, but non-identical when polluted by different Cu2+ concentrations, so that the Cu2+ content in soil can be also retrieved based on the IMF components. Accordingly, the HHT analysis in time-frequency provides a new method and idea for information extracting, spectral diagnosing and Cu2+ content retrieving of soil spectra.

Keyword: Heavy metal copper pollution; Soil spectra; Empirical mode decomposition; Hilbert-Huang transform; Weak information detection; Time-frequency analysis

引 言

土壤和农作物中重金属含量微弱, 但重金属在土壤和农作物中的积累对土壤有机质、 农作物生长和生态环境可持续发展等方面有极其严重影响[1]。 传统环境监测技术虽然精度高, 但是由于耗时费力周期长等缺点而难以满足当代精准农业快速、 大量监测环境污染信息的需求。 高光谱遥感因光谱分辨率高、 无损无污染等[2, 3, 4]特点在环境监测中优势显而易见。 近年来对土壤和农作物的相关研究表明[5, 6], 对高光谱数据进行微分、 分数阶微分、 吸收峰深度、 倒数、 包络线去除等变换, 可以增强高光谱数据之间的相关性, 筛选出敏感的特征波段和确定红边参数、 红边“ 蓝移” 指数等一些相关参数[7], 以建立反演含量模型进行污染监测等。 总结出此类研究均是在光谱域和空间域中分析, 可使高光谱数据的研究思路拓展到频率域中。

电磁遥感信号获取过程中不可避免地会受到外部通讯系统干扰和仪器电子热运动[8]的影响而产生噪声, 噪声极大干扰重金属在高光谱上污染信息的解读, 同样导致诊断高光谱重金属污染信息不准确。 传统电磁信号降噪的中值滤波、 平滑去噪和算术平均去噪等方法虽然在高信噪比小时有一定抑制噪声的效果, 但针对含有微弱重金属污染信息的高光谱降噪过程, 很容易淹没此类微弱的有用信息, 因而很难实现有效降噪。 伴随数学方法连续小波变换[9, 10](continue wavelet transform, CWT)、 短时傅里叶变换[11](short-time Fourier transform, STFT)、 傅里叶变换(Fourier transform, FT)和谐波分析[12](harmonic analysis, HA)等在电磁信号中应用, 可将时频分析方法引入高光谱遥感中的研究效果显著, 但对土壤或农作物重金属污染信号产生机理的研究尚少, 未能从根本上挖掘重金属污染信息, 从而寻求一种有效提取土壤或农作物重金属污染微弱信息的方法对高光谱遥感重金属污染监测具有意义重大。 本工作从频率域研究, 将希尔伯特-黄变换(Hilbert-Huang transform, HHT)[13]引入到高光谱信号降噪和微弱的重金属污染信息挖掘中, 以实现高光谱遥感对土壤或农作物重金属污染机理等根本因素的分析。

1 理论与方法
1.1 经验模态分解

经验模态分解法(empirical mode decomposition, EMD)是美籍华人黄锷于1996年提出的信号分解算法, 又称Huang算法。 EMD可从复杂信号“ 筛分” 处理过程中实现信号分解, 从而分离出固有模态函数(intrinsic mode function, IMF)。

对任意给定的光谱信号矢量f(x), 采用三次样条曲线连接所有极大值形成上包络线, 连接所有极小值形成下包络线, 并计算光谱信号f(x)与上、 下包络线的均值之差得到疑似IMF分量。 依据IMF分量的两个重要条件对疑似IMF分量进行判断, 直到第n阶残余信号是单调函数而不能继续分解出IMF时终止筛分。 筛分条件如下[14]:

(1) 均值线的平均值趋近于0, 一般与0相差小于0.1即可认为符合条件。

(2) 原始信号极值点个数与原始信号同横轴交点个数之差不能大于1。

数学中, 通过公式简化表示该分解过程是n个IMFj(x)分量和一个残余项R(x)之和, 见式(1)

f(x)=j=1nIMFj(x)+R(x)(1)

1.2 Hilbert变换算法

希尔伯特变换(Hilbert transform, HT)由“ 无冕之王” David Hilbert提出, 当信号经过Hilbert变换后会出现90° 相移[15]。 定义高光谱矢量由f(x)表示, 则f(x)的Hilbert变换是H[f(x)]

H[f(x)]=1π-+f(t)x-tdt(2)

根据信号理论可知, 式(2)中H[f(x)]是f(x)和1/π x的卷积, x表示波段号, t为积分变量。 式(2)查表并结合欧拉公式可得

θ(x)=-π2, x> 0+π2, x< 0(3)

由式(3)可看出, 频率大于0时, 相位左移90° ; 频率小于0时, 相位右移90° 。 通常将Hilbert变换结果用复光谱信号表示, 即

Z(x)=f(x)+iH[f(x)]=a(x)eiθ(x)(4)a(x)=f2(x)+H2[f(x)](5)θ(x)=arctanH[f(x)]f(x)(6)r(x)=12π×dθ(x)dx(7)

式(4)— 式(7)中, Z(x)被称为f(x)的解析信号, f(x)被称为复信号Z(x)的实部, H[f(x)]被称为复信号Z(x)的虚部, a(x)是瞬时幅值, θ (x)是瞬时相位, r(x)是瞬时频率, i表示虚数。

1.3 Hilbert-Huang变换

Hilbert变换在信号解析应用中需要信号同时满足“ 线性” 和“ 稳态” 这两个苛刻的要求, 自然界中绝大多数信号不是“ 线性非稳态” 就是属于“ 非线性稳态” , 甚至很多信号属于“ 非线性非稳态” 。 然而黄锷的EMD算法很好地解决了Hilbert变换存在的这一问题, 因此, Hilbert-Huang变换也就应运而生, Hilbert-Huang变换(Hilbert-Huang transform, HHT)的算法流程如图1所示。

图1 Hilbert-Huang变换算法流程Fig.1 The algorithm flow of Hilbert-Huang transform

2 实验部分
2.1 方法

(1)实验材料与设备。 培育农作物的土壤是直径小于1 mm沙土, 用CuSO4· 5H2O配置重金属污染溶液, 培育辅助设备: 电热板、 电子分析天平。 土壤光谱测量采用350~2 500 nm光谱范围的SVC HR-1024I高性能地物光谱仪。

(2)污染梯度设置。 实验设计中唯一变量是土壤中铜离子浓度, 设置铜离子浓度梯度为100 μ g· g-1, 浓度范围在0~1 000 μ g· g-1, 即0, 100, …, 900, 1 000 μ g· g-1共11个梯度, 同一梯度有3组平行实验, 共33组实验对象。

(3)土壤污染设置步骤。 采用CuSO4溶液与沙土相拌匀的方法, 以构成重金属胁迫农作物生长的土壤环境。 每个污染梯度沙土重3 kg, 计算出3 kg土最大持水量, 将对应梯度所需CuSO4· 5H2O溶解于800 mL(80%持水量)的蒸馏水中, 均匀浇入土壤中, 再加入20%蒸馏水至土壤最大持水量, 放置24 h后再种植农作物。

2.2 土壤光谱数据

采用SVC HR-1024I地物光谱仪测定不同Cu2+浓度污染下土壤光谱数据。 在采集光谱数据时, 为了防止太阳光对土壤光谱的影响, 土壤光谱测定在不透光的室内进行, 采用了功率为50 W的卤素灯光源和垂直于土壤表面40 cm的4° 视场角光纤探头。 对每个梯度土壤中3组光谱取平均得到11个梯度土壤均值光谱数据如图2所示。

图2 不同浓度Cu2+污染下的土壤光谱Fig.2 Soil spectra polluted by Cu2+ under different concentrations

3 结果与讨论
3.1 时频分析土壤光谱

根据信号理论可知, 土壤光谱属于“ 非线性非稳定” 信号, 不能满足Hilbert变换的要求, 因而需对光谱信号进行经验模态分解成IMF后再做Hilbert变换和频谱分析。 采集原始光谱存在一定的噪声, 在进行时频分析时, 应当剔除光谱中噪声的干扰, 然而经过EMD得到的IMF可以有效消除原始信号中噪声的影响。 随机选择无污染的土壤光谱中一组进行HHT分析得到9个IMF和对应每个IMF及原始光谱的频谱图如图3所示。

图3 土壤光谱平行组HHT时频分析结果Fig.3 Results of HHT time frequency analysis on parallel group soil spectra

3.2 时频分析不同污染梯度土壤光谱

土壤中Cu2+含量会导致光谱局部产生突变, 在经过去噪处理后的土壤光谱中找出导致土壤光谱突变的点, 是解决重金属污染微弱信息难挖掘的关键。 研究表明, 不同污染梯度的土壤光谱存在细微的差异, 但并不能充分表达土壤污染光谱关键所在。 实验对11个梯度土壤光谱进行HHT时频分析后, 不同的污染梯度得到IMF数量不同, 存在一定差异, 每个梯度得出的最高频率IMF所属第几IMF也不尽相同。 但从土壤HHT时频分析结果可以看出IMF瞬时频率、 调制信号和IMF包络线等均存在很大的差异, 计算各对应曲线间相关系数都很小。 因实验中Cu2+污染梯度和每个梯度光谱EMD分解出的IMF多, 故任意选取梯度0, 300, 500和800 μ g· g-1中的IMF5分量的结果进行绘图对比分析, 如图4。

图4 不同Cu2+浓度污染时土壤光谱的IMF5分析结果Fig.4 The IMF5 analysis results on spectral differences of soil polluted by different Cu2+concentrations

3.3 时频分析相同污染梯度土壤光谱

不同浓度Cu2+污染土壤时, 反映在土壤光谱上信息均存在微小的差异, 验证该方法能否有效正确估算土壤中Cu2+含量, 通过实验分析11个梯度中的每组平行样本, 对相同浓度Cu2+污染土壤下不同样本光谱的时频分析结果进行对比检验。 以Cu2+浓度含量0 μ g· g-1的土壤光谱的IMF5分量为例, 采集平行实验中不同花盆中土壤样本进行HHT变换和频谱分析, 研究结果表明同一梯度中不同花盆的土壤光谱调制信号、 瞬时频率等曲线十分相似, 计算各同类曲线相关系数均达到0.9以上, 表明同类曲线之间十分相似, 具体如图5。

图5 相同Cu2+浓度污染时土壤光谱的IMF5分析结果Fig.5 IMF5 analysis results on spectral differences of soil polluted by identical Cu2+concentrations

4 结 论

从频率域思考, EMD可以有效去除信号中噪声分解成若干IMF分量, 而结合EMD算法的HHT变换是时频分析的研究方法之一。 故将时频分析方法HHT尝试性地运用到高光谱数据信息提取中, 最大化地消除噪声的影响, 在微弱的光谱差异中挖掘有利于诊断光谱差异和识别污染光谱信息以便后续研究。 采用HHT对不同Cu2+浓度污染时土壤高光谱数据进行时频分析, 重复随机地选取光谱数据对比分析平行实验组内、 外光谱数据的HHT分析结果表明, 组内时频分析各IMF分量评价曲线相似度高于0.9, 组外则相似度很低。 因此, HHT算法中EMD去噪所得各IMF分量的调制信号、 包络线和频谱等曲线能够准确诊断不同浓度Cu2+污染土壤时光谱之间的差异, 同时, 根据分析结果建立Cu标准曲线反演土壤含量Cu2+含量。

The authors have declared that no competing interests exist.

参考文献
[1] Conforti M, Buttafuoco G, Leone A P, et al. Catena, 2013, 110(11): 44. [本文引用:1]
[2] Koppnen S M, Brezonik P L, Olmanson L G, et al. Remote Sensing of Environment, 2002, 82: 38. [本文引用:1]
[3] Koppnen S, Pulliainen J. Remote Sensing of Environment, 2002, 79: 51. [本文引用:1]
[4] Pablo H Rosso, James C Pushnik, Mui Lay, et al. Environmental Pollution, 2005, 137: 241. [本文引用:1]
[5] LIU Mei-ling, LIU Xiang-nan, LI Ting, et al(刘美玲, 刘湘南, 李婷, ). Transactions of the Chinese Society of Agricultural Engineering(农业工程学报), 2010, 26(3): 191. [本文引用:1]
[6] ZHANG Long, PAN Jia-rong, ZHU Cheng(张龙, 潘家荣, 朱诚). Journal of Zhejiang University·Agric. & Life Sci. (浙江大学学报·农业与生命科学版), 2013, 39(1): 50. [本文引用:1]
[7] ZHU Ye-qing, QU Yong-hua, LIU Su-hong, et al(朱叶青, 屈永华, 刘素红, ). Journal of Remote Sensing(遥感学报), 2014, 18(2): 335. [本文引用:1]
[8] LI Cheng-wu, DONG Li-hui, WANG Qi-fei, et al(李成武, 董利辉, 王启飞, ). Journal of China Coal Society(煤炭学报), 2016, 41(8): 1933. [本文引用:1]
[9] YU Lei, HONG Yong-sheng, ZHOU Yong, et al(于雷, 洪永胜, 周勇, ). Spectroscopy and Spectral Analysis(光谱学与光谱分析), 2016, 36(5): 1428. [本文引用:1]
[10] Liao Q H, Wang J H, Yang G J, et al. Journal of Applied Remote Sensing, 2013, 7(1): 1. [本文引用:1]
[11] PANG Cun-suo, LIU Lei, SHAN Tao(庞存锁, 刘磊, 单涛). Acta Electronica Sinica(电子学报), 2014, 42(2): 347. [本文引用:1]
[12] Jakubauskas M E, Legates D R, Kastens J H. Photogrammetric Engineering and Remote Sensing, 2001, 4: 461. [本文引用:1]
[13] Huang N E, Shen Z, Long S R, et al. Proc. R. Soc. Lond, 1998, A: 903. [本文引用:1]
[14] Xue Yajuan, Cao Junxing, Tian Renfei. [J]. Appl. Geophys. , 2013, 89: 108. [本文引用:1]
[15] Xu Chuanlong, Liang Cai, Zhou Bin, et al. Chemical Engineering Science, 2010, 65(4): 1334. [本文引用:1]