递推最小二乘辨识
遗忘因子递推最小二乘法辨识rc电路
遗忘因子递推最小二乘法辨识rc电路遗忘因子递推最小二乘法（Recursive Least Squares，简称RLS）是一种用于参数辨识的方法，主要针对线性时变系统的辨识问题。
在本文中，我们将探讨如何使用RLS方法辨识RC电路。
首先，我们需要了解什么是RC电路。
RC电路由一个电阻和一个电容组成，其工作原理是通过电阻阻碍电流的流动，而电容则储存电荷。
RC电路可以用于滤波器、积分器等应用中。
在辨识RC电路时，我们可以将其视为一个一阶系统，其传递函数可以表示为：H(s) = 1 / (RCs + 1)其中，s是频率，RC是电阻和电容的乘积。
我们的目标是辨识出RC的值。
接下来，我们将介绍如何使用RLS方法辨识RC电路的参数。
首先，我们需要采集一些输入-输出数据对。
我们可以通过施加不同的输入信号，并记录输出信号来获取这些数据对。
这些输入信号可以是脉冲信号、阶跃信号或正弦信号等。
设输入信号为u(t)，输出信号为y(t)，则我们可以得到如下的离散模型：y(k) = wT(k) · x(k)其中，y(k)是第k个采样点的输出，w(k)是待辨识的参数向量，x(k)是输入-输出数据对构成的特征向量。
在训练过程中，我们需要定义一个递归的更新规则来更新参数向量w(k)。
这里，我们引入遗忘因子λ，用于控制历史数据对参数的影响程度。
遗忘因子的取值范围为0到1之间，较大的值表示较快的遗忘，较小的值表示较慢的遗忘。
参数向量的更新公式如下：w(k) = w(k-1) + P(k-1) · x(k) · e(k)其中，P(k-1)是协方差矩阵，它与遗忘因子有关，可以通过递推方式得到：P(k) = 1/λ · (P(k-1) - P(k-1) · x(k) · x(k)T · P(k-1) / (λ + x(k)T · P(k-1) · x(k)))e(k)是输出误差，可以通过以下方式计算：e(k) = y(k) - wT(k-1) · x(k)以上就是使用遗忘因子递推最小二乘法辨识RC电路的基本原理和公式。
现代控制理论_第14章_最小二乘法辨识
y n 2 ai y n 2 i bi u n 2 i n 2
i 1 i 0
n
y n N ai y n N i bi u n N i n N
即
y k ai y k i bi u k i v k ai v k i
i 1 i 0 i 1
n
n
n
(14-3)
假设v k k 1,2,, n 是均值为零的独立分布的平稳随机序列，且与 序列u k k 1,2,, n 相互独立。设
ˆ 表示 y 的最优估值，则有 设ˆ 表示 的最优估值， y
ˆ ˆ y
(14-12)
式中
ˆ n 1 y a ˆ ˆ y n 2 ˆ ˆ y , b ˆ ˆ y n N
T 的展开式如下所示：
y n 1 y n y n y n 1 y 1 y 2 T u n 1 u n 2 u n 1 u n u 2 u 1 y n N 1 n 1 y n N 2 n 2 yN u n N u n N 1 n N uN
1
因为ˆ 有解与 T 正定等价，所以可以保证 T 正定来确定对输 入 u k 序列的要求。由式(14-9)可知
Y U
(14-20)
则
YT U YT U Y T Y U T T T U Y U U U
递推阻尼最小二乘法辨识算法公式的详细推导与说明
控制理论与控制工程学位课程《系统辨识》考试报告递推阻尼最小二乘法公式详细推导专业：控制理论与控制工程班级：2011双控（研）学生姓名：江南学号：20110201016任课教师：蔡启仲老师2012年06月29 日摘要在参数辨识中，递推最小二乘法是用得最多的一种算法。
但是，最小二乘法存在一些缺点，如随着协方差矩阵的减小，易产生参数爆发现象；参数向量和协方差矩阵的处置选择不当会使得辨识过程在参数收敛之前结束；在存在随机噪声的情况下，参数易产生漂移，出现不稳定等。
为了防止参数爆发现象，Levenberg 提出在参数优化算法中增加一个阻尼项，以增加算法的稳定性。
本文在一般的最小二乘法中增加了阻尼因子，构成了阻尼最小二乘法。
又根据实时控制的要求，详细推到了递推阻尼最小二乘公式，实现在线辨识。
关键字：系统辨识，最小二乘法，递推算法正文1.题目的基本要求已知单入单出系统的差分方程以及噪声，在应用最小二乘法进行辨识的时候，在性能指标中加入阻尼因子，详细推导阻尼最小二乘法的递推公式。
2.输入辨识信号和系统噪声的产生方法和理论依据 2.1系统辩识信号输入选择准则（1）输入信号的功率或副度不宜过大，以免使系统工作在非线性区，但也不应过小，以致信噪比太小，直接影响辩识精度；（2）输入信号对系统的“净扰动”要小，即应使正负向扰动机会几乎均等； （3）工程上要便于实现，成本低。
2.2白噪声及其产生方法 （1） 白噪声过程（2）白噪声是一种均值为0、谱密度为非0常数的平稳随机过程。
（3）白噪声过程定义：如果随机过程()t ω的均值为0，自相关函数为()()2R t t ωσδ= （2.2.1）式中()t δ 为狄拉克（Dirac) 分布函数，即(){(),00,01t t t dt δδ∞∞=≠∞==⎰-且t （2.2.2）则称该随机过程为白燥声过程。
2.3白噪声序列 (1) 定义 如果随机序列{()}w t 均值为0，并且是两两不相关的，对应的自相关函数为()2,0,1,2w l R l l σδ==±± 式中{1,00,0l l l δ=≠=则称这种随机序列{()}w t 为白噪声序列。
递推最小二乘估计及模型阶次辨识
实验二 递推最小二乘估计(RLS)及模型阶次辨识(F-Test ）1 实验方案设计1.1 生成输入数据和噪声用M 序列作为辨识的输入信号，噪声采用标准正态分布的白噪声。
生成白噪声时，首先利用乘同余法生成U[0,1]均匀分布的随机数，再利用U[0,1]均匀分布的随机数生成标准正态分布的白噪声。
1.2 过程仿真辨识模型的形式取)()()()()(11k e k u z B k z z A +=--，为方便起见，取n n n b a == 即nn n n zb a b z b z B z a a a z a z A ------++++=++++=...1)(...1)(22112211用M 序列作为辨识的输入信号。
1.3 递推遗忘因子法数据长度L 取534，初值⎪⎪⎪⎩⎪⎪⎪⎨⎧⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡==1000010000100001)0(001.0)0(P θ 1.4 计算损失函数、噪声标准差损失函数⎥⎦⎤⎢⎣⎡+---+-=μθμττ)()1()()]1(ˆ)()([)1()(2k h k P k h k k h k z k J k J噪声标准差θλdim )(ˆ-=L L J1.6 F-Test 定阶法计算模型阶次统计量t)22,2(~222)1()1()()1,(----++-=+n L F n L n J n J n J n n t其中，)(∙J 为相应阶次下的损失函数值，L 为所用的数据长度，n 为模型的估计阶次。
若a t n n t >+)1,(，拒绝00:n n H >，若a t n n t <+)1,(，接受00:n n H >，其中αt 为风险水平α下的阀值。
这时模型的阶次估计值可取1+n 。
1.6 计算噪信比和性能指标噪信比22ye σση= 参数估计平方相对偏差i i i ni i i θθθθθδˆ~,~1221-=⎪⎪⎭⎫ ⎝⎛=∑= 参数估计平方根偏差ii i n i ini iθθθθθδˆ~,)()~(2122122-==∑∑== 2 编程说明M 序列中，M 序列循环周期取15124=-=p N ，时钟节拍t ∆=1Sec ，幅度1=a ，特征多项式为1)(56⊕⊕=s s s F 。
应用最小二乘一次完成法和递推最小二乘法算法地系统辨识
1最小二乘法的理论基础1.1最小二乘法设单输入单输出线性定长系统的差分方程表示为：其中δ(k)为服从N(0,1)的随机噪声,现分别测出n+N 个输出输入值y(1),y(2),…,y(n+N),u(1),u(2),…,u(n+N),则可写出N 个方程，写成向量－矩阵形式()()()()()()()()1201121n n y k a y k a y k a y k n b u k b u k b u k n k ξ=-------++-++-+L L ()()()()()()()()()()()()()()()()()()10111121222112n n y n y n y u n u y n y n y u n u y n N y n N y N u n N u N a n a n b n N b ξξξ+--+⎡⎤⎡⎤⎢⎥⎢⎥+-+-+⎢⎥⎢⎥=⨯⎢⎥⎢⎥⎢⎥⎢⎥+-+--+⎢⎥⎢⎥⎣⎦⎣⎦⎡⎤⎢⎥+⎡⎤⎢⎥⎢⎥⎢⎥+⎢⎥+⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥+⎢⎥⎣⎦⎢⎥⎢⎥⎣⎦LLL L M M MMLL M M M(4.1.1)则式(1.1.1)可写为 （4.1.2）式中：y 为N 维输出向量；ξ为N 为维噪声向量；θ为（2n+1）维参数向量；Φ为N ×（2n+1）测量矩阵。
因此，式（4.1.1）是一个含有（2n+1）个未知参数，由N 个方程组成的联立方程组。
11y θφφξ--=-在给定输出向量y 和测量矩阵Φ的条件下求参数θ的估计，这就是系统辨识问题。
设 表示 θ 的估计值，ŷ表示y 的最优估计，则有 （4.1.3） 式中：()()()10ˆˆ1ˆˆ2ˆˆ,ˆˆˆn n ay n a y n y b y n N b θ⎡⎤⎢⎥+⎡⎤⎢⎥⎢⎥⎢⎥+⎢⎥⎢⎥==⎢⎥⎢⎥⎢⎥⎢⎥+⎢⎥⎣⎦⎢⎥⎢⎥⎣⎦M M M ()()()()()()()()()()()()()()()()()()101122,,11112221n n a y n n y n a n y b y n N n N b y n y u n u y n y u n u y n N y N u n N u N ξξθξξφ⎡⎤⎢⎥++⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥++⎢⎥⎢⎥===⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥++⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥⎢⎥⎣⎦--+⎡⎤⎢⎥-+-+⎢⎥=⎢⎥⎢⎥-+--+⎢⎥⎣⎦M M M M L LL L M MMLL y φθξ=+ˆθˆˆyθ=Φ设e(k)=y(k)- ŷ(k), e(k)称为残差，则有e=y- ŷ=y-Φθ 最小二乘估计要求残差的平方和最小，即按照指数函数（4.1.4）求Ｊ对 的偏导数并令其等于0可得：（4.1.5）由式（4.1.5）可得的 θ 最小二乘估计：（4.1.6）J 为极小值的充分条件是：即矩阵ΦT Φ为正定矩阵，或者说是非奇异的。
第四章-最小二乘参数辨识方法及原理
脉冲传递函数描述：
Y ( z ) b0 b1 z 1 bn z n B( z 1 ) G( z) = 1 n 1 U ( z ) 1 a1 z an z A( z )
2.随机模型
v(k )
u(k )
G (k )
x(k )
y(k )
观测值可表示为：
本章的学习目的
1、掌握最小二乘参数辨识方法的基本原理
2、掌握常用的最小二乘辨识方法 3、熟练应用最小二乘参数辨识方法进行模型参数辨识 4、能够编程实现最小二乘参数辨识
回
顾
辨识目的：根据过程所提供的测量信息，在某种准则意 义下，估计模型的未知参数。 Input
Process 目 的
Output
工程实践
88 1010
95.7 1032
表 1 热敏电阻的测量值
t (C ) R ()
20.5 765
26 790
32.7 826
40 850
51 873
61 910
73 942
80 980
88 1010
95.7 1032
R a bt
N N N N ˆ Ri 2 R a 702ti.762i ti ti i 1 i 1 a i 1 i 1 ˆ 2 N N ˆ N. t i2 t i b 3 4344 i 1 i 1 N N N N Ri t i Ri t i i 1 i 1 b i 1 ˆ R 943 .1682 N N 2 N ti ti i 1 i 1
上述N个方程可写成下列向量－矩阵形式
y (n) y (1) u (n 1) y (n 1) y (n 2) y (n 1) y (2) u ( n 2) y (n N ) y (n N 1) y ( N ) u (n N )
递归最小二乘法辨识参数
递归最小二乘法辨识参数递归最小二乘法（Recursive Least Squares, RLS）是一种参数辨识方法，它使用递归算法来求解最小二乘法中的参数。
在许多领域中，例如系统辨识、自适应控制、信号处理等，递归最小二乘法都是一个广泛使用的方法。
递归最小二乘法的基本思想是：通过递归迭代来更新参数估计值，使其逼近最优解。
在递归过程中，每一次迭代时，都会通过当前的测量值来更新参数的估计值，同时保留历史测量值的影响，从而获得更精确的估计值。
具体地说，在递归过程中，首先需要定义一个初始参数向量，然后通过观测数据序列来递归更新参数向量。
假设有一个如下所示的线性关系：y(k) = Φ(k) * θ + v(k)其中，y(k)是被观测到的输出值，Φ(k)是与该输出值相关的输入向量，θ是待辨识的参数向量，v(k)是误差项。
递归最小二乘法的目标就是通过观测数据来估计θ的值。
在递归最小二乘法中，首先需要定义一个初始的参数向量θ0，然后通过数据序列递归地更新θ的值。
每一次迭代时，都会用最新的观测数据来更新参数向量，使得估计值更接近真实值。
具体来说，每次观测到新的数据之后，都会根据当前参数估计值和新的观测值来计算估计误差，并更新参数向量。
具体的迭代步骤如下：1.从数据序列中读取观测值y(k)和输入向量Φ(k)；2.计算估计值y(k)hat和估计误差e(k)：y(k)hat = Φ(k) * θ(k-1)e(k) = y(k) - y(k)hat3.计算卡尔曼增益K(k)和参数估计值θ(k)：K(k) = P(k-1) * Φ(k) / (λ + Φ(k)' * P(k-1) * Φ(k))θ(k) = θ(k-1) + K(k) * e(k)其中，P(k-1)是先前迭代步骤中的误差协方差矩阵，λ是一个小的正数，用于确保逆矩阵的存在性。
需要注意的是，递归最小二乘法的计算量相对较大，因此通常需要对算法进行优化，以提高计算效率和精度。
广义递推最小二乘辨识
广义递推最小二乘辨识一、实验目的1 通过实验掌握广义最小二乘辨识算法；2 运用MATLAB编程，掌握算法实现方法。
二、实验原理广义最小二乘法的基本思想是基于对数据先进行一次滤波预处理，然后利用普通最小二乘法对滤波后的数据进行辨识。
如果滤波模型选择得合适，对数据进行了较好的白色化处理，那么直接利用普通最小二乘法就能获得无偏一致估计。
广义最小二乘法所用的滤波模型实际上就是一种动态模型，在整个迭代过程中不断靠偏差信息来调整这个滤波模型，使它逐渐逼近于一个较好的滤波模型，以便对数据进行较好的白色化处理，使模型参数估计称为无偏一致估计。
理论上说，广义最小二乘法所用的动态模型经过几次迭代调整后，便可对数据进行较好的白化处理，但是，当过程的输出噪信比比较大或模型参数比较多时，这种数据白色化处理的可靠性就会下降。
此时，准则函数可能出现多个局部收敛点，因而辨识结果可能使准则函数收敛于局部极小点上而不是全局极小点上。
这样，最终的辨识结果往往也会是有偏的。
其收敛速度比较慢，需要经过多次迭代计算，才能得到较准确的参数估计值。
一般情况下，经过多次迭代后，估计值便会收敛到稳态值。
但在某些情况下（如噪声比较低时）存在局部极小值，估计值不一定收敛到准则函数的全局极小值上。
为了防止参数估计值收敛到局部极小值，最好选定初值接近最优解，一般可以用最小二乘法的批处理估计值作为初值。
如果系统是时变的，或为了克服数据饱和现象，可以在两次RLS算法中分别引进遗忘因子。
三、实验内容<1> 数据获取：实验数据按照表9-1，为二阶线性离散系统的输入输出数据<2> 数据处理：为了提高辨识精度，实验者必须对原始数据进行剔除坏数据、零均值化、工频滤波等处理。
实验进行了白化滤波处理。
<3> 辨识算法：利用处理过的数据（取适当的数据长度），选择某种辨识方法（如RLS递推最小二乘法、RELS、RIV或RML等参数估计算法及F-检验或AIC定阶法），估计出模型参数和阶次，同时分析辨识结果。
递推最小二乘辨识
• 令
P(k ) (ΦΦk )-1 k
(2)
将Φk展开,故有
P(k ) ([Φ-1 (k -1)][Φ-1 (k -1)] )-1 k k [Φ-1Φk -1 (k -1) (k -1)]-1 k
[ P -1 (k -1) (k -1) (k -1)]-1 (3)
– 遗忘因子法
• 通过对不同时刻的数据赋予一定的加权系数,使得对旧 数据逐渐遗忘,加强新数据对当前辨识结果的影响,从而 避免协方差矩阵P(k)与增益矩阵K(k)急剧衰减而失去对 参数估计的修正能力,使算法始终保持较快的收敛速度. • 遗忘因子法在后面将重点讨论.
L=[(0), (1), ..., (L-1)]T, L×(na+nb) (k - 1) [ y (k - 1), , y (k - na ) u (k - 1), , u (k - nb )]
θ [-a1 , , - ana
b1 , , bnb ]
• 设在k-1时刻和k时刻,系统的参数估计结果为
在第k次递推时,我们已获得新的观测数据向量 (k-1)和 y(k),则记 Φk-1=[(0), (1), ..., (k-2)]T
ˆ θ(k -1) (Φ1Φk 1 )1 Φ1Yk 1 k k 1 ˆ θ(k ) (Φk Φk ) Φk Yk
使用最小二乘法法进行系统辨识的两种方法
递推最小二乘法辨识与仿真现在有如下的辨识仿真对象：图中， )(k v 是服从N )1,0(分布的不相关随机噪声。
且)(1-zG )()(11--=z A z B ，)(1-z N )()(11--=zC zD ， （1）⎪⎪⎩⎪⎪⎨⎧=+==+-=--------1)(5.00.1)()(7.05.11)(121112111z D z zz B z C z z a z A选择上图所示的辨识模型。
仿真对象选择如下的模型结构：)()2()1()2()1()(2121k w k u b k u b k y a k y a k y +-+-=-+-+可得系统模型为：)()2(5.0)1()2(7.0)1(5.1)(k w k u k u k y k y k y +-+-=-+-- 递推最小二乘法的推导公式如下：)]1(ˆ)()()[()1(ˆ)(ˆ--+-=k k k z k k k θh K θθτ )1()]()([)]()()1([)(11--=+-=--k k k k k k k τP h K I h h PP τ1]1)()1()()[()1()(-+--=k k k k k k h P h h P K τ相关程序如下：% exp053 %%递推最小二乘法程序%clear%清理工作间变量L=55;% M序列的周期y1=1;y2=1;y3=1;y4=0;%四个移位寄存器的输出初始值for i=1:L;%开始循环，长度为Lx1=xor(y3,y4);%第一个移位积存器的输入是第3个与第4个移位积存器的输出的“或”x2=y1;%第二个移位积存器的输入是第3个移位积存器的输出x3=y2;%第三个移位积存器的输入是第2个移位积存器的输出x4=y3;%第四个移位积存器的输入是第3个移位积存器的输出y(i)=y4;%取出第四个移位积存器幅值为"0"和"1"的输出信号，if y(i)>0.5,u(i)=-0.03;%如果M序列的值为"1"时,辨识的输入信号取“-0.03”else u(i)=0.03;%当M序列的值为"0"时,辨识的输入信号取“0.03”end%小循环结束y1=x1;y2=x2;y3=x3;y4=x4;%为下一次的输入信号做准备end%大循环结束，产生输入信号uw=normrnd(0, sqrt(0.1), 1, 55);%加入白噪声figure(1);%第1个图形，伪随机序列stem(u),grid on%以径的形式显示出输入信号并给图形加上网格z(2)=0;z(1)=0;%取z的前两个初始值为零for k=3:55;%循环变量从3到55z(k)=1.5*z(k-1)-0.7*z(k-2)+u(k-1)+0.5*u(k-2)+w(k);%给出理想的辨识输出采样信号endc0=[0.001 0.001 0.001 0.001,0.001]';%直接给出被辨识参数的初始值,即一个充分小的实向量p0=10^6*eye(5,5);%直接给出初始状态P0，即一个充分大的实数单位矩阵E=0.000000005;%相对误差E=0.000000005c=[c0,zeros(5,54)];%被辨识参数矩阵的初始值及大小e=zeros(5,55);%相对误差的初始值及大小for k=3:55; %开始求Kh1=[-z(k-1),-z(k-2),u(k-1),u(k-2),w(k)]'; x=h1'*p0*h1+1; x1=inv(x); %开始求K(k)k1=p0*h1*x1;%求出K的值d1=z(k)-h1'*c0; c1=c0+k1*d1;%求被辨识参数ce1=c1-c0;%求参数当前值与上一次的值的差值e2=e1./c0;%求参数的相对变化e(:,k)=e2; %把当前相对变化的列向量加入误差矩阵的最后一列c0=c1;%新获得的参数作为下一次递推的旧参数c(:,k)=c1;%把辨识参数c 列向量加入辨识参数矩阵的最后一列p1=p0-k1*k1'*[h1'*p0*h1+1];%求出 p(k)的值p0=p1;%给下次用if e2<=E break;%若参数收敛满足要求，终止计算end%小循环结束end%大循环结束c;%显示被辨识参数e;%显示辨识结果的收敛情况%分离参数a1=c(1,:); a2=c(2,:); b1=c(3,:); b2=c(4,:);d1=c(5,:);ea1=e(1,:); ea2=e(2,:); eb1=e(3,:); eb2=e(4,:);figure(2);%第2个图形i=1:55;%横坐标从1到55plot(i,a1,'r',i,a2,':',i,b1,'g',i,b2,':',i,b1,'k') %画出a1，a2，b1，b2的各次辨识结果title('系统辨识结果')%图形标题如果系统中的参数为：⎪⎪⎩⎪⎪⎨⎧+-=+==+-=----------211211121112.01)(5.00.1)()(7.05.11)(z z z D z zz B z C z z a z A 那么系统模型机构为：)2(2.0)1()()2(5.0)1()2(7.0)1(5.1)(-+--+-+-=-+--k w k w k w k u k u k y k y k y 相关程序如下： % exp054.m % %递推最小二乘法编程%clear%清理工作间变量 L=55;% M 序列的周期y1=1;y2=1;y3=1;y4=0;%四个移位寄存器的输出初始值 for i=1:L;%开始循环，长度为Lx1=xor(y3,y4);%第一个移位积存器的输入是第3个与第4个移位积存器的输出的“或” x2=y1;%第二个移位积存器的输入是第3个移位积存器的输出 x3=y2;%第三个移位积存器的输入是第2个移位积存器的输出 x4=y3;%第四个移位积存器的输入是第3个移位积存器的输出 y(i)=y4;%取出第四个移位积存器幅值为"0"和"1"的输出信号，if y(i)>0.5,u(i)=-0.03;%如果M序列的值为"1"时,辨识的输入信号取“-0.03”else u(i)=0.03;%当M序列的值为"0"时,辨识的输入信号取“0.03”end%小循环结束y1=x1;y2=x2;y3=x3;y4=x4;%为下一次的输入信号做准备end%大循环结束，产生输入信号uw=normrnd(0, sqrt(0.1), 1, 55);figure(1);%第1个图形，伪随机序列stem(u),grid on%以径的形式显示出输入信号并给图形加上网格z(2)=0;z(1)=0;%取z的前两个初始值为零for k=3:55;%循环变量从3到55z(k)=1.5*z(k-1)-0.7*z(k-2)+u(k-1)+0.5*u(k-2)+w(k)-w(k-1)+0.2*w(k-2);%给出理想的辨识输出采样信号endc0=[0.001 0.001 0.001 0.001,0.001]';%直接给出被辨识参数的初始值,即一个充分小的实向量p0=10^6*eye(5,5);%直接给出初始状态P0，即一个充分大的实数单位矩阵E=0.000000005;%相对误差E=0.000000005c=[c0,zeros(5,54)];%被辨识参数矩阵的初始值及大小e=zeros(5,55);%相对误差的初始值及大小for k=3:55; %开始求Kh1=[-z(k-1),-z(k-2),u(k-1),u(k-2),w(k),w(k-1),w(k-2)]';x=h1'*p0*h1+1; x1=inv(x); %开始求K(k)k1=p0*h1*x1;%求出K的值d1=z(k)-h1'*c0; c1=c0+k1*d1;%求被辨识参数ce1=c1-c0;%求参数当前值与上一次的值的差值e2=e1./c0;%求参数的相对变化e(:,k)=e2; %把当前相对变化的列向量加入误差矩阵的最后一列c0=c1;%新获得的参数作为下一次递推的旧参数c(:,k)=c1;%把辨识参数c 列向量加入辨识参数矩阵的最后一列p1=p0-k1*k1'*[h1'*p0*h1+1];%求出 p(k)的值p0=p1;%给下次用if e2<=E break;%若参数收敛满足要求，终止计算end%小循环结束end%大循环结束c;%显示被辨识参数e;%显示辨识结果的收敛情况%分离参数a1=c(1,:); a2=c(2,:); b1=c(3,:); b2=c(4,:);d1=c(5,:);d2=c(6,:);d3=c(7,:); figure(2);%第2个图形i=1:55;%横坐标从1到55plot(i,a1,'r',i,a2,':',i,b1,'g',i,b2,':',i,d1,'k',i,d2,':',i,d3,'*') %画出a1，a2，b1，b2的各次辨识结果title('系统辨识结果')%图形标题。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
实验4 递推最小二乘法的实现
实验报告
哈尔滨工业大学
航天学院控制科学与工程系
专业：自动化
班级：1040101
姓名：
日期：2013 年10 月23 日
1．实验题目： 递推最小二乘法的实现 2．实验目的：
(n a y k -()(n b u k n ξ+-+,
,n a ，0,n b b 为待辨识的未知参数，(k ξ为系统的输出，为系统的输入。
分别测出n N +个输出、(3),
(),(1),(2),
()y y n N u u u n N ++，则可写出矩阵形式，有
0()
(1)(1)(1)((1)(2)
(2)
(2)((1)
()()
()(n n y n y u n u n a y n y u n u n b y n N y N u n N u N n b ξξξ⎢⎥
--+⎤⎡⎤⎡⎢⎥⎢⎥⎢⎢⎥-+-+⎢⎥⎢+⎢⎥⎢⎥⎢⎢⎥⎢
⎥⎢⎢⎥-+--+⎣⎦⎣⎢⎥⎢⎥⎣⎦
（2）
101)(1),,()()n n a n a b y n N n N b ξξθξξ⎡⎤⎢⎥
++⎤⎡⎤⎢⎥⎥⎢⎥⎢⎥
⎢
⎥⎢⎥==⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥++⎣⎦⎣⎦⎢⎥
⎢⎥⎣⎦，()
(1)(1)(1)(1)(2)
(2)
(2)1)
()()
()y n y u n u n y u n u N y N u n N u N -+⎤⎥+-+⎥⎥⎥+--+⎦
y θξ=Φ+ （维输出向量；ξ为N 维噪声向量；θ为21n +维参数向量；Φ测量矩阵。
为了尽量减小噪声ξ对θ估值的影响，应取2N n >+
()T θ=ΦΦ为了实现实时控制，必须采用递推算法，这种辨识方法主要用于在线辨设已获得的观测数据长度为N ，将式N N Y θ=ΦN θ表示θ的最小二乘估计，则
(T
N N
θ=ΦΦ(T N N P =ΦΦT N N N P θ=Φ如果再获得一组新的观测值(u n N +11T
N ψθ++=+(1),n N ++](1)(1)
(1)y N u n N u N -++++
）合并，并写成分块矩阵形式，可得
T 111N N N N N Y y ξθξψ+++⎡⎤
Φ⎡⎤⎤⎢⎥⎢⎥⎢⎥=+⎢
⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦
（9）于是，类似地可得到新的参数估值
1
T T
1T T T 1111N N
N N N N N N N Y y θψψψ-+++++⎧⎫⎡⎤⎡⎤⎡⎤ΦΦΦ⎡⎪⎪⎢⎥⎢⎥⎢⎥⎪⎪⎢⎥=⎨⎬
⎢⎥⎢⎥⎢⎥⎢⎥⎪⎪⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦
⎣⎦⎣⎦⎣⎦
⎪⎪⎩⎭ )
T
1T 1111N N N N N T N N N N Y P y Y y ψψ+++++⎡⎤Φ⎡⎤
⎢⎥⎢⎥⎢⎥
⎢⎥⎢⎥⎢⎥⎣⎦
⎣⎦
Φ+
)
T T 111T 11
N N
N N N N N ψψψψ
++-++⎪⎢⎥⎢⎥⎪⎬
⎢⎥⎢⎥⎪⎢⎥⎢⎥⎣⎦⎣⎦⎪⎭+ )
1T 11N N ψ
-++ 与P 的递推关系式出发，导，最终可求得递推最小二乘法辨识公式：
()T
1111N N N N N N K y θθ
ψθ++++=+- ()1
T 11111N N N N N N K P P ψψψ-++++=+ ()1
T T 11111N N N N N N N N P P P P ψψψψ-++++=-+ 为了进行递推计算，需要给出N P 和N θ的初值0P 和0θ。
2
000,P c I θ==是充分大单位矩阵，则经过若干次递推之后能得到较好的参数估计。
4．实验对象或参数
5．程序框图
6．程序代码
7．实验结果及分析
本实验给定系统为
12012()(1)(2)()(1)(2)()y k a y k a y k b u k b u k b u k k ξ=----++-+-+
要求运用递推最小二乘法对所给系统参数进行在线识别，其主要目的是设计状态反馈控制律，因此首先采用状态空间分析法对系统进行分析。
闭环系统原理图如下图所示。
根据系统辨识模型，对于该二阶系统，不妨设状态量1
()x y k =，
21(1)x x y k ==+，所以该系统状态空间表达式为：
12(1)0(2)1x y k X AX BV x y k ξ+⎡⎤⎡⎤⎡⎤
===++⎢⎥⎢⎥
⎢⎥+⎣⎦⎣⎦⎣⎦ Y CX =
其中：210
1A a a ⎡⎤=⎢
⎥--⎣⎦，01B ⎡⎤
=⎢⎥⎣⎦
，[]10C =，[]0
1
2(2)(1)()u k V b b b u k u k +⎡⎤
⎢⎥=+⎢⎥⎢⎥⎣⎦。
暂不考虑干扰信号，设状态反馈控制律为：
()U k KX V =+
[]1212K k k ⨯=为反馈增益系数矩阵，11V ⨯为系统参考输入，则状态反馈闭环
系统的状态空间表达式为
()X A BK X BV
Y CX
=++=
闭环系统传递函数：
[]1
2
1221()()1
()()
k W s C sI A BK B
s a k s a k -=-+=+-+-
设期望闭环极点1,2
2s j =-±1216k a =-，214k a =-。
通过最小二乘法辨识出1ˆa 、2ˆa 、0ˆb 、1ˆb 、2ˆb ，即可以得到实际控制律： [][]120
12ˆˆˆˆˆ()164(2)(1)()T U k a
a X
b b b u k u k u k ⎡⎤=--+++⎣⎦
在实验给定的系统辨识模型以及参数下，按照实验要求的正弦信号参考输入与白噪以及理论计算获得的控制律基础上进行MATLAB 仿真，具体得到以下结果分析。
输出data ：
下面为循环次数为25，c 为1000时的参数收敛情况：
由上图可知，在第14步时，参数基本稳定。
当c取不同值的时候，探索对系统辨识的影响。
如下表：
从表中可以看出，无论 c 取多大的值，最终都能得到系统参数真实值，参数收敛的快慢不同，c 较小时所需步数较多，c 较大时所需步数较少。
当 c 达到一定程度之后，收敛递推步数逐渐稳定，但是不可能小于N=2n+1=5 步。
方程个数 N只需要满足 N=2n+1 的条件即可以获得方程的解；同样的，在递推最小二乘法中，递推步数同样要大于 2n+1才能获得较为准确的结果。
8．结论
本实验中，递推最小二乘法的在线算法需要设计状态反馈控制规律。
经过状态反馈后，无论c取值多少，都能在一定的步数后得到真实的辨识参数。
c值越大，收敛越迅速。
但是收敛步数不可能小于N=2n+1=5 步。
相对最小二乘法，递推最小二乘法可以实现实时控制，所需要的数据量较少。
经过这次实验，不仅最小二乘法的递推过程更加熟悉，而且认识了状态反馈控制规律在辨识中的重大作用。