To extract the frequency domain characteristics from your EMG signal in MATLAB, follow the steps below:
1. Apply a bandpass filter (e.g., 20–450 Hz) to your raw EMG signal using the “butter” and “filtfilt” functions to remove noise and artifacts.
[b, a] = butter(4, [20 450]/(fs/2), 'bandpass'); %Example code
emg_filtered = filtfilt(b, a, emg_raw);
2. Compute the Fast Fourier Transform (FFT) of the filtered signal using the “fft” function to obtain the frequency and power spectrum.
3. Extract the positive frequencies by indexing the first half of the FFT output.
N = length(emg_filtered); %Example Code
emg_fft = fft(emg_filtered);
psd = abs(emg_fft).^2 / N;
f = (0:N-1) * (fs/N);
% Keep only the positive frequencies
half_N = floor(N/2) + 1;
f = f(1:half_N);
psd = psd(1:half_N);
4. Calculate the mean frequency using basic arithmetic operations:
sum(f .* psd) / sum(psd)
5. Find the median frequency by computing the cumulative sum with “cumsum” and using “interp1” to find the frequency at which the cumulative power reaches half the total power.
6. Identify the peak frequency by finding the maximum value in the power spectrum with the “max” function.
7. Compute the total power by summing the power spectrum using the sum function.
% Example Code
% Mean Frequency
mean_freq = sum(f .* psd) / sum(psd);
% Median Frequency
cumulative_power = cumsum(psd);
median_freq = interp1(cumulative_power, f, cumulative_power(end)/2);
% Peak Frequency
[~, idx_max] = max(psd);
peak_freq = f(idx_max);
% Total Power
total_power = sum(psd);
For more information on the functions used, please refer to the below documentation links:
- “butter”: https://www.mathworks.com/help/signal/ref/butter.html
- “fft”: https://www.mathworks.com/help/matlab/ref/fft.html
- ”cumsum”: https://www.mathworks.com/help/matlab/ref/double.cumsum.html
- “max”: https://www.mathworks.com/help/matlab/ref/double.max.html
Hope this is helpful!
