社区首页 >专栏 >自适应滤波器（二）NLMS自适应滤波器

自适应滤波器（二）NLMS自适应滤波器

猫叔Rex

发布于 2022-01-24 05:42:18

1.8K03

代码可运行

文章被收录于专栏：科学计算科学计算

运行总次数：3

代码可运行

前一篇文章我们讲了LMS自适应滤波器，我们先回顾一下LMS算法流程：

影响LMS性能的因素，也就是最后一个公式的三个因素：

步长，它是由我们事先指定
输入向量
估计误差

如果过大，那么

的结果中，占的权重就太大了，而是带噪信号，这样梯度噪声就被放大了。为了克服这个问题，可使用归一化LMS滤波器。在迭代时，对输入向量欧式范数（就是模值）的平方进行归一化（Normalized LMS)。

归一化LMS滤波器是最小化干扰原理的一种表现形式，这个原理可以表述如下：

从一次迭代到下一次中，自适应滤波器的权向量应当以最小方式改变，而且受到更新的滤波器输出所施加的约束。

用

表示第n次迭代滤波器的权向量，

表示第n+1次迭代滤波器的权向量，那么NLMS设计准则可表述为约束优化问题：给定输入向量和目标响应，确定更新抽头向量

，以使如下增量

的欧式范数最小化，并受制于以下约束条件

使用拉格朗日乘子法来解决这个约束问题，那么代价函数为：

其中，为复数拉格朗日乘子，∗表示复共轭，表示取实部运算，约束对代价函数的贡献是实值的；

表示欧式范数的平方运算，其结果也是实数。因此，代价函数是实值的二次函数，且可表示为：

采用如下步骤寻找最小的最优更新权向量：

代价函数对

求导，可得：

令其为零，即得最优解为：

将上式带入约束条件，求解未知乘子：

可求得：

其中，

是误差信号。

结合前两步的结果，可得：

为了对一次迭代到下一次迭代抽头权向量的增量变化进行控制而不改变向量的方向，引入一个正的实数标度因子，该增量可以写为：

等价的，我们可以写出：

这个公式就是归一化LMS算法抽头权向量的递归公式，为什么叫归一化呢？因为公式中对输入向量欧式范数的平方就行了归一化。

当输入向量较小时，

的值过小，可能导致数值计算困难的情况，为了克服这个情况，将上面的表达式改为：

其中，

我们总结NLMS算法的步骤如下：

% 输入参数:
%   xn   输入的信号，列向量
%   dn   所期望的响应
%   M    滤波器的阶数
%   mu   收敛因子(步长)
% 输出参数:
%   W    滤波器系数矩阵  
%   en   误差序列 
%   yn   滤波器输出        
function [yn, W, en]=nlmsFunc(xn, dn, M, mu, delta)
itr = length(xn);
en = zeros(itr,1);            
W  = zeros(M,itr);    % 每一列代表-次迭代,初始为0
% 迭代计算
for k = M:itr                  % 第k次迭代
    x = xn(k:-1:k-M+1);        % 滤波器M个抽头的输入
    y = W(:,k-1).' * x;        % 滤波器的输出
    en(k) = dn(k) - y ;        % 第k次迭代的误差
    % 滤波器权值计算的迭代式
    W(:,k) = W(:,k-1) + mu*en(k)*x/(delta+x.' * x);
end

yn = inf * ones(size(xn)); % 初值为无穷大是绘图使用，无穷大处不会绘图
for k = M:length(xn)
    x = xn(k:-1:k-M+1);
    yn(k) = W(:,end).'* x;  % 最终输出结果
end

产生测试信号：

fs = 1;
f0 = 0.02;
n = 1000;
t = (0:n-1)'/fs;
xs = cos(2*pi*f0*t);
ws = awgn(xs, 20, 'measured');

M  = 20 ;   % 滤波器的阶数
xn = ws;
dn = xs;
% rho_max = max(eig(ws*ws.'));   % 输入信号相关矩阵的最大特征值
% mu = (1/rho_max) ;    % 收敛因子 0 < mu < 1/rho
mu1 = 0.01;
mu2 = 1;
delta = 1e-3;
[yn,W,en] = lmsFunc(xn,dn,M,mu1);
[yn2,W2,en2] = nlmsFunc(xn, dn, M, mu2, delta);