LSM和RLS算法实现系统辨识和逆辨识

用MATLAB脚本实现，使用LMS和RLS两种算法分别实现系统辨识和逆辨识，通过脚本实现过程可比对LMS和RLS两种算法收敛速度、稳态误差变化趋势等特征，模仿本段代码可作为理解自适应滤波算法的入门练习。

%%
%Adaptive Filter Test
%Copyright By popwaic
%2019.11
%In ZhengZhou
%%
%
clc;
clear;

d=randn（110000）‘;


h=[0.2 0.66 0.31 0.77 0.91 0.65 0 0 0 0 0 0 0];
N=11;
h1=[0 0 0 0.3 0.6 0.3 0.7 0.9  0 0 0];

w=zeros（113）;  %init
w1=zeros（12*N+3）;  %init
w1（N+2）=1;

r=zeros（1length（d））;

u=0.2;
u1=0.015;

 Lmin=7000;

%%
%System Identification
    for i=1:length（d）-13
        xx=d（i+12:-1:i）;  
        dr=h*xx;       
        r（i）=dr;    
    end

    %LMS
    [w_est_Lee_L]=LMS_t（rdwu1）;
    
     %use RLS
     R=eye（length（w））;
     t=0.0001;
     R=t*R;
     R=inv（R）;
     j=0.995;

     [w_est_Ree_R]=RLS_t（d（1:Lmin）rwRj）;
    
%drow picture
    subplot（221）;
    stem（h）; 
    title（‘h（n）‘）;
    subplot（222）;
    stem（w_est_L）; 
    hold on;
    stem（w_est_R）;
    hold off;
    title（‘w（n）‘）;
    subplot（223）;
    stem（h-w_est_L）;
    hold on;
    stem（h-w_est_R‘）;
    hold off;
    title（‘h（n）-w（n）‘）;
    subplot（224）;
    plot（ee_L（1:Lmin-100））; 
    hold on 
     plot（ee_R（1:Lmin-100））; 
    title（‘MSE of e*e‘）;
    

%%
%System inverse identification

for i=1:length（d）-N
    xx=d（i+N-1:-1:i）;  
    dr=h1*xx;       
    r（i）=dr;    
end

r=r*0.7359;         %Normalization

eee=zeros（1length（r））;
 
 %use LMS
 [w_est_Lee_L]=LMS_t（d（17:Lmin）rw1u10）;
 
 %use RLS
 R=eye（length（w1））;
 t=0.0001;
 R=t*R;
 R=inv（R）;
 j=0.995;
 
 [w_est_Ree_R]=RLS_t（rd（17:Lmin）w1Rj）;
 
figure;

subplot（221）;
stem（h1）; 
title（‘h1（n）‘）;
subplot（222）;
%stem（w1）
stem（w_est_L）; 
hold on;
stem（w_est_R）;
hold off;

title（‘w1（n）‘）;

subplot（223）;
stem（conv（h1w_est_L（:）））; 
hold on;
stem（conv（h1w_est_R（:）））; 
hold off;

%stem（conv（h1w1（:）））; 
title（‘h1（n）*w1（n）‘）;
subplot（224）;
plot（ee_L（1:Lmin-100））; 
hold on;
plot（ee_R（1:Lmin-100））; 
hold off;
title（‘MSE of e*e‘）;


%%noise 


%%
%functrion of LMS
% r:y
% d:d
% w:estimate of h
% u: step length
% mod : 0:basic LMS   1:normalzial LMS

function [w_estee]=LMS_t（rdwumod）
    N=length（w）;
    n=30;
    for i=1:length（r）-2*N
        xx=d（i+N-1:-1:i）; %
        
        %
        [m h]=size（xx）;        
        if m==1
            xx=xx‘;
        end
        
        y=w*xx;        
        e=r（i）-y;
       
        if mod==1
            w=w+（u/（sum（xx.^2）+0.01）*xx*e）‘; 
        elseif mod==0
            w=w+（u*xx*e）‘;    
        end

    %calculate  MSE
    Le=i;
    eee（i）=e.^2;  
    if i     ee（i）=10*log10（sum（eee（1:Le）.^2）/Le）;   
    else
        ee（i）=10*log10（sum（eee（Le-n+1:Le）.^2）/Le）;   
    end
       
    end
    
    w_est=w;   
end

%%function
%RLS 
% x: x           N*1 D
% d: d           N*1 D
% w: w（n-1）      N*1 D
% R:  the inverse of square matrix  Rxx（n-1）    N*N D

% t:             a very small param
% I:             N*N D  R（0）=t*I
% j:             forget param

%init: w（0）=0 R0 jt

f

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----
     文件        3855  2019-12-08 15:05  Adaptive_Filter_LSM_RLS_SysIdent_SysEqu.m