Chapter 14: when doing bandpass filter, why fft() plus ifft() gives different results from filtfilt()

[Chao Han]

Hi Mike,

Chapter 14 showed how to use filtfilt() to filter the data. I thought filtfilt() in theory works in the same way as convolution, so I tried fft() plus ifft() to get the filtered results instead of using filtfilt(). But it looks like they give slightly different results, with the results from fft() plus ifft() giving larger amplitude. I was wondering, is the difference due to some mistake in my code, or the difference is real? If it's real, what caused the difference?

Here is the code adapted from your codes to compare the different results:

% load data

load sampleEEGdata

 

% specify Nyquist freuqency

nyquist = EEG.srate/2;

 

% filter frequency band

filtbound = [4 10]; % Hz

 

% transition width

trans_width = 0.2; % fraction of 1, thus 20%

 

% filter order

filt_order = round(3*(EEG.srate/filtbound(1))); % A typical filter order: 3 times the ratio of the sampling rate to the lower edge of the frequency band.

 

% frequency vector (as fraction of Nyquist

ffrequencies  = [ 0 (1-trans_width)*filtbound(1) filtbound (1+trans_width)*filtbound(2) nyquist ]/nyquist;

 

% shape of filter (must be the same number of elements as frequency vector

idealresponse = [ 0 0 1 1 0 0 ];

 

% get filter weights

filterweights = firls(filt_order,ffrequencies,idealresponse);

 

% apply filter to signal

signal = EEG.data(47,:,1);

 

% use filtfilt()

filtfilt_function = filtfilt(filterweights,1,signal);

 

% use fft() and ifft()

n_data = EEG.pnts;

n_kernel = length(filterweights);

half_kernel = floor(n_kernel/2);

n_conv = n_data + n_kernel - 1;

n_conv_pow2 = pow2(nextpow2(n_conv));

fft_kernel = fft(filterweights, n_conv_pow2);

fft_kernel = fft_kernel / max(fft_kernel);

fft_data = fft(signal, n_conv_pow2);

conv_result = ifft(fft_kernel .* fft_data);

conv_result = conv_result(1:n_conv);

conv_result = conv_result(half_kernel+1 : end-half_kernel);

fft_function = real(conv_result)*2;

 

% plot to compare

figure(1),clf

plot(EEG.times, signal), hold on

plot(EEG.times,filtfilt_function, 'linew', 2)

plot(EEG.times, fft_function, 'k', 'linew', 2)

legend({'raw data';'result from filtfilt()';'result from fft() and ifft()'})

Thank you!

Chao

[Mike X Cohen]

Hi Chao. Indeed, the results will be similar but not necessarily identical. MATLAB's filtfilt function calls the filter function, which is causal and does some normalization and padding. Based on the plot your code generates, I wouldn't worry about it -- some changes to filter characteristics would produce similar kinds of minor differences.

Mike

--
You received this message because you are subscribed to the Google Groups "AnalyzingNeuralTimeSeriesData" group.
To unsubscribe from this group and stop receiving emails from it, send an email to analyzingneuraltimes...@googlegroups.com.
Visit this group at https://groups.google.com/group/analyzingneuraltimeseriesdata.
To view this discussion on the web visit https://groups.google.com/d/msgid/analyzingneuraltimeseriesdata/5b9d2da4-4a0a-46a3-ab62-527091673294%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

--

Mike X Cohen, PhD
Fresh look: mikexcohen.com

[Chao Han]

Hi Mike,

The MATLAB's documentation of filtfilt() functions says "Do not use filtfilt with differentiator and Hilbert FIR filters, because the operation of these filters depends heavily on their phase response." Is this a concern for the Filter Hilbert method?

Chao

To unsubscribe from this group and stop receiving emails from it, send an email to analyzingneuraltimeseriesdata+unsub...@googlegroups.com.

Visit this group at https://groups.google.com/group/analyzingneuraltimeseriesdata.
To view this discussion on the web visit https://groups.google.com/d/msgid/analyzingneuraltimeseriesdata/5b9d2da4-4a0a-46a3-ab62-527091673294%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Mike X Cohen]

Hi Chao. That refers to the symmetric shape of the filter frequency response. It's unrelated to the Hilbert transform, just named after the same person (I guess...).

Chao

To unsubscribe from this group and stop receiving emails from it, send an email to analyzingneuraltimes...@googlegroups.com.

Visit this group at https://groups.google.com/group/analyzingneuraltimeseriesdata.
To view this discussion on the web visit https://groups.google.com/d/msgid/analyzingneuraltimeseriesdata/5b9d2da4-4a0a-46a3-ab62-527091673294%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

--

Mike X Cohen, PhD
Fresh look: mikexcohen.com

--
You received this message because you are subscribed to the Google Groups "AnalyzingNeuralTimeSeriesData" group.

To unsubscribe from this group and stop receiving emails from it, send an email to analyzingneuraltimes...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/analyzingneuraltimeseriesdata/47f98fe6-1f35-43f5-be48-0b67557288bf%40googlegroups.com.