使用Matlab2019b实现最小均方误差算法（LMS）

可以与之前的dsp.LMSFilter作为对比，从原理来实现自适应滤波器中的最小均方误差算法（Least Mean Square）

注意：本程序仅兼容2018a以后版本

这个程序有一个地方目前还存在一些问题，就是求取自相关函数R时所用到的方法，我通过尝试自行编写自相关函数求取方程和matlab自带的corrmtx函数作为对比，发现数值有所不同，应该是matlab自带的corrmtx函数和书本上的求解方式并不相同，所以使用时需要自行判断使用哪种方式。

clc
clear
%%构建信号%%
N=256;%信号长度
for j=1:50  %进行50次lms运算
t=1:N;
d=sin(0.1*pi*t);
figure(1)
plot(t,d);
axis([0,N,-2,2]);
%%加噪声%%
x=awgn(d,20);%awgn 加性高斯白噪声  信噪比20DB

% 一种简单的计算自相关矩阵方法
% R=x'*x;%计算自相关矩阵

%%直接计算自相关矩阵 设自相关矩阵长度为L自己写
% L=256;%自相关长度
% Rxx=zeros(L,L);
% for m=1:L
%     for n=1:L
%         mm=abs(m-n);
%         temp=0;
%         for nn=1:N-mm
%             temp=temp+x(nn)*x(nn+mm);
%         end
%         Rxx(m,n)=temp/(N-mm);
%     end
% end
% [V,D]=eig(Rxx);%特征值分解求u
% ma=max(max(D));
% u=1/ma;

% 一维实值信号x的自相关矩阵Rxx应为实对称的toeplitz矩阵，
%而一维实值信号x,y的互相关矩阵Rxy为非对称的toeplitz阵，
%matlab提供的corrmtx函数产生的并非通常意义下的autocorrelation matrix

m= 255;% dfine the time lag m+1,and m+1<=n
[rx,r]=corrmtx(x,m);%corrmtx为自相关矩阵
R=rx'*rx;
[V,D]=eig(R);%提取自相关矩阵的特征值对角阵D
maxlamuda=max(D,[],'all');%提取自相关矩阵最大特征值
u=1/maxlamuda-0.01;%保证LMS算法收敛并保持稳定的条件

%%LMS迭代滤波%%
T=50;
k=25;
p=zeros(T,N-k);%50次 误差矩阵
e=zeros(1,N);%误差矩阵
w=zeros(1,k);%滤波器长度系数
y=zeros(1,N);%滤波输出
y(1:k)=x(1:k);%假设滤波器前25点为输入值
    for i=(k+1):N%开始迭代更新参数
        Xn=x(i:-1:(i-k+1));%输入
        y(i)=w*Xn';%输出
        e(i)=d(i)-y(i);%误差
        w=w+2*u*e(i)*Xn;%权值更新
    end
    p(j,:)=(e(k+1:N)).^2;%误差信号功率
end
for i=1:N-k
    E(i)=sum(p(:,i))/T;%50次误差取平均
end
figure(2)
t=1:N-k;
plot(t,E);

Copyright © 2020 by RichardYang. All rights reserved.
仅供参考，严禁转载，感谢。

算法工作模式

参数：

• 抽头权系数向量w(n)
• 输入向量x(n)
• 期望输出d(n)
• 滤波器输出y(n)
• 滤波器抽头权更新参数w(n+1)

流程：

• 滤波：$y(n)={\omega }^T(n)x(n)$
• 误差估计：$e(n)=d(n)-y(n)$
• 抽头权向量更新：$\omega (n+1)=\omega (n)+2\mu e(n)x(n)$

注意：

• 为了保证算法的收敛和稳定性必须有：$0<\mu <\frac{1}{{\lambda }_{\max }}$
• ${\lambda }_{\max }$为输入自相关矩阵R的最大特征值（公式证明，但实际使用时可能会出现问题，通常不这么求最大步长，但测试可以使用 ¹）

---

1. 自适应滤波器原理及Matlab仿真应用（原书第2版） ↩︎