how to do digital filter on this EEG data

Illustration
Delilah - 2021-03-05T10:35:30+00:00
Question: how to do digital filter on this EEG data

Hi! I recently recorded an EEG  signal at 50K Hz sampling rate. The signal was band passed at 1-500Hz before it was digitalized. Now I wanted to do a digital fitlering to this data, which is a low pass filter to cut off bands with frequency larger than 10Hz. But I found it hard to achieve. My code was like this:     filter_n=20*sample_rate; cutoff=10; %cutoff frequency lowpass=fir1(filter_n,cutoff*2/sample_rate,'low'); y_lowpass=filter(lowpass,1,y); But I found the low-passed signal still had much power in bands larger than 10Hz. Can someone help me out? Thanks!  

Expert Answer

Profile picture of John Williams John Williams answered . 2025-11-20

Your issue may be due to insufficient filter order or improper design of the FIR filter. To apply an effective low-pass filter to your EEG signal, here are a few key suggestions:

Improvements to Your Code

  1. Increase Filter Order:

    • The filter order should be sufficiently high to achieve the desired attenuation for frequencies above the cutoff. However, an overly high order can make the filter unstable or cause excessive computation time.
  2. Check the Normalized Cutoff Frequency:

    • In your code, the normalized cutoff frequency is calculated as cutoff×2/sample_rate\text{cutoff} \times 2 / \text{sample\_rate}. This is correct, but ensure the resulting value is ≤1\leq 1.
  3. Use a Window Method or Different Design Approach:

    • The fir1 function uses a window-based design method, but you might achieve better performance with an alternative design approach, like using designfilt.
  4. Analyze the Filter's Frequency Response:

    • Plot the frequency response of the filter to confirm that it meets your requirements.

Corrected Code with Suggestions

Here’s an updated and more robust way to design and apply the filter:

 

% Given parameters
sample_rate = 50000; % 50 kHz
cutoff = 10; % Cutoff frequency in Hz
filter_order = 500; % Increase filter order for sharper cutoff

% Normalize cutoff frequency
normalized_cutoff = cutoff / (sample_rate / 2);

% Design the low-pass FIR filter
lowpass = fir1(filter_order, normalized_cutoff, 'low');

% Analyze the frequency response of the filter
fvtool(lowpass, 1);

% Apply the filter to your EEG data
y_lowpass = filter(lowpass, 1, y);

% Plot the original and filtered signals
t = (0:length(y)-1) / sample_rate; % Time vector
figure;
subplot(2, 1, 1);
plot(t, y);
title('Original Signal');
xlabel('Time (s)');
ylabel('Amplitude');

subplot(2, 1, 2);
plot(t, y_lowpass);
title('Low-pass Filtered Signal');
xlabel('Time (s)');
ylabel('Amplitude');

Key Points:

  1. Filter Order:

    • A higher order (e.g., 500) provides a sharper transition band but increases computation time. You can tune this value based on your needs.
  2. Frequency Normalization:

    • Ensure that the cutoff frequency is normalized with respect to the Nyquist frequency (sample_rate/2\text{sample\_rate}/2).
  3. Filter Verification:

    • Use fvtool to visualize the frequency response of the filter and verify that it attenuates frequencies above 10 Hz effectively.

Alternative: Using designfilt

Instead of fir1, you can use designfilt for more control over the filter design:

 

% Design a low-pass FIR filter with designfilt
d = designfilt('lowpassfir', ...
               'FilterOrder', filter_order, ...
               'CutoffFrequency', cutoff, ...
               'SampleRate', sample_rate);

% Analyze the frequency response
fvtool(d);

% Apply the filter
y_lowpass = filter(d, y);

Verify the Filtered Signal

To ensure that frequencies above 10 Hz are removed:

  1. Perform a Fourier Transform of the original and filtered signals to compare the frequency spectra.
  2. Use fft and plot to visualize the power in each frequency band.

 

% FFT of the original signal
Y_orig = fft(y);
f = linspace(0, sample_rate / 2, floor(length(Y_orig) / 2) + 1);

% FFT of the filtered signal
Y_filtered = fft(y_lowpass);

% Plot the magnitude spectra
figure;
plot(f, abs(Y_orig(1:length(f))), 'b');
hold on;
plot(f, abs(Y_filtered(1:length(f))), 'r');
title('Frequency Spectrum');
xlabel('Frequency (Hz)');
ylabel('Magnitude');
legend('Original', 'Filtered');

Troubleshooting

  1. If the filtered signal still contains power above 10 Hz, increase the filter order or consider using an IIR filter (e.g., Butterworth) for a sharper cutoff.
  2. Ensure that your signal is not being corrupted by noise or artifacts during recording or preprocessing.

Let me know if you need additional guidance!


Not satisfied with the answer ?? ASK NOW

Get a Free Consultation or a Sample Assignment Review!