广义最小二乘

2024-11-23

广义最小二乘(精选7篇)

广义最小二乘 篇1

0引言

在系统理论、信号处理和矩阵方程求解等许多应用领域中,往往需要给出形如:

的矩阵的某种类型的广义逆,其中表示样本数目的k通常都很大,而且随时间增长而增大,因而计算量也迅速增加,这在实时应用中是应该避免的。

本文给出了某些类型的广义逆的算法,并用所得结果,给出了求解由多个矩阵方程所构成的系统的最小二乘估计,从而推广了[1]的工作,并简化了[2]中有关求多个方程公共解的算法,文中使用广义逆的记号均引自[3]。

1广义逆的计算

设k≥1为一整数,为具有相同行数的矩阵且满足

再设的某种广义逆,我们将证明,存在着满足

的矩阵Mk和Gk,使得

为的同种广义逆。

1°,如果(·)-表{1}-逆;

2°,如果(·)-表{1,2}-逆;

1°的证明:注意到Qk-1为对称幂等矩阵且Ak=(AkTAk)(1)AkT=(AkTAk)(1)HkTQk-1[3],我们有

2°的证明:注意到对任意A(1,3)∈A{1,3},总存在着矩阵Z,使得A(1,3)=A-+(I-A-A)Z,因此AA-=AA+,于是由(1),(2)和(3)可知

为对称阵,另外由1°可知:

3°的证明:由2°,我们只需再进一步证明

引理3[5]如果Gk如下确定:

从引理1到引理3知,为获得的广义逆,必须求出

容易看出,以上两式的计算量将随着k的增大而迅速增大,为避免这一困难,注意到

或等价地,

显然,(5)~(8)的计算量与k无关.

综上所述,我们证明了以下两定理。

定理1设Mk和Zk由(1),(2)和(3)所定义且对k=0,1,…,令

Q-1=I,Ak=Qk-1Hk,Qk=Qk-1-AkGk

定理2设Mk和Zk由(1),(2)和(3)所定义且对k=0,1,…,令

下面,我们将定理1与定理2应用于

其中和Hk为具有相同列数的矩阵,且注意到如下事实[6]

我们有

定理3设由(9)所定义,且令

定理4设和Zk由(9)和(10)所定义且对k=0,1,…,令

2最小二乘估计

考虑观测系统

其中Yk∈Rm×p,Hk∈Rm×n为已知,Vk∈Rm×p为零均值随机矩阵,X为待估未知矩阵,对任何k,称X赞k为X的最小二乘估计(LSE),如果X赞k为最优化问题

的最小范数解,而

及加权矩阵wi∈Rm×n皆正定,i=0,1,…,k.

首先考虑无加权情形,即的情形。

为使估计误差的方差最小[6,7],我们有

其中的Moore-Penrose逆.

为得到(14)的递推式,注意到定理4以及(9)、(10)和(13),我们有

定理5系统(11)的X的无加权的极小范数最小二乘估计可如下递推地给出

其中Gk由定理4递推而得。

证明:

如果m=p=1,则定理5便退化为[1]的结果。定理5还有如下简单的推论:

若(11)中的Vk=0,坌k,则定理5给出多个方程所构成的系统

的公共解,亦即我们有

则(16)的公共解可按

递推地求得。

该推论中要求计算Ak+,考虑广义逆的性质和定理3,如果Gk=Qk-1Ak(1),则(17)给出(16)的公共解,从而只需计算A(1)∈A{1}而使计算量大为减少。类似地,如果Gk=Ak(1,4),则(17)给出(16)的具有极小范数的公共解。

现在再考虑一般加权最小二乘估计,因

故可通过无加权问题(18)的右端来处理加权问题。其中

或等价地,

应用定理5于无加权问题(18),便有

定理6令P-1=0,Q-1=I,且对k=0,1,…令

H′kTCk,则系统(11)的加权极小范数二乘估计可如下递推地求得:

为使(11)的加权最小二乘估计是最优估计[8],其中加权矩阵应如何选择?答案如下:

设而其方差设为

则由[4]中的推论6.4.4知,当时,为(11)在给定Yo,…,Yk时的最优最小二乘估计。

关于定理3和定理4的更多应用,可进一步参考[6,9,10]等。

3结论

应用定理1、2、3的结果,将[1]中关于单变量情形下的最小二乘估计推广到了多变量的情形,给出了后者的最优最小二乘估计的算法,同时还简化了[2]中求多个相容矩阵方程公共解的算法。

摘要:给出了某些广义逆的计算方法,将这些结果应用于观测系统的参数估计,给出了系统无加权和加权最小二乘的最新估计结果,推广了[1]的工作,同时简化了[2]中求多个相容矩阵方程公共解的算法。

关键词:算法,最小二乘估计,广义逆,方差,加权

参考文献

[1]Albert A,Sittle R.W.A method for computing leastsquares estimators that keep up with data SIAM[J].Control,1965,3:384-417.

[2]Morris G.L,Odell P.L.Common solutions for a matrixequations with applications[J].Assoc comput.math,1968(15):272-274.

[3]A.B.Israel,T.N.E.Greville.Generalized inverses:Theory and Applications[M].New York:John Wiley,1974.

[4]Campbell S.L,Meyer,C.D.Jr.Generalized Inverses of lin-ear transformations[M].London:Pitman publishing Ltd,1979:104-115.

[5]Cline R.E.Representations for the generalized inverses ofa partitioned matrix[J].Soc.Indust.Appl.Math.1964(12):588-600.

[6]X.Liu(刘轩黄).Recursive computation of generalizedinverses with applications to optimal state estimation[J].Control theory and advanced technology.1995(10):1485-1497.

[7]Rao G.R.Linear statistical inference 2nd its applications(2nd edtion)[M].New York:John Wiley&Sons,1973.

[8]刘轩黄.无初始状态统计知识时多变量系统的最优平滑与预报[J].海南大学学报,1998,16(3):210-214.

[9]X.Liu(刘轩黄).Optimal and minimum energy optimaltracking problems[J].苏州科技学院学报,2005,22(1):1-9;22(2):9-18.

[10]X.Liu(刘轩黄).Optimal state estimation without the requirement of a priori statistic information of the ini-tial state[J].IEEE.Trans on auto.control.1994,39:2087-2091.

广义最小二乘 篇2

在实验数据的处理、分析时, 数据拟合是经常采用的一种方法。数据拟合的目标是找到能反映变量之间关系的一种表达式, 使其在某种准则下最佳地接近已知数据。其原理有最小二乘法、契比雪夫法等, 且以最小二乘法最为常见和重要。

随着人类认识能力的不断进步以及计算技术的快速发展, 对于变量之间的未知关系, 应用曲线拟合的方法揭示其内在规律具有重要的理论和现实意义。针对最小二乘曲线拟合的有关问题以及相应的MATLAB实现进行探讨。

1 最小二乘曲线拟合

曲线拟合是指:已知n个数据点 (xi, yi) , i=1, 2, ……, n, 其中xi不全相同, 寻求函数f (x;a1, a2, …, am) 的待定参数a1, a2, …, am的一组取值, 使得在这组取值之下, 函数f (x;a1, a2, …, am) 与已知n个数据点整体上最为接近。

最小二乘曲线拟合方法根据已知数据, 首先构造出能够反映含有待定参数的函数f (x;a1, a2, …, am) 与n个数据点 (xi, yi) , i=1, 2……, n偏离程度的函数:

然后应用数学方法求函数J (a1, a2, …, am) 的最小值a1, am2, i…n, am=J (a1, a2, …, am) , 此时a1, a2, …, am的取值就是所求的待定值。这样一组取值使得函数f (x;a1, a2, …, am) 与n个数据点在二次平方和意义下最为接近。

2 最小二乘曲线拟合的MATLAB实现

2.1 线性最小二乘曲线拟合的MATLAB实现

线性最小二乘曲线拟合的含有待定参数的函数形式为:

其中rk (x) 为事先选定的一组已知函数, ak为待定系数, k=1, 2, …, n, m

该方程组可简化为RTRA=RTY, 其中:

当r1 (x) , …, rm (x) 线性无关时, RTR可逆, 方程组RTRA=RTY有唯一解。

根据线性最小二乘拟合理论可得关于待求参数构成的矩阵A的解, 利用MATLAB下命令A=RY, 可以直接求得待求参数, 下面举例说明。

例1下表给出了一组实测数据

已知数据 (xi, yi) 的函数原型为:

试用已知数据进行曲线拟合, 求出待定参数的值。

在MATLAB中输入以下程序:

运行程序, 输出待定参数结果为:

下面画出拟合得到的曲线及已知的数据散点图:

由图1可见, 曲线拟合效果非常好。事实上, 所给数据就是由已知曲线:

2.2 非线性最小二乘曲线拟合

在最小二乘曲线拟合中, 如果待定参数不能全部以线性形式出现, 如指数拟合函数f (x) =a0+a1e-a2·x, 这便是非线性最小二乘曲线拟合。

非线性最小二乘曲线拟合与线性最小二乘曲线拟合的原理没有什么区别, 但是最小二乘的解常常难以通过人的手算实现, 从而制约了该方法的应用。随着计算机技术的进步、专业软件的不断涌现, 这一问题的求解已不再困难。但是, 非线性曲线拟合中初值的选取是一个重要的问题, 目前为止还没有固定的理论或方法给出一般性的结论。

在MATLAB中实现非线性最小二乘曲线拟合有三个常用的命令。

(1) lsqcurvefit () 命令, 其使用格式为x=lsqcurvefit (fun, x0, xdata, ydata) , 其中fun是要拟合的非线性函数, x0是初始参数, xdata, ydata是拟合点的数据, 该函数最终返回系数矩阵。

(2) nlinfit () 命令, 其应用格式为beta=nlinfit (x, y, fun, beta0) , 其中x和y是拟合点数据, fun是被回归 (拟合) 的函数, beat0是初始参数。

(3) lsqnonlin () 命令, 其应用格式为x=lsqnonlin (fun, x0) , 其中fun的定义与lsqnonlin () 函数中fun的定义有差别, x0仍为初始参数向量, 将输出的系数结果放在变量x中。

例2假设已知

并已知该函数原型为y (x) =a1exp (a2x) +0.54exp (a3x) .sin (a4x) , 其中ai为待定系数。

在MATLAB中输入以下命令:

%建立函数原型, 则可以根据他来进行下面的求取系数的计算

运行结果为:

所求得的系数与原式中的系数相近。

如果想进一步提高精度, 则需修改最优化的选项, 函数的调用格式也将随之改变。

在MATLAB中输入以下命令:

%修改精度, 即是修改其限制条件

%两个空矩阵表示系数向量的上下限a', res运行结果为:

下面画出由拟合得到的曲线及已知的数据散点图, 如图2。

3 结论

首先介绍了最小二乘法的原理, 其次对线性最小二乘曲线拟合和非线性最小二乘曲线拟合应用MATLAB进行了具体实现, 使得相关曲线拟合理论更加生动易懂, 这使我们可以在今后的研究和工作中建立曲线模型对相应对象进行曲线拟合, 并应用MATLAB来实现, 从而找到更好更形象的反映变量之间关系的曲线, 也就找到最合适的模型了。

参考文献

[1]薛定宇, 陈阳泉.高等应用数学问题的MATLAB求解[M].北京:清华大学出版社, 2004.

[2]刘琼荪, 何中市, 傅鹂等.数学实验[M].北京:高等教育出版社, 2000.

[3]梁国业, 廖健平.数学建模[M].北京:冶金工业出版社, 2004.

广义最小二乘 篇3

如何确定视频的拍摄地点和拍摄日期是视频数据分析的重要方面。其中一种简单的情形是:根据相关影长数据确定参照物位置, 这一技术也称为太阳影子定位技术。在本文中, 根据2015 年全国大学生数学建模竞赛A题附件1 的数据[1], 确定出该地点的位置。

2 影长与时间的函数关系

2.1 影响影长的各个参数

太阳高度角是指太阳光的入射方向和地平面之间的夹角, 可以用太阳高度角来估计影子场电影。 太阳高度角h与实测地的地理纬度θ, 太阳直射角纬度 φ, 时角 Ω 满足以下函数关系式[2]:

其中, 关于各个位置参数确定如下:

太阳直射点每时都在向西移动, 每小时移过15 度经度, 在地理题的计算中可粗略取每天移动0.25 度纬度[3]。 据此求解出太阳直射角纬度与日期之间的关系, 进而通过MATLAB编程, 得到了每天太阳直射点的纬度, 部分表见表1, 在4 月18 日的太阳直射点纬度为 φ=7。

一个天体的时角表示该天体是否通过了当地的子午圈, 其数值则表示了该天体与当地子午圈的角距离, 以小时来计量 (1 时角=15 度) 。根据以上理论, 建立时角的关系表达式[4], 即:

其中Cr表示北京时;LC表示经度校正 (4min/度) , 经度大于120 的经度校正为正, 小于120 的经度校正为负;EQ为时差。

2.2 建立影长随时间变化函数

根据式 (1) , 我们得到了太阳高度角与时角、太阳直射点纬度以及实测地的地理纬度的函数关系式, 进而我们可以通过太阳高度角与影子长度之间的关系, 建立影长随时间变化的函数关系式, 其步骤如下:

根据太阳高度角的定义, 得到太阳高度角与影子长度之间满足以下关系:

其中L为影子长度;H为旗杆高度3m;h为太阳高度角;

根据式 (1) , 将各个参数值代入公式进行整理, 得到太阳高度角与各个参数之间的函数关系:

再根据式 (4) , 将 (5) 转化为:

其中:

3 利用最小二乘拟合确定直杆高度以及所处位置

3.1 模型的建立

已知对应时刻的影子长度 (ti, Li) , 影长与时间的函数关系, 见公式 (6) , 利用最小二乘原理, 建立如下模型:

其中c1, c2, c3的初始值设为1, 0.1, 10

3.2 参数及经度的确定

根据参数的初始值, 将21 组 (ti, Li) 代入拟合方程中进行拟合。 在拟合之前, 本文根据c1, c2, c3的表达式做了以下猜想:发现c1, c2关于其自身的自变量的表达式十分简单, 而c3表达式较为复杂, 并且c3在拟合方程中的表达式也较为复杂, 所以关于c3的拟合效果相对于c1, c2来说, 拟合效果较差, 可能得不到精确的解。

而在实际拟合过程也证明了我们的猜想。 在实际拟合过程中, 随着初始值的不断变化, c1, c2变化不大或者只在固定的值的范围内小幅度变化, 而c3的变化幅度较大, 且不能找到唯一确定的值, 因此通过求出的纬度以及直杆高度来确定经度。经求解:c1=2.1316m, c2=0.313或0.1287, 即H=2.1316m, φ=18.24N或7.3969N。

将H, φ 中的数据代入 (5) 中, 用MATLAB编程计算, 我们得到21个时间点所对应的经度位置, 我们发现, 这21 个经度变化都十分小, 故我们经度该位置的经度就选取求出的21 个经度的平均值, 得到4 个经度, 进而得到四个固定直杆可能所处的地理位置: (18.24N, 105.9528E) , (18.24N, 37.4972E) , (7.3969N, 107.1988E) , (7.3969N, 36.2512E) 。

4 地理位置的检验

根据测定的关于x, y的坐标, 我们发现随着时间的增加, 影子长度也在不断增加, 说明当地时间已经过了正午时间12:00, 而此时的北京时间为14:42-15:42, 结合北京时间以东八区120 度为标准, 每向西15 度, 当地时间往前1 小时的规律, 我们计算出在79.5 E时, 北京时间为14:52, 当地时间恰好为12:00, 因此我们所确定的经度必须大于79.5E。 进而确定可能的地理位置为: (18.24N, 105.9528E) 和 (7.3969N, 107.1988 E) 。

根据地图确定地理位置为 (18.24N, 105.9528E) 的地点为越南, 地理位置为 (7.3969N, 107.1988E) 的地点为泰国湾 (为靠近马来西亚的一个海峡) , 故这个位置不符合我们的标准, 故我们确定的直杆所处位置大致在越南。

5 模型的检验

我们将所确定的地理位置代入关于太阳高度角的式子 (5) 中, 画出14:42-15:42的影子长度的变化曲线l1, 同时拟合出附件1中影子长度随时间变化曲线l2, 两者之间的误差, 可以这21个时刻的误差分别为:0.0004, 0.0005, 0.0005, 0.0006, 0.0006, 0.0006, 0.0006, 0.0006, 0.0006, 0.0006, 0.0005, 0.0005, 0.0005, 0.0005, 0.0004, 0.0004, 0.0004, 0.0004, 0.0004, 0.0004, 0.0004。两条影子长度随时间变化的曲线 (见图1) 。

根据图1, 结合21 个点的误差, 可以发现两条曲线误差最大的为0.0006, 误差非常小, 认为该模型合理。

参考文献

[1]http://www.shumo.com/2015cumcm.html[OL].

[2]http://baike.baidu.com/view/86609.htm[OL].

[3]http://baike.baidu.com/link?url=_a Uu Zt IGSo3DBZgj Wi SGC2Nq H9NNBj Ec XX37JYu5MWs F6M22qio Ndju Y81e Nje P31Llpv JKTBhf_e PWF9g5CXa#1[OL].

[4]王国安, 米鸿涛, 邓天宏, 李亚男, 李兰霞.太阳高度角和日出日落时刻太阳方位角一年变化范围的计算[J].气象与环境科学, 2007, 30:161-163.

广义最小二乘 篇4

卫星导航定位在现代军事和经济领域扮演着不可或缺的重要角色, 卫星导航系统除了提供导航功能外, 还必须具有在系统不能使用时及时向用户发出告警的能力, 这种能力叫做系统的完好性。完好性和定位精度是卫星导航中两个重要指标。当卫星导航用于航空、武器系统等生命相关场合时, 完好性成为决定性指标。因此, 卫星导航系统的完好性问题已成为卫星导航领域重点解决的问题, 并将存在于卫星导航整个生命周期。

引起接收机伪距观测量偏差的主要因素有较大的卫星钟漂、导航电文的不正确上传、各种欺骗与干扰以及卫星组成部分的故障等, 这里统称为卫星故障, 在GPS系统中这种故障发生的概率大约为10~4/小时。对这种故障, 导航系统本身做出告警反应的时间较长, 如GPS系统的反应时间为15分钟到数小时, 对于生命相关场合显然太慢, 不能满足需求。因此, 卫星故障的快速检测只有在用户端进行, 因而出现了接收机自主完好性检测技术 (Receiver Autonomous Integrity Monitoring, 英文简写RAIM) 。本文针对导航接收机常用的位置解算方法, 即最小二乘法, 设计了基于最小二乘残差法的RAIM检测方案, 并对方案进行了仿真分析。

1 基于最小二乘残差法的RAIM检测方案

导航接收机位置解算常用最小二乘法, 针对最小二乘位置解算方法, 设计了相应的RAIM检测方案, 即基于最小二乘残差法的RAIM检测方案。基于最小二乘残差法的RAIM方案除了对当前可视卫星的几何结构有很强的依赖性, 还对可视卫星数目提出了很高要求。如果接收机要完成故障检测, 最少需要5颗可视卫星;如果接收机要完成故障排除, 则最少需要6颗卫星。

1.1 最小二乘位置解算法[5,6]

对于伪距定位方程:

undefined. (1)

在X0= (x0, y0, z0, B0) ′ (B=cΔt, 看作一个未知数, 以米为单位) 处进行泰勒级数展开[4,5], 略去二次项及高次项得到:

undefined

令undefined

则:undefined. (3)

undefined. (4)

undefined. (5)

将式 (3) (4) (5) 代入式 (2) , 方程可整理为:

ρi-ρ0=αi (x-x0) +βi (y-y0) +γi (z-z0) + (B-B0)

=αiΔx+βiΔy+γiΔz+ΔB. (6)

测量四颗卫星的伪距, 可获得四个类似式 (6) 的方程, 联立这四个方程得到一个四元线性方程组:

undefined

ΔX=ΔX ΔY ΔZ ΔBT.

Y=T.

则式 (7) 变为: GΔX=Y. (8)

解此方程得: ΔX=G-1Y. (9)

若测量伪距的卫星数n大于4, 则G为n×4的矩阵, 解此方程得:

ΔX= (GTG) -1GTY. (10)

所以 X=X0+ΔX. (11)

再以新的X作为X0, 继续上述解算, 迭代数次后, 当满足设定的误差时, 迭代结束, 得到所求的解, 实践经验表明, 迭代次数为5左右。

1.2 故障检测与识别

基于最小二乘位置解算原理, 利用最小二乘迭代解算过程的中间结果, 就可以完成RAIM故障检测与识别。引入RAIM故障检测与识别的最小二乘位置解算过程如图1所示, 其中完好性处理部分的详细流程如图2所示。

完好性处理部分首先对观测量Y进行故障检测后, 观测量Y才能进入最小二乘迭代运算过程。完好性处理主要完成可视卫星数目的判断、水平定位误差保护级HPL、检验统计量及其限值的计算, 并对故障进行检测、识别等。具体过程如下:在RAIM检测之前, 有两个重要的限值需要计算, 它们是水平定位误差保护限值HAL和故障检测门限T。水平定位误差保护限值和故障检测门限可通过文献中给出的算法计算得到, 或者可通过文献提出的基于最小风险准则来确定。

在RAIM检测开始时, 首先根据获取的观测信息, 判断可视卫星数目是否达到故障检测与识别要求, 判断可视卫星数目的流程如图3所示。判断可视卫星的主要依据是计算过境卫星的仰角, 如果仰角大于5度, 认为可视。理论上卫星导航系统中最多可视卫星为12颗, 所以判断可视卫星时最多进行12次循环判断。

如果前面判断的可视卫星数目符合要求, 对水平定位误差保护级HPL进行计算, 并与HAL比较, 看是否超限, 以保证RAIM检测的有效性。当HPL小于HAL时, 对检验统计量进行计算。如果检验统计量小于检测门限值T, 则认为无故障, 否则, 进入故障识别与隔离阶段。如果可视卫星数目只有5颗, 那么只能完成故障的检测, 并对检测到的故障向用户发出告警;如果可视卫星数目大于等于6颗, 则可以对检测到的故障进行隔离。最后, 输出完好性处理结果。

2 方案仿真

为了验证基于最小二乘残差法的RAIM方案对故障检测与识别的性能, 利用STK5.0软件产生星座数据, 采用MATLAB软件对该方案进行仿真。由于GPS系统已有较为准确的误差统计数据, 因此仿真采用的卫星星座是GPS卫星星座。利用了在GPS系统时间1 Jan 2010 13:15:00.00至1 Jan 2010 13:17:00.00之间, 在某地区所有的可视卫星, 它们分别是GPS星座系统第1、2、3、6、9、11、15、20颗卫星。仿真条件如下:故障检测误警率为0.005, 由此计算得到故障检测门限值为10.343 m。设观测噪声均值为0, 标准差为5 m。在观测时段内第30秒时刻, 对其中一颗卫星对应的伪距加入一个100 m突变, 模拟此时卫星出现故障, 相应的伪距出现异常。

在上面仿真条件下, 仿真结果如图4所示。横轴为观测时间, 单位为秒。左边纵轴表示检验统计量, 单位为米, 右边纵轴为观测伪距残差, 单位为米。从仿真结果图上可以看出, 在30秒时刻观测伪距发生异常时, 相应的伪距残差出现突变, 从而导致了检验统计量的突变。图中所示的检验统计量在30秒时刻已经超过了故障检测门限值, 因此, 此时的故障完全能够被检测出来。仿真结果表明, 本文所设计的基于最小二乘残差法的RAIM故障检测方案是可行的, 有效的。

3 结束语

完好的系统故障检测性能、快速的故障反馈能力是现代卫星导航定位系统接收机自主完好性检测的发展趋势。本文提出的基于接收机自主完好性检测方案具有告警时问短、运算量小的特点。接收机位置解算如果采用最小二乘法, 就可以很容易地将该方案嵌入其中。我国二代卫星导航定位系统正在建设当中, 建议在系统建设中考虑好系统的完好性检测扩展能力, 在接收机的设计中加入自主完好性检测功能。

摘要:为了使导航接收机的自主完好性检测算法与最小二乘位置解算法结合起来, 提供准确可靠的定位定时服务, 对最小二乘位置解算法和基于最小二乘残差法的接收机自主完好性检测算法进行了研究。针对导航接收机最小二乘位置解算法, 提出了基于最小二乘残差法的卫星导航接收机自主完好性检测方案, 并对方案进行了仿真分析, 仿真结果表明方案是可行的, 有效的。

关键词:导航接收机,自主完好性,最小二乘法

参考文献

[1]Brown R G.Global Positioning System Theory and Appli-cations[M].American:American Institute of Aeronauticsand Astronantics, 1996:143-164.

[2]Xie G, Pullen S, Luo M, et al.Integrity Design and Up-dated Test Results for the Stanford LAAS Integrity MonitorTestbed[C].Albuquerque, NM:Proceedings of ION 57thAnnual Meeting and CIGTF 20th Biennial Guidance TestSymposium, 2001:681-693.

[3]冯庆堂, 沈林成, 常文森.一种监测GPS完好性的新方法[J].导航与控制, 2004, 3 (1) :4-5.

[4]Bradford W, Parkinson, Penina Axelrad.AutonomousGPS Integrity Monitoring Using the Pseudorange Residual[J].Journal of the Institute of Navigation, 1988, 35 (2) :456-460.

[5]王惠南.GPS导航原理与应用[M].北京:科学出版社, 2003.

[6]刘基余.GPS卫星导航定位原理与方法[M].北京:科学出版社, 2003.

[7]黄晓瑞, 田巍, 李波.GPS接收机的自主完善性监测算法研究[J].遥测遥控, 2003, 24 (1) :1-3.

广义最小二乘 篇5

考虑下列二阶抛物问题

{ut(x,t)-(a(x)u(x,t))=f(x,t),(x,t)Ω×J,u(x,t)=0,(x,t)Γ×J,u(x,0)=u0(x),xΩ(1)

(1)式中J=[0,T],ΩRn为有界凸域,边界Γ是Lipschitz连续的。假定a(x)有界,即存在常数α*和α*,使得

0≤α*≤a(x)≤α*,∀xΩ (2)

对上述问题人们提出了许多数值模拟格式,如有限元格式[1],混合有限元格式[2],扩展混合有限元格式[3],最小二乘混合有限元格式[4]等。现将最小二乘思想和扩展混合有限元方法相结合,提出并分析了最小二乘扩展混合有限元格式。该格式继承了最小二乘和扩展混合元方法的优点,即格式能同时高精度逼近未知函数,未知函数的梯度以及流体的通量,所形成的有限元方程组是对称正定的,有限元空间无须满足LBB相容性条件。数值算例说明了该格式的有效性。

1 最小二乘扩展混合元格式

定义空间Λ=(L2(Ω))n,V=H01(Ω),W=H(div;Ω)={q∈(L2(Ω))n,∇·qL2(Ω)},相应范数为μ0,Ω=(μ,μ)12,v1,Ω=(v0,Ω2+v0,Ω2)12,qΗ(div;Ω)=(q0,Ω2+q0,Ω2)12,其中(·,·)表示空间(L2(Ω))n中的内积。作时间剖分τ=T/N,N为正的常数,在时间tn=,n=1,2,…,N处逼近方程的解。Ω=ΚτhΚΩ的一个剖分,h=max{diam(K):K∈Th}。设Vh,ΛhWh分别是V,ΛW的有限维子空间,具有通常逼近性质。Vh,ΛhWh的一个标准的选法是k,rs次的分片多项式空间。

λ=-∇u,σ=,并用一阶差分近似ut,则式(1)可变成离散形式I

{¯tun=σn=fn+R1n,xΩ,λn+un=0,xΩ,aλn-σn=0,xΩ,un(x)=0,xΩ,u0(x)=u0(x)xΩ(3)

(3)式中¯tun=(un-un-1)/τ,R1n=¯tun-utn

对于(u,λ,σ)∈V×Λ×W和(v,u,q)∈V×Λ×W,定义双线性形式

aΝ(u,λ,σ;υ,μ,q)=(u+τσ,υ+τq)+τ(λ+u,μ+υ)+τ(a-1(aλ-σ),aμ-q)(4)

根据离散形式式(3),省略时间截断误差定义最小二乘扩展混合元格式

格式Ⅰ 给定初始值uh0Vh,对于n=1,2,…,N,求(u,λ,σ)∈Vh×Λh×Wh,使得

aΝ(unn,λnn,σnn;υh,μh,qh)=(uhn-1+τfn,υh+τqh),(υh,μh,qh)Vh×Λh×Wh(5)

2 解的存在唯一性

为证明格式Ⅰ解的唯一性,首先说明双线性形式aN(·;·)是强制的。

引理1 对任意的(v,μq)∈V×A×W存在c>0,使得aN(v,μ,q;v,μq)≥c(‖v1,n2+‖μ0,n2+‖qH(div n))。 (6)

证明 引入适当选取的常数β>0,式(4)可整理为

aΝ(υ,μ,q;υ,μ,q)=(τq+(1-β)υ,τq+(1-β)υ)+τ(q-a(μ+βυ),a-1(q-a(μ+βυ)))+2τ(μ+(1-aβ)υ,μ+(1-aβ)υ)+(2β-β2)(υ,υ)(7)

由于υ(x)=0,xΓ,所以由Poincaré-Friedrichs不等式知,存在正常数CF,使得

v0,n2CF‖∇v0,n2 (8)

结合式(2)有

aΝ(υ,μ,q;υ,μ,q)β(2(τα*+CF)-β(τα*2+τα*+CF))υ0,Ω2(9)

β=τα*/(τα*2+τα*+CF),则

aΝ(υ,μ,q;υ,μ,q)Cυ0,Ω2Cυ1,Ω2(10)

从而由式(10)和式(4),知

μ0,Ω22μ+υ0,Ω2+2υ0,Ω2CaΝ(υ,μ,q;υ,μ,q),(11)q0,Ω22aμ-q0,Ω2+2aμ0,Ω2CaΝ(υ,μ,q;υ,μ,q),(12)q0,Ω22τυ+τq0,Ω2+2τυ0,Ω2CaΝ(υ,μ,q;υ,μ,q),(13)

结合式(10)—式(13)便证得式(6)成立。

定理1 设fL2(Ω),则格式Ⅰ存在唯一解(μhn,λhn,σhn)∈Vh×Λh×Wh

证明 考虑空间V×Λ×W,其范数定义为‖υ‖1,Ω+‖μ‖0,Ω+‖QH(div;Ω)。由引理1知双线民生形式αN(·;·)在空间Vh×Λh×Wh中是强制的;而且由Cauchy不等式易证αN(·;·)在Vh×Λh×Wh中是有界的。因此由Lax-Milgram引理知格式Ⅰ解存在唯一。

3 数值算列

为验证所提格式的有效性,在区域[0,1]×[0,1]上考虑方程

{ut-(aux)x=f,(x,t)(0,1)×(0,1],u(0,t)=u(1,t)=0,t[0,1],u(x,0)=sinπx,x[0,1](14)

(14)式中a=1100,f=(1+π2100)etsinπx

容易验证方程的真解为u=etsinπx,λ=-πetcosπx,σ=-1100πetcosπx。取Λh为分片常数多项式空间,VhWh为分片线性多项式空间,用MATLAB编程计算近似解uh,λhσh,绘出真解和最小二乘扩展混合元解的比较图,见图1—图3。

摘要:对二阶抛物问题提出了最小二乘扩展混合元数值模拟格式。数值算例验证了格式的有效性。

关键词:二阶抛物问题,最小二乘,扩展混合有限元方法

参考文献

[1]Thome V.Galerkin finite element methods for parabolic problems.Lecture Notes in Mathematics,Springer,1984

[2]Johnson C V.Thome,Error estimates for some mixed finite element methods for parabolic type problem.RAIRO Numer Anal,1981;15:41—78

[3] Guo Ling,Chen Huanzhen.An expanded characteristic mixed finite element method for a convection-dominated transport problem.Journal of Coputational Mthematics,2005;23:479—490

广义最小二乘 篇6

支持向量机

支持向量机(Support Vectormachine,SVM)作为一类新型机器学习方法,由Vapnik等人提出的是,这种方法对小样本、非线性及高维等模式识别问题有更好的解决办法。该方法具有良好的泛化能力,因而在模式识别中得到了广泛应用。

1最优分类面

SVM算法是假设有一个两类样本的分类问题,给定训练样本{xi,yi},{xi,yi},i=,1,2∧n,x∈Rn,yi∈{-,1+1},存在一个超平面可以将它划分。支持向量机的基本思路是寻找一个最优超平面,使它的分类间隙最大。

2广义最优分类面

如果一个超平面不能把两类点彻底分开时,可引入松弛变量ξi(ξi≥,0i=,1∧,n),,使超平面wTx+b=0满足:

当0<ξi<1时,可以对样本点进行正确分类。针对包含噪音的数据,导致训练集达到零误差,会影响模型发生过拟合和较差的泛化能力。根据实际情况分析,为解决该问题并允许存在一些样本的错分,目标函数变为:

i1式中c>0为规定常数,为样本的错分上界,也称为损失函数。

在式(2)约束下,求式(1)目标函数的极小值,在线性不可分条件下可以得到最优超平面。线性不可分情况下最优超平面的对偶问题与线性可分情况大致相同,只在条件上有所变化。

3支持向量机的核函数

目前,主要有以下四种常用的核函数形式。

(1)线性核函数

该核函数没有待定参数。

(2)多项式核函数(Polynomial)

d为其待定系数。

(3)径向基形式核函数

其优越性体现在将原有空间中的非线性问题转换为其他特征空间中的线性问题,而且实际涉及的计算又全部在原空间进行。

最小二乘支持向量机

介于SVM在其训练时总能找到全局最优解和良好泛化能力和,其广泛应用在经验建模领域。但其有约束的二次规划问题,导致了训练时间较长,从而不被接受。为提高训练效率,1999年J.A.KSuyken等人提出一种新的最小二乘支持向量机(简称LSSVM)。该方法的训练只需求解一个线性方程组,使SVM容易实现,大幅度提高SVM的训练效率,因此在模式识别领域得到广泛应用。

1构造多类分类最小二乘支持向量机

经典支持向量机和最小二乘支持向量机的基本思想相同,但后者为每个数据点加入一个改正量ei,即把求最优分类超平面问题就转化为求解凸优化问题。

本文把车型分为轿车、货车和客车三类,因此首先要构造多类分类的LS-SVM。

当LS-SVM在多类分类问题应用时,假设给出的c-分类问题训练样本为{yim,xi}i=,1,2∧,n;m=,1∧,c,n为训练样本数,yim为第i个样本,其识别可以表示为求解下面的问题

约束条件:

定义拉格朗日函数,求解得:就可以得到多元分类最小二乘支持向量机的决策函数为:

2分类方法的选取

本文采用的是基于LS-SVM的多类识别。在识别时,常用的有三种方法,本文对这三种方法进行了研究和对比。

(1)一对多方法(OAA法)

该方法基本思想是如果数据有类,则需要构造个k分类器,将第i个分类的样本数据记为正k类,不属于类别i的样本数据记为负类。测试中,对测试数据计算各子分类器的判别函数值,将最大判别函数值所对应的类别做为测试数据的类别。该算法缺点是每个分类函数都要所有的样本参与,训练时间与类别数量成正比,扩展能力差,存在有不可区分的区域。

(2)一对一方法(OAO法)

该方法基本思想是选取2个不同类别构成一个子分类器,总共有k(k-)1/2个子分类器。在构造类别i和类别j的子分类器时,样本数据集选取属于类别i、类别j的样本数据作为训练样本数据,并将属于类别i的数据标记为正,将属于类别i的数据标记为负。该算法可导致过学习问题,而且训练时间随和类别数呈超线性关系。

(3)多对多方法(ECOC法)

该方法的基本思想是按规则构造n个两分类器,要求对每个类的样本两分类器都判断正确,从而得到一个k行n列的编码矩阵。

实验结果及分析

1实验环境

本文采用海康4004HC/B视频监控,Windows XP SP3操作系统,Matlab 2010软件。

2实验步骤

本文将从监控系统中得到的视频图像作为实验数据,采集客车30辆、货车29辆、轿车35辆。

(1)视频图像的背景差分及运算,消除图像中的孤立噪声,确定识别车辆的边界,将识别车辆目标从图像序列中分割出来,得到车辆几何特征。

(2)在提取车辆尺寸和长宽比的基础上引入伪Zernike矩车型特征。

(3)确定LS-SVM的惩罚因子c和径向基核函数参数σ2等参数,进行测试,确定最优的LSSVM分类模型。

3车型分类结果及分析

由实验结果得出,ECOC方法具有较高的分类精度,完成分类的时间最短。

结论

广义最小二乘 篇7

总体最小二乘[1]认为平差模型中观测向量和系数矩阵同时含有误差, 并对其进行改正。针对系数矩阵含有误差的平差问题, 在理论上总体最小二乘较最小二乘更为严密。总体最小二乘的经典解法是奇异值分解法 (SVD) , 但由于其计算复杂且不利于编程, 测量学者提出了一些基于平差模型的迭代解法[2~7], 并对其实际应用作了一些研究[8~10]。在测量数据处理中, 线性拟合是总体最小二乘一个最重要的应用, 但总体最小二乘的常规解法 (即奇异值分解法和迭代解法) 都无法顾及线性回归中系数矩阵的常数列而是对其进行了改正, 这对线性回归模型来讲是不合理的。文献[5]以直线必须通过数据点的中心, 才能使其偏差最小为约束条件, 推导了一种求解线性模型的总体最小二乘算法, 虽然该算法提高了线性回归的逼近精度, 但仍然存在偏差。文献[7]通过将数据中心化后再采用奇异值分解法对回归参数进行解算, 则分离了系数矩阵的常数列, 得到的参数估值合理可靠。但其缺陷是计算复杂难以编程实现。鉴于此, 本文在基于平差模型的基础上, 推导了一种求解线性回归参数的总体最小二乘算法, 该算法易于编程实现。算例分析表明其可靠、合理。

1 总体最小二乘常规解法

1.1 SVD解法

线性回归模型为[1,2,6]:

式中:Δ观测量y的真误差, ak为回归参数, 当同时考虑到自变量xk也含有误差, 此时系数矩阵A中的元素含有误差, 则式 (2) 变换为总体最小二乘的平差模型:

采用奇异值分解 (SVD) 法将式 (2) 变换为:

求得参数估值为[1,2,6]:

式 (4) 中, n为非固定列对应参数的个数, σ2n+1为增广矩阵[A L]的最小特征值。其单位权中误差计算公式为σ02=σ2n+1/ (m-n-1) , 其中m为观测个数, n+1为参数的个数。

SVD解是将系数矩阵中所有元素看成是含有误差的, 这对线性回归模型来讲是不合理的。而且SVD解法较为复杂, 不利于测量数据处理的编程实现。

1.2 常规迭代解法

总体最小二乘的迭代解法为[2]:

进行迭代解算时, 首先按最小二乘法计算参数的初值X (0) , 然后根据式 (5) 反复迭代直到‖X (i+1) -X (i) ‖<ε则停止计算输出参数值, 单位权中误差按σ0=V/ (m-n) 计算。迭代解法最大的优点是充分运用到测量平差模型, 其迭代格式简单易于编程。但常规的总体最小二乘迭代解法都是将系数矩阵所有元素看成是含误差的, 这对线性回归模型也是不合理的。

2 线性回归参数的总体最小二乘迭代解法

2.1 算法推导

根据前文所述, 有必要建立基于总体最小二乘的线性回归模型和解法, 线性回归数学模型为:

可将式 (6) 进行等价转换为:

表示为矩阵形式为:

设v=vec (EA) , vec表示向量拉直运算, vec (EA) 即将矩阵从左至右逐列拉直成一列向量。在不考虑权重时, 根据总体最小二乘原理, 相应的误差期望和方差为:

式 (9) 中, 为矩阵的克罗内克积, Im和In分别为m和n阶的单位矩阵, 相应的其平差准则为:

根据式 (8) 与式 (10) 构造拉格朗日目标函数:

式 (11) 中, k为 (n+1) ×1的拉格朗日乘数;, In+1为n+1阶的单位矩阵。根据拉格朗日函数求极值的必要条件, 将式 (11) 分别对v、x求导并令其等于0, 化简整理后可得:

由式 (12) 的第一式可得:

式 (13) 中, vec-1表示vec的逆运算, 即将矩阵重新构造为m× (n+1) 的矩阵。

将式 (12) 的第一式代入到式 (8) 中, 并顾及, 整理后可得:

根据式 (11) 第二式和式 (14) 则可得:

根据式 (15) 即可得到参数^x的表达式:

根据上述推导, 其参数解算步骤如下:

(1) 对回归参数附初值^x (0)

(2) 按式 (17) 计算k值和新的回归参数值

(3) 重复步骤 (2) , 直到‖x (i+1) -x (i) ‖<ε, 则停止迭代

按上述迭代方法即可求得式 (7) 的方程, 相应的线性回归方程的参数解:

2.2 单位权中误差评定

根据线性回归模型式 (7) , 当同时考虑自变量误差即采用总体最小二乘求取回归参数值时, 其单位权中误差的计算公式为:

根据上文的迭代解再结合式 (12) 可得:

则按本文的推导, 其单位权中误差计算公式为:

在计算线性回归单位权中误差时, 按常规的总体最小二乘解法算得的结果比按式 (19) 计算的结果要小, 这与实践是不不相符的。因为, 常规的总体最小二乘解法没有顾及线性回归模型中系数矩阵常数列。而本文推导的解法从理论上比常规总体最小二乘解法更严密, 究其原因, 是因为按常规总体最小二乘解法对线性回归模型的单位权中误差评定有误, 下面则具体进行讨论。

按常规总体最小二乘法将线性回归系数矩阵中的常数列也进行了改正并参与精度评定, 其改正后的模型可表示为:

其单位权中误差评定公式为:

从式 (22) 可以看出, 其将常数列的改正值看成是自变量的改正值, 这显然是不对的。从改正模型来看, 应该将常数列的改正值归入到因变量的改正值中, 则变量的改正值为:

按式 (24) 进行改正值修正后, 再按式 (19) 进行单位权中误差的评定。如此得到的结果才与线性回归模型相符。这也解释了前文提出的疑问, 同时也表明不能单一以单位权中误差来评定一种平差模型和算法的优劣。因为它们的单位权中误差评定方式不一定都相同。

3 算例分析

3.1 算例1

运用文献[6]例5~10的数据, 用三个点 (1, 2) 、 (2, 6) 、 (6, 1) 来拟合一元线性回归方程。运用本文的解法得到的参数估值为6、1。这与文献[5]中将数据中心化后再采用奇异值分解法得到的结果相同, 而按常规总体最小二乘的奇异值分解法和迭代解法得到的参数估值皆为6.74、0.99。这是因为总体最小二乘的常规解法将系数矩阵的常数列也进行了改正, 得到的结果有偏差。而本文的解法则将系数矩阵的常数列分离开不加改正, 计算得到的参数估值较之可靠。

3.2 算例2

为进一步验证本文方法对线性回归模型的统一性和可靠性, 运用MATLAB模拟一平面方程。设平面方程为:z=1.5+x+2y, 在x和y没有误差时求得z的值上加上均值为0, 方差为0.005的随机误差, 组成观测值。然后分别对x和y添加均值为0, 方差为0.03的随机误差, 组成新的观测值, 如表1所示。分别采用最小二乘法 (LS) 、常规总最小二乘法、本文解法解算线性回归方程的参数值, 并按文中所述的单位权中误差评定公式计算其单位权中误差, 结果如表2所示。

从表2中的参数估值结果不难看出, 采用本文的总体最小二乘解法解算的回归参数值与真值最为接近, 而常规的总体最小二乘解法由于将线性回归模型中的系数矩阵中的常数列也进行了改正导致其结果与真实值有偏差, 最小二乘法由于没有考虑到线性回归自变量的误差即系数矩阵元素的误差, 计算的结果与真值偏差最大。这也可以从其单位权中误差中得出, 本文的总体最小二乘解法所得单位权中误差最小, 常规总体最小二乘解法次之, 最小二乘法最小。需要指出的是, 如果按常规总体最小二乘的单位权中误差计算公式得到的单位权中误差为0.0227, 比本文方法得到的值要小, 这就会误以为按常规总体最小二乘解法解算线性回归参数得到的精度比本文要高。其实不然, 这是因为其针对线性回归模型的单位权中误差计算公式有误, 通过本文的纠正得出的值才符合其实际。这也可以通过参数的估值结果与真值的比较中看出。

4 结束语

本文给出了一种线性回归参数估计的总体最小二乘算法, 该算法即能同时考虑线性回归自变量的误差即平差模型系数矩阵元素的误差, 又能顾及线性回归中系数矩阵常数列即对常数项不加改正。这弥补了常规总体最小二乘解法针对线性回归模型解算的不足, 且对基于总体最小二乘的线性回归单位权中误差进行了探讨, 纠正了针对线性回归模型的常规总体最小二乘单位权中误差评定的偏颇。

摘要:根据线性回归模型的特点, 推导了一种解线性回归参数的总体最小二乘算法。并对其单位权中误差的评定进行了探讨, 通过理论推导表明按常规的总体最小二乘解法求得的线性回归单位权中误差与实际不符, 并对其进行了纠正。通过算例分析且与常规的总体最小二乘解法进行比较, 结果表明了算法的正确性和可靠性。

关键词:总体最小二乘,线性回归,平差模型,迭代解法,奇异值分解

参考文献

[1]GOLUB G H, VAN L C F.An Analysis of the Total Least Squares Problem[J].SIAM J Numer.Anal, 1980, 17:883~893.

[2]鲁铁定, 周世健.总体最小二乘的迭代解法[J].武汉大学学报 (信息科学版) , 2010, 35 (11) :1351~1354.

[3]孔建, 姚宜斌, 吴寒.整体最小二乘的迭代解法[J].武汉大学学报 (信息科学版) , 2010, 35 (6) :711~714.

[4]许超钤, 姚宜斌, 张豹等.基于整体最小二乘的参数估计新方法及精度评定[J].测绘通报, 2011, (10) :1~4.

[5]邱卫宁, 齐公玉, 田丰瑞.整体最小二乘求解线性模型的改进算法[J].武汉大学学报 (信息科学版) , 2010, 35 (6) 708~710.

[6]丁士俊, 姜卫平, 杨颜梅.整体最小二乘线性回归模型与算法[J].测绘通报, 2012, (12) :8~10.

[7]邱卫宁, 陶本藻, 姚宜斌等.测量数据处理理论与方法[M].武汉:武汉大学出版社, 2008.

[8]陈义, 陆珏, 郑波.总体最小二乘方法在空间后方交会中的应用[J].武汉大学学报 (信息科学版) , 2008, 33 (12) :1271~1274.

[9]陆珏, 陈义, 郑波.总体最小二乘方法在三维坐标转换中的应用[J].大地测量与地球动力学, 2008, 28 (5) :77~81.

上一篇:铜绿假单胞菌治疗下一篇:音乐载体