递推最小二乘法算法

合集下载
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

题目:

(递推最小二乘法)

考虑如下系统:

)()4(5.0)3()2(7.0)1(5.1)(k k u k u k y k y k y ξ+-+-=-+--

式中,)(k ξ为方差为0.1的白噪声。

取初值I P 610)0(=、00=∧

)(θ。选择方差为1的白噪声作为输入信号)(k u ,采用PLS 法进行参数估计。

Matlab 代码如下:

clear all

close all

L=400; %仿真长度

uk=zeros(4,1); %输入初值:uk(i)表示u(k -i)

yk=zeros(2,1); %输出初值

u=randn(L,1); %输入采用白噪声序列

xi=sqrt(0.1)*randn(L,1); %方差为0.1的白噪声序列

theta=[-1.5;0.7;1.0;0.5]; %对象参数真值

thetae_1=zeros(4,1); %()θ初值

P=10^6*eye(4); %题目要求的初值

for k=1:L

phi=[-yk;uk(3:4)]; %400×4矩阵phi 第k 行对应的y(k-1),y(k-2),u(k-3), u(k-4) y(k)=phi'*theta+xi(k); %采集输出数据

%递推最小二乘法的递推公式

K=P*phi/(1+phi'*P*phi);

thetae(:,k)=thetae_1+K*(y(k)-phi'*thetae_1);

P=(eye(4)-K*phi')*P;

%更新数据

thetae_1=thetae(:,k);

for i=4:-1:2

uk(i)=uk(i -1);

end

uk(1)=u(k);

for i=2:-1:2

yk(i)=yk(i -1);

end

yk(1)=y(k);

end

plot([1:L],thetae); %line([1,L],[theta,theta]);

xlabel('k'); ylabel('参数估计a 、b');

legend('a_1','a_2','b_0','b_1');

%axis([0 L -2 2]);

结果

分析如下

系统方程为

)()4(5.0)3()2(7.0)1(5.1)(k k u k u k y k y k y ξ+-+-=-+--

由CAR 模型:

)()4(5.0)3()2(7.0)1(5.1)(k k u k u k y k y k y ξ+-+-+---=

可以得到:

2117.05.11)(----+=z z z A

)5.01()(131---+=z z z B

我们可以知道纯迟延为3T ,na=2,nb=1,d=3,

phi(k,:)=[-yk;uk(3:4)]'; 400×4矩阵phi第k行对应的y(k-1),y(k-2),u(k-3), u(k-4)

y(k)=phi*theta+xi(k);输出等于矩阵phi与真值矩阵相乘加上白噪声。在这里phi是列向量不需要用phi(k,:)形式

这是循环的主体,接下来

一开始我们赋给初值给输入以及输出

初始输入u(0)=u(-1)=u(-2)=u(-3)=0;初始输出y(0)=y(-1)=0;

所以

=u

-

+

u

+

y

y

y

+

)2-(

5.0

)3-(

)1-(

7.0

)1(

5.1

)0(

)1(ξ

经过白噪声的更新代码后,以及输入输出值得更新后我们有,也就是把计算出来的y(1)给到下一组y(2)的计算中区去:

-

+

+

=u

y

+

y

y

u

)2(

)1-(

5.0

)2-(

)2(ξ

)0(

5.1

)1(

7.0

于是我让程序中显示phi以便只观的表现递推最小二乘法和批处理最小二乘法的不同:

>> phi

phi =

-1.9333

-3.3952

-0.4448

-0.4710

递推最小二乘法phi截图如下:

相比于第一次作业的批处理最小二乘法的phi:

可以发现递推最小二乘法和批处理最小二乘法的运算结果在L数据长度是一样的情况下是相同的,两者的实质原理相同,仅仅是数学原理上的差别。但是优点在于后一种办法,取一组数据用一组,避免了矩阵求逆的过程,减轻了计算机软件的负担,因为在matlab中求逆是应该避免的。

相关文档
最新文档