MATLAB中基于CEEMDAN与FastICA算法的单通道信号规律噪声的滤除方法

本文章思想主要来源于下面的一篇文献，将其用matlab实现。
[1]罗志增,严志华,傅炜东.基于CEEMDAN-ICA的单通道脑电信号眼电伪迹滤除方法[J].传感技术学报,2018,31(08):1211-1216.

1.CEEMDAN

基于自适应噪声的完备经验模态分解CEEMDAN ( Complete Ensemble Empirical Mode Decomposition with Adaptive Noise) 是对经验模态分解EMD,集合经验模态分解EEMD算法的改进方法。

  • EMD：EMD 方法可以将非线性、非平稳信号分解为有限个本征模态函数 IMF( Intrinsic Mode Function) ，每个 IMF 分量代表了不同时间尺度上的样本特征，其存在“模态混叠”问题。

  • EEMD： EEMD 方法是对 EMD 方法的改进，它利用白噪声频率均匀分布的特性，在原始信号中加入正态分布的白噪声，使信号在整个频带中极值点间隔分布均匀，在不同尺度上具有连续性，从而降低了模态混叠效应。

  • CEEMDAN： CEEMDAN 方法则在 EEMD 的基础上，通过自适应加入白噪声，克服了 EEMD 在加入白噪声后分解失去完备性、产生重构误差的问题。

2.FastICA

首先先了解ICA（独立成分分析）

独立成分分析法(ICA) 是一种常见的盲源分离方法，在多个源信号和传输通道参数未知的情况下， 通过观测信号来估计源信号，从而将源信号恢复。 它的分离模型如图所示。 观测信号 X(t) = { x1(t) ， x2(t) ， …， xn(t) }经源信号 S(t) = { s1(t) ， s2(t) ，…， sm(t)}通过某一未知混合矩阵A得到， 即X=AS。ICA 就是在S和A均为未知的情况下，求解解混矩阵W，通过W从观测信号中得出源信号S的估计 Y，ICA 的前提条件是观测信号的数量不少于源信号的数量。

更具体的请参考
独立成分分析ICA与FastICA： https://blog.csdn.net/sinat_27763257/article/details/80568792

在诸多ICA算法中，固定点算法 (也称FastlCA)以其收敛速度快、分离效果好被广泛应用于信号处理领域。该算法能很好地从观测信号中估计出相互统计独立的、被未知因素混合的原始信号。

3.该算法主要思路

针对单通道信号，首先通过CEEMDAN算法将该信号分解为有限个本征模态函数 IMF，然后将其作为观测信号 X(t) = { x1(t) ， x2(t) ， …， xn(t) }通过FastICA算法得出源信号估计S(t) = { s1(t) ， s2(t) ，…， sm(t)}，而规律噪声信号就包含在其中。那如何寻找该规律噪声呢？这里我们用到了模糊熵。

  • 模糊熵
    信息论中的熵可以用来衡量数据的规律性，通过熵值判断我们可以寻找出相对于有价值信号熵值更低的规律噪声。而其中的模糊熵比样本熵和近似熵抗噪能力更强，稳定性好，所以我们选用模糊熵。

通过模糊熵算法，我们可以得出各个源信号的模糊熵。利用文献中提及的GomezHerrero提出的阈值判别式，我们可以将熵值较小的规律噪声找出并将其置零，即达到了去除规律噪声的目的。

最后，我们需要通过ICA逆变换得出新的IMF，再将IMF直接累加达到所需结果。

4.程序代码

%% 导入单条时间序列信号x（本例是脑电信号）%% CEEMDAN得出IMF[modes its]=ceemdan(x,0.2,100,20); %modes就包含了各个IMF%% 快速独立成分分析[icasig,A,W]=fastica(modes) ;% icasig就是源信号估计%% 求出各源信号估计对应模糊熵值[n,m]=size(icasig);for i=1:nfuzzyen(i)=Fuzzy_Entropy(4,0.2,icasig(i,:));end%% 利用判别法求出规律噪声熵值new_fuzzyen=sort(fuzzyen)k=2; while new_fuzzyen(k+1)-new_fuzzyen(k)<new_fuzzyen(k)-new_fuzzyen(k-1)            if k<n/2            k=k+1;            else k=fix(n/2);break;            end        endk=k-1;%1：k都为伪迹成分%% 将对应规律噪声置零icasigCopy=icasig;for i=1:k    for j=1:n        if fuzzyen(j)==new_fuzzyen(i)            icasigCopy(j,:)=0        end    endend%% ICA逆变换并叠加new_modes=A*icasigCopy;for j=1:mnew_x(j)=sum(new_modes(:,j));endnew_x=round(new_x)';% 四舍五入取整%% 画图比较观察plot(x)hold onplot(new_x)12345678910111213141516171819202122232425262728293031323334353637383940414243444546

  • 资源
    包含论文，CEEMDAN、FastICA、模糊熵的函数代码，以及练手用的脑电信号。

赞 (0)

相关推荐