系统辨识大作业最小二乘法及其相关估值方法应用
学号：**********
姓名：***
基于最小二乘法的多种系统辨识方法研究
1.最小二乘法的引出
在系统辨识中用得最广泛的估计方法是最小二乘法(LS)。
设单输入-单输出线性定长系统的差分方程为
x(k)+a1x(k−1)+⋯+a n(k−n)=b0u(k)+⋯+b n u(k−n)，k=1，2，3，⋯
(5.1.1)
式中: u(k)为随机干扰；x(k)为理论上的输出值。x(k)只有通过观测才能得到，在观测过程中往往附加有随机干扰。x(k)的观测值y(k)可表示为
y(k)=x(k)+n(k)
(5.1.2)
式中: n(k)为随机干扰。由式(5.1.2)得
x(k)=y(k)−n(k)
(5.1.3)
将式(5.1.3)带入式(5.1.1)得
y(k)+a1y(k−1)+⋯+a n y(k−n)
n
=b0u(k)+b1u(k−1)+⋯+b n u(k−n)+n(k)+∑a i(k−i)
i=1
(5.1.4)
我们可能不知道n(k)的统计特性，在这种情况下，往往把n(k)看做均值为0的白噪声。
设
n
ξ(k)=n(k)+∑a i(k−i)
i=1
(5.1.5)
则式(5.1.4)可写成
y(k)=−a1y(k−1)−a2y(k−2)−⋯−a n y(k−n)+b0u(k)+b1u(k−1)+⋯
+b n u(k−n)+ξ(k)
(5.1.6)
在观测u(k)时也有测量误差，系统内部也可能有噪声，应当考虑它们的影响。因此假定ξ(k)不仅包含了x(k)的测量误差，而且包含了u(k)的测量误差和系统内部噪声。假定ξ(k)是不相关随机序列(实际上ξ(k)是相关随机序列)。
现分别测出n+N个随机输入值y(1)，y(2)，⋯，y(n+N)，u(1)，u(2)，⋯，u(n+ N)，则可写成N个方程，即
y(n+1)=−a1y(n)−a2y(n−1)−⋯−a n y(1)+b0u(n+1)+b1u(n)+⋯+b n u(1) +ξ(n+1)
y(n+2)=−a1y(n+1)−a2y(n)−⋯−a n y(2)+b0u(n+2)+b1u(n+1)+⋯+b n u(2) +ξ(n+2)
⋮
y(n+N)=−a1y(n+N−1)−a2y(n+N−2)−⋯−a n y(N)+b0u(n+N)
+b1u(n+N−1)+⋯+b n u(N)+ξ(n+N)
上述N个方程可写成向量-矩阵形式
[y(n+1) y(n+2)
⋮
y(n+N)]=[
−y(n)
−y(n+1)
⋯
⋯
−y(1)u(n+1)⋯u(1)
−y(2)u(n+2)⋯u(2)
⋮⋮ ⋮ ⋮⋮
−y(n+N−1)⋯−y(N)u(n+N)⋯u(N)
]×
[
a1
⋮
a n
b0
⋮
b n]
+[
ξ(n+1)
ξ(n+2)
⋮
ξ(n+3)
]
(5.1.7)
设
y=[y(n+1) y(n+2)
⋮
y(n+N)]，θ=
[
a1
⋮
a n
b0
⋮
b n]
，ξ=[
ξ(n+1)
ξ(n+2)
⋮
ξ(n+N)
]
Φ=[
−y(n)
−y(n+1)
⋯
⋯
−y(1)u(n+1)⋯u(1)
−y(2)u(n+2)⋯u(2)⋮⋮ ⋮ ⋮⋮
−y(n+N−1)⋯−y(N)u(n+N)⋯u(N)
]
则式(5.1.7)可写为
y=Φθ+ξ
(5.1.8)
式中：y为N维输出向量；ξ为N维噪声向量；θ为(2n+1)维参数向量；Φ为N×(2n+1)测量矩阵。因此式(5.1.8)是一个含有(2n+1)个未知参数，由N个方程组成的联立方程组。如果N<2n+1，方程数少于未知数数目，则方程组的解是不定的，不能唯一地确定参数向量。如果N=2n+1，方程组正好与未知数数目相等，当噪声ξ=0时，就能准确地解出
θ=Φ−1y
(5.1.9)
如果噪声ξ≠0，则
θ=Φ−1y−Φ−1ξ
(5.1.10)
从上式可以看出噪声ξ对参数估计是有影响的，为了尽量较小噪声ξ对θ估值的影响。在给定输出向量y和测量矩阵Φ的条件下求系统参数θ的估值，这就是系统辨识问题。可用最小二乘法来求θ的估值，以下讨论最小二乘法估计。
2.最小二乘法估计算法
设θ̂表示 θ的最优估值，ŷ表示 y的最优估值，则有
ŷ=Φθ̂
(5.1.11)
ŷ=[ŷ(n+1) ŷ(n+2)
⋮
ŷ(n+N)
]，θ̂=
[a1̂⋮
a n̂
b0̂
⋮
b n̂]写出式(5.1.11)的某一行，则有
ŷ(k)=−a1̂y(k−1)−a2̂y(k−2)−⋯−a n̂y(k−n)+b0̂u(k)+b1̂u(k−1)+⋯
+b n̂u(k−n)+ξ(k)=−∑a îy(k−i)
n
i=1+∑b îu(k−i)
n
i=0
，k
=n+1，n+2，⋯，n+N
(5.1.12) 设e(k)表示y(k)与ŷ(k)之差，即
e(k)= y(k)-ŷ(k)=
y(k)−ŷ(k)=y(k)—∑a îy(k−i)
n
i=1+∑b îu(k−i)
n
i=0
=
(1+a1̂z−1+⋯+a n̂z−n)y(k)−(b0̂+b1̂z−1+⋯+b n̂z−n)u(k)
=â(z−1)y(k)−b̂(z−1)u(k)，k=n+1，n+2，⋯，n+N
(5.1.13)式中
â(z−1)=1+a1̂z−1+⋯+a n̂z−n
b̂(z−1)=b0̂+b1̂z−1+⋯+b n̂z−n
e(k)成为残差。把k=n+1，n+2，⋯，n+N分别代入式(5.1.13)可得残差e(n+
1)，e(n+2)，⋯，e(n+N)。设
e=[e(n+1)e(n+2)⋯e(n+N)]T
则有
e(k)= y−ŷ=y−Φθ̂
(5.1.14) 最小二乘估计要求残差的平方和为最小，即按照指数函数
J=e T e=(y−Φθ̂)T(y−Φθ̂)
(5.1.15) 为最小来确定估值θ̂。求J对θ̂的偏导数并令其等于0可得
∂J
∂θ̂
=−2ΦT(y−Φθ̂)=0
(5.1.16)
ΦTΦθ̂=ΦT y
(5.1.17) 由式(5.1.17)可得θ的最小二乘估计
θ̂=(ΦTΦ)−1ΦT y
(5.1.18) 3.递推最小二乘法
为了实现实时控制，必须采用递推算法，这种辨识方法主要用于在线辨识。
设已获得的观测数据长度为N，将式(5.1.8)中的y和ξ分别用Y N，ΦN，ξ̅N来代替，
即
Y N=ΦNθ+ξ̅N
(5.3.1) 用θ̂N表示 θ的最小二乘估计，则
θ̂N=(ΦN TΦN)−1ΦN T Y N
(5.3.2)
