【心电信号】心电信号去除基线漂移【含Matlab源码 955期】

最新推荐文章于 2024-07-01 21:26:23 发布

Matlab领域

最新推荐文章于 2024-07-01 21:26:23 发布

阅读量1.2k

点赞数 27

分类专栏： Matlab信号处理（初级版）文章标签： matlab

本文链接：https://blog.csdn.net/weixin_63266434/article/details/137014735

版权

Matlab信号处理（初级版）专栏收录该内容

342 篇文章 16 订阅

订阅专栏

在这里插入图片描述

⛄一、获取代码方式

获取代码方式1：
完整代码已上传我的资源：【心电信号】基于matlab心电信号去除基线漂移【含Matlab源码 955期】
点击上面蓝色字体，直接付费下载，即可。

获取代码方式2：
付费专栏Matlab信号处理（初级版）

备注：
点击上面蓝色字体付费专栏Matlab信号处理（初级版），扫描上面二维码，付费29.9元订阅海神之光博客付费专栏Matlab信号处理（初级版），凭支付凭证，私信博主，可免费获得1份本博客上传CSDN资源代码（有效期为订阅日起，三天内有效）；
点击CSDN资源下载链接：1份本博客上传CSDN资源代码

⛄二、心电信号中基线漂移的去除方法比较及算法实现简介

0 引言
心电信号是生物医学信号中一种低幅低频的微弱信号。基线漂移一般是由于人体呼吸、电极运动等引起, 其一般小于1 Hz。常用的去除方法有中值滤波法、小波变换法、形态学滤波法等。本文简单介绍了几种常用方法并用Matlab实现。以下用于处理的原始心电信号是截取自MIT-BIH数据库的序号为103和114两种数据中的一段数据, 我们通过对比基线漂移严重和微弱两种情况, 比较各种方法的优劣之处, 其中103为较严重漂移的心电数据, 114为漂移较为微弱的心电数据。

1 方法
1.1 中值滤波法
中值滤波的定义就是取信号中的某一点前的N (N为奇数) 个数据, 将这N个数据进行排序, 取数据的中值, 将数据的中值作为该点的数据。将整段信号都做以上处理, 提取出漂移信号, 将原始信号减去漂移信号, 便得到滤波后的效果。测试效果如图1、图2所示。在这里插入图片描述
图1 中值滤波去除严重基线漂移

1.2 高通滤波法
基线漂移的频率成分较低, 所以用高通的方法对基线漂移的去除也有一定的效果。IIR滤波器能用较低的阶数获得较高的频率选择性, 但不具有线性相位。由于心电信号的处理对线性相位的要求不是很高, 所以可以选择IIR滤波器。这里我们选用butterworth滤波器, 截止频率为1 Hz。测试效果如图3、图4所示。
在这里插入图片描述
图2 中值滤波去除微弱基线漂移

图3 高通滤波去除严重基线漂移

图4 高通滤波去除微弱基线漂移
1.3 小波变换法
根据小波变换多分辨分析理论, 对含噪信号进行多尺度小波分解后, 可以利用信号和噪声在频谱和能量分布上的不同, 直接将噪声所对应的小波分解尺度上的细节分量去除, 然后利用小波逆变换进行信号重构, 就可以有效去除信号中的噪声成份。经过试验得到, 采用小波函数db8对原始信号进行8个尺度的分解, 第8层分解的信号频率与基线漂移的频率接近, 因此将前7层信号重构, 去除第8层信号, 即可得到去除基线漂移的信号。测试效果如图5、图6所示。
在这里插入图片描述
图5 小波变换去除严重基线漂移

图6 小波变换去除微弱基线漂移

⛄三、部分源代码

clc;
clear;
close all;

%% 提取信号
M = importdata(‘3.txt’);
fsample=1000;%采样率为1KHz

[mx,my]=size(M);
Signal=M(:,2);%M的第一列为时间，第二列为信号

length=floor(mx/2);%取原始信号的一半。
S=Signal(1:length);

%% 高通滤波，去除基线漂移的影响
disp(‘-------------------------------------------’);
disp(‘1：工具箱巴特沃斯高通滤波器’);
disp(‘2：IIR高通滤波’);
disp(‘3：FIR高通滤波’);
disp(‘4：中值滤波’);
disp(‘5：稀疏小波滤波’);
disp(‘6：中值+小波滤波’);
disp(‘-------------------------------------------’);
choose=input(‘选择滤波方式choose=’);
function [s, err_mse, iter_time]=sparseOMP(x,A,m,varargin)
[n1 n2]=size(x);
if n2 == 1
n=n1;
elseif n1 == 1
x=x’;
n=n2;
else
error(‘x must be a vector.’);
end

sigsize = x’*x/n;
initial_given=0;
err_mse = [];
iter_time = [];
STOPCRIT = ‘M’;
STOPTOL = ceil(n/4);
MAXITER = n;
verbose = false;
s_initial = zeros(m,1);
GradSteps = ‘auto’;
alpha = 1;
weakness = 1;

if verbose
display(‘Initialising…’)
end

switch nargout
case 3
comp_err=true;
comp_time=true;
case 2
comp_err=true;
comp_time=false;
case 1
comp_err=false;
comp_time=false;
case 0
error(‘Please assign output variable.’)
otherwise
error(‘Too many output arguments specified’)
end

% Put option into nice format
Options={};
OS=nargin-3;
c=1;
for i=1:OS
if isa(varargin{i},‘cell’)
CellSize=length(varargin{i});
ThisCell=varargin{i};
for j=1:CellSize
Options{c}=ThisCell{j};
c=c+1;
end
else
Options{c}=varargin{i};
c=c+1;
end
end
OS=length(Options);
if rem(OS,2)
error(‘Something is wrong with argument name and argument value pairs.’)
end
for i=1:2:OS
switch Options{i}
case {‘stopCrit’}
if (strmatch(Options{i+1},{‘M’; ‘corr’; ‘mse’; ‘mse_change’},‘exact’));
STOPCRIT = Options{i+1};
else error(‘stopCrit must be char string [M, corr, mse, mse_change]. Exiting.’); end
case {‘stopTol’}
if isa(Options{i+1},‘numeric’) ; STOPTOL = Options{i+1};
else error(‘stopTol must be number. Exiting.’); end
case {‘P_trans’}
if isa(Options{i+1},‘function_handle’); Pt = Options{i+1};
else error(‘P_trans must be function _handle. Exiting.’); end
case {‘maxIter’}
if isa(Options{i+1},‘numeric’); MAXITER = Options{i+1};
else error(‘maxIter must be a number. Exiting.’); end
case {‘wf’}
if isa(Options{i+1},‘numeric’); alpha = Options{i+1};
if alpha <1 weakness =0; else alpha =1; weakness = 1; end
else error(‘wf must be a number. Exiting.’); end
case {‘verbose’}
if isa(Options{i+1},‘logical’); verbose = Options{i+1};
else error(‘verbose must be a logical. Exiting.’); end
case {‘start_val’}
if isa(Options{i+1},‘numeric’) && length(Options{i+1}) == m ;
s_initial = Options{i+1};
initial_given=1;
else error(‘start_val must be a vector of length m. Exiting.’); end
case {‘GradSteps’}
if isa(Options{i+1},‘numeric’) || strcmp(Options{i+1},‘auto’) ;
GradSteps = Options{i+1};
else error(‘start_val must be a vector of length m. Exiting.’); end
otherwise
error(‘Unrecognised option. Exiting.’)
end
end

if strcmp(STOPCRIT,‘M’)
maxM=STOPTOL;
else
maxM=MAXITER;
end

if nargout >=2
err_mse = zeros(maxM,1);
end
if nargout ==3
iter_time = zeros(maxM,1);
end

if isa(A,‘float’) P =@(z) A*z; Pt =@(z) A’z;
elseif isobject(A) P =@(z) Az; Pt =@(z) A’*z;
elseif isa(A,‘function_handle’)
try
if isa(Pt,‘function_handle’); P=A;
else error(‘If P is a function handle, Pt also needs to be a function handle. Exiting.’); end
catch error(‘If P is a function handle, Pt needs to be specified. Exiting.’); end
else error(‘P is of unsupported type. Use matrix, function_handle or object. Exiting.’); end

if initial_given ==1;
IN = find(s_initial);
Residual = x-P(s_initial);
s = s_initial;
oldERR = Residual’*Residual/n;
else
IN = [];
Residual = x;
s = s_initial;
sigsize = x’*x/n;
oldERR = sigsize;
end

    mask=zeros(m,1);
    mask(ceil(rand*m))=1;
    nP=norm(P(mask));
    if abs(1-nP)>1e-3;
        display('Dictionary appears not to have unit norm columns.')
    end

if verbose
display(‘Main iterations…’)
end
tic
t=0;
normA=(sum(A.^2,1)).0.5;
NI = 1:size(A,2);
p=zeros(m,1);
DR=Pt(Residual);
[v I]=max(abs(DR(NI))./normA(NI)');
INI = NI(I);
IN=[IN INI];
NI(I) = [];
if weakness ~= 1
[vals inds] = sort(abs(DR),‘descend’);
I=inds( find( vals >= alpha * v ) );
end

IN = union(IN,I);
if strcmp(STOPCRIT,‘M’) & length(IN) >= STOPTOL
IN=IN(1:STOPTOL);
end
MASK=zeros(size(DR));
pDDp=1;
done = 0;
iter=1;

while ~done

% Select new element
if isa(GradSteps,'char')
    if strcmp(GradSteps,'auto')

% finished=0;
% while ~finished
% Update direction
if iter==1
p(IN)=DR(IN);
Dp=P§;
else
MASK(IN)=1;
PDR=P(DR.MASK);
b=-Dp’PDR/pDDp;
p(IN)=DR(IN) +bp(IN);
Dp=PDR +b Dp;
end
% Step size
% Dp=P§; % =P(DR(IN)) +b P(p(IN));
pDDp=Dp’Dp;
a=Residual’Dp/(pDDp);
% Update coefficients
s=s+ap;
% New Residual and inner products
Residual=Residual-aDp;
DR=Pt(Residual);
% select new element
[v I]=max(abs(DR(NI))./normA(NI)');
INI = NI(I);
if weakness ~= 1
[vals inds] = sort(abs(DR),‘descend’);
I=inds( find( vals >= alpha * v ) );
end
IN = union(IN,INI);
NI(I) = [];
if strcmp(STOPCRIT,‘M’) & length(IN) >= STOPTOL
IN=IN(1:STOPTOL);
end

    else
        
        
        error('Undefined option for GradSteps, use ''auto'' or an integer.')
    end
elseif isa(GradSteps,'numeric') 

            
    % Do GradSteps gradient steps
    count=1;
    while count<=GradSteps
        % Update direction    
             if iter==1
                 p(IN)=DR(IN);
                 Dp=P(p);
             else
                 MASK(IN)=1;
                 PDR=P(DR.*MASK);
                 b=-Dp'*PDR/pDDp;
                 p(IN)=DR(IN) +b*p(IN);
                 Dp=PDR +b* Dp;
             end
         % Step size

% Dp=P§;
pDDp=Dp’Dp;
a=Residual’Dp/(pDDp);
% Update coefficients
s=s+ap;
% New Residual and inner products
Residual=Residual-aDp;
DR=Pt(Residual);
count=count+1;
end
% select new element
[v I]=max(abs(DR(NI))./normA(NI)');
INI = NI(I);
if weakness ~= 1
[vals inds] = sort(abs(DR),‘descend’);
I=inds( find( vals >= alpha * v ) );
end
IN = union(IN,INI);
NI(I) = [];
if strcmp(STOPCRIT,‘M’) & length(IN) >= STOPTOL
IN=IN(1:STOPTOL);
end

 else
      error('Undefined option for GradSteps, use ''auto'' or an integer.')
 end


 ERR=Residual'*Residual/n;
 if comp_err
     err_mse(iter)=ERR;
 end
 
 if comp_time
     iter_time(iter)=toc;
 end


 if strcmp(STOPCRIT,'M')
     if length(IN) >= STOPTOL
         done =1;
     elseif verbose && toc-t>10
        display(sprintf('Iteration %i. --- %i selected elements',iter ,length(IN))) 
        t=toc;
     end
elseif strcmp(STOPCRIT,'mse')
     if comp_err
        if err_mse(iter)<STOPTOL;
            done = 1; 
        elseif verbose && toc-t>10
            display(sprintf('Iteration %i. --- %i mse',iter ,err_mse(iter))) 
            t=toc;
        end
     else
         if ERR<STOPTOL;
            done = 1; 
         elseif verbose && toc-t>10
            display(sprintf('Iteration %i. --- %i mse',iter ,ERR)) 
            t=toc;
         end
     end
 elseif strcmp(STOPCRIT,'mse_change') && iter >=2
     if comp_err && iter >=2
          if ((err_mse(iter-1)-err_mse(iter))/sigsize <STOPTOL);
            done = 1; 
         elseif verbose && toc-t>10
            display(sprintf('Iteration %i. --- %i mse change',iter ,(err_mse(iter-1)-err_mse(iter))/sigsize )) 
            t=toc;
         end
     else
         if ((oldERR - ERR)/sigsize < STOPTOL);
            done = 1; 
         elseif verbose && toc-t>10
            display(sprintf('Iteration %i. --- %i mse change',iter ,(oldERR - ERR)/sigsize)) 
            t=toc;
         end
     end
 elseif strcmp(STOPCRIT,'corr') 
      if max(abs(DR)) < STOPTOL;
         done = 1; 
      elseif verbose && toc-t>10
            display(sprintf('Iteration %i. --- %i corr',iter ,max(abs(DR)))) 
            t=toc;
      end
 end
 
% Also stop if residual gets too small or maxIter reached
 if comp_err
     if err_mse(iter)<1e-16
         display('Stopping. Exact signal representation found!')
         done=1;
     end
 else


     if iter>1
         if ERR<1e-16
             display('Stopping. Exact signal representation found!')
             done=1;
         end
     end
 end

 if iter >= MAXITER
     display('Stopping. Maximum number of iterations reached!')
     done = 1; 
 end
 

 if ~done
    iter=iter+1;
    oldERR=ERR;
 end

end

if nargout >=2
err_mse = err_mse(1:iter);
end
if nargout ==3
iter_time = iter_time(1:iter);
end
if verbose
display(‘Done’)
end