最小二乘法的基本原理和多项式拟合

数值分析2013年12月最小二乘法的基本原理和多项式拟合摘要：随着科技发展和社会进步，尤其是计算机大范围的普及，计算机应用逐渐由大规模科学计算的海量数据处理转向日常学习中遇到的小型数据的处理与计算，这就产生了解决各种实际问题需要的各种应用程序，特别是在电磁检验实验中的应用更为广泛。最小二乘法（又称最小平方法）是一种数学优化技术，是利用最小化误差的平方和寻找数据的最佳函数匹配的一种计算方法[1]，本文对最小二乘法进行了深入细致的研究，利用Matlab语言编制程序和数学知识相结合，通过实验数据的输入，实现多项式和曲线拟合的输出，并利用设计的程序实现了一些实际问题的求解和处理。关键词：最小二乘法，曲线拟合，多项式一、前言曲线拟合又称作函数逼近，是求近似函数的一类数值方法．它不要求近似函数在每个节点处与函数值相同，只要求其尽可能的反映给定数据点的基本趋势以及某种意义上的无限“逼近”．在需要对一组数据进行处理、筛选时，我们往往会选择合理的数值方法，而曲线拟合在实际应用中也倍受青睐．采用曲线拟合处理数据时，一般会考虑到误差的影响，于是我们往往基于残差的平方和最小的准则选取拟合曲线的方法，这便是经常所说的曲线拟合的最小二乘法．通过对一些文献的分析和整理，不仅有数据处理（数据采集），还有模型建立，可以了解到曲线拟合的最小二乘法的应用领域较为广泛。二、最小二乘法法的介绍1），最小二乘法的基本原理从整体上考虑近似函数同所给数据点(i=0,1,…,m)误差(i=0,1,…,m)(i=0,1,…,m)绝对值的最大值，即误差向量的∞—范数；二是误差绝对值的和，即误差向量r的1—范数；三是误差平方和的算术平方根，即误差向量r的2—范数；前两种方法简单、自然，但不便于微分运算，后一种方法相当于考虑2—范数的平方，因此在曲线拟合中常采用误差平方和来度量误差(i=0，1，…，m)的整体大小。数据拟合的具体作法是：对给定数据(i=0,1,…，m)，在取定的函数类中,求,使误差(i=0,1,…,m)的平方和最小，即=从几何意义上讲，就是寻求与给定点(i=0,1,…,m)的距离平方和为最小的曲线（图6-1）。函数称为拟合函数或最小二乘解，求拟合函数的方法称为曲线拟合的最小二乘法。可有不同的选取方法.6—12），多项式拟合假设给定数据点(i=0,1,…,m)，为所有次数不超过的多项式构成的函数类，现求一,使得(1)当拟合函数为多项式时，称为多项式拟合，满足式（1）的称为最小二乘拟合多项式。特别地，当n=1时，称为线性拟合或直线拟合。［显然为的多元函数，因此上述问题即为求的极值问题。由多元函数求极值的必要条件，得(2)即(3)（3）是关于的线性方程组，用矩阵表示为(4)式（3）或式（4）称为正规方程组或法方程组。可以证明，方程组（4）的系数矩阵是一个对称正定矩阵，故存在唯一解。从式（4）中解出(k=0,1,…，n)，从而可得多项式(5)可以证明，式（5）中的满足式（1），即为所求的拟合多项式。我们把称为最小二乘拟合多项式的平方误差，记作由式(2)可得(6)多项式拟合的一般方法可归纳为以下几步：(1)由已知数据画出函数粗略的图形——散点图，确定拟合多项式的次数n；(2)列表计算和；(3)写出正规方程组，求出；(4)写出拟合多项式。在实际应用中，或；当时所得的拟合多项式就是拉格朗日或牛顿插值多项式。3），最小二乘拟合多项式的存在唯一性定理1设节点互异，则法方程组（4）的解存在唯一。证由克莱姆法则，只需证明方程组（4）的系数矩阵非奇异即可。用反证法，设方程组（4）的系数矩阵奇异，则其所对应的齐次方程组（7）有非零解。式(7)可写为（8）将式（8）中第j个方程乘以(j=0,1,…，n)，然后将新得到的n+1个方程左右两端分别相加，得因为其中所以(i=0,1,…,m)是次数不超过n的多项式，它有m+1＞n个相异零点，由代数基本定理，必须有，与齐次方程组有非零解的假设矛盾。因此正规方程组（4）必有唯一解。定理2设是正规方程组（4）的解，则是满足式（1）的最小二乘拟合多项式。证只需证明，对任意一组数组成的多项式，恒有即可。因为(k=0,1,…，n)是正规方程组（4）的解，所以满足式（2），因此有故为最小二乘拟合多项式。2算法实现在化学反应中，由实验测得分解物浓度与时间的关系如下表1所示表1浓度（y）与时间(t)的关系实验数据表t0510152025y01.272.162.863.443.87t303540455055y4.154.374.514.584.624.643思路解答根据题设要求即根据所给数据求一未知的多项式拟合曲线，由于曲线未知，只能通过已知数据去确定待定的多项式参数，进而使数据尽可能的逼近所求曲线．最后再通过总误差比较，取最小误差对应多项式即可．不妨设所求多项式为1010mjmjmjyaaxaxax[fiv]，式中m表示多项式的最高次数．说明：表格(1)中iy是题设给定的函数值，f0，f1，f2，f5，f9分别表示拟合多项式最高次数m=1,2,3,6,10时的拟合值．由于数据不是很多，可以直接用最小二乘法找出最好的拟合曲线．具体就是比较总误差的大小，然后取误差最小的那条拟合曲线即可．当然，可根据均方差2121]))(([iniitfyQ和偏差G=yf。计算求得．同时根据法方程组或者多项式公式，可求出其多项式系数表a据表a可得一次多项式系数a1=0.0765，a0=1.2667，故2677.10765.00fx可得二次多项式系数a2=—0.0022，a1=0.1992，a0=0.2453，故三次多项式系数a3=0.000032499，a2=-0.0049，a1=0.2557，a0=0.0122，故同理六次和十次多项式拟合曲线分别为6-2344-55-6-67-88-109-8101410*8.76060.24900.01810.0058-10*6.519210*4.0014-10*1.480910*3.3893-10*4.697110*6140.310*1852.19xxxxxxxxxxfa.用matlab绘出多项式拟合的直观图：itiyiti^2tiyiti^3ti^4ti^5ti^6ti^2yiti^3yi00000000000151.27256.3512562531251562531.75158.752102.1610021.610010000100000100000021621603152.8622542.933755062575937511390625643.59652.54203.4440068.880001600003200000640000001376275205253.8762596.751562539062597656252441406252418.7560468.756304.15900124.5270008100002430000072900000037351120507354.371225152.954287515006255252187518382656255353.25187363.88404.511600180.4640002560000102400000409600000066402656009454.582025206.191125410062518452812583037656259274.5417352.510504.622500231125006250000312500000156250000001155057750011554.643025255.216637591506255032843752768064062514036771980和33040.47126501386.554311002498375011933625005859321875055274.752431806图一4总误差比较与分析：用matlab程序[tw]计算可得最大偏差和均方差分别为Gw1=1.2677，Gw2=0.2453，Gw3=0.0661，Gw6=0.0188，Gw10=0.00400074.0))((100342.0))((60933.0))((34471.0))((20858.2))((1211219211215211212121211211210iiiiiiiiiifyEfyEfyEfyEfyE根据图一可知一次多项式拟合和二次多项式拟合偏差太大，三次、六次、十次多项式比较贴近真实值，再结合多项式拟合的最大偏差和均方差，可清晰知道十次多项式拟合的数据最小故十次多项式拟合更接近实际值，综合可知10次多项式拟合曲线比1、2、3和6次拟合更加接近题设所给数据．三、参考文献[1]李庆扬，王能超，易大义.数值分析[M].清华大学出版社:2009年：第三章[2]姜健飞，胡良剑，唐俭.数值分析及其MATLAB实验[M].科学出版社：2004年：第四章[3]罗建军，杨琦，冯博琴.MATLAB数值分析[M].机械工业出版社：2012年：第三、七章[4]张铁，阎家斌．数值分析［M］．北京:冶金工业出版社，2005．第二、三章[5]李岳生，黄友谦.数值逼近.人民教育出版社,1978，（3）：36-41[6]陈杰，张增强，于锋.MATLAB宝典[第3版].北京：电子工业出版社，2011.1.第六章第二三节四、附录t=[0510152025303540455055]y=[01.272.162.863.443.874.154.374.514.584.624.64]figure(1)plot(t,y,'bo')gridonholdonxs0=polyfit(t,y,1);%插值点，其中1表示一次多项式，即直线a=xs0(1)%插值求得一次项系数b=xs0(2)%插值求得常数项项系数f0=polyval(xs0,t)%求得拟合出来的多项式在输入值为t下的y值xs1=polyfit(t,y,2);%插值点，其中1表示一次多项式，即直线a=xs1(1)%插值求得二次项系数b=xs1(2)%插值求得一次项系数c=xs1(3)%插值求得常数项系数f1=polyval(xs1,t)%求得拟合出来的多项式在输入值为t下的y值xs2=polyfit(t,y,3);a=xs2(1)%插值求得三次项系数b=xs2(2)%插值求得二次项系数c=xs2(3)%插值求得一次项系数d=xs2(4)%插值求得常数项系数f2=polyval(xs2,t)%求得拟合出来的多项式在输入值为t下的y值xs5=polyfit(t,y,6);a=xs5(1)%插值求得六次项系数b=xs5(2)%插值求得五次项系数c=xs5(3)%插值求得四次项系数d=xs5(4)%插值求得三次项系数e=xs5(5)%插值求得二次项系数f=xs5(6)%插值求得一次项系数g=xs5(7)%插值求得常数项系数f5=polyval(xs5,t)%求得拟合出来的多项式在输入值为t下的y值xs9=polyfit(t,y,10);a=xs9(1)%插值求得十次项系数b=xs9(2)%插值求得九次项系数c=xs9(3)%插值求得八次项系数d=xs9(4)%插值求得七次项系数e=xs9(5)%插值求得六次项系数f=xs9(6)%插值求得五次项系数g=xs9(7)%插值求得四次项系数h=xs9(8)%插值求得三次项系数i=xs9(9)%插值求得二次项系数j=xs9(10)%插值求得一次项系数l=xs9(11)%插值求得常数项系数f9=polyval(xs9,t)%求得拟合出来的多项式在输入值为t下的y值plot(t,f0,'r-')%用红色画出通过数据点而拟合出来的曲线holdonplot(t,f1,'g')%用绿色画出通过数据点而拟合出来的曲线holdonplot(t,f2,'b-')%用蓝色画出通过数据点而拟合出来的曲线holdonplot(t,f5,'k-')%用黑色画出通过数据点而拟合出来的曲线holdonplot(t,f9,'m-