我想在 Python 中应用自适应滤波器，但在网上找不到任何关于如何实现这种算法的文档或示例。我熟悉使用scipy.signal工具箱设计“静态”滤波器，但我不知道如何设计自适应滤波器。
澄清一下：我有一个S包含噪音的记录信号。在这段录音中有一个我想访问的“真实”函数，称之为T. 我也有一个估计T。我想设计一个过滤器，使过滤器S和过滤器之间的误差T最小化。请注意，在这种情况下，静态滤波器没有用，因为我正在尝试过滤非平稳信号。
问问题
20003 次
2 回答
18
这是一个使用 Numpy 的 Python 中的基本 LMS 自适应滤波器。
欢迎评论，欢迎测试用例。
""" lms.py: a simple python class for Least mean squares adaptive filter """
from __future__ import division
import numpy as np
__version__ = "2013-08-29 aug denis"
#...............................................................................
class LMS:
""" lms = LMS( Wt, damp=.5 ) Least mean squares adaptive filter
in:
Wt: initial weights, e.g. np.zeros( 33 )
damp: a damping factor for swings in Wt
# for t in range(1000):
yest = lms.est( X, y [verbose=] )
in: X: a vector of the same length as Wt
y: signal + noise, a scalar
optional verbose > 0: prints a line like "LMS: yest y c"
out: yest = Wt.dot( X )
lms.Wt updated
How it works:
on each call of est( X, y ) / each timestep,
increment Wt with a multiple of this X:
Wt += c X
What c would give error 0 for *this* X, y ?
y = (Wt + c X) . X
=>
c = (y - Wt . X)
--------------
X . X
Swings in Wt are damped a bit with a damping factor a.k.a. mu in 0 .. 1:
Wt += damp * c * X
Notes:
X s are often cut from a long sequence of scalars, but can be anything:
samples at different time scales, seconds minutes hours,
or for images, cones in 2d or 3d x time.
"""
# See also:
# http://en.wikipedia.org/wiki/Least_mean_squares_filter
# Mahmood et al. Tuning-free step-size adaptation, 2012, 4p
# todo: y vec, X (Wtlen,ylen)
#...............................................................................
def __init__( self, Wt, damp=.5 ):
self.Wt = np.squeeze( getattr( Wt, "A", Wt )) # matrix -> array
self.damp = damp
def est( self, X, y, verbose=0 ):
X = np.squeeze( getattr( X, "A", X ))
yest = self.Wt.dot(X)
c = (y - yest) / X.dot(X)
# clip to cmax ?
self.Wt += self.damp * c * X
if verbose:
print "LMS: yest %-6.3g y %-6.3g err %-5.2g c %.2g" % (
yest, y, yest - y, c )
return yest
#...............................................................................
if __name__ == "__main__":
import sys
filterlen = 10
damp = .1
nx = 500
f1 = 40 # chirp
noise = .05 * 2 # * swing
plot = 0
seed = 0
exec( "\n".join( sys.argv[1:] )) # run this.py n= ... from sh or ipython
np.set_printoptions( 2, threshold=100, edgeitems=10, linewidth=80, suppress=True )
np.random.seed(seed)
def chirp( n, f0=2, f1=40, t1=1 ): # <-- your test function here
# from $scipy/signal/waveforms.py
t = np.arange( n + 0. ) / n * t1
return np.sin( 2*np.pi * f0 * (f1/f0)**t )
Xlong = chirp( nx, f1=f1 )
# Xlong = np.cos( 2*np.pi * freq * np.arange(nx) )
if noise:
Xlong += np.random.normal( scale=noise, size=nx ) # laplace ...
Xlong *= 10
print 80 * "-"
title = "LMS chirp filterlen %d nx %d noise %.2g damp %.2g " % (
filterlen, nx, noise, damp )
print title
ys = []
yests = []
#...............................................................................
lms = LMS( np.zeros(filterlen), damp=damp )
for t in xrange( nx - filterlen ):
X = Xlong[t:t+filterlen]
y = Xlong[t+filterlen] # predict
yest = lms.est( X, y, verbose = (t % 10 == 0) )
ys += [y]
yests += [yest]
y = np.array(ys)
yest = np.array(yests)
err = yest - y
averr = "av %.2g += %.2g" % (err.mean(), err.std())
print "LMS yest - y:", averr
print "LMS weights:", lms.Wt
if plot:
import pylab as pl
fig, ax = pl.subplots( nrows=2 )
fig.set_size_inches( 12, 8 )
fig.suptitle( title, fontsize=12 )
ax[0].plot( y, color="orangered", label="y" )
ax[0].plot( yest, label="yest" )
ax[0].legend()
ax[1].plot( err, label=averr )
ax[1].legend()
if plot >= 2:
pl.savefig( "tmp.png" )
pl.show()
于 2013-08-20T09:08:32.243 回答
10
您应该检查 Padasip 库https://pypi.python.org/pypi/padasip
在这个库中实现了多个自适应滤波器。您可以直接使用它，也可以调查那里的过滤器是如何创建的（源代码在 github https://github.com/matousc89/padasip上）
如果您只想针对您的案例测试任何已实现的自适应滤波器的适用性，请参阅本教程Padasip Adaptive Filters Basics - Noise Cancelation, System Identification and Signal Prediction
于 2016-05-12T19:58:26.370 回答