系统辨识之最小二乘
系统辨识之最小二乘
周耀龙 144173248
摘 要 最小二乘法是是一种经典的数据处理方法。在系统辨识领域中，最小二乘法是一种得到广泛应用的估计方
法, 可用于动态系统、静态系统、线性系统、非线性系统。在随机的环境下，利用最小二乘法时，并不要求观测数据提供其概率统计方面的信息，而其估计结果，却有相当好的统计特性。
关键词 最小二乘法；系统辨识法
在科学研究中，为了揭示某些相关量之间的关系，找出其规律，往往需要求解其函数解析式。一种方法是采用插值逼近法，即所构造的近似函数f (x ) 在已知节点x i 上必须满足f (x i ) =y i 要求逼近函数f (x i ) 与被逼近函数f (x ) 在各已知点x i 处的误差为零，即要求f (x ) 的曲线必须通过所有的点，常用的插值法有拉格朗日(Lagrange)插值，牛顿(Newton)插值，埃尔米特(Hermite)插值等。
另一方面，由于观测数据较多，一般不用插值法，而是用拟合的方法。即只要找到一条曲线，即能反映给定数据的一般趋势，又不出现局部较大的波动即可，只要f (x ) 与f (x ) 的偏差满足某种要求就行了。这种数据间的非确定关系需要统计方法来描述，最常用的方法就是数据拟合。数据拟合就是找一种函数的解析表达式或近似表达式来描述这组数据间的函数关系，通常用到最小二乘法。数据拟合的最小二乘法通过最小化误差的平方和，寻找数据的最佳函数。利用最小二乘法可以简便地求得未知的数据，并使得这些求得的数据与实际数据之间误差的平方和为最小。
最小二乘理论是有高斯（K.F.Gauss ）在1795年提出：“未知量的最大可能值是这样一个数值，它使各次实际观测值和计算值之间的差值的平方乘以度量其精度的数值以后的和最小。”这就是最小二乘法的最早思想。
最小二乘辨识方法提供一个估算方法，使之能得到一个在最小方差意义上与实验数据最好拟合的数学模型。递推最小二乘法是在最小二乘法得到的观测数据的基础上，用新引入的数据对上一次估计的结果进行修正递推出下一个参数估计值，直到估计值达到满意的精确度为止。
对工程实践中测得的数据进行理论分析，用恰当的函数去模拟数据原型是一类十分重要的问题，最常用的逼近原则是让实测数据和估计数据之间的距离平方和最小，这即是最小二乘法。最小二乘法是一种经典的数据处理方法。在随机的环境下，利用最小二乘法时，并不要求观测数据提供其概率统计方面的信息，而其估计结果，却有相当好的统计特性。
图 1曲线拟合流程图
图 2最佳平方逼近流程图
一、最小二乘法的系统辨识
设单输入—单输出线性定常系统的差分方程为：
x (k ) +a 1x (k -1) +×××+a n x (k -n ) =b 0u (k ) +×××+b n u (k -n ), k =1, 2, 3, ××× （1）
上式中：u (k ) 为输入信号；x (k ) 为理论上的输出值。x (k ) 只有通过观测才能得到，在观测过程中往往附加有随机干扰。x (k ) 的观测值y (k ) 可表示为
y (k ) =x (k ) +n (k ) （2）
将式（2）代入式（1）得
y (k ) +a 1y (k -1) +! +a n y (k -n ) =b 0u (k ) +b 1u (k -1) +! +
b n u (k -n ) +n (k ) +
åa n (k -i ) (3)
i i =1
n
我们可能不知道n (k ) 的统计特性，在这种情况下，往往把n (k ) 看做均值为0的白噪声。 设
x(k ) =n (k ) +åa i n (k -i ) （4）
i =1
n
则式（3）可以写成
y (k ) =-a 1y (k -1) -! -a n y (k -n ) +b 0u (k ) +b 1u (k -1) +! +
b n u (k -n ) +x(k ) (5)
在测量u (k ) 时也有测量误差，系统内部也可能有噪声，应当考虑它们的影响。因此假定而且还包含了u (k ) 的测量误差和系统内部噪声。假定x(k ) x (k ) 不仅包含了x (k ) 的测量误差，
是不相关舱机序列(实际上x(k ) 是相关随机序列) 。
分别测出n+N个输出输入值y (1), y (2), ! y (n +N ), u (1), u (2), ! , u (n +N ) , 则可以写出N 个方程，即
y (n +1) =-a 1y (n ) -a 2y (n -1) -! -a n y (1) + b 0u (n +1) +b 1u (n ) +! +b n u (1) +x (n +1) y (n +2) =-a 1y (n +1) -a 2y (n ) -! -a n y (2) +
b 0u (n +2) +b 1u (n +1) +! +b n u (2) +x (n +2)
!
y (n +N ) =-a 1y (n +N -1) -a 2y (n +N -2) -! -a n y (N ) +
b 0u (n +N ) +b 1u (n +N -1) +! +b n u (N ) +x(n +N )
上述N 个方程可写成向量—矩阵形式
éy (n +1) ùêy (n +2) ú
ú = ê
êú! êúy (n +N ) ëû-y (n ) ! -y (1) u (n +1) u (1) ùé
ê-y (n +1) ú! -y (2) u (n +2) u (2) êú ê" " " " úêú-y (n +N -1) ! -y (N ) u (n +N ) u (N ) ëû
´
éa 1ù
ê! ú
êúéx(n +1) ùêb n úêx(n +2) ú
ú （6） êú+ê
ú! êb 0úê
úê! úêx(n +N ) ëû
êúêëb n úû
设
éa 1ù
ê! ú
éy (n +1) ùéx(n +1) ùêúêy (n +2) úêx(n +2) úêb n ú
ú , q = êú, x=êú y =ê
êúêú! ! êb 0úêúêúê! úëy (n +N ) ûëx(n +N ) û
êúêb n ûúë
-y (n ) ! -y (1) u (n +1) u (1) ùé
ê-y (n +1) ú! -y (2) u (n +2) u (2) ú F=êê" " " " úêú-y (n +N -1) ! -y (N ) u (n +N ) u (N ) ëû
则式（6）可以写成 y =Fq+x (7)
上式中：y 为n 维输出向量；x为N 维噪声向量；q 为（2n+1）维参数向量；F为N ´(2n +1) 测量矩阵。因此N
q=F-1y
如果噪声x¹0，则
q=F-1y -F-1x
从上式可以看出噪声x对参数估计有影响，为了尽量减小噪声x对q估值的影响，应取N>(2n+1),即方程数大于未知数数目。在这种情况下，不能用解方程的办法来求q，而要采用数理统计的办法，以便减小噪声对q估值的影响。在给定输出向量y 和测量矩阵F 的条件下求系统参数q的估值，这就是系统辨识问题。可用最小二乘法来就q的估值。
二、最小二乘法的原理
以一个SISO 系统为例说明最小二乘法的原理。 假设一个SISO 系统如下图所示：
其离散传递函数为：
b n z -1 B ( z ) b 1 z + b 2 z + … + 3.1
G (z ) =
输入输出的关系为：
-1-1
-1-2-n
A (z )
=
1+a 1z -1+a 2z -2+…+a n z -n
u (k ) •G (z -1) +e (k ) =y (k ) 3.2
进一步，我们可以得到：
y (k ) ×A (z -1) =u (k ) ×B (z -1) +e (k ) 3.3
其中，扰动量e (k ) 为均值为0，不相关的白噪声。 将式3.3写成差分方程的形式：
y (k ) =-a 1y (k -1) -a 2y (k -2) -…-a n y (k -n ) +b 1u (k -1) +b 2u (k -2) +…b n u (k -n ) +e (k )
3.4 令j(k ) =[-
y (k -1) -y (k -2) …-y (k -n ) u (k -1) u (k -2) …u (k -n ) ]T
q =[a 1a 2…a n b 1b 2…b n ]
则式3.4可以写为：
y (k ) =jT (k ) q+e (k ) 3.5
将上述式子扩展到N 个输入、输出观测值{u (k ), y (k ) }，k=1,2,…，N+n。将其代入到式3.5中，写成矩阵的形式为：
Y =Fq+e 3.6
取泛函J (q) 为
-y (n -1) é-y (n )
ê-y (n -1) -y (n ) F=ê
ê" " ê
ë-y (n +N ) -y (n +N -1)
u (1) ùu (2) úú
…" " #" ú
ú
! -y (N ) u (n +N ) ! u (N ) û……
……
-y (1) -y (2)
u (n ) u (n +1)
J (q) =å(Y -Fq) 2=åe 2(n +i ) =e T •e =(Y -Fq) T (Y -Fq)
i =1
N
最小二乘法原理既是使J (q) 最小，对其求极值得：
由此可得系统的最小二乘法估计值为：
¶J ¶=[(Y -Fq) T (Y -Fq)]=0¶q¶q
q=(FT F) T FT Y
这样，我们就得到了系统的最小二乘估计值。
三、最小二乘法系统辨识的应用举例
系统辨识是通过建立动态系统模型, 在模型输入输出数据的基础上, 运用辨识方法对模型参数进行辨识，从而得到一个与所观测的系统在实际特性上等价的系统。应用最小二乘法对系统模型参数进行辨识的方法有离线辨识和在线辨识两种。离线辨识是在采集到系统模型所需全部输入输出数据后，用最小二乘法对数据进行集中处理，从而获得模型参数的估计值；而在线辨识是一种在系统运行过程中进行的递推辨识方法，所应用的数据是实时采集的系统输入输出数据，应用递推算法对参数估计值进行不断修正, 以取得更为准确的参数估计值。由于在线辨识方法具有实时采集系统输入输出数据，实时辨识模型参数, 且占据计算机存储量小的优点，因此与离线辨识相比, 在线辨识方法得到了更为广泛的应用。在线辨识的参数估计的最小二乘递推算法如下： ^θ(k+1) = ^θ(k)+K(k+1)[y(k+1)-xT(k+1)^θ(k)] K(k+1) = P(k)x(k+1)[1+xT(k+1)P(k)x(k+1)]-1 P(k+1) = P(k)-K(k+1)xT(k+1)P(k)
2
递推初值:^θ(0) =任意值；P(0) =αI, α取计算机容许的最大值。 式中x 与y 分别为系统的输入输出，θ为参数估计值,K 为增益矩阵， P(m) = (xm x m )
T
-1
其最优性准则函数为: J =
åe
i =1
m
2
(i )
其中m 为数据采集的次数，e 为残差向量。由于上述递推算法无法反映参数随时间变化的特点，新数据被大量的老数据所淹没，对于慢时变参数的辨识来说，这必然得不到跟踪参数变化的实时估计，因此又进一步有了改进的最小二乘递推算法，即带遗忘因子的渐消记忆的递推算法，该算法贬低老数据的作用，强调新数据的作用，选取遗忘因子λ，得到渐消记忆的最小二乘递推算法如下： ^θ(k+1) = ^θ(k)+K(k+1)[y(k+1)-xT(k+1)^θ(k)] K(k+1) = P(k)x(k+1)[λ+xT(k+1)P(k)x(k+1)]-1 P(k+1) =1λ[P(k)-K(k+1)xT(k+1)P(k)]
递推初值：^θ(0) =任意值；P(0) =αI ，α取计算机容许的最大值。 其最优性准则函数为： J =
2
ål
i =1
m
m -i
e 2(i )
其中加权系数0
四、最小二乘法例题分析
(一) 题目概述
求0-100o C 之间水的定压比热变化的数学模型问题。
(二) 原理
实验数据的平滑处理—插值与逼近
所谓“插值”，就是求取两测量点之间“函数值”的计算方法，常用的有“线性插值”和“三次
样条插值”。
线性插值所建立的数学描述/模型在插值点上是“非光滑的” 。三次样条插值可以较完美地逼近理想的数学描述/模型，其代价是计算量与存储空间的增加。
实验数据的统计处理—最小二乘法
对于随机型系统，其数据处理需要依据“数理统计”的理论与方法来处理，常用的方法是“最
小二乘法”。 目标
(x i , y i ), i =1,2,..., n
y =j (x )
要求是某给定函数类H 中的一个函数，并要求 能使 与 的差的平方和相对于同一函数类中的其他函数而言是最小的，即
j (x )
y i
j (x i )
å[y i -j(x i )]2=Min å[y i -j
(x i )]2
i =1
n n
jÎ
H i =1
(三) 辨识方法
“系统辨识”的基本原理与三要素
“数据、假设模型、准则”是系统辨识建模过程中的“三要素”。
(四) Matlab 程序
clear;clc; format long
%键入数据点 T=0:5:100;
C=[1.00762 1.00392 1.00153 1 0.99907 0.99852 0.99826 0.99818 0.99828 0.99849 0.99878 0.99919 0.99967 1.00024 1.00091 1.00167 1.00253 1.00351 1.00461 1.00586 1.00721];
a=polyfit(T,C,3);%曲线拟合 b=polyval(a,T);%真实数据点
plot(T,C,'-' ,T,b, '*'); %绘制拟合曲线以及数据点 %绘制误差曲线图 diff=abs(C-b); plot(T,diff);
(五) 结果
(六) 结果分析
最小二乘法特点： 上图为最小二乘法拟合误差图。由图可知，最大误差出现在T=0℃处，误差值约为0.001174。
l
l
l
原理易于理解（不需要数理统计方面知识）； 应用广泛（动态/静态系统，线性/非线性系统的辨识）； 所得的“估计值”具有条件最优的统计特性（即具有一致/无偏/有效的性质）。
最小二乘法缺陷：
l
l
当模型噪声是有色噪声时，最小二乘估计不是无偏、一致估计； 随着数据的增长，将出现所谓的“数据饱和”现象。
缺陷的解决方案：
遗忘因子法、限定记忆法、偏差补偿法、增广最小二乘、广义最小二乘、辅助变量法、二步法及多级最小二乘法等。
(七) 总结
当已经若干离散的工程数据时，我们可以通过应用最小二乘原理，对杂乱无章的数据进行离线拟合，最终画出平滑曲线。该曲线在一定程度上，可以预测更多工况情况下的数据输出，也能反映多数据的相关特性。书上给出的数学递推式相对复杂，我们可以应用 matlab 对数据进行简单快速地拟合操作，并方便地绘制平滑曲线图。通过拟合曲线与真实采样数据点的绘制，我们可以发现，拟合曲线与采样数据偏差不大，最大误差出现在温度值为零摄氏度时。故最小二乘法确实能够实现对离散数据的较好拟合，且 matlab 实现方便，该原理也因此得到广泛的认知与应用。
综上所述，最小二乘法就是以最小二乘原理为依据，同时解出一组未知参量的最佳值，最后确
定函数解析式的方法。同时，它所做出的曲线拟合能够清晰地表现出变量间的函数关系，还能通过数据点与曲线的偏离程度大致估计观测误差。特别适用于通过实验求解未知形式的函数关系式或者用简单的解析式近似复杂的表达式，是科学实验中对数据处理的一种重要方法。
五、参考文献
[1] 王青建. 数学史简编[M].北京. 科学出版社.1998.
[2] 丁天彪. 数值计算方法[M].河南. 黄河水利出版社,2003.
[3] 张池平. 计算方法[M].北京. 科学出版社,2006.
[4] 张韵华. 数值计算方法与算法[M].北京. 科学出版社,2006.