递归：带有`scipy.lfilter`的IIR过滤器答案

【问题标题】：Recursion: IIR filter with `scipy.lfilter`递归：带有`scipy.lfilter`的IIR过滤器
【发布时间】：2018-02-08 17:52:16
【问题描述】：

给定一些数据x：

from pandas_datareader.data import DataReader as dr
x = np.squeeze(dr('DTWEXB', 'fred').dropna().values)

我想计算另一个向量y如下：

在这种情况下，alpha 等于 0.03。

我可以用scipy.lfilter? 来做这件事吗？类似的问题here，但在这种情况下，结果的起始值为 0，这会导致一些问题。

我的尝试：

from scipy.signal import lfilter
a = 0.03
b = 1 - a
y0 = x[0]
y = lfilter([a], [y0, -b], x)

结果应该是：

true_y = np.empty(len(x))
for k in range(0, len(true_y)):
    if k == 0:
        true_y[k] = x[0]
    else:
        true_y[k] = a*x[k] + b*true_y[k-1]
print(true_y)
[ 101.1818      101.176862    101.16819314 ...,  120.9813121   120.92484874
  120.85786628]

【问题讨论】：

标签： python python-3.x numpy recursion scipy

【解决方案1】：

对过滤器进行 z 变换，得到分子 b 和分母 a 的值：

       alpha
y(z) = ------------------ x(z)
       1 - (1-alpha) z^-1

所以你跑了

import scipy.signal
x[0] /= alpha
y = scipy.signal.lfilter([alpha], [1, - 1 + alpha], x)

产量

array([ 101.1818    ,  101.176862  ,  101.16819314, ...,  120.9813121 ,
        120.92484874,  120.85786628])

请注意，我已将 x[0] 缩放以考虑您想要的初始条件。

【讨论】：

似乎在正确的轨道上，但有些东西不对劲。结果y 应该等于上面的true_y。我也从设置a[0]=0 即[0, 1 - alpha] (ValueError: BUG: filter coefficient a[0] == 0 not supported yet) 得到一个错误
递推关系可以写成y[k] - (1 - a)*y[k-1] = a*x[k]。取that的z变换，求解Y(z)得到传递函数。
我似乎需要提高我的 z 变换技能：-D
除了@WarrenWeckesser 的回答之外，它还有效。我只是建议lfilter([a], [1., -1. + a], np.insert(x[1:], 0, x[0] / a)) 避免就地修改x。

【解决方案2】：

传递函数系数的正确参数是[a] 和[1, -b]。

要处理您想要的初始条件，您可以使用scipy.signal.lfiltic 为过滤器创建正确的初始状态：

zi = lfiltic([a], [1, -b], y=[x[0]])

然后使用zi 参数调用lfilter：

y, zo = lfilter([a], [1, -b], x, zi=zi)

这是x、y（使用lfilter 和zi 计算）和您的true_y：

In [37]: x
Out[37]: array([ 3.,  1.,  2.,  0., -1.,  2.])

In [38]: y
Out[38]: 
array([ 3.        ,  2.94      ,  2.9118    ,  2.824446  ,  2.70971262,
        2.68842124])

In [39]: true_y
Out[39]: 
array([ 3.        ,  2.94      ,  2.9118    ,  2.824446  ,  2.70971262,
        2.68842124])

【讨论】：