how to do digital filter on this EEG data

% 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 ($\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!