SciPy 沿 N 维数组的任意轴应用初始条件的 lfilter

Question

根据 SciPy 文档 lfilter：

zi : array_like, 可选滤波器延迟的初始条件。它是一个长度为 max(len(a),len(b))-1 的向量（或 N 维输入的向量数组）。如果 zi 是 None 或未给出，则假定初始休息。有关详细信息，请参阅 lfiltic。

以下代码调用 lfilter，并使用 lfilter_zi 传递 zi，使得 zi 最后一维的长度为 max(len(a),len(b))-1。但是它会引发错误，具体取决于应用轴：

import numpy as np
import scipy.signal as sig

def apply_filter(B, A, signal, axis=-1):
   # apply filter, setting proper initial state (doesn't assume rest)
   filtered, zf = sig.lfilter(B, A, signal, 
             zi=sig.lfilter_zi(B, A) * np.take(signal, 0, axis=axis)[..., np.newaxis], axis=axis)
   return filtered

B, A = sig.butter(1, 0.5)
x = np.random.randn(12, 50)
apply_filter(B, A, x, axis=1)    # works without error
apply_filter(B, A, x, axis=0)    # raises ValueError

ValueError: The number of initial conditions must be max([len(a),len(b)]) - 1

如何避免错误，并在不假设初始静止的情况下沿任何轴应用过滤器？

Answer 1

zi中的初始条件必须与给定lfilter的轴在同一个轴上。改变这个：

np.take(signal, 0, axis=axis)[..., np.newaxis]

到

np.take(signal, [0], axis=axis)

np.take(signal, 0, axis=axis)和np.take(signal, [0], axis=axis)的区别在于后者保留了维数。例如

In [105]: signal
Out[105]: 
array([[ 0,  1,  2,  3,  4],
       [ 5,  6,  7,  8,  9],
       [10, 11, 12, 13, 14]])

In [106]: signal.shape
Out[106]: (3, 5)

如果我们 take 来自轴 1 的第一个索引，我们得到一个形状为 (3,) 的一维数组：

In [107]: a = np.take(signal, 0, axis=1)

In [108]: a.shape
Out[108]: (3,)

In [109]: a
Out[109]: array([ 0,  5, 10])

如果我们在 indices 参数中使用列表 [0]，我们将得到一个形状为 (3, 1) 的数组：

In [110]: b = np.take(signal, [0], axis=1)

In [111]: b.shape
Out[111]: (3, 1)

In [112]: b
Out[112]: 
array([[ 0],
       [ 5],
       [10]])

SciPy 沿 N 维数组的任意轴应用初始条件的 lfilter

SciPy lfilter with initial conditions applied along any axis of N-D array

python

filtering

numpy

scipy