Python: Analysing EMG signals – Part 4
We have seen how Python can be used to process and analyse EMG signals in lessons 1, 2 and 3. When EMG signals are filtered, how does changing filter settings change the appearance of the filtered EMG signal?
A low pass filter allows frequencies below the cut-off frequency to pass through (ie. higher frequencies are removed). The filtered EMG signal can be used to interpret different physiological properties. For example, scientists investigating muscle force and muscle activity often use a low pass filter to capture the shape or “envelope” of the EMG signal as this is thought to better reflect force generated by a muscle. Let’s see how changing the low pass cut-off frequency changes the rectified signal.
We want to write a bit of code to filter and rectify an EMG signal and plot graphs comparing the unfiltered and filtered signal, for different low pass cut-off frequencies. Since this involves doing things repeatedly (eg. EMG processing, plotting) for different values (low pass cut-offs), we will write a function to do the EMG processing etc. and send the function different values in a
for loop. For quick refreshers, see this series on functions and this post on loops in Python.
Here is our function and the
for loop at the end. Key bits of code are discussed below:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72
def filteremg(time, emg, low_pass=10, sfreq=1000, high_band=20, low_band=450): """ time: Time data emg: EMG data high: high-pass cut off frequency low: low-pass cut off frequency sfreq: sampling frequency """ # normalise cut-off frequencies to sampling frequency high_band = high_band/(sfreq/2) low_band = low_band/(sfreq/2) # create bandpass filter for EMG b1, a1 = sp.signal.butter(4, [high_band,low_band], btype='bandpass') # process EMG signal: filter EMG emg_filtered = sp.signal.filtfilt(b1, a1, emg) # process EMG signal: rectify emg_rectified = abs(emg_filtered) # create lowpass filter and apply to rectified signal to get EMG envelope low_pass = low_pass/(sfreq/2) b2, a2 = sp.signal.butter(4, low_pass, btype='lowpass') emg_envelope = sp.signal.filtfilt(b2, a2, emg_rectified) # plot graphs fig = plt.figure() plt.subplot(1, 4, 1) plt.subplot(1, 4, 1).set_title('Unfiltered,' + '\n' + 'unrectified EMG') plt.plot(time, emg) plt.locator_params(axis='x', nbins=4) plt.locator_params(axis='y', nbins=4) plt.ylim(-1.5, 1.5) plt.xlabel('Time (sec)') plt.ylabel('EMG (a.u.)') plt.subplot(1, 4, 2) plt.subplot(1, 4, 2).set_title('Filtered,' + '\n' + 'rectified EMG: ' + str(int(high_band*sfreq)) + '-' + str(int(low_band*sfreq)) + 'Hz') plt.plot(time, emg_rectified) plt.locator_params(axis='x', nbins=4) plt.locator_params(axis='y', nbins=4) plt.ylim(-1.5, 1.5) plt.plot([0.9, 1.0], [1.0, 1.0], 'r-', lw=5) plt.xlabel('Time (sec)') plt.subplot(1, 4, 3) plt.subplot(1, 4, 3).set_title('Filtered, rectified ' + '\n' + 'EMG envelope: ' + str(int(low_pass*sfreq)) + ' Hz') plt.plot(time, emg_envelope) plt.locator_params(axis='x', nbins=4) plt.locator_params(axis='y', nbins=4) plt.ylim(-1.5, 1.5) plt.plot([0.9, 1.0], [1.0, 1.0], 'r-', lw=5) plt.xlabel('Time (sec)') plt.subplot(1, 4, 4) plt.subplot(1, 4, 4).set_title('Focussed region') plt.plot(time[int(0.9*1000):int(1.0*1000)], emg_envelope[int(0.9*1000):int(1.0*1000)]) plt.locator_params(axis='x', nbins=4) plt.locator_params(axis='y', nbins=4) plt.xlim(0.9, 1.0) plt.ylim(-1.5, 1.5) plt.xlabel('Time (sec)') fig_name = 'fig_' + str(int(low_pass*sfreq)) + '.png' fig.set_size_inches(w=11,h=7) fig.savefig(fig_name) # show what different low pass filter cut-offs do for i in [3, 10, 40]: filteremg(time, emg_correctmean, low_pass=i)
Line 1. Define a function called
filteremg to accept
emg values (for plotting on the x- and y-axes) with default values for a low pass filter (for the EMG envelope), sampling frequency, and high and low pass filter frequencies. The default values mean that when we call the function, if we don’t specify different values, the function uses the defaults during code execution.
Line 28 onwards. Plot 3 subplots to see (1) the unfiltered, unrectified EMG signal, (2) the filtered, rectified signal, (3) the rectified signal with a low pass filter to get the EMG envelope and (4) a zoomed-in section of the signal from (3) over the time period indicated by the red line to see the underlying shape of the final signal. This line says to index values from 0.9 to 1.0 seconds, convert these numerical float values to integers, and plot time on the x- and EMG signal on the y-axis.
Line 70. Write a
for loop that says: for each low pass cut-off value in the list
[3, 10, 40], when
i is 3, call the function
filteremg and assign the value of
i to the local variable
low. When i is 10, call the function and do the same. Etc.
For each of the cut-off frequencies supplied, we see in Figure 1 below that low cut-off frequencies produce a “smoothed” EMG signal which has a more noticable shape (envelope) compared to when higher cut-off frequencies are used.
There is still discussion around which lowpass cut-off frequencies are appropriate to get the EMG envelope under different conditions (eg. sustained contractions, gait). The biomechanist David Winter argues that choosing a lowpass cut-off value should be based on how the biomechanics of production of muscle force is modelled; see here for a discussion on the linear EMG envelope.
Changing the cut-off frequencies changes the shape and appearance of the filtered EMG signal. We now have a function that can be used to see how changing these frequencies changes the appearance of the signal. As a learning activity, stick the function into an IPython notebook (or similar) and play around with the filter settings to see just how much the shape can change.