数值计算方法第三章函数逼近与曲线拟合当函数只在有限点集上给定函数值,要在包含该点集的区间上用公式给出函数的简单表达式,这些都涉及到在区间[a,b]上用简单函数逼近已知复杂函数的问题,这就是函数逼近问题。插值法就是函数逼近问题的一种仍然是已知x1…xm;y1…ym,求一个简单易算的近似函数P(x)f(x)。但是①m很大;②yi本身是测量值,不准确,即yif(xi)这时没必要取P(xi)=yi,而要使P(xi)yi总体上尽可能小。使误差在某种度量意义下最小实例:考察某种纤维的强度与其拉伸倍数的关系,下表是实际测定的24个纤维样品的强度与相应的拉伸倍数是记录:编号拉伸倍数强度编号拉伸倍数强度11.91.41355.5221.3145.2532.11.81565.542.52.5166.36.452.72.8176.5662.72.5187.15.373.531986.583.52.72087944218.98.51043.52298114.54.2239.58.1124.63.524108.1iiyxiiyx1234567891012345678912345678910123456789纤维强度随拉伸倍数增加而增加系要关系应是线性关的主与拉伸倍数因此可以认为强度xy并且24个点大致分布在一条直线附近xxy10)(为待定参数其中10,越接近越好样本点与所有的数据点我们希望),)(()(10iiyxxxy必须找到一种度量标准来衡量什么曲线最接近所有数据点(1)最小二乘法的基本概念iiiyxy)(令一般使用mii0222在回归分析中称为残差miiiyxy02))((准偏离程度大小的度量标与数据点作为衡量),()(iiyxxy称为平方误差在回归分析中称为残差平方和从而确定(1)中的待定系数mii0222miiiyxy02))((注意(1)式是一条直线关系的关系并不一定是线性但yx,因此将问题一般化)(,xSyyx的关系为设来自函数类其中)(xS来自线性函数类中如)()1(xy为给定的一组数据设),,1,0)(,(miyxii),,1,0)((nixi的基函数为设函数类mn一般要求即生成的函数集是由也称,),,1,0)((nixi)}(,),(),({10xxxspannmii0222miiiyxS02))((仍然定义平方误差njjjxaxS0)()(我们选取的度量标准是)(*xS中选取一个函数在函数类njjjxaxS0*)()(*)(*)(*)(*1100xaxaxann22*miiiyxS02))(*(miiixSyxS02)())((min22)(minxS中的任意函数为其中mjjjxaxS0)()(---------(2)---------(3)数据拟合的最小二乘法的方法为的求函数称满足条件njjjxaxS0*)()(*)3(为最小二乘解njjjxaxS0*)()(*为拟合系数为拟合函数),,1,0(,)()(0njaxaxSjnjjj(),(0,1,,)jSxajn在确定了拟合函数后如何求拟合系数呢?满足拟合条件使得)3()()(*0*njjjxaxS误差称为最小二乘解的平方22*miinjijjyxa020))((miiiyxS02))((法方程组22njjjxaxS0)()(由的函数为拟合系数),,1,0(njaj可知因此可假设),,,(10naaamiinjijjyxa020))((因此求最小二乘解转化为二次函数的问题点极小值的最小值求*,*,*,)(),,,(1010nnaaaaaa由多元函数取极值的必要条件0),,,(10knaaaank,,1,0)]())((2[00ikmiinjijjxyxaka0得即miikimiiknjijjxyxxa000)()()(0)]()()([00ikmiinjikijjxyxxa),,,(10naaamiinjijjyxa020))((miikimiiknjijjxyxxa000)()()(miikinjjikmiijxyaxx000)()]()([nk,,1,0---------(4)miikiikmiinnikmiiikmiixyxxaxxaxxa00011000)()()()()()()(nk,,1,0即元线性方程组的是一个关于显然1,,,)4(10naaan引入记号))(,),(),((10mrrrxxxr),,,(10myyyf)()(),(0ijmiikjkxx则由内积的概念可知imiikkyxf0)(),(---------(5)---------(6)),(jk),(kj显然内积满足交换律方程组(4)便可化为),(),(),(),(1100faaaknknkknk,,1,0---------(7)的线性方程组常数项为这是一个系数为),(),,(fkjk将其表示成矩阵形式naaa10),(),(),(10fffn),(),(),(01000n),(),(),(11101n),(),(),(10nnnn-----(8)上的法方程组在点式为函数序列称mnxxxxxx,,,)(,),(),()8(1010的基为函数类由于)(,),(),(10xxxn必然线性无关因此)(,),(),(10xxxn并且其系数矩阵为对称阵所以法方程组的系数矩阵非奇异,即0])),det[((nnji根据Cramer法则,法方程组有唯一解*,*,*,1100nnaaaaaa*),*,*,(10naaamiinjijjyxa020))((),,,(10naaa即是的最小值22*miiiyxS02))(*(miiixSyxS02)())((min22)(minxS所以miinjijjyxa020))(*(miinjijjxSyxa020)())((minmiinjijjyxa020))(*(为最小二乘解njjjxaxS0*)()(*因此的拟合函数作为常使用多项式),,1,0)(,()()(miyxxPxSiin作为一种简单的情况,的基函数为拟合函数)()(xPxSn,1)(0x,)(1xx,,)(,kkxxnnxx)(基函数之间的内积为)()(),(0ijmiikjkxxmijikixx0mijkix0imiikkyxf0)(),(miikiyx022*平方误差miiiyxS02))(*(njjjfaff0),(*),(例1.回到本节开始的实例,从散点图可以看出纤维强度和拉伸倍数之间近似与线性关系xaaxy10)(故可选取线性函数为拟合函数,其基函数为1)(0xxx)(1建立法方程组根据内积公式,可得24),(005.127),(1061.829),(111.113),(0f6.731),(1f法方程组为61.8295.1275.1272410aa6.7311.1131505.00a即为所求的最小二乘解xxy8587.01505.0)(*8587.01a解得6615.5*22平方误差为1234567891012345678912345678910123456789拟合曲线与散点的关系如右图:例2.求拟合下列数据的最小二乘解x=.24.65.951.241.732.012.232.522.772.99y=.23-.26-1.10-.45.27.10-.29.24.561解:从数据的散点图可以看出xxycos之间具有三角函数关系与xexy系之间还具有指数函数关与xxyln系之间还具有对数函数关与因此假设拟合函数与基函数分别为xcexbxaxScosln)(xex)(2xxln)(0xxcos)(100.511.522.53-1.5-1-0.500.51xy6.7941-5.347563.2589-5.34755.1084-49.008663.2589-49.00861002.51.6163-2.382726.7728通过计算,得法方程组的系数矩阵及常数项矩阵为00.511.522.53-1.5-1-0.500.51xyGo!用Gauss列主元消去法,得cba-1.0410-1.26130.030735xexxxS030735.0cos2613.1ln0410.1)(*的最小二乘解是关于xy22*20))(*(miiiyxS20)030735.0cos2613.1ln0410.1(miixiiyexxi92557.0拟合的平方误差为图象如图例3.在某化学反应里,测得生成物浓度y%与时间t的数据如下,试建立y关于t的经验公式x=1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16y=4.00,6.40,8.00,8.80,9.22,9.50,9.70,9.86,10.00,10.20,10.32,10.42,10.50,10.55,10.58,10.60解:的散点图与浓度画出时间yx例:xy(xi,yi),i=1,2,…,m方案一:设baxxxPy)(求a和b使得最小。miiiiybaxxba12)(),(线性化/*linearization*/:令,则xXyY1,1bXaY就是个线性问题将化为后易解a和b。),(iiYX),(iiyxTakeiteasy!Wejusthavetolinearizeit…方案二:设xbeaxPy/)((a0,b0)线性化:由可做变换xbaylnlnbBaAxXyY,ln,1,lnBXAY就是个线性问题将化为后易解A和B),(iiYX),(iiyxxbAeaxPBbea/)(,,定义权函数:①离散型/*discretetype*/根据一系列离散点拟合时,在每一误差前乘一正数wi,即误差函数,这个wi就称作权/*weight*/,反映该点的重要程度。),...,1(),(niyxiiniiiiyxPw12])([②连续型/*continuoustype*/在[a,b]上用广义多项式P(x)拟合连续函数f(x)时,定义权函数(x)C[a,b],即误差函数=。权函数(x)必须满足:非负、可积,且在[a,b]的任何子区间上(x)0。dxxyxPxba2)]()([)(加权最小二乘法),,1,0)(,(miyxii对于一组给定的数据点中在拟合的数据点),,1,0)(,(miyxii各点的重要性可能是不一样的的重度表示数据点假设),(iiiyx重度:即权重或者密度,统称为权系数mk,,1,0定义加权平方误差为miii0222miiiiyxy02))((-----(9)例:用来拟合,w12210xaxaayx1234y4101826解:0(x)=1,1(x)=x,