[SciPy-User] how to use signal.lfiltic?

Skipper Seabold jsseabold@gmail....
Wed Nov 24 15:37:22 CST 2010


On Wed, Nov 24, 2010 at 11:09 AM, Skipper Seabold <jsseabold@gmail.com> wrote:
> On Wed, Nov 24, 2010 at 10:39 AM, Skipper Seabold <jsseabold@gmail.com> wrote:
>> This is mainly for my own understanding of going back and forth
>> between signal processing language and time series econometrics.
>>
>> I don't see how to use lfiltic.
>>
>> Say I have a (known) output vector of errors and an input y vector.  y
>> follows a mean zero ARMA(p,q) process given by b and a below where p =
>> q = 2.  If I want to use lfilter to recreate the errors, forcing the
>> first p outputs (errors) to be zero, then I need to solve the
>> difference equations for the zi that do this which is given by zi
>> below.  But if I try to recreate this using lfiltic, it doesn't work.
>> Am I missing the intention of lfiltic?  Matlab's documentation, which
>> a lot of the signal stuff seems to be taken from, suggests that y and
>> x in lfiltic need to be reversed, but this also doesn't give the zi I
>> want.
>>
>> from scipy import signal
>> import numpy as np
>>
>> errors = np.array([ 0.,  0.,  0.00903417,  0.89064639,  1.51665674])
>> y = np. array([-0.60177354, -1.60410646, -1.16619292, 0.44003132, 2.36214611])
>>
>> b = np.array([ 1.        , -0.8622494 ,  0.34549996])
>> a = np.array([ 1.        ,  0.07918344, -0.81594865])
>>
>> # zi I want to produce errors = 0,0,...
>>
>> zi = np.zeros(2)
>> zi[0] = -b[0] * y[0]
>> zi[1] = -b[1] * y[0] - b[0] * y[1]
>>
>> zi
>> # array([ 0.60177354,  1.08522758])
>>
>> e = signal.lfilter(b, a, y, zi=zi)
>> e[0]
>> # array([ 0.        ,  0.        ,  0.00903417,  0.89064639,  1.51665674])
>>
>> zi2 = signal.lfiltic(b,a, errors[:2], y[:2])
>>
>> zi2
>> # array([-0.03533984, -0.20791273])
>>
>> e2 = signal.lfilter(b, a, y, zi=zi2)
>>
>> e2[0]
>> # array([-0.63711338, -1.24269149, -0.41241704, -0.0899541 ,  1.25042151])
>>
>> Skipper
>>
>
> Basically, I think this line in signal.lfiltic
>
> for m in range(M):
>    zi[m] = sum(b[m+1:]*x[:M-m],axis=0)
>
> should be
>
>  for m in range(M):
>    zi[m] = sum(-b[:m+1][::-1]*x[:m+1],axis=0)
>
> I'm not sure about the next loop, since my output are zero, I didn't
> have to solve for it.
>
> Is this a bug or am I misunderstanding?
>

Misunderstanding.  Matlab gives the same thing, but I still don't
understand why lfiltic doesn't give the same zi that I used to
actually create the output.

Skipper


More information about the SciPy-User mailing list