递推最小二乘法算法
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 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中求逆是应该避免的。