一、理论基础

肌电信号(EMG)是一种非常微弱的生物电信号，它与神经肌肉活动密切相关，其中蕴涵着很多与肢体运动相关联的信息。肌电信号可通过表面肌电拾取电极或针式肌电拾取电极加以引导、记录，通过表面电极拾取的肌电信号称之为表面肌电信号，而通过针电极拾取的肌电信号称之为针肌电信号。目前，肌电信号已经深入应用到临床医学、运动医学、生物医学工程等领域。特征提取是肌电信号分析的基础，肌电分析被广泛地应用于肌肉疲劳、重症肌无力、肌强直和肌萎缩等各种肌肉疾病的临床诊断；肌电图也是运动员训练中动作和强度分析的依据；另外，还可利用人体表面肌电的某些特征进行模式分类，进而驱动假肢或其它运动机械(如机器人)的各种动作，实现仿生控制。由此可见，提取有效的肌电信号特征对肌电信号分析意义重大。随着检测技术、信息处理技术及计算机技术的发展，使肌电信号的定量分析及更深入的研究成为可能。

二、案例背景

本项目的目的是开发利用小波分析处理针肌电信号的信号处理方法，这将有助于神经肌肉疾病患者的诊断。

使用小波分析的方法分析肌电图信号，从而实现病人数据的分析。

临床肌电图是评估神经肌肉疾病患者的标准方法。目前没有自动诊断患者的方法，因此诊断必须由临床医生通过目视检查信号或在扬声器上收听其属性来确定。我们已经研究了小波变换作为一种从信号中获得比仅时域频域（FFT）分析更多信息的手段。组合的时频信息为创建患者诊断方法提供了真正的可能性。我们以前发表过几篇关于这个问题的文章。

此段说明，采用传统的快速傅立叶变换方法已经无法达到所要的目的。

有许多不同的小波可用。该项目将首先使用B-sline小波，因为它在时间局部化方面具有理想的特性，但许多包括其他小波。到目前为止，我们一直在使用复杂的Labview程序在PC上收集肌电图数据。新的EMG机器将允许在机器本身上收集数据，但我们需要一个带有医生界面的Labview程序来控制信号向患者的呈现，以便在最大自主收缩（MVC）和大约50%的MVC时收集信号。我们希望对两种信号的比较分析将揭示患者和正常受试者之间的差异。

这段的要求，就是使用B样条小波分析/或其余小波分析的方法来分析肌电图数据；

该项目将从为医生设置控制程序开始。一旦在阿伯丁收集到数据，这些数据将被匿名并发送到邓迪，在那里将使用小波分析方法进行分析。学生将开发分析方法。分析程序将用Matlab编写。该项目将适合具有编程经验和信号处理知识的人员。小波理论需要良好的数学技能才能理解使用MATLAB编写相关代码来实现分析。

三、matlab程序

clc;
clear;
close all;
warning off;[signals,Fs]=wavread('data\triceps1.wav');%由于数据采样率太高，对于数据较大的数据，只选择其中一段
if length(signals) > 1000000signal = signals(1:1000000);
elsesignal = signals;
end
clear signals;figure;
plot(signal);title('Initial data');
save mat\signal.mat signal
save mat\Fs.mat Fs
clear signal%Step1:利用小波变换做预处理，去噪
%Step1:利用小波变换做预处理，去噪
%Step1:利用小波变换做预处理，去噪
%画出原始信号
load mat\signal.mat
s = signal';
figure;
subplot(221); plot(s); title('原始信号');grid;
%用db1小波对原始信号进行3层分解并提取系数
[c,l]=wavedec(s,3,'db1');
ca3=appcoef(c,l,'db1',3);
cd3=detcoef(c,l,3);
cd2=detcoef(c,l,2);
cd1=detcoef(c,l,1);
%对信号进行强制性消噪处理并图示结果
cdd3=zeros(1,length(cd3));
cdd2=zeros(1,length(cd2));
cdd1=zeros(1,length(cd1));
c1=[ca3 cdd3 cdd2 cdd1];
s1=waverec(c1,l,'db1');
subplot(2,2,2); plot(s1); title('强制消噪后的信号');grid;
%用ddencmp函数获得信号的默认阈值,使用wdencmp命令函数实现消噪过程
[thr,sorh,keepapp]=ddencmp('den','wv',s);
s2=wdencmp('gbl',c,l,'db1',3,thr,sorh,keepapp);
subplot(2,2,3); plot(s2);title('默认阈值消噪后的信号');grid;
%用给定的软阈值进行消噪处理
cd1soft=wthresh(cd1,'s',1.465);
cd2soft=wthresh(cd2,'s',1.823);
cd3soft=wthresh(cd3,'s',2.768);
c2=[ca3 cd3soft cd2soft cd1soft];
s3=waverec(c2,l,'db1');
subplot(2,2,4); plot(s3); title('给定软阈值消噪后的信号');grid
set(gcf,'color','w');save mat\delete.mat s s1 s2 s3
clear all%Step2:直接对信号进行功率谱的计算分析
%Step2:直接对信号进行功率谱的计算分析
%Step2:直接对信号进行功率谱的计算分析
load mat\delete.mat
load mat\Fs.mat
signals = s2;
figure;
plotspec(signals,1/Fs)
set(gcf,'color','w');
save mat\signals.mat signals
clear all%Step3:对该信号进行小波变换
%Step3:对该信号进行小波变换
%Step3:对该信号进行小波变换
%db10小波进行4层分解
%一维小波分解
load mat\delete.mat
load mat\Fs.mat
load mat\signals.mattype    = 'db5';
n       = 7;
signals = s2;
[c,l]   = wavedec(signals,n,type);
%分解
a7 = wrcoef('a',c,l,type,7);
a6 = wrcoef('a',c,l,type,6);
a5 = wrcoef('a',c,l,type,5);
a4 = wrcoef('a',c,l,type,4);
a3 = wrcoef('a',c,l,type,3);
a2 = wrcoef('a',c,l,type,2);
a1 = wrcoef('a',c,l,type,1);
%重构第1～4层细节系数
d7 = wrcoef('d',c,l,type,7);
d6 = wrcoef('d',c,l,type,6);
d5 = wrcoef('d',c,l,type,5);
d4 = wrcoef('d',c,l,type,4);
d3 = wrcoef('d',c,l,type,3);
d2 = wrcoef('d',c,l,type,2);
d1 = wrcoef('d',c,l,type,1);
%显示细节信号
figure
subplot(5,3,1);plot(d7);title('d7');
subplot(5,3,2);plot(d6);title('d6');
subplot(5,3,3);plot(d5);title('d5');
subplot(5,3,4);plot(d4);title('d4');
subplot(5,3,5);plot(d3);title('d3');
subplot(5,3,6);plot(d2);title('d2');
subplot(5,3,7);plot(d1);title('d1');
subplot(5,3,9) ;plot(a7);title('a7');
subplot(5,3,10);plot(a6);title('a6');
subplot(5,3,11);plot(a5);title('a5');
subplot(5,3,12);plot(a4);title('a4');
subplot(5,3,13);plot(a3);title('a3');
subplot(5,3,14);plot(a2);title('a2');
subplot(5,3,15);plot(a1);title('a1');
set(gcf,'color','w');
save mat\d.mat d1 d2 d3 d4 d5 d6 d7
clear allload mat\delete.mat
load mat\signals.mat
load mat\d.mat
load mat\Fs.mat
type    = 'db5';
n       = 7;%7层小波包分解
T=wpdec(s2,n,type);
%重构低频信号
y1=wprcoef(T,[7,1]);
figure;
plotspec(y1,1/Fs)
set(gcf,'color','w');
%第1层的频谱图
figure;
plotspec(d1,1/Fs);title('d1');
%第2层的频谱图
figure;
plotspec(d2,1/Fs);title('d2');
set(gcf,'color','w');
%第3层的频谱图
figure;
plotspec(d3,1/Fs);title('d3');
set(gcf,'color','w');
%第4层的频谱图
figure;
plotspec(d4,1/Fs);title('d4');
set(gcf,'color','w');
%第5层的频谱图
figure;
plotspec(d5,1/Fs);title('d5');
set(gcf,'color','w');
%第6层的频谱图
figure;
plotspec(d6,1/Fs);title('d6');
set(gcf,'color','w');
%第7层的频谱图
figure;
plotspec(d7,1/Fs);title('d7');
set(gcf,'color','w');clear alltype    = 'db5';
n       = 7;
load mat\delete.mat%计算各层能量分布
for i=1:nwpt=wpdec(s2,i,type);   e=wenergy(wpt);E=zeros(1,length(e));for j=1:2^iE(j)=sum(abs(wprcoef(wpt,[i,j-1])).^2);endfigure(20);subplot(7,1,i);bar(e);axis([0 length(e) 0 130]);title(['第 ',num2str(i), ' 层']);
endclear all;

clc;
clear;
close all;
warning off;[signals,Fs]=wavread('data\EDC3.wav');%由于数据采样率太高，对于数据较大的数据，只选择其中一段
if length(signals) > 1000000signal = signals(1:1000000);
elsesignal = signals;
end
clear signals;
x=signal;l0=length(x);%计算B样条小波函数的尺度特性曲线
m     = 3;%B样条阶数
N     = func_bspline(2*m,2*m+1);%求样条函数在整数点的值
w     = 0:0.1*pi:10*pi;
%求尺度函数
mea_w = fun_bspline_wavelet(m,w,N);
figure;
subplot(131);
plot(w/pi,abs(mea_w));title('尺度函数的幅频图');mea_ang=atan(imag(mea_w)./(real(mea_w)+eps));%计算相位角
subplot(132),plot(w/pi,mea_ang);title('尺度函数的相频图');%求反傅立叶变换，求出尺度函数的时域函数。
mea_x=ifft(mea_w);
mea_x=ifftshift(mea_x);
subplot(1,3,3);
plot(w/pi,real(mea_x));title('尺度函数的时域图');
set(gcf,'color','w');x=x';
%三次样条小波的滤波器系数
ld=[0,0.0625,0.25,0.375,0.25,0.0625,0,0];
hd=[-0.00008,-0.01643,-0.10872,-0.59261,0.59261,0.10872,0.01643,0.00008];   %第1层小波分解
ld=rot90(ld,2);
hd=rot90(hd,2);
a1=wconv1(x,ld);
d1=wconv1(x,hd);
l1=length(a1);
a1=a1([l1-l0+1:l1]);
d1=d1([l1-l0+1:l1]);
a1=dyaddown(a1);
d1=dyaddown(d1);
l1=length(a1);
%第2层小波分解
a2=wconv1(a1,ld);
d2=wconv1(d1,hd);
l2=length(a2);
a2=a2([l2-l1+1:l2]);
d2=d2([l2-l1+1:l2]);
a2=dyaddown(a2);
d2=dyaddown(d2);
l2=length(a2);
%第3层小波分解
a3=wconv1(a2,ld);
d3=wconv1(d2,hd);
l3=length(a3);
a3=a3([l3-l2+1:l3]);
d3=d3([l3-l2+1:l3]);
a3=dyaddown(a3);
d3=dyaddown(d3);
l3=length(a3);
%第4层小波分解
a4=wconv1(a3,ld);
d4=wconv1(d3,hd);
l4=length(a4);
a4=a4([l4-l3+1:l4]);
d4=d4([l4-l3+1:l4]);
a4=dyaddown(a4);
d4=dyaddown(d4);
l4=length(a4);%显示细节信号
figure
subplot(421);plot(a4);title('a4');
subplot(422);plot(a3);title('a3');
subplot(423);plot(a2);title('a2');
subplot(424);plot(a1);title('a1');
subplot(425);plot(d4);title('d4');
subplot(426);plot(d3);title('d3');
subplot(427);plot(d2);title('d2');
subplot(428);plot(d1);title('d1');set(gcf,'color','w');%第1层的频谱图
figure;
plotspec(d1,1/Fs);title('d1');
set(gcf,'color','w');
%第2层的频谱图
figure;
plotspec(d2,1/Fs);title('d2');
set(gcf,'color','w');
%第3层的频谱图
figure;
plotspec(d3,1/Fs);title('d3');
set(gcf,'color','w');
%第4层的频谱图
figure;
plotspec(d4,1/Fs);title('d4');
set(gcf,'color','w');

四、仿真结论

·输入的信号：

·通过小波变换的方法去噪声，其仿真结果如下所示：

这里给出了三种不同方法进行滤波了，通过上面的仿真结果可知：

从图可以看出，经全局阈值和分层阈值方法降噪后的信号都很好的保留了信号发展初期的高频特性。分层阈值法虽然损失了部分性能（与原信号的相似性：经计算，分层阈值降噪信号与原信号的标准差较全局阈值法较大），但降噪信号更光滑，而其信号发展初期的高频特性也几乎不受影响，最大限度地反映了原信号本身的特性。

·首先对信号直接进行频谱分析

正如题目要求可知，如果对该信号直接采用fft，无法判断出信号的具体特征，下面的仿真通过小波变换进行分析。

·小波变换

d1~d7为1~7层的细节信号，a1~a7为1~7层的逼近信号。

·对1~7层分别进行频谱分析

得到了如下的仿真结果：

·小波变换

分别对d1~d4进行频谱分析

所用到的MATLAB内部函数简介：

[C,L]=wavedec(X,N,'wname')

利用小波'wname'对信号X进行多层分解

A=appcoef(C,L,'wname',N)

利用小波'wname'从分解系数[C,L]中提取第N层近似系数

[C,L]=wavedec(X,1,’wname’)中返回的近似和细节都存放在C中，即C=[cA,cD]，L存放是近似和各阶细节系数对应的长度。

A09-02

基于小波变换的EMG信号病人数据matlab仿真分析相关推荐

【小波滤波】基于小波变换的噪声信号滤波处理matlab仿真
1.软件版本 MATLAB2021a 2.核心代码 % 小波分解与程序,Xk0是要分解的原始信号,step是表示要分解的层数 function [Xh,D]=decomposition(Xk0,ste ...
MATLAB基于小波变换的语音信号去噪算法改进
MATLAB基于小波变换的语音信号去噪算法改进概述 0. 需要调用的子函数 0.1 Gnoisegen函数 0.2 snrr函数 1. 语音信号输入和加噪 1.1 语音信号输入 1.2 语音信号加噪 ...
基于GRNN网络和小波变换的ECG信号睡眠监测matlab仿真
目录一.理论基础二.核心程序三.仿真测试结果作者ID :fpga和matlab CSDN主页:https://blog.csdn.net/ccsss22?type=blog 擅长技术: 1.无 ...
基于小波变换的脉搏信号滤波matlab仿真
目录 1.算法仿真效果 2.MATLAB源码 3.算法概述 4.部分参考文献 1.算法仿真效果 matlab2022a仿真结果如下:
【信号处理】基于小波变换的音频水印嵌入提取matlab源码
较早利用分块DCT的水印技术,他们的水印方案是用一个密钥随机的选择图像的一些分块,在频域的中频上稍稍改变一个三元组来隐藏二进制序列信息.这种方法对有损压缩和低通滤波是稳健的.Cox等[提出了着名的基于 ...
【智能优化算法】基于遗传算法实现城市交通信号优化附matlab代码
1 简介本文设计实时优化的配置方案对道路畅通的应急决策管理具有重要意义.本文在分析交通控制基本理论的基础上,根据交叉口的实际情况并考虑信号灯的转换与车辆的启动损失时间,采用四相位对称式放行方案,以车 ...
【故障诊断分析】基于小波变换实现外圈轴承故障诊断含Matlab源码
1 简介在滚动轴承的故障诊断时,传统的频谱分析法通常采用共振解调技术 , 具有良好的效果 ,但当内圈 . 滚动体或多点故障时,解调谱线却很难分辨故障类型小波包是小波理论在信号处理应用领域的又一重大发 ...
【信号去噪】基于NLM时间序列心电信号去噪附matlab代码
1 简介作为一种信号预处理手段,信号去噪在众多信号处理应用中发挥着重要的作用.到目前为止,信号去噪问题被大量研究,并取得了许多重要成果,涌现出了包括非局部均值(NLM)去噪算法在内的一批优秀的去噪方 ...
码元速率 matlab,[转载]基于小波变换的移相键控信号符号速率估计（matlab仿真）...
西安电子科技大学西电大宝在电子技术迅猛发展的当代社会,空间中充满了各种各样不同频率.不同调制类型的通信信号.正常通信条件下,发送方和接受方进行的是合作通信,即接受方预先知道发送信号的频率.调制类 ...

基于小波变换的EMG信号病人数据matlab仿真分析

一、理论基础

二、案例背景

三、matlab程序

四、仿真结论

基于小波变换的EMG信号病人数据matlab仿真分析相关推荐

最新文章

热门文章