正如我在评论中指出的那样,您可以使用
scipy.signal.lfilter。在这种情况下(假设
A是一维numpy数组),您需要做的是:
B = lfilter([a], [1.0, -b], A)
这是一个完整的脚本:
import numpy as npfrom scipy.signal import lfilternp.random.seed(123)A = np.random.randn(10)a = 2.0b = 3.0# Compute the recursion using lfilter.# [a] and [1, -b] are the coefficients of the numerator and# denominator, resp., of the filter's transfer function.B = lfilter([a], [1, -b], A)print B# Compare to a simple loop.B2 = np.empty(len(A))for k in range(0, len(B2)): if k == 0: B2[k] = a*A[k] else: B2[k] = a*A[k] + b*B2[k-1]print B2print "max difference:", np.max(np.abs(B2 - B))
该脚本的输出为:
[ -2.17126121e+00 -4.51909273e+00 -1.29913212e+01 -4.19865530e+01 -1.27116859e+02 -3.78047705e+02 -1.13899647e+03 -3.41784725e+03 -1.02510099e+04 -3.07547631e+04][ -2.17126121e+00 -4.51909273e+00 -1.29913212e+01 -4.19865530e+01 -1.27116859e+02 -3.78047705e+02 -1.13899647e+03 -3.41784725e+03 -1.02510099e+04 -3.07547631e+04]max difference: 0.0
另一个示例,在IPython中,使用pandas Dataframe而不是numpy数组:
如果你有
In [12]: df = pd.Dataframe([1, 7, 9, 5], columns=['A'])In [13]: dfOut[13]: A0 11 72 93 5
并且您想要创建一个新列,
B这样
B[k] = A[k] + 2*B[k-1](
B[k] == 0对于k <0而言),您可以编写
In [14]: df['B'] = lfilter([1], [1, -2], df['A'].astype(float))In [15]: dfOut[15]: A B0 1 11 7 92 9 273 5 59



