数值建模课件-拟合.ppt

1、1 这样，最小二乘问题就转化为求多元函数),(10naaaIminjiijjixfxax002)()()(（4）的极小点 问题.),(*1*0naaa 用最小二乘法求拟合曲线的问题，就是在形如(2)的 中求一函数 ，)(xS)(*xSy 由求多元函数极值的必要条件，有 minjikiijjikxxfxaxaI000)()()()(2).,1,0(nk使（3）取得最小.2曲线拟合的最小二乘法曲线拟合的最小二乘法 1 最小二乘法及其计算最小二乘法及其计算 在函数的最佳平方逼近中 如果 只在一组离散点集 上给定，这就是科学实验中经常见到的实验数据 的曲线拟合.,)(baCxf)(xf,1,0,mix

2、i,1,0),(miyxii3 记误差,1,0,)(*miyxSiii则 的各分量分别为 个数据点上的误差.T10),(mm 问题为利用 求出一个函数,1,0),(mixfyii)(*xSy 与所给数据 拟合.,1,0),(miyxii4 设 是 上线性无关函数族，)(,),(),(10 xxxn,baC在 中找一函数 ，)(,),(),(10 xxxspann)(*xS使误差平方和miiimiiyxS02*0222)(,)(min02)(miiixSyxS（1）这里)()()()(1100 xaxaxaxSnn).(mn（2）5 这个问题称为最小二乘逼近，几何上称为曲线拟合的最小二乘法.用最

3、小二乘求拟合曲线时,首先要确定 的形式.)(xS 确定 的形式问题不仅是数学问题,还与问题的实际背景有关.)(xS 通常要用问题的运动规律及给定的数据进行数据描图,确定 的形式,然后通过实际计算选出较好的结果.)(xS6 为了使问题的提法更有一般性，通常在最小二乘法中考虑加权平方和)()()()(1100 xaxaxaxSnn)(mn（2）,)()(0222miiiiyxSx（3）这里 是 上的权函数，它表示不同点 处的数据比重不同.0)(x,ba)(,(iixfx就是 次多项式.)(xSn 若 是 次多项式，)(xkk 的一般表达式为（2）表示的线性形式.)(xS7其中,),(,),(T10

4、T10nndddaaada.),(),(),(),(),(),(),(),(),(101110101000nnnnnnG（7）,dGa 要使法方程(6)有惟一解，就要求矩阵 非奇异，G而 在 上线性无关不能推出)(,),(),(10 xxxn,ba矩阵 非奇异，必须加上另外的条件.Gknjjjkda0),().,1,0(nk（6）8若记,)()()(),(0miikijikjxxx（5）kmiikiikdxxfxf0)()()(),().,1,0(nk上式可改写为 knjjjkda0),().,1,0(nk（6）这个方程称为法方程法方程，可写成矩阵形式9一般可取 ，但这样做当 时，,1nxxs

5、pan3n通常对 的简单情形都可通过求法方程（6）得到 1n).(*xS 给定 的离散数据 ，,1,0),(miyxii)(xf 例如，bxaxSe)(,ln)(lnbxaxS求解法方程（6）将出现系数矩阵 为病态的问题，G 有时根据给定数据图形，其拟合函数 表面上)(xfy 不是（2）的形式，但通过变换仍可化为线性模型.若两边取对数得knjjjkda0),().,1,0(nk（6）)()()()(1100 xaxaxaxSnn)(mn（2）10 显然 在任意 个点上满足哈尔条件.nxx,1)(nmm哈尔条件，则法方程(6)的系数矩阵(7)非奇异，如果 在 上满足,)(,),(),(10bax

6、xxnmix0函数 的最小二乘解为)(xf 定义定义1 1设 的任意线,)(,),(),(10baxxxn性组合在点集 上至多只有 个)(,1,0,nmmixin不同的零点，则称 在点集 )(,),(),(10 xxxn,1,0,mixi上满足哈尔哈尔(Haar)条件条件.,1,0,*nkaakk方程(6)存在惟一的解从而得到于是knjjjkda0),().,1,0(nk（6）11,)()()()()()(0202*miiiimiiiixfxSxxfxSx这样得到的 ，)(*xS对任何形如(2)的 ，)(xS).()()()(*1*10*0*xaxaxaxSnn都有故 确是所求最小二乘解.)(

7、*xS)()()()(1100 xaxaxaxSnn)(mn（2）12 例例7 7113125.8865.4454321iiifx这样就变成了形如（2）的线性模型.此时，若令 则,ln),(ln)(bBaAxSxS,)(BxAxS已知一组实验数据如下，求它的拟合曲线.13 解解 从图中看到各点在一条直线附近，故可选择线性函数作拟合曲线，将所给数据在坐标纸上标出，见图1.图114 令,)(101xaaxS,8),(4000ii,22),(),(400110iiix,74),(40211iiix,47),(400iiiff.5.145),(401iiiifxf,1)(,1,40 xnm这里故,)(

8、1xx15.5.1457422,472281010aaaa解得.13.1,77.210aa.13.177.2)(*1xxS由（6）得方程组 于是所求拟合曲线为knjjjkda0),().,1,0(nk（6）16 关于多项式拟合，Matlab中有现成的程序 ),(polyfitmyxa 其中输入参数 为要拟合的数据，为拟合多项式的次数，yx,m输出参数 为拟合多项式的系数.a 利用下面的程序，可在Matlab中完成上例的多项式拟合.17x=1 1 2 3 3 3 4 5;f=4 4 4.5 6 6 6 8 8.5;aa=poly(x,f,1);y=polyval(aa,x);plot(x,f,r

9、+,x,y,k)xlabel(x);ylabel(y);gtext(y=s1(x)）18结果如下：19 例例1 1设数据 由表3-1给出，)4,3,2,1,0)(,(iyxii,ebxay 用最小二乘法确定 及 .ab46.845.753.679.510.500.275.150.125.100.143210iiyxi1表3 解解,lniiyy表中第4行为通过描点可以看出数学模型为它不是线性形式.,ebxay 用给定数据描图可确定拟合曲线方程为两边取对数得20 若令,ln,lnaAyy先将 转化为),(iiyx),(iiyx为确定 ，bA,根据最小二乘法，取,1)(,)(,1)(10 xxxx.

10、lnlnbxay.,1,xbxAy则得数据表见表3-1.得,5),(00,5.7),(4010iix,875.11),(40211iix135.2008.2876.1756.1629.146.845.753.679.510.500.275.150.125.100.143210iiiyyxi1表3 21.422.14),(401iiiyxy,404.9),(400iiyy故有法方程.422.14875.1150.7,404.950.75bAbA解得.071.3e,505.0,122.1AabA.e071.3505.0 xy 于是得最小二乘拟合曲线为 22miikimiikiikkkkxxxxfx

11、fa020*)()()()()(),(),().,1,0(nk（9）则方程（6）的解为 且平方误差为.)(02*2222nkkkaAfknjjjkda0),().,1,0(nk（6）23 接下来根据给定节点 及权函数 mxxx,10,0)(x构造带权 正交的多项式 .)(x)(xPn 注意 ，用递推公式表示 ，即mn)(xPk)()()()(),()()(,1)(1110110 xPxPxxPxPxxPxPkkkkk).1,2,1(nk（10）这里 是首项系数为1的 次多项式，)(xPkk根据 的)(xPk正交性，得24 利用下面的程序，可在Matlab中完成曲线拟合.x=1.00 1.25

12、1.50 1.75 2.00;y=5.10 5.79 6.53 7.45 8.46;y1=log(y);aa=poly(x,y1,1);a=aa(1);b=exp(aa(2);y2=b*exp(a*x);plot(x,y,r+,x,y2,k)xlabel(x);ylabel(y);gtext(y=a*exp(bx);25结果如下：26 而 ，11ks,0),(),(skskxPPPxP于是由（5.12），当 时，2 ks.0),(1skPP 另外，是首项系数为1的 次多项式，它可由)(xxPs1s由归纳法假定，当 时20ks,0),(slPP.0),(1skPP110,sPPP的线性组合表示.

13、由归纳法假定又有),(),)(),(111skkskkskPPPPxPP（12）).,(),(),(11skkskkskPPPPPxP272 用正交多项式做最小二乘拟合用正交多项式做最小二乘拟合 如果 是关于点集)(,),(),(10 xxxn,0,0)()()(),(0kmiikijikjAxxx,kj,kj（8）用最小二乘法得到的法方程组（6），其系数矩阵 是病态的.G带权 正交的),1,0(mixi),1,0()(mixi函数族，即knjjjkda0),().,1,0(nk（6）28),(),(),(),(111kkkkkkkkkkPPPPPxPPP至此已证明了由（10）及（11）确定的

14、多项式 )(xPk组成一个关于点集 的正交系.ix 用正交多项式 的线性组合作最小二乘曲线拟合，)(xPk只要根据公式（10）及（11）逐步求 的同时，)(xPk相应计算出系数*ka.0),(),(),(),(kkkkkkkkPPPPPxPPxP最后，由 和 的表达式（11）有 kk29),(),(11kkkkPPPP（11）).1,2,1(nkmiikimiikikxPxxPx02102)()()()(miikimiikiikxPxxPxx02021)()()()(下面用归纳法证明这样给出的 是正交的.)(xPk)(),()(),(xPxPxPxxPkkkk),(),(kkkkPPPxP30

15、),(),(),(0010010PPxPPPP 假定 对 及)(0),(slPPsl1,1,0ls,1,0kl要证 对 均成立.0),(1skPPks,1,0由（10）有),(),)(),(111skkskkskPPPPxPP 由（10）第二式及（11）中 的表达式，有 1),(),(),(),(00000000PPPPPxPxPP.0nk 均成立，（5.12）).,(),(),(11skkskkskPPPPPxP)()()()(),()()(,1)(1110110 xPxPxxPxPxxPxPkkkkk).1,2,1(nk（10）)()()()(),()()(,1)(1110110 xPxP

16、xxPxPxxPxPkkkkk).1,2,1(nk（10）31由假定有),(),(11kkkkxPPPxP 再考虑（13）),(),(),(),(1111111kkkkkkkkkkPPPPPxPPP),(10kjjjkkPcPP).,(kkPP利用（11）中 表达式及以上结果，得 k),(),(),(11111kkkkkkkPPPxPPP.0),(),(kkkkPPPP),(),(11kkkkkPPPP32),(),(*kkkkPPPfa),1,0(nk并逐步把 累加到 中去，最后就可得到所求的)(*xPakk)(xS).()()()(*1*10*0 xPaxPaxPaxSynn 用这种方法编

17、程序不用解方程组，只用递推公式，并且当逼近次数增加一次时，只要把程序中循环数加1，其余不用改变.这里 可事先给定或在计算过程中根据误差确定.n)()()()()(200ikmiimiikiixPxxPxfx拟合曲线33 只要 在 上分段连续，则级数（3）一致收敛到 .)(xf 2,0)(xf 对于最佳平方逼近多项式（6.1）有.)()()()(222222xSxfxSxfnn由此可以得到相应于（4.11）的贝塞尔不等式.d)(1)(2120212220 xxfbaankkk因为右边不依赖于 ，左边单调有界，所以级数 n10)sincos(21kkkkxbkxaa（3）nxbnxaxbxaaxSnnnsincossincos21)(110（6.1）.)()(22122*xfxankkk（11）