《参数估计理论与应用(第三章)课件.ppt》由会员分享，可在线阅读，更多相关《参数估计理论与应用(第三章)课件.ppt（77页珍藏版）》请在三一办公上搜索。
1、2022/11/26,第三章 参数估计理论与应用,第三章 参数估计理论与应用,3.1 参数估计的评价准则3.2 基于统计分布的参数估计方法3.3 基于模型的参数最小二乘估计 本章小结,2022/11/26,第三章 参数估计理论与应用,在许多情况下，观测数据所服从的概率模型已知的，而模型的未知部分是以未知参数形式出现的。 参数估计的基础是优化理论，即被估计的参数应该在某种准则下是最优的，以及任何获得最优的估计。 非参数估计方法不假定观测数据服从某种特定的概率模型。例如，频域上的谱估计与谱线拟合就是典型的非参数估计方法。,2022/11/26,第三章 参数估计理论与应用,3.1 参数估计的评价准则
2、,参数估计是通过样本去估计总体的某些数字特征或统计量。任何一个统计量都可作为参数的估计量，但其效果的优劣有所差别。3.1.1 无偏性、有效性与相容性 （1）无偏性 设样本的总体分布密度函数为 p(x;), 是未知参数。从总体中抽取容量为 N 的样本 x=x1, , xN , 用样本的估计量 来估计，如果希望多次估计中，平均的估计值没有偏差，即则称 是的无偏估计量。,2022/11/26,第三章 参数估计理论与应用,例3-1 样本均值是总体数学期望的无偏估计。 设x1, , xN 是随机过程 xk 的N个独立观测样本，如果参数是总体的数学期望Ex，即用样本的均值作为的估计量，对该估计量取期望值，
3、有 一个无偏估计量在多次估计中将不会产生系统偏差，但并不意味着有偏估计就不好。如果一个有偏估计是渐进无偏的，即,2022/11/26,第三章 参数估计理论与应用,那么它仍然有可能是一个好的估计。 考虑实随机过程xk的相关函数的两种估计量： 假定数据xk是独立观测的，容易验证 式中，Rx()=Exk+ xk 是随机数据xk的相关函数。 以上二式表明，估计量 1() 是无偏的，而 2()则是有偏的。但是， 2()是渐进无偏的，即,2022/11/26,第三章 参数估计理论与应用,渐进无偏估计量 2()是半正定的，而无偏估计量 1()却不一定是半正定的，故 2()的使用场合较多。 （2）有效性 如果
4、 1 和 2 是两个根据N个独立观测样本得到的无偏估计量，无疑地，对 的平均偏差较小是选择的标准之一。例如，如果则 1的值比 2 的值更密集地聚集在真值的附近。通常将方差（或协方差阵）在所有的无偏估计量中达到最小的 称为有效估计量。 例3-2 设x1,xN 是N个独立观测样本，若被估计参数,2022/11/26,第三章 参数估计理论与应用,=Ex，则对任何满足都是的无偏估计量。利用不等式 可得在估计总体的数学期望时，简单的算术平均比加权平均好。 （3）一致性 估计量的精度是与样本的容量 N 有关系的。一般说来，总是认为N 越大估计的效果应该越好。如果记依赖样本容量 N 的估计为 N ，当满足,
5、2022/11/26,第三章 参数估计理论与应用,则称 N 是的一致性估计量，或相容估计。 例3-3 设总体 x 具有均匀分布，分布密度为其中，1 和2 是未知参数。 总体样本的均值和二阶矩分别为（严格按定义计算）解得,2022/11/26,第三章 参数估计理论与应用,按矩的估计方法，用独立样本的均值和独立样本的二阶矩，分别作为总体均值和二阶矩的估计量，就有 下面说明 1 和 2 分别是1 和2 的相容估计。 设 y1,yN 是具有同分布的独立观测样本，根据大数定律，有令y=x2, 就有,2022/11/26,第三章 参数估计理论与应用,于是3.1.2 Fisher信息和Cramer-Rao不
6、等式 通常希望获得有效的参数估计量。但是，由于不存在导致最小方差无偏估计量的最佳算法，所以通常采用参数无偏估计的Cramer-Rao下限(或CR下界), 作为评价参数估计性能的测度。为了简洁叙述这一的评价测度，先定义一个重要的概念。 Fisher 信息 Fisher 信息用J（）表示，定义为（3.1.1）,2022/11/26,第三章 参数估计理论与应用,当考虑 N 个观测样本 X= x1,xN , 此时，联合条件分布密度函数可表示为 将式（3.1.1）中的p(x|)改为p(X|)就可给出N个样本变量X的Fisher信息的表达式。 定理（Cramer-Rao不等式） 设观测样本X= x1,xN
7、 , 若参数估计 是真实参数 的无偏估计，并且条件分布密度函数的p(X|) 对参数 的一、二阶偏导数存在，则有（3.1.2） 参数 的方差所能达到的下限（称为CR下限），即上式等号成立的充要条件是,2022/11/26,第三章 参数估计理论与应用,其中, 函数K()0，并与样本向量 X 无关。 当 为有偏估计量时，Cramer-Rao 不等式为 （3.1.3） 式中()为估计偏差，即()=E -，并假定b()是可微分的。 对于多个参数的情况，记=1，p，则用矩阵J() 表示Fisher信息，其元素Jij() 定义为（3.1.4）,2022/11/26,第三章 参数估计理论与应用,且Cramer
8、-Rao不等式变为矩阵不等式：（3.1.5） 上式表示无偏估计量的协方差矩阵cov( )与逆Fisher信息阵之差是一半正定矩阵。 Fisher信息是描述从观测数据中得到的 的 “信息” 测度，它给出利用观测数据估计参数的方差下界。但是，满足这一下界的估计量有的时候可能不存在。3.2 基于统计分布的参数估计方法 参数估计量的优劣取决于所采用的评价准则（或代价函数）和估计算法。现在介绍已知总体统计分布的两种最有效的参数估计方法：Bayes 估计和最大似然估计。,2022/11/26,第三章 参数估计理论与应用,3.2.1 Bayes 估计 在参数估计中，估计误差 - 通常不为零。因此，除了采用前
9、面介绍的无偏、有效和相容估计作为评价准则外，还可以利用估计误差的变化范围作为参数估计的测度，这种测度叫做代价函数，用符号C( ,)表示。常用的代价函数有绝对型、二次型和均匀型三种。,2022/11/26,第三章 参数估计理论与应用,本节仅介绍最常用的二次型代价函数，即 当总体的分布密度函数p(X|)已知时，利用X= x1,xN 进行参数估计，通常是采用代价函数的期望值作为评价参数估计量效果的测度，并称之为风险函数。使风险函数最小的参数估计叫做 Bayes 估计；基于二次型风险函数最小的估计称为最小均方误差（minimum mean square error, MMSE）估计。二次型风险函数定义
10、为（3.2.1） 根据条件概率公式，有,2022/11/26,第三章 参数估计理论与应用,其中，p( | x1,xN )是给定N个观测样本X= x1,xN 条件下 的后验分布密度函数。于是，式（3.2.1）可以写成（3.2.2） 为使风险函数RM M S E 最小，对上式取 的偏导，并令其结果为零，便得到由于p(x1,xN ) 是非负的，因此，RM M S E / =0, 等价于上式中=0。故有,2022/11/26,第三章 参数估计理论与应用,（3.2.3） 注意，在式（3.2.3）中，利用了以下事实： 由此可得出重要的结论：未知参数 的MMSE估计是给定样本X条件下的条件均值。 例3-4
11、某一随机参量x 服从高斯N(mx,Cx)分布，用仪器可测量其线性组合y ，即（1）式中，yN 维，kNM 维，x M维，e N 维。,2022/11/26,第三章 参数估计理论与应用,其中，测量误差 e 服从高斯N(0,Ce)分布；k 为给定的常数阵。假设 () e 与 x 独立； () e 与 x 相关，互协方差函数为Cxe 。 试分别求出两种情况下的MMSE估计x(y)和估计误差x (y)的协方差R x(y)。 解 如果将 x 看作未知参数，那么，根据上面讨论, x 的MMSE估计是给定观测样本y1,yN 时 x 的条件均值。因此，可利用公式（1.4.16）和（1.4.17）pp.29 （
12、2）（3）来求解。,2022/11/26,第三章 参数估计理论与应用,对式（1）两边取均值，得到 （4） 将式（1）和（3）代入有关定义式，得（5）（6）（7）,2022/11/26,第三章 参数估计理论与应用,（i）当 e与 x 互相独立，Cxe=0。将式（4）（7）代入式（2）和（3），得到x(y)的估计及协方差R x(y)（ii）当 e 与 x 相关，只需注意Cxe 0即可。 这个问题留给读者解决。请构造一组数据，在Matlab 平台上仿真这两种的估计结果。3.2.2 最大似然估计 最大似然估计（maximum likelihood estimate, ML估计）的基本思路是：在给定参数
13、条件下，将观测样本 x,K,2022/11/26,第三章 参数估计理论与应用,联合条件概率密度函数p(x|)视为真实参数 的函数，即似然函数L(x，) （包含未知参数的可能性函数），然后利用容量为 N 的观测样本x= x1,xN ，求出使L(x，)达到最大化的参数 作为=1，p的估计值。在数学上，通常将未知参数 的最大似然估计量记为式中是参数 的值域。故ML估计量 ML就是p(x|)的全局极大点。 由于对数函数是严格单调的，故 L(x，) 的极大点与ln L(x，)的极大点是一致的。通常，将ln L(x，)称为对数似然函数。于是，ML估计量 ML可由（3.2.4）,2022/11/26,第三章
14、 参数估计理论与应用,确定。如果 x1,xN 是N个独立的观测样本，则对数似然函数可写作（3.2.5） ML估计量 ML只要能够求出来，总是比较好的估计，它具有以下性质：最大似然估计是有效和一致估计；对于大的N，ML估计量 ML服从高斯分布，并且是无偏的，方差可达CR下界。 例3-5 设样本x= x1,xN 服从高斯分布N（m，），则其对数似然函数为,2022/11/26,第三章 参数估计理论与应用,分别求 lnL 关于 m 和2 的偏导，并令它们等于零，得到解得显然有 可见，均值的ML估计 ML 是无偏的，而方差的ML估计 ML是有偏的。但若将 ML N / (N-1)作为新的估计量，则该估
15、计是无偏的。,2022/11/26,第三章 参数估计理论与应用,计算L(x，)的相对于m 的二阶偏导数，有 由式（3.1.1）得Fisher 信息：Cramer-Rao不等式为等号成立的充要条件是事实上，我们有,2022/11/26,第三章 参数估计理论与应用,因此，只要取K(m)=N/2，ML估计 ML就可达CR下界2/N。这表明ML估计 ML是一有效估计量。 例3-6 （二元阵最大似然测向系统 ） 设二元阵布置在 x轴上，两个基元坐标分别为x1 和x2 ，如图3-2所示。如果取x1=0，则 x2=d，d为两传感器的位置间隔。假设信号为平面波，入射角为，则传感器1相对于传感器2的信号时延为
16、（3.2.6）式中，c 为声速。我们的问题是如何利用二元阵中两个输入过程的时差来测定目标的方位角。,2022/11/26,第三章 参数估计理论与应用,解 设两传感器的零均值接收过程可分别表示为其中，si 为单频平面波信号，wi (i=1，2)为零均值高斯噪声，二者互相独立。 如果采用图3-3所示的时延补偿方法，则单频平面波信号的归一化声程补偿（或指向）向量 v 在所考虑的二元阵中可表示为 下面，我们来推导信号的协方差矩阵和噪声的协方差矩阵，以便于求出观测样本的似,2022/11/26,第三章 参数估计理论与应用,函数。记输入信号和输入噪声的傅立叶系数为设信号和噪声的功率谱分别为S (n) 和N
17、 (n) ,那么，由公式（1.4.6） pp.26, n= 2n / T） 信号和噪声的协方差矩阵可分别表示为（3.2.7）,2022/11/26,第三章 参数估计理论与应用,于是，观测样本的似然函数可表示为（3.2.8）式中，X (1)= X1 (1)，X2(1)T ， ，X (TW) =X1 (TW)，X2(TW)T 是传感器的接收过程x=x1，x2T的傅立叶系数阵； T 是过程的持续时间（采样数据的长度），W 是接收过程的带宽。 容易验证，行列式| Cw + Cs | 与时延无关。于是，ML估计就是选择，使ln p(X |)最大,也即使式（3.2.8）的指数函数（3.2.9）,2022/
18、11/26,第三章 参数估计理论与应用,最大。下面，我们从式（3.2.9）出发，推导时延参数 的最大似然估计的等效形式。为此，首先引进下列求逆公式（3.2.10）式中，A为nn非奇异矩阵；g为n1列向量。证明留给请读者课外练习【利用恒等式 g(1+gHA-1g)=(A+ gHg)A-1g）】。 利用求逆公式，可知 ,2022/11/26,第三章 参数估计理论与应用,将上式代入式（3.2.9），略去与无关的量T/N (n)。因此，选择使式（3.2.9）最大，等价于使下式（3.2.11） 最大。现引入记号在此将 X(n，) 视为某时间函数 x(t, )在时间（t-T，t）内的傅立叶系数。将上述替换
19、量代入式（3.2.11）后，再应用周期函数的 Parseval 公式，就有,2022/11/26,第三章 参数估计理论与应用,略去无关紧要的常数项1/2，计算 z(x, ) 的结构如图3-4所示。调节时延，使输出 z(x, )达到最大，相应的时延就是真实时延的ML估计 ML。 根据ML估计的传递性，由式（3.2.6）可得真实方位的ML估计 （3.2.12）,xH0(t),z(x,),x1(t),x2(t),H0(),()2,图3-4 二元阵最大似然测向系统,exp(-j),2022/11/26,第三章 参数估计理论与应用,二元阵最大似然测向系统与二元阵似然比检测系统具有完全相同的结构。这是因为
20、：在 H1 情况下， p(X | ) 等价于 p(X | H1)，后者也可看作是时延参数的函数；而在 H0 情况下，p(X | H0) 与无关。因此，选取 使似然函数最大，也就是使似然比p(X | H1)/ p(X | H0)最大。由此可见，检测问题与参数估计问题是密切相关。 顺便指出，可用测向测距近似公式（3.2.13） 构成最大似然联合测向测距系统。其中，di 表示第i 个传感器与“基准” 传感器位置的间距；D 表示目标与“基准” 传感器位置之间的距离。,2022/11/26,第三章 参数估计理论与应用,3.3 基于模型的参数最小二乘估计,最小二乘法（Least square method
21、，LS）是一种不需要任何先验知识的参数估计方法。在被测系统的静态（稳态）模型和动态模型的参数辨识中，最小二乘法是最常用的参数估计方法，在测控技术领域获得了广泛的应用。3.3.1 最小二乘估计器及其统计特性 在一般的最小二乘问题中，线性系统的参数化模型可以表示为(3.3.1) 其中，u=u1,upT 是模型的输入向量，f1,fn 是u 的已知函数，也可以是未知输入的观测数据； 1, ,n 是待估计,2022/11/26,第三章 参数估计理论与应用,的参数，又称为回归系数； y 是系统的输出。 当 f1,fn 是u 的稳态响应状态或是实测的确定性变量，且 y 是系统的稳态输出，则称式(3.3.1)
22、是描述线性系统的静态模型；当 y 是u 的动态响应或瞬态观测数据，那末式(3.3.1)就是描述线性系统的动态模型。 为了估计未知参数i, 必须做实验来获得数据对u i yi 或 fk (u i) yi, i=1,2,N, k=1,2,n；N n 以构成训练数据。将数据对代入方程 (3.3.1)，可以获得一组线性方程： 用矩阵表示方法，将上式写成更简洁的形式，即,2022/11/26,第三章 参数估计理论与应用,(3.3.2)其中 为了唯一地识别出未知参数，通常要求 N n，即数据对的数目多于拟合参数的数目。满足所有 N 个方程的精确解是不可能的，因为观测数据难免受到噪声的污染，或者描述系统的参
23、数化数学模型不够精确。故必须考虑随机噪声或建模误差，在方程(3.3.2)中引入随机误差向量e，得到 (3.3.3),2022/11/26,第三章 参数估计理论与应用,参数的最小二乘估计 LS ，就是使目标函数 (3.3.4) 达到最小值的参数估计。为此，通常都采用求极值的方法。 将式(3.3.4)展开后，得到 对 求导数，有 J 极小化的条件是一般均假设T非奇异，于是，LS 有唯一的解：,2022/11/26,第三章 参数估计理论与应用,(3.3.5) 式中+表示的伪逆。 上述表示误差向量对整体平方误差有相同权重。可以进一步扩展，令每个误差项有不同的权重。设W 为所需的权值矩阵，它是对称和正定
24、的，则加权的目标函数为(3.3.6) 按上述求极小值的方法，可得加权的最小二乘估计量：(3.3.7)显然，当W 选为单位矩阵时， WLS = LS。 例3-7 考虑最简单的一维线性模型（静态的），即只有一个控制变量 u 的情形，这时模型的形式是,2022/11/26,第三章 参数估计理论与应用,求未知参数0 和1的LS估计量。 解 实际过程输出是模型的输出加上一随机误差项，即观测数据对ui, yi的结构应为式中，ei 称为模型的残差或观测噪声，一般认为是零均值、相互独立的随机序列，并具有相同的方差2。将上式写成矩阵形式： 根据式(3.3.5)，可得LS估计量：,2022/11/26,第三章 参
25、数估计理论与应用,如果进一步假定 ei 的分布是正态的，则容易验证，方差2 的ML估计量是 作为练习，请读者在Matlab平台上输入以下数据和函数：x=1 2 3 4 5; y=1.3 1.8 2.2 2.9 3.5; p, s=polyfit(x, y, 1)% 生成拟合一次多项式运行结果是：p=0.55 0.69，s=0.1643。即y=0.55x+0.69标准差为0.1643。,2022/11/26,第三章 参数估计理论与应用,例5-8 （可线性化的非线性静态模型曲线回归）假设有一个非线性模型的输出为 其中，x1，x2 为确定性输入变量，a，b和 c为待估计参数。 解 上式两边经简单的代
26、数运算，再同时取自然对数，可转化为一个线性模型：这说明变换后的输出 ln (y-1-1) 可以显式地表达为以 lnx1和 x2 为输入、以 lna，b 和 c 为参数的线性模型。因此，就可以按变换后的线性模型用最小二乘法来估计变换后的未知参数，然后，再根据变换后的估计参数计算出原参数。,2022/11/26,第三章 参数估计理论与应用,判定输入x-输出y之间的关系能否用一个线性模型来描述的标准，通常用互相关系数的大小来衡量：(3.3.8)xy 的绝对值越大，表示变量之间的线性关系越密切，因而线性回归的效果就越好。 例3-9 设某一结构参数 n，m 和 d 已知的离散线性系统，其差分方程的形式为
27、： (3.3.9),2022/11/26,第三章 参数估计理论与应用,其中，e(k)为噪声，(k)为输入-输出观测向量，为未知参数向量，且要求根据 N 次数据对 y(i), u(i), i= 1, 2 , N；N n+m+1来估计对未知参数。 解 将式(3.3.9)改写成矩阵形式，得到将数据写成下标形式，就有这样，未知参数向量 可按式(3.3.5)进行估计。,2022/11/26,第三章 参数估计理论与应用,考虑如下单输入-单输出系统：用Matlab中rarx函数进行系统辨识，程序如下：A=1 -1.5 0.7; % a0=1, a1=-1.5, a2=0.7B=0 0.3 0.2 0.5;
28、% b0=0, b1=0.3, b2=0.2, b3=0.5th0=arx2th(A,B,1,1); % 实际系统的ARX模型e=randn(200,1);u=idinput(200，prbs); % 高斯噪声和伪随机信号y=idsim(u e,th0);z=y u; % 模型仿真；输入-输出信号zna=2; nb=3; nk=1 % ARX模型的阶次thm,yhat=rarx(z,na nb nk,ng,0.1); % 根据z进行ARX模型参数辨识plot(y,-);grid % 作图，实际系统的输出曲线hold onplot(yhat,:) % 作图，辨识系统的输出曲线参数辨识结果thm：
29、 a1= -1.3798，a2=0.7039 ，b1=0.3007 ， b2=0.1170， b3=0.4243。,2022/11/26,第三章 参数估计理论与应用,应当指出，要求观测数据容量N n+m+1是为了保证T非奇异，降低过程噪声序列e(k)的影响，从而提高参数估计的精度。不论 e(k)是何种形式的噪声序列，式(3.3.5)总是成立的。换言之，噪声性质仅影响LS估计的统计特性。 下面介绍LS估计的统计特性。如果观测噪声或建模误差序列e(k)具有零均值和相同的方差，即则LS估计量 LS是无偏、有效和相容的，并具估计误差的协方差为 2(T) -1。 对于动态控制系统的辨识，输入信号u(t)
30、 必须满足持续激励条件，也即输入信号u(t)的频谱必须包含足够丰富（Sufficient rich）的频率成分，以保证充分激励受控对象的,2022/11/26,第三章 参数估计理论与应用,所有振型，从而使观测数据载有动态系统的主要信息。LS估计在满足持续激励条件时，是渐进无偏的，也称为估计的一致性。 式(3.3.9)，也称为CAR模型（即受控的AR模型），可以写成更简洁的形式 (3.3.10)式中 q 表示时间算子，d为整数，表示系统的滞后量；A(),B()分别为 q1 的降次幂多项式。 CAR模型满足一致估计（或相容估计）的条件为：（1）e(k)是白噪声序列；（2）u(k)的均值和协方差有界
31、；且满足（m+1）阶持续激励条件（或正定条件）：,2022/11/26,第三章 参数估计理论与应用,（3）u(k) 和 e(k) 相互独立。通常，u(k) 都采用伪随机二元序列。 只要选择恰当的模型阶次或最小二乘多项式阶次(参见taylor.m, Matlab)， 最小二乘法总是可以很好地拟合数据，但是，如果观测数据波动较大，将严重影响参数估计的准确性。对此，可采用数据预处理和数字滤波的方法加以解决。 检验模型准确性的最简单方法是准备另外一组输入-输出数据对，称为检验数据集，在参数估计时不用，待模型建立后，用这组数据对来验证所得模型的普适性或泛化能力。,2022/11/26,第三章 参数估计理
32、论与应用,上述讨论，均假设数学模型的阶次是已知的。实际上，对于动态系统，模型的阶次很少是预知的。检验模型阶次是否合适的一种简单而有效的方法是：评估不同阶次的理论模型对观测数据的拟合度，用拟合误差函数来描述。通常，当 n 或（n, m） 增大时 J(n) 或 J(n, m) 就会减小；而当 n (或n, m)大于模型的真实阶次 n0 (或n0, m0) 时，J 的减小就不显著了。由此，可以很方便地用多次实验的方法来确定模型的恰当阶次。 注意，对于J(n, m) 形式的拟合误差函数，一般应按正交实验法来确定模型的恰当阶次，以减少实验的次数。,2022/11/26,第三章 参数估计理论与应用,3.3
33、.2 递推最小二乘估计,在测控系统中，被测对象通常可以不断提供新的输入-输出数据。如果希望利用新的信息来改善估计精度，那么就应当采用递推估计算法，这不仅可避免观测数据矩阵的行数的不断“膨胀” ，而且可减少参数估计的计算量。 在推导最小二乘递推算法前，先引入一个与式(3.2.10)类似的矩阵求逆定理。设A和I+CA-1B 均是非奇异方阵，则 (3.3.11) 下面介绍最小二乘递推算法。为了简化符号，以下推导均用 代替 LS。 设N是时刻N为止的观测数据，N+1时刻的LS估计为,2022/11/26,第三章 参数估计理论与应用,式中参见式(3.3.2)和(3.3.9)于是有(3.3.12) 令 P
34、N =N TN-1 ，由求逆公式(3.3.11)知 (3.3.13)定义增益向量为,2022/11/26,第三章 参数估计理论与应用,将式(3.3.13)和(3.3.14)代入(3.3.12)，得到(3.3.15) 上式表明，新的估计量 N+1 等于前一时刻的估计量 N 与修正项 KN+1(yN+1 - TN+1 N) 之和，这是一切递推公式的共同特征。如果令代表基于前一时刻的估计量 N 对N+1时刻的预测。那么，递推估计提供的新息,2022/11/26,第三章 参数估计理论与应用,就是预测误差或拟合误差。因此，修正量的大小与新息成正比，而各校正分量的权，由增益向量决定。 在启动上述递推算法时
35、，必须知道初值 0 和P0 ，通常令其中，2 1。然后从得到的一组数据，按式(3.3.16)开始递推运算。 从物理上看，这种初值选取虽然初始误差较大，但校正的作用也大，因此这种递推算法是有效的。此外，还可以先取得N m+n+1组数据，算出,2022/11/26,第三章 参数估计理论与应用,作为初值，然后，按式(3.3.16)进行递推运算。 增益向量KN+1在递推运算过程是怎样变化的？先考察PN =N TN-1如果 N 大于估计参数的数目，而且，输入输出数据对含有足够的“信息”（满足充分激励条件），则NTN 通常是正定的。显然，当 N 趋于无穷大，NTN / N 接近于非奇异的常数阵。于是，有
36、可见，式(3.3.16)中自适应增益向量KN+1 随着每次迭代而递减，这意味着递推运算过程将逐渐收敛于参数空间的最优点。事实上，在白噪声或低噪声条件下，递推最小二乘,2022/11/26,第三章 参数估计理论与应用,估计是一种简便而又有效的算法。这种递推算法在递推过程中虽然没有保存全部先前的数据，但所有先前数据的影响却一直在起作用，故称为无限增长记忆的递推最小二乘估计。3.3.3 卡尔曼滤波器的递推算法（状态估计） 与参数估计器不同，卡尔曼滤波器主要是解决如何从被噪声污染的观测数据中估计出已知动态系统模型的状态，而不是动态系统模型的未知参数。然而，仅从算法上看，这二者是非常相似的。为了便于比较
37、二者的异同，我们在此不加证明地列出卡尔曼滤波器算法。 设离散定常系统的状态方程和观测方程分别为,2022/11/26,第三章 参数估计理论与应用,式中，X(k)是 n 维状态向量， u(k)是 m 维控制向量，w(k)是 p 维过程噪声，y(k)是 r 维观测向量，v(k) 是 r 维观测噪声； A，B，C，分别是相应维数的系统矩阵、控制矩阵、观测矩阵和过程噪声权矩阵。噪声的统计特性满足 假定到 k 时刻为止，观测数据为y(1),y(k)，要估计l时刻的状态，就有三种情况：l k，称为预测问题；l=k，称为滤波问题；lk，称为平滑问题。 下面主要介绍滤波问题。,2022/11/26,第三章 参
38、数估计理论与应用,滤波算法：预测算法：滤波增益：滤波误差的协方差：预测误差的协方差：初始条件：,2022/11/26,第三章 参数估计理论与应用,例3-10 考虑某一动态系统表示的空间导航问题，其中加速度为白噪声。被白噪声污染之后观测其位置。因此，过程的状态方程为（1）式中，X(t)=x1(t), x2(t)T, 表示过程的位置和速度；而 w(t) 是一均值为 0、方差为 1 的高斯白噪声。观测方程为（2）其中 v(t) 是一均值为0、方差为 10 的高斯白噪声。 解 首先把连续的状态方程离散化（参见Matlab中c2d函数）,2022/11/26,第三章 参数估计理论与应用,不妨设采样间隔T
39、=1，则可写出状态方程（1）和观测方程（2）的离散化表达式进一步假设则可按前面介绍的滤波算法进行递推计算：X0X(1|0), P(1|0)K(1)X(1|1); X(2|1), P(1|1)P(2|1)K(2)X(2|2)。 在Matlab中用Kalman函数仿真卡尔曼滤波器的设计。请读者在Matlab平台上完成例3-10的仿真计算。,2022/11/26,第三章 参数估计理论与应用,3.3.4 限定记忆的递推最小二乘估计 以上介绍的递推最小二乘法适用于定常系统的参数估计。对于时变系统，由于参数时变的信息显然更多地蕴藏在当前的观测数据中，而与先前观测数据的关系将逐渐减弱。因此，利用一切观测数据
40、对来决定的自适应增益向量KN+1，显然会削弱递推过程跟踪变化参数的能力。解决这一问题的一种简单方法是：当我们怀疑观测数据发生显著的变化时，就将当前的PN 设置为P0，重新进行参数估计，这是因为LS估计能快速地收敛到当前的最优参数。处理这一问题的另一种方法是对过去的数据引入带遗忘因子，逐渐削弱它们在参数估计中的作用。为此，采用加权的目标函数：,2022/11/26,第三章 参数估计理论与应用,其中 (3.3.17) 而且0 1，当=1，就退化为基本的递推最小二乘法。由式(3.3.7)知，在时刻N，使上述加权的目标函数最小化的估计为 (3.3.18) 每当取得一个新的数据后，就对沟渠以前的加权矩阵
41、乘以，于是N+1时刻的估计量为,2022/11/26,第三章 参数估计理论与应用,与前面略有不同，在此令 PN =N T WNN-1 ，则由求逆公式(3.3.11)，可知于是，按前面推导出式(3.3.15)思路，可得式中 遗忘因子将老的数据逐渐从“记忆”中去掉，因此这种使用数据信息的方式也叫做“渐消记忆”法，相应的算法称为带遗忘因子的递推最小二乘法，即,2022/11/26,第三章 参数估计理论与应用,关于的选取通常由经验或实验确定，一般范围为0.95 0.99 取得愈小，最新数据的权重就愈大，也就更适合于跟踪大的时变参数，但与此同时，估计器也可能会发生较大的波动从而加大估计的误差。 例5-1
42、0 考虑模型采样300次后变为,2022/11/26,第三章 参数估计理论与应用,试用两组模拟数据，一组不考虑噪声，一组是带观测噪声的数据，分别用不同的遗忘因子，对时变模型进行参数估计，并讨论估计结果。 解 用MATLAB rarx函数 进行带遗忘因子的系统辨识算法。程序如下：e=randn(300,1);u=idinput(300,prbs); % 产生高斯噪声和伪随机信号A1=1 0.8; B1=0 0.5; th0=arx2th(A1,B1,1,1); % a1=1, 0.8 ;b2=0, 0.5y1=idsim(u e,th0); % 初始模型仿真；A2=1 0.6; B2=0 0.3
43、; th0=arx2th(A2,B2,1,1); % a2=1, 0.8 ;b2=0, 0.5y2=idsim(u e,th0); % 采样300次后模型仿真；for k=1:300 y(k,1)=(300-k)/300)*y1(k,1)+(k/300)*y2(k,1); % 时变模型 end,2022/11/26,第三章 参数估计理论与应用,z=y u；na=1; nb=1; nk=1; % 产生输入-输出信号z和ARX模型的阶次thm,yhat=rarx(z,na nb nk,ff,0.97); % 带遗忘因子0.97的ARX参数辨识plot(thm(:,1),-);grid % 作图，时
44、变系统的a(k)输出曲线hold onplot(thm(:,2), :) % 作图，时变系统的b(k)输出曲线,图3-5 带遗忘因子的时变系统参数辨识结果,2022/11/26,第三章 参数估计理论与应用,从仿真结果（图3-5）可以看出，采用带遗忘因子的系统辨识算法能较好地跟踪系统时变参数的变化。3.3.5 广义最小二乘估计 上述推导中，假定模型(3.3.10)中的噪声项 e(t) 为白噪声。如果 e(t) 为有色噪声，则被测对象就应当的用CARMA 模型（受控的ARMA）来表示，即与CAR 模型（基本最小二乘估计的数学模型）比较，此处C(q-1)1，而有色噪声 e(t) 可以认为是由白噪声(
45、t) 通过成型滤波器而产生的输出（以下写成离散数据的形式）：,2022/11/26,第三章 参数估计理论与应用,下面介绍广义最小二乘法，用以改善CARMA 模型参数估计的统计特性。 将CARMA 模型写成常规的最小二乘结构，即(3.3.20) 因为(k) 是不可测量的，所以在初始估计时只能用(k)估计值代替(k)中有关(k) 的分量。按基本最小二乘法估计出，然后，计算(k) 的估计量(3.3.21)再构造,2022/11/26,第三章 参数估计理论与应用,这样，就可以按前面介绍的递推公式(3.3.22) 对CARMA 模型进行参数估计。估计的初始值可按下式选取其中，2 1。在实际应用中，(k)
46、的估计有两种形式，一种是预测的形式：,2022/11/26,第三章 参数估计理论与应用,为了避免这种问题的产生，最常采用的方法是对 PN 采用U-D分解法, 下面就介绍这种方法 。3.3.6 改进数值稳定性的U-D分解法 为了保证在递推过程中始终保持PN 的非负定性，可把矩阵PN 进行U-D分解，即【在此，将PN 写成P(N)】其中，U是对角元全为1的上三角矩阵，D为对角矩阵。这种分解只需要实时修正U，而无需修正P，因而能保证P的非负定性。下面，我们不加证明地将U-D分解基本算法的计算步骤总结如下：,2022/11/26,第三章 参数估计理论与应用,(1) 设置初值(0)，U(0)和D(0)，
47、输入初始数据； (2) 采样当前的输入和输出（假设数据的行数为p）； (3) 按下式计算 f(N) 和 g(N) (4) 定义由此可得,2022/11/26,第三章 参数估计理论与应用,按上式计算i(N)；（5）计算（6）计算参数估计（7）按下式计算di(N)式中，di是D的对角元素；（8）引入上三角阵V，它的元素为,2022/11/26,第三章 参数估计理论与应用,由此可得其中，uij是U的元素。按上式计算vii(N+1) 和 vij(N+1)；（9）按下式计算uij(N)(10) 返回（2），直到算法收敛或满足精度要求为止。 为了简化书写，在上述各个符号中，不加区分地直接写出矩阵、向量和标
48、量。请读者注意区分。下面用一个例子说明U-D分解算法的计算步骤。 例3-11 设置初值,2022/11/26,第三章 参数估计理论与应用, 采样当前的输入-输出数据对T（N+1）=1/3，1; 计算f(N+1) 和g(N+1)；令=1经计算得 计算(N+1)：,2022/11/26,第三章 参数估计理论与应用, 计算 K(N+1)： 计算估计参数 计算di(N+1)(下次迭代用)： 计算vii(N+1) 和 vij(N+1)：,2022/11/26,第三章 参数估计理论与应用, 计算 uii(N+1) 和 uij(N+1) (下次迭代用)： 返回步骤，直到算法收敛（被估计参数的两次迭代计算没有
49、显著的差异）或满足精度要求为止。 这种分解不仅能够保证P非负定性，而且，在每步的递推过程中还同时生成低阶模型的估计参数（参见。但在一般情况下，前面介绍的最小二乘法已能满足要求。在计算精度要求不高时，当检验到数据异常时，也可以人为地重新设置PN ，使之满足正定条件。,2022/11/26,第三章 参数估计理论与应用,本章小结,本章的核心是参数估计的基本理论与算法，以及这些理论方法在检测、系统辨识等专题中的应用。首先，我们介绍了参数估计量的选择标准：无偏性、有效性和一致性（或相容性）以及估计方差可能达到的CR下界。在随后的几节中，依次介绍了已知样本总体统计分布的Bayes估计、最大似然估计和样本总
50、体分布未知的最小二乘估计。 在工程实际中，必须具备用数学语言描述物理系统的能力，建立适合于参数估计的理论模型。此外，还要考虑训练数据的获取与安排，即实验设计问题。关于这方面问题很多，至今尚无系统的教科书。因此，建议读者多做访真练习，以加深对上述算法理解及其应用范围。,2022/11/26,第三章 参数估计理论与应用,习题,3-1 某一随机过程由x(n)=a+bn+v(n)描述，其中，v(n)N(0,2)；a 和b 是未知参数。求估计量a和b的估计方差的CR下界。 3-2 令 y 为一测量数据向量，其观测方程为其中C为观测矩阵，x是不可观测的状态向量，而 v 是加性观测向量。假定观测噪声向量服从
