2动态过程数学模型参数估计的最小二乘方法

第二章 参数估计的最小二乘方法Least Squares§2—１静态线性模型参数的最小二乘估计（多元线性回归）一、 什么是最小二乘估计系统辨识三要素：模型，数据，准则。

例： y = ax + ε其中：y 、x 可测；ε — 不可测的干扰项；a —未知参数。

通过 N 次实验，得到测量数据 y k 和x k k = 1、2、3 …，确定未知参数 a 称“参数估计”。

使准则 J 为最小 ：令：∂ J / ∂ a = 0 , 导出 a = ?称为“最小二乘估计”，即残差平方总和为最小的估计，Gauss 于 1792年提出。

min)(21=-=∑=k N k k ax y J 0)(21=--=∂∂∑=k k N k k ax y x a J二、多元线性回归线性模型 y = a 0+ a 1x 1+ + a n x n + ε 式(2 - 1- 1)引入参数向量: θ = [ a 0,a 1, a n ]T (n+1)*1进行 N 次试验，得出N 个方程：y k = ϕk T θ + εk ; k=1、2…、N 式(2 -1- 2)其中：ϕk = [ 1,x 1,x 2, ,x N ] T (n+1) *1方程组可用矩阵表示为y = Φ θ + ε 式（2 -1- 3）其中：y = [ y 1,y 2, 。

,y N ]T (N *1) ε = [ ε1, ε2, 。

,ε N ]T (N *1) N *(n+1) 估计准则有：= （y — Φ θ）T ( y — Φ θ)(1*N) ( N *1)⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=T N T T nN N n n x x x x x x ϕϕϕφ....1...........1...121121211121)(θϕT k N k k y J -=∑=[]⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡----=)(..)(*)(...)(1111θϕθϕθϕθϕT N N T T N N T y y y y JJ = y T y + θT ΦT Φ θ -y T Φ θ - θT ΦT y= y T y + θT ΦT Φ θ - 2 θT ΦT y 式（2 -1- 4）假设：（ΦT Φ）(n+1)(n+1) 满秩，由 利用线性代数的以下两个矩阵对向量求偏导数的公式:A x A x T =∂∂)( 和 Ax xAx x T 2)(=∂∂(必须A 为对称阵) 由ΦΦT 为对称阵 有: y y T T T ΦΦ=∂∂θθ)( 和 θθθθΦΦ=∂ΦΦ∂T T T2)( 所以:y y y y J T T T T T T T ΦΦΦΦΦΦ22)2(-=-+∂∂=∂∂θθθθθθ 令上式等于零，解出参数估计向量：θ Ls =（ΦT Φ）-1 ΦT y 式（2 -1- 5）令：P = （ΦT Φ）-1 则参数估计向量 θ Ls = P ΦT y参数估计向量 θ Ls 被视为以下“正则方程”的解：（ΦT Φ）θ = ΦT y 式（2 -1- 6）注：为了便于区别, 我们用红体字符表示估计量或计算值，而用黑体表示为参数真值或实际测量值。

又 02)22()(22>ΦΦ=Φ-ΦΦ∂∂=∂∂T T T y J θθθθ 所以θ Ls 使J 为最小。

0=∂∂θJ三、关于参数最小二乘估计θLs性质的讨论以上求解参数最小二乘估计θLs时并为对{ εk }的统计特性做任何规定，这是最小二乘估计的优点。

当{ εk }为平稳零均值白噪声时，则θLs有如下良好的估计性质：a) 参数最小二乘估计θLs是y的线性估计θLs= PΦT y 是y 的线性表出；b) 参数最小二乘估计θLs是无偏估计，即E θLs= θ（参数真值）[ 证明]：E θLs= E[ PΦT y ]= PΦT E( y ) = PΦT E ( Φθ + ε ) =PΦTΦθ + E( ε ) = θ + 0 = θc) 最小二乘估计θLs 的估计误差协方差阵是σ2P (n+1)(n+1)即：E [ ( θLs- θ ) (θLs-θ )T ] = σ2P[ 证明]：E [ ( θLs- θ ) (θLs - θ )T ] = E [ PΦT ( y -Φθ) ( y-Φθ)TΦP] = E [ PΦT εεT ΦP ] = PΦT E ( εεT) ΦP =PΦTσ2 I N*NΦP = PΦTΦσ2P=σ2Pd) 若{ εk }为正态分布零均值白噪声时，则θLs是线性无偏最小方差估计（证明从略）。

如若{ εk }是有色噪声，则θLs不具有上述性质，即为有偏估计。

四、最小二乘估计θLs的几何意义和计算问题1. 最小二乘估计的几何意义最小二乘估计的模型输出值为y k= ϕk T θLs k = 1,2,…N输出实际测量值与模型输出值之差叫残差：εk = y k –y k模型输出向量为y = ΦθLs ，而残差向量为：ε= y –y= y–ΦθLsΦT ε= ΦT y –ΦTΦ（ΦTΦ）-1 ΦT y=ΦT y –ΦT y= 0即残差向量ε与由测量数据矩阵Φ的各个向量：Φ1, Φ2 ,…, Φ张成的超平面（估计空间）正交，而最小二乘模型输出向量y 为N实际输出向量y 在估计空间上的正交投影，这就是最小二乘估计的几何意义。

最小二乘估计的几何意义2. 关于最小二乘估计计算中的病态问题估计参数向量θLs一般是求解正则方程：（ΦTΦ）θ = ΦT y 式（2 -1- 6）得出。

可以利用消元法等一系列求解多元线性一次方程组的方法，计算得出，其有解的条件是（ΦTΦ）=P –1 矩阵非奇异（行列式数值大于零）。

但有时在求解式（2 -1- 6）方程组是会出现矩阵接近于奇异（行列式数值接近于零），即所谓“病态”的情况。

由此导致参数估计的结果不稳定，不可信。

出现上述情况的原因可能是由于①被辨识的过程受到的外加激励不够；②采样间隔太密；③A/D转换的位数太短，计算舍入误差累计所致。

为解决最小二乘计算中可能出现的病态问题，提出了不少改进算法，例如：Householder变换法、改进的平方根法和U—D分解算法。

后者是Bierman 1977提出的改善P阵计算性质（对称性、正定性和稳定性）而又不增加计算量的算法。

正定P阵可以分解成一个上三角阵U（其对角线元素都为1）和一个对角阵DP = U D U T由此可解决最小二乘计算中可能出现的病态问题，具体可参阅关于《计算方法》的文献。

总之，我们在使用最小二乘的辨识方法时，应该注意避免出现和克服病态问题。

应用举例在建立生产过程的静态模型时，特别是在机理不清之时常用多元线性回归方法，例如：水泥凝固放热量与水泥成分的关系模型y = a0+ a1x1+ a2x2+ a3x3+ a4x 4 +y 水泥凝固时的放热量（卡/克）; x1 ~ x1 水泥的几种成分。

五、非线性最小二乘法（Nonlinear Least Square）以“误差平方总和为最小”的估计准则，估计非线性模型参数的方法。

假设非线性静态系统模型为y = f ( x, θ )+ ε非线性模型 f 的形式是已知的，参数 θ 未知。

经过N 次实验，取得N 组数据（x 1, y 1） （x N , y N ）。

准则：需要用优化算法求解，常用的算法有两类——搜索法和迭代法。

前者如单纯型法；后者如梯度法、高斯法、牛顿—拉夫森法、变尺度法等等。

该类方法也还可应用于动态模型和时间连续模型的参数估计。

单纯型法是先给出参数空间的几个猜想点，构成正多面体，计算各点的目标函数值，比较各值后舍去最差的点，按照反射、开拓、收缩等步骤确定新的估计点，直到预定的精度要求后停止搜索。

迭代法是先猜想一个估计的初值，确定一个向量为可接受的方向和步长，进行迭代计算 θk+1 –θ k =μ•v … 。

具体内容可阅有关计算方法的参考书籍。

近年来发展出一系列基于生物进化论的优化新算法，如遗传算法（Genetic Algorithms ）(基于改进遗传算法的系统辨识方法. 北方工业大学学报. 第10卷, 第1期,P62-67.)和免疫算法(Immune Algorithms)，使得优化算法的性能得到了很大的改善。

∑=-=N k k k x f y J 12)]([θ，§2—2 动态过程参数估计的线性最小二乘法一、模型类：考虑CAR 模型式（2 -2- 1）其中{ y(k) }和{ u(k) }为可测的输出和输入，{ ε(k) }为不可测的随机干扰。

上式还可表示成：A(z -1)y(k)= B(z -1)u(k)+ ε(k) 式（2 -2- 2） 其中：A(z -1) = 1+ a 1z - 1+ a 2z - 2+ 。

+a n z - nB(z -1) = b 1z - 1+ b 2z - 2+ +b n z - n 还可表示为式（2 -2- 3） 其中：当进行了 k = 1-n ,2-n ,..,0,1,2,…,N 共计（N+n ）次采样，得到N 个方程：)()(k k y T kεθϕ+=)](),...,1(),(),...,1([],...,,,,....,,[2121n k u k u n k y k y b b b a a a T kn n T k ------==ϕθ∑∑==+-=-+n i n i i i k i k u b i k y a k y 11)()()()(εy(1)= -a 1y(0) -… -a n y(1-n)+b 1u(0)+ … +b n u(1-n)+ ε(1)y(2)= -a 1y(1) -… -a n y(2-n)+b 1u(1)+ … +b n u(2-n)+ ε(2)………………………………………………y(N)=-a 1y(N-1)-…-a n y(N -n)+b 1u(N-1)+…+b n u(N-n)+ε(N)用矩阵表示成y N = ΦN θ + ε N 式（2 -2- 4） 其中：y N = [ y(1), y(2), … , y(N) ]Tε N = [ ε (1), ε (2), … , ε (N) ] T二、 参数最小二乘估计 θ Ls 的导出估计准则为式（2 -2- 5）由 解出： ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡--------------=)(...)1(.........)2(...)1()1(...)0(,)(..)1(......,)2(..)1(,)1(..)0(n N u N u n u u n u u n N y N y n y y n y y N φ0=∂∂θJ∑=-=N k T k k y J 12])([θϕθ N =（Φ N T Φ N ）-1 Φ N T y N 式（2 -2- 6） 上式可视为以下正则方程的解（Φ N T Φ N ）θ N = Φ N T y N 式（2 -2- 7）称为最小二乘的“一次完成算法”，是离线算法，有解的条件是（Φ N T Φ N ）2n *2n 满秩。