前一篇文章我们讲了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);