EMG data comparison problem
19 views (last 30 days)
Show older comments
I have some sEMG data sets that I want to compare. How can I do this with Matlab?
%% Analysis of EMG data
%An EMG signal has a frequency between 0 to 500 Hz; but the usable
%frequency is between 50 and 150 Hz.
% EMG signal is produced using an instrument called the electromyograph that records electric activity
% produced by muscles.
% The amplitude for surface EMG signals is within the range of 0 to 10 mV.
% sEMG (surace EMG) signals tend to have a DC bias, therefore a high pass
% filter is applied to a cutoff frequency of 50Hz. This high pass filter
% removes any DC bias.
% The root mean-square (RMS): calculated by finding the square root of: squaring each data point,
%summing the squares, dividing the sum by the number of observations.
% RMS value of a signal is a measure of the power of the signal and it
% reflects on the activity in the motor unit during contraction.
% This makes RMS useful in producing waveforms.
% Amplitude of an EMG signal is random with Gaussian Distribution/
% RMS and MA are two parameters used to measure the amplitude.
% MA short for: Mean Absolute value: is a mean of the absolute value of
% EMG signal
% MA is a measure of the area under the rectified EMG.
% In order to calculate MA remove all of the negative phases of the raw EMG
% This is known as: half-wave rectification.
% The MA and RMS provide an estimate of the amplitude of the raw EMG signal.
clear,clc
load('A1_Ball_high_t2.mat')
t = t2(:,1);
y1 = t2(:,2);
N = length(y1);% find the length of the data per second
ls = size(y1); %% size
f = 1/N;% find the sampling rate or frequency
%Note: sampling frequency should be at least twice the incoming signal.
fs = 2000;
T = 1/fs % period between each sample
t1 = (0 : N-1) *T;
t = (0:1:length(y1)-1)/fs; % sampling period
% nyquist frequency is half of the sampling rate of a signal. Also known as:
% folding frequency.
Nyquist = fs/2;
figure;
%
subplot (3,1,1), plot(t,y1,'b');
title ('EMG signal of single muscle 25 years old patient ');
xlabel ('time (sec)');
ylabel ('Amplitute (V)');
grid on;
Y= abs(fft(y1));
Y(1) = [];
power = abs(Y(1:N/2)).^2;
nyquist = 350;
freq = (1:N/2)/(N/2)*nyquist;
subplot(212), plot(freq,power), grid on
xlabel('Sample number (in Frequency)')
ylabel('Power spectrumen');
title({'Single-sided Power spectrum' ...
' (Frequency in shown on a log scale)'});
axis tight
%%% RMS of the signal
rms_y1 = sqrt(mean(y1.^2));
msgbox(strcat('RMS of EMG signal is = ',mat2str(rms_y1), ''));
rms_emg = rms (y1);
%%%%% MA of the signal
ma_y1 = abs(mean(y1));
msgbox(strcat('MA of EMG signal is = ',mat2str(ma_y1), ''));
%remove any DC offset of the signal
% use the detrend function.
% It subtracts the mean or a best-fit line from data.
y2 = detrend(y1); %% y2 is the singal after DC offset removed.
figure;
rec_y = abs(y2);
plot (rec_y);
xlabel('Sample number (in Frequency)')
ylabel('Rectified EMG signal');
title({'Rectified EMG signal' ...
' (Frequency shown on a log scale)'});
figure;
xdft = fft(y1);
xdft = xdft(1:N/2+1);
psdx = (1/(fs*N)) * abs(xdft).^2;
psdx(2:end-1) = 2*psdx(2:end-1);
freq = 0:fs/length(y1):fs/2;
plot(freq,10*log10(psdx))
grid on
title(' Power spectrum FFT')
xlabel('Frequency (Hz)')
ylabel('Power/Frequency (dB/Hz)')
%amplitude of the entire signal
amplitude = max(y2(:)) - min(y2(:))
%%% DFT to describe the signal in the frequency
NFFT = 2 ^ nextpow2(N);
% The reason I used this function above:
%The execution time for fft depends on the length of the transform.
% It is fastest for powers of two. It is almost as fast for lengths that
%have only small prime factors. It is typically several times slower for
%lengths that are prime or which have large prime factors.
Y = fft(y1, NFFT) / N;
f = (fs / 2 * linspace(0, 1, NFFT / 2+1))'; % Vector containing frequencies in Hz
amp = ( 2 * abs(Y(1: NFFT / 2+1))) % Vector containing corresponding amplitudes
figure;
plot (f, amp);
title ('plot single-sided amplitude spectrum of the EMG signal')
xlabel ('frequency (Hz)')
ylabel ('|y(f)|')
grid on;
References:
Ali H. Al-Timemy, Rami N. Khushaba, Guido Bugmann, and Javier Escudero, " Improving the Performance against Force Variation of EMG Controlled Multifunctional Upper-Limb Prostheses for Transradial Amputees", IEEE Transactions on Neural Systems and Rehabilitation Engineering, Accepted for Publication, 2015.
http://biomedicalsignalandimage.blogspot.com/2016/02/matlab-code-to-study-emg-signal.html
0 Comments
Answers (1)
Francisco Romero
on 11 Feb 2022
Hi José,
I saw to this question by looking for other things.
You may find this link useful:
0 Comments
See Also
Categories
Find more on Spectral Estimation in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!