search

nilearn / nilearn

Machine learning for NeuroImaging in Python
http://nilearn.github.io
Other
1.15k stars 581 forks source link

Filtering looks broken #482

Closed GaelVaroquaux closed 9 years ago

GaelVaroquaux commented 9 years ago

Our usage of butterworth looks like it is not doing what is should.

GaelVaroquaux commented 9 years ago

Related to #374 ?

AlexandreAbraham commented 9 years ago

Can you be more specific?

GaelVaroquaux commented 9 years ago

Sorry, I wrote this issue live during the course at Leipzig, in front of the audience. I dug a bit more, and the following seems really wrong: import numpy as np from nilearn.signal import butterworth from matplotlib import pyplot as plt

t = np.linspace(0, 10, 2000) y = np.cos(t)

plt.clf() plt.plot(y) plt.plot(butterworth(y[:, np.newaxis], sampling_rate=1., low_pass=10000., ).squeeze())

plt.show()

Even though I am putting a very high cut off frequency for the low_pass filter, I am getting the signal entirely cut. Even stranger, if I change the low_pass to a high_pass, I am getting the opposite: no cut at all. It seems to me that there is something deeply broken here.

bthirion commented 9 years ago

On 04/03/2015 21:39, Gael Varoquaux wrote:

plt.plot(butterworth(y[:, np.newaxis], sampling_rate=1., low_pass=10000., ).squeeze()) On numpy 1.9.1, this throws an error ("wn should be btw 0 and 1").

AlexandreAbraham commented 9 years ago

On numpy 1.9.1, this throws an error ("wn should be btw 0 and 1").

Not on my box.

GaelVaroquaux commented 9 years ago

On 04/03/2015 21:39, Gael Varoquaux wrote:

plt.plot(butterworth(y[:, np.newaxis], sampling_rate=1., low_pass=10000., ).squeeze()) On numpy 1.9.1, this throws an error ("wn should be btw 0 and 1").

Yes, we have had this in the workshop also. It clearly seems to me that this code is broken.

AlexandreAbraham commented 9 years ago

OK, I think that I got it. In fact, Wn (see scipy.signal.butter) must be in [0, 1] for a digital filter. The problem is that when the value is outside of the bounds, the filter seems to have an unexpected behavior.

You can see it easily: with low_pass=0.5 you get the expected behavior (signal is the same, Wn=1.). But with low_pass=0.501 you can see that the filter goes completely crazy.

I still don't know why I don't get the warning associated with Wn value.

AlexandreAbraham commented 9 years ago

Related to #374 ?

No. This issue has been clearly identified and we only need to agree on the behavior to adopt if a high pass filter is applied on data with a linear trend.

GaelVaroquaux commented 9 years ago

OK, I think that I got it. In fact, Wn (see scipy.signal.butter) must be in [0, 1] for a digital filter.

It seems that we should clip it (and maybe warn) for values too high or too low.

AlexandreAbraham commented 9 years ago

:+1:

banilo commented 9 years ago

On numpy 1.9.1, this throws an error ("wn should be btw 0 and 1").

I also get this error on Mac + numpy '1.9.2'