5
$\begingroup$

I'm trying to port some MATLAB code to Python and am running into some strange behaviour. I am implementing a 5th order Butterworth bandpass filter. The sampling rate is 30 Hz.

Running MATLAB R2012b, Spyder 2.2.0 with Python 2.7, SciPy 0.12.0 on Windows 7 x64.

In MATLAB:

[b,a] = butter(5, [0.75*2/30, 5.0*2/30], 'bandpass');
y = filtfilt(b, a, input_signal)

This is the raw signal:

raw signal

and the filtered signal:

MATLAB filtered signal

and its power spectrum:

MATLAB power spectrum

which makes sense since the normalized bandpass frequencies are 0.05 - 0.33.

I found before that SciPy's butter function does not give the same coefficients as MATLAB so I exported the filter coefficients from MATLAB to Python using hdf5write (see here: https://stackoverflow.com/questions/7117797/export-matlab-variable-to-text-for-python-usage)

In Python:

y = signal.filtfilt(b, a, input_signal, padtype = None)

and the output is:

Python filtered signal

and its power spectrum:

Python power spectrum

I used padtype = None because by default it is padtype = 'odd'. However, I've tried all the different padding options and they all look more or less the same.

I'm not entirely sure what's going wrong...any help would be greatly appreciated.

EDIT: Added graph for padtype = "odd" and b and a filter coefficients used

Signal with padtype = "odd":

padtype = odd signal

Seems like there is some underlying signal so I don't think it's the impulse response although that large transient at the beginning is strange.

Python filter coefficients (5th order Butterworth filter):

passband = [0.75*2/30, 5.0*2/30]
b, a = scipy.signal.butter(5, passband, 'bandpass')

b, a are arrays of type float64.

b = array([  5.49209388e-03,   0.00000000e+00,  -2.74604694e-02,
    -1.97776791e-17,   5.49209388e-02,   2.47220989e-17,
    -5.49209388e-02,  -1.97776791e-17,   2.74604694e-02,
     0.00000000e+00,  -5.49209388e-03])
a = array([  1.        ,  -6.52098852,  19.50384534, -35.47189804,
    43.65758795, -38.07760914,  23.83047021, -10.55670367,
     3.16710078,  -0.58122912,   0.04945954])

MATLAB filter coefficients

b =   0.0055, 0, -0.0275, 0, 0.0549, 0, -0.0549, 0, 0.0275, 0, -0.0055
a =   1.0000, -6.5210, 19.5038, -35.4719, 43.6576, -38.0776, 23.8305, 
      -10.5567, 3.1671, -0.5812, 0.0495

Turns out the coefficients are about the same, I believe the issues I had before were due to the frequency being much higher (coefficients must be generated differently in MATLAB and Python).

The input signal is a list of float32, although I've tried converting to array with numpy.array and the result is the same.

EDIT: More information about padtypes and Python vs. MATLAB

Python SciPy's filtfilt function includes a parameter called padtype which indicates the type of padding extended on both sides of the signal. This padding serves to reduce transients. Odd and even are descriptors of the type of symmetry these extensions have with the endpoints (http://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.filtfilt.html).

Some illustrations of the difference in padtypes:

padtype = "even"

even padding

padtype = "odd"

odd padding

padtype = None

no padding

Based on the results, it appears that MATLAB uses odd padding (which is also the default for Python):

MATLAB

$\endgroup$
  • $\begingroup$ Your python output looks like the impulse response. Do you have the input_signal in the right orientation? (i.e. should it be transposed?). Try setting axis=0 or axis=1 in the call to filtfilt. $\endgroup$ – Peter K. Nov 4 '13 at 19:28
  • $\begingroup$ It appears to be the correct orientation - the output is the same length as the input. If I put axis = 0, it's the same as before. Setting axis = 1 gives me an error that tuple index is out of range. $\endgroup$ – limi44 Nov 4 '13 at 20:57
  • $\begingroup$ Did you try generating a and b like this. $\endgroup$ – Peter K. Nov 5 '13 at 12:40
  • $\begingroup$ Please include b and a in your question. Are they formatted correctly? I suspect padtype='odd' is the same as Matlab, since it's the default and Matlab docs say "filtfilt minimizes start-up and ending transients by matching initial conditions", but I'm not sure. $\endgroup$ – endolith Nov 5 '13 at 14:29
  • 2
    $\begingroup$ I uploaded the data to see if anyone else can replicate my problem. Data file. Thank you endolith and Peter for all your help so far. $\endgroup$ – limi44 Nov 5 '13 at 15:27
3
$\begingroup$

It works for me:

from scipy.io import loadmat
from scipy.signal import butter, filtfilt
from matplotlib.pyplot import plot

signaldata = loadmat('signaldata.mat')

input_signal = signaldata['input_signal'][0]

passband = [0.75*2/30, 5.0*2/30]
b, a = butter(5, passband, 'bandpass')

y = filtfilt(b, a, input_signal)
plot(y)

enter image description here

I don't know what your issue is but maybe you can figure it out from here.

$\endgroup$
  • $\begingroup$ I ran the code you posted and I'm still getting the same thing as before...are you running a different version of SciPy? $\endgroup$ – limi44 Nov 5 '13 at 16:25
  • $\begingroup$ @limi44: Perhaps? Python 2.7.2, NumPy 1.7.1, SciPy 0.12.0 $\endgroup$ – endolith Nov 5 '13 at 16:26
  • 1
    $\begingroup$ I installed the newest SciPy package 0.13.0 and now it works! I originally installed SciPy as part of Python(x,y) so perhaps part of the installation got botched. Either way, thanks for confirming I wasn't going crazy. $\endgroup$ – limi44 Nov 5 '13 at 16:55
  • $\begingroup$ @limi44: I wouldn't expect that sort of thing to differ from one version to the next, but if you can run the latest version and that fixes it, I guess that's ok. Can you confirm that the default padtype matches Matlab? What effect do the other types have? $\endgroup$ – endolith Nov 5 '13 at 19:36
  • $\begingroup$ I don't think it was a version issue, probably just some sort of issue with my installation. I can confirm that MATLAB matches the default padtype (odd). $\endgroup$ – limi44 Nov 5 '13 at 21:02
1
$\begingroup$

According to the newsgroup https://mail.python.org/pipermail/scipy-user/2014-April/035648.html,

The difference is in the default padding length. In matlab's filtfilt, it is 3*(max(len(a), len(b)) - 1), and in scipy's filtfilt, it is 3*max(len(a), len(b)).

Therefore, I modified the padlen when invoking filtfilt in Python. And the results between MATLAB and Python will be the same.

MATLAB code

input_signal = [1,2,3,4,5,1,2,3,4,5,1,2,3,4,1,2,3,4,5,1,2,3,4,5,1,2,3,4,5,1,2,3,4,5]
[b,a] = butter(5, [0.75*2/30, 5.0*2/30], 'bandpass');
y = filtfilt(b, a, input_signal)

Python with SciPy code

from scipy.signal import butter, filtfilt
input_signal = [1,2,3,4,5,1,2,3,4,5,1,2,3,4,1,2,3,4,5,1,2,3,4,5,1,2,3,4,5,1,2,3,4,5]
[b,a] = butter(5, [0.75*2/30, 5.0*2/30], 'bandpass');
y = filtfilt(b, a, input_signal, padtype = 'odd', padlen=3*(max(len(b),len(a))-1))
$\endgroup$
1
$\begingroup$

I got the same results in Paco's example, but changing the input signal to a vector of ones does not yield same results in MATLAB and Python (Scipy 1.0.0):

MATLAB:

input_signal = ones(10000,1);
[b,a] = butter(5, [0.75*2/30, 5.0*2/30], 'bandpass');
y = filtfilt(b, a, input_signal);
y(1:5)'

ans =

    1.0e-13 *

    -0.306211857054885  -0.468000282200120  -0.612592864567431  -0.736885525141862  -0.838459113089290

Python:

from numpy import ones
from scipy.signal import butter, filtfilt
input_signal = ones((10000,))
[b,a] = butter(5, [0.75*2/30, 5.0*2/30], 'bandpass');
y = filtfilt(b, a, input_signal, padtype = 'odd', padlen=3*(max(len(b),len(a))-1))
y[:5]

Out[49]: 
array([-2.52233927e-14, -1.79580673e-14, -1.07857124e-14, -3.97094772e-15, 2.30161154e-15])
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.