系统辨识报告.docx
《系统辨识报告.docx》由会员分享，可在线阅读，更多相关《系统辨识报告.docx（19页珍藏版）》请在冰豆网上搜索。
系统辨识报告
一、
题目要求：
（1）一钢球以零初速由电视天线上自由下落，实验测得下落时间和距离的数据如下：
t（s）
1
2
3
4
5
6
l（m）
8.94
20.05
50.65
72.19
129.85
171.56
假设时间t无量测误差，而距离l有量测误差。
试用最小二乘法确定重力加速度g。
分析解答：
分析可得，我们已知重力加速度乘以下落时间t的平方除以二等于其降落距离，因此，此处使用最小二乘法的基本算法，即可求得g的估计值，同时经画图可在坐标上描出对应的点，通过其曲线对结果进行相应的验证。
（1）Matlab程序如下：
t=[123456]';
T=0.5*t.^2;
l=[8.9420.0550.6572.19129.85171.56]';
H=[ones（6,1）,T];
Z=inv（H'*H）*H'*l;
l0=Z
（1）,g=Z
（2）
L=l0+g*T;
plot（T,L,'r',T,l,'b+'）;holdon
xlabel（'t/s'）;ylabel（'l/m'）;
title（'下落时间和距离关系'）;
holdoff
结果：
a=
9.7729
（2）做出曲线如下图所示：
图一：
小球下落时间与距离关系曲线图
（2）
题目要求：
设两国军备竞赛模型为
x（k）=ax（k-1）+by（k-1）+f
y（k）=cx（k-1）+dy（k-1）+g
式中x（k）和y（k）为两国军事费用（单位：
百万美元），已知数据如下：
K
伊朗
伊拉克
北约
华约
1972
2891
909
216478
112893
1973
3982
1123
211146
115020
1974
8801
2210
212267
117169
1975
11230
2247
210525
119612
1976
12178
2204
205717
121461
1977
9867
2303
212009
123561
1978
9165
2179
215988
125498
1979
5080
2675
218561
127185
1980
4040
225411
129000
1981
233957
131595
试用最小二乘法确定模型参数a,b,c,d,f和g.
分析解答：
仍然可直接使用最小二乘法进行估计计算.
1.伊朗和伊拉克：
Matlab程序如下：
y=[9091123221022472204230321792675];
x=[2891398288011123012178986791655080];
fork=2:
8
h（k,:
）=[x（k-1）y（k-1）1];
z（k,:
）=[x（k）y（k）];
end
guji=inv（h'*h）*h'*z;%算出a、b、f、c、d、g估计值，为三行两列的矩阵a=guji（1,1）,b=guji（2,1）,f=guji（3,1）,c=guji（1,2）,d=guji（2,2）,g=guji（3,2）%分别取到a、b、f、c、d、g
plot（[1:
8],x,'r'）;holdon%画出实际与各自的估计曲线
hg=h*guji;
abf=（hg（:
[1]））';
plot（[1:
8],abf,'b'）;holdon
plot（[1:
8],y,'r'）;holdon
hg=h*guji;
cdg=（hg（:
[2]））';
plot（[1:
8],cdg,'g'）;holdon
xlabel（'k（1-8代表1972到1980）'）;ylabel（'军事费用/百万美元'）;%写出坐标表头
legend（'伊朗实际值','伊朗估计值','伊拉克实际值','伊拉克估计值'）;
title（'伊朗与伊拉克军事费用对比'）;holdoff
a=
0.5481
b=
-0.2034
f=
4.4475e+003
c=
-0.0306
d=
0.7443
g=
987.1565
：
图二：
伊朗与伊拉克军事费用比较图
2、北约和华约：
Matlab程序如下：
x=[216478211146212267210525205717212009215988218561255411233957];
y=[112893115020117169119612121461123561125498127185129000131595];
fork=2:
10
h（k,:
）=[x（k-1）y（k-1）1];
z（k,:
）=[x（k）y（k）];
end
guji=inv（h'*h）*h'*z;%用最小二乘法求出估计值
a=guji（1,1）,b=guji（2,1）,f=guji（3,1）,c=guji（1,2）,d=guji（2,2）,g=guji（3,2）
plot（[1:
10],x,'r'）;holdon
hg=h*guji;
abf=（hg（:
[1]））';%取guji矩阵的第一列即a，b,f的值
plot（[1:
10],abf,'b'）;holdon
plot（[1:
10],y,'g'）;holdon
cdg=（hg（:
[2]））';
plot（[1:
10],cdg,'b'）;holdon
xlabel（'k（1-8代表1972到1980）'）;ylabel（'军事费用/百万美元'）;%写出坐标表头
legend（'北约实际值','北约估计值','华约克实际值','华约克估计值'）;
title（'北约与华约事费用对比'）;holdoff
a=
0.1540
b=
1.7438
f=
-2.5468e+004
c=
0.0181
d=
0.9675
g=
2.0756e+003
图三：
华约与北约军事费用对比图
最小二乘法辨识结果
参数
A
b
f
c
d
g
伊朗与伊拉克
0.5481
-0.2034
4447.5
-0.0306
0.7443
987.1565
北约与华约
0.1540
1.7438
-70770
0.0181
0.9675
2075.6
分析：
由以上对比曲线可以看出，经过最小二乘法估计得到的数据与实际数据之间虽然存在区别，但是
二：
题目要求：
考虑理想数学模型为
式中，
是服从正态分布的白噪声
。
输入信号采用4阶
序列，其幅值为1。
选择如下的辨识模型进行增广递推最小二乘参数辨识。
给出各参数的辨识曲线和辨识误差曲线。
分析解答：
分析：
应该使用增广递推最小二乘算法
n=60;
%4位移位寄存器产生的M序列的周期
y1=1;y2=1;y3=1;y4=0;
fori=1:
n;
x1=xor（y3,y4）;
x2=y1;
x3=y2;
x4=y3;
y（i）=y4;
ify（i）>0.5,u（i）=-1;
elseu（i）=1;
end
y1=x1;y2=x2;y3=x3;y4=x4;
end
figure
（1）;subplot（2,1,1）;stem（u）,gridon
title（'M序列'）%M序列
%白噪声的产生
A=19;x0=12;M=1024;
fork=1:
n
x=A*x0;
x1=mod（x,M）;
v（k）=x1/512;
x0=x1;
end
%辨识主程序
z=zeros（7,n）;zs=zeros（7,n）;zm=zeros（7,n）;zmd=zeros（7,n）;
z
（1）=0;z
（2）=0;
zs
（1）=0;zs
（2）=0;
zm
（1）=0;zm
（2）=0;
zmd
（1）=0;zmd
（2）=0;
theta=zeros（7,1）;
p0=10^6*eye（7,7）;
the=[theta,zeros（7,14）];
e=zeros（7,15）;
fork=3:
n;
z（k）=1.5*z（k-1）-0.7*z（k-2）+u（k-1）+0.5*u（k-2）+v（k）-v（k-1）+0.2*v（k-2）;
h=[-z（k-1）,-z（k-2）,u（k-1）,u（k-2）,v（k）,v（k-1）,v（k-2）]';
x=inv（h'*p0*h+1）;
q=p0*h*x;
d1=z（k）-h'*theta;
theta1=theta+q*d1;
zs（k）=-1.5*z（k-1）+0.7*z（k-2）+u（k-1）+0.5*u（k-2）;
zm（k）=[-z（k-1）,-z（k-2）,u（k-1）,u（k-2）]*[theta1
（1）;theta1
（2）;theta1（3）;theta1（4）];
zmd（k）=h'*theta1;
e（:
k）=theta1-theta;
theta=theta1;
the（:
k）=theta1;
p1=p0-q*q'*[h'*p0*h+1];
p0=p1;
end
figure
（1）;
i=1:
n;
plot（i,the（1,:
）,'r',i,the（2,:
）,'r:
',i,the（3,:
）,'b',i,the（4,:
）,'b:
',i,the（5,:
）,'g',i,the（6,:
）,'g:
',i,the（7,:
）,'g+'）
title（'递推增广最小二乘辨识方法'）
figure
（2）;
i=1:
n;
plot（i,e（1,:
）,'r',i,e（2,:
）,'r:
',i,e（3,:
）,'b',i,e（4,:
）,'b:
',i,e（5,:
）,'g',i,e（6,:
）,'g:
',i,e（7,:
）,'r+'）
title（'辨识误差'）
运行程序，得到如下曲线：
产生的M序列信号
递推增广最小二乘辨识方法
辨识误差
增广最小二乘参数辨识结果
参数
a1
a2
b1
b2
d1
d2
d3
真值
1.5
0.7
1.0
0.5
1.2
-1.0
0.2
估计值
1.5
0.7
1.0
0.5
1.2
-1.0
0.2
三：
题目要求：
完善试验2，给出由相关分析法给出的脉冲响应函数辨识曲线，在试验2的基础上，给出拟合二阶及三阶系统的最小二乘辨识参数曲线，及由脉冲响应函数反求传函
。
1、由实验二已经求得，相关分析法得到的脉冲响应函数辨识曲线
相关分析法所得脉冲响应函数辨识曲线
2、二阶脉冲响应辨识：
Matlab程序：
n=2;Ts=2;np=15;%n阶次，Ts采样时间
U=UY（61:
300,1）;%输入初值U
Y=UY（61:
300,2）;%输入初值Y
f=zeros（120,4）;
f（:
1）=-1*Y（61:
180）;
f（:
2）=-1*Y（60:
179）;
f（:
3）=U（61:
180）;
f（:
4）=U（60:
179）;
yy=Y（62:
181）;
q=inv（f'*f）*f'*yy;
a=q（1:
2）;
b=q（3:
4）;
G=tf（b',[1a'],2）%传递函数GLS（z）
yc=zeros（n+np,1）;%脉冲响应
uc=zeros（n+np,1）;
uc（n+1）=1/Ts;
fori=1:
np-1
yc（n+i）=-1*[yc（i+n-1:
-1:
i）]'*a+[uc（i+n-1:
-1:
i）]'*b;
end
yimp=yc（n+1:
n+np）;
figure
（1）;
stem（0:
Ts:
（np-1）*Ts,2*yimp,'r'）
holdon
plot（0:
Ts:
（np-1）*Ts,2*yimp,'r^-'）
impulse（G）
得到的传递函数为：
Transferfunction:
0.3882z+0.2645
-----------------------
z^2-0.9641z+0.2137
Samplingtime:
2
得到辨识曲线如图所示：
改为三阶系统的传递函数
运行结果：
0.3888z^2+0.2581z+0.002418
G（z）=-------------------------------------
z^3-0.9819z^2+0.2452z-0.01913
三阶系统的脉冲响应函数
四：
比较最小二乘法及其改进算法的性能、优缺点。
答：
从最小二乘的定义中得知，最小二乘的思想就是寻找一个
的估计值
，使得各次测量的
与由估计
确定的量测估计
之差的平方和最小，由于此方法兼顾了所有方程的近似程度，使整体误差达到最小，因而对抑制误差是有利的。
但是。
最小二乘一般方法的估计精度不够高，这是由于对各个测量数据同等对待，而各次测量数据一般不会在相同的条件下获得，造成测量数据的置信度不变较大。
因此，最小二乘法有很多改进算法，虽然没有一个是完美的，但是能够适应不同的情况，条件，对应选择不同的算法，其各自的性能及优缺点如下：
1）基本的最小二乘法（LS，RLS）
对低噪声水平最有效，参数估计可很快收敛到真值，所需计算量相对较少。
对于有色噪声的情况只能得到有偏估计。
几乎不需要特殊的先验假设，在递推算法中只要求给定和，RLS具有可靠的收敛性。
2）广义最小二乘法（GLS）
一般可给出好的结果，但计算量较大，信噪比较小时收敛不到真值。
可能出现局部最小值的情况。
3）多级最小二乘法（MSLS）
分三步获得系统模型参数和噪声模型参数的估计，计算量小，速度快，但精度难提高，工程应用不便。
4）增广最小二乘法（RELS）
算法简单，可同时得到参数和噪声模型的估计，工程应用效果很好。
5）加权最小二乘法可对不同置信度的测量值采用加权的办法分别对待，置信度高的，权重取得大些；置信度低的，权重取的小些。
但加权最小二乘法仅能用于事先能估计方程误差对参数估计的影响。
6）递推最小二乘法的可以对计算量大、存储大的数据进行批处理，在线辨识，可以减少数据存储量，避免矩阵求逆，从而减少计算量。
矩形窗（限定记忆）RLS方法需要保留一定的数据存储量，此存储量大小取决于矩形窗宽度，因而在应用范围上有一定程度的限制。
五：
题目要求：
时变系统的参数估计问题。
已知系统y（k）+a（k）y（k-1）=b（k）u（k-1）+v（k）
其中a（k）、b（k）具有以下数值：
a（k）＝0.8、b（k）＝0.5当
a（k）＝0.6、b（k）＝0.3当
v（k）为零均值白噪声，输入
采用四级逆M序列，系统用数字计算机模拟，利用渐消记忆法估计时变系统参数：
分别取0.9、0.95、0.99，试说明取值不同对参数估计的影响。
分析解答如下：
分析题意可知，此题应该使用带遗忘因子的递推最小二乘法
L=100;%设置仿真长度
x1=1,x2=1,x3=1,x4=0;%四位移位寄存器产生M序列
s=1;
fork=1:
L
im=xor（s,x4）;
if（im==0）
u（k）=1;
else
u（k）=1;
end
s=not（s）;
m（k）=xor（x3,x4）;
x4=x3;x3=x2;x2=x1;x1=m（k）;
end%M序列产生结束
stairs（u）;grid;
v=randn（L）;
num=v;
k1=k;
y
（1）=0;
y=zeros（1,L）;
fork=2:
L%输出响应
ifk<300
y（k）=-0.8*y（k-1）+0.5*u（k-1）+v（k-1）;
else
y（k）=-0.6*y（k-1）+0.3*u（k-1）+v（k-1）;
end
end
L=1000;
Z=zeros（3,1）;
P=10^6*eye（3）;%设置初值
Q=zeros（3,n）;
r=0.98;
fork=1:
L
K=P*Q/（r+Q'*P*Q）;
Q=[-y（（k-1）:
-1:
（k-na））u（（k-1）:
-1:
（k-nb））];
C（:
k）=Z+K*（y（k）-Q*Z）;
P=（eye（3）-K*Q）*P/r;
Z=C（:
k）;
end
k=1:
L
plot（k,C）;
holdon
1）遗忘因子
=0.9
2）遗忘因子
=0.95
2）遗忘因子
=0.99
分析：
由上述辨识曲线可以看出，当遗忘因子增大时，其估计所得参数曲线变化并不明显，但是能够看出估计所得参数的变化则越平滑；λ越小，则跟踪越好，但参数会出现一定的跳跃。
但遗忘因子越大，辨识精度及速度都会降低。
但是由图1、图2、图3可以得到相应的结论。
六、广义最小二乘法和辅助变量法的区别
★辅助变量法（工具变量法）：
某一个变量与模型中随机解释变量高度相关，但却不与随机误差项相关，那么就可以用此变量与模型中相应回归系数的一个一致估计量，这个变量就称为工具变量，这种估计方法就叫工具变量法。
辅助变量法的辨识精度叫最小二乘法高，尤其适当造声势有色噪声，且噪声结构不好确定时，辅助变量法更显示出它的优势，因为这种方法无需知道噪声结构能获得系统的无偏一致估计。
缺点：
工具变量法的关键是选择一个有效的工具变量，由于工具变量选择中的困难，工具变量法本身存在两方面不足：
一是由于工具变量不是惟一的，因而工具变量估计量有一定的任意性；
二是由于误差项实际上是不可观测的，因而要寻找严格意义上与误差项无关而与所替代的随机解释变量高度相关的变量事实上是困难的。
★广义最小二乘法：
其基本思想在于对数据先进行一次白化滤波处理，然后利用基本的最小二乘法对滤波后的数据进行辨识。
如滤波模型选的合适，对数据进行了较好的白化处理，则直接基本的最小二乘法就能得到无偏一致估计。
所用的滤波模型实际上是一种动态模型，经过几次迭代调整后，便可对数据进行较好的白化处理。
其收敛速度比较慢，需要经过多次迭代计算，才能得到较准确的参数估计值。
一般情况下，经过多次迭代后，估计值便会收敛到稳态值。
但在某些情况下（如噪声比较低时）存在局部极小值，估计值不一定收敛到准则函数的全局极小值上。
为了防止参数估计值收敛到局部极小值，最好选定初值接近最优解，一般可以用最小二乘法的批处理估计值作为初值。
如果系统是时变的，或为了克服数据饱和现象，可以在两次RLS算法中分别引进遗忘因子。
升级会员 