RCLS.rar_MATLAB偏差补偿_偏差估计_偏差补偿_最小二乘估计_模型偏差补偿

clc clear L=400; %系统输入信号为M序列 u=[1 0 1 0]; for i=5:15+L; u(i)=xor(u(i-1),u(i-4)); end u=1-2*u; v=0.1*randn(1,15+L); %产生N(0,1)分布的随机输入噪声 %输出观测值的长度 z=zeros(1,15+L); for j=3:15+L; z(j)=1.5*z(j-1)-0.7*z(j-2)+u(j-1)+0.5*u(j-2)+v(j); %系统在M序列输入下的输出值 end n=2; %二阶系统 %求H H=[]; for i=1:14+L for j=1:n Z(j)=z(n+i-j); U(j)=u(n+i-j); end h=[-Z U]; H(i,:)=h; end H=H'; I=eye(n+n); p=10^6*I; %初始化p ThtaLs=zeros(1,2*n)';%最小二乘初值 ThtaC=zeros(1,2*n)';%偏差补偿初值 D=[eye(n) zeros(n,n);zeros(n,n) zeros(n)]; J=0; %按6.3.20式计算ThtaC for i=1:L; k=p*H(:,i)/(1+(H(:,i))'*p*H(:,i)); ThtaLs=ThtaLs+k*(z(n+i)-(H(:,i))'*ThtaLs); p=(I-k*H(:,i)')*p; J=J+(z(n+i)-H(:,i)'*ThtaLs)^2/(1+(H(:,i))'*p*H(:,i)); s2=J/(i*(1+ThtaC'*D*ThtaLs)); ThtaC=ThtaLs+i*s2*p*D*ThtaC; end a=ThtaC(1:n,:)' b=ThtaC(n+1:2*n,:)'