Ch09：数值计算方法之最小二乘法与曲线拟合市公开课特等奖市赛课微课一等奖课件

第9章最小二乘法与曲线拟合最小二乘法是一个非常有用数学方法，它直接应用是求矛盾线性方程组最小二乘解.在工程中，它可用来求经验公式，对试验数据进行曲线拟合;在统计学中，它可用来求最小二乘预计，多元回归等;另外，在数值分析领域，它可用来进行误差分析。能够说，伴随我们未来知识增多，最小二乘法应用领域将愈来愈大。正是因为这个原因，我们把这个词出现在章标题中。第1页9.1问题提出为了让大家对最小二乘法适用实际问题有一个感性认识，我们首先构想一个简单工程问题：怎样尽可能准确地测算出一个弹簧长度与所受到外界拉力之间函数关系。我们知道，在弹性程度以内，弹簧长度改变量与所受到外界拉力成正比，亦即弹簧长度Y与所受到外界拉力X满足一个线性关系:Y=a0+a1x我们要作就是综合物理和数学方法求出a0和a1。为了处理这个问题，我们可按下面步骤来进行。第2页1.设计试验装置首先我们要设计一个试验装置，方便能方便地测量出在不一样拉力下弹簧长度来。我们能够简单地在工作台上固定一个支架，弹簧一端固定在支架上，在工作台边缘固定一个定滑轮，一根连线下面吊一个法码盘，连线另一端绕过滑轮与弹簧相连。这么就完成了如图所表示试验装置安装。第3页2.统计试验数据有了这个装置，我们能够设计一个试验方案，比如控制外界拉力分别为x0,x1,…,xn,对应地测量出弹簧长度分别为y0,y1,…,yn，从而得到n+1对试验数据列表：弹簧性能测试数据列表

这么，我们就完成了需要物理方法完成工作，接下来就是用数学方法来处理试验数据。Xx0x1…xnYy0y1…yn第4页3.组织线性方程组所谓分析弹簧性能，就是要尽可能准确地确定出a0和a1。为此我们能够把上面n+1对观察值代入到y=a0+a1x中，从而得到两个变元线性方程组

记我们得到方程组（1）简写形式Xa=y。第5页4.写出误差平方和在普通情况下，方程组（1）无解。这是因为y0,y1,…,yn都是观察值，它们不可防止地带有观察误差，但我们有理由寻找尽可能准确解。记εi=a0+a1xi-yi，i=0,1,…,n,显然，它就是yi误差。记（2）它就是n个误差平方和。注意在这里，x0,x1,…,xn和y0,y1,…,yn都是常数，而a0，a1都是待确定系数，所以应把J看作是（a0,a1）函数。提醒：在普通情况下，J作为（a0,a1）函数还是正定二次函数，它有全局唯一极小点。第6页5.能够转向求二次函数极小值即使我们不能求出（1）解，但我们能够求（2）极小解，也就是寻找（a0,a1）,使J(.,.)取极小值，并把这么解称为是（1）最小二乘解。（2）极小解可表为下面方程组解：（3）第7页6.正规方程利用（2）式对两个变元求偏导数，整理，我们能够得到（3）详细形式为

（4）结论：上面（4）解就是（1）最小二乘解。术语：上面方程组（4）普通称为正规方程。在不一样应用领域内，也有些人称之为法方程、回归方程。第8页7.求解正规方程我们完全能够利用上一章学到方法求解上面给出正规方程，从而基本完成问题求解。第9页8.检验解有效性假如我们得到了（1）最小二乘解（a0,a1）,那么我们也就找到了y与x近似函数关系y=a0+a1x.在普通情况下这个结果是很有用，但有时候依然会产生较大误差。为此，我们可作深入检验，这即使超出本课程范围，但也不算困难。记则e就是用a0+a1x表示y平均误差。在实际工作中，我们能够要求出一个临界值e0,实际也就是误差限，当e小于这个值时就认为结果有效，不然作深入考虑。第10页9.2最小二乘法设有线性方程组Xa=Y(1)其中a是方程组未知向量,且在普通情况下，n>m,X是满列秩矩阵，方程组无解。最小二乘法简单说来就是求（1）尽可能准确解数学方法。第11页1．化为二次型形式定义n+1维向量ε=(ε0,ε1,…,εn)T,令

ε=Xa-Y（2）我们能够把ε看作是用y近似表示Xa误差向量。思绪：假如能找到a,使得ε为零向量，那么a就是（1）准确解；假如a使得ε非常靠近零向量，那么a就是（1）很好近似解。为了度量这种靠近程度，我们将（2）两边左乘以各自转置向量得

εT·ε=(Xa-Y)T·(Xa-Y)（3）整理后即得εT·ε=aT(XTX)a-2aTXTY+yTY注意到εT·ε就是ε各项平方和，我们可用它来评价ε靠近零程度。第12页2.正规方程矩阵形式记（3）

利用数学分析中多变量函数取极值必要条件（在我们这里也是充分条件）能够推出（3）最优解也就是正规方程（4）解，我们称之为（1）最小二乘解。我们把求一个线性方程组最小二乘解方法称为最小二乘法。第13页3.算法说明假如我们兴趣仅在于求（1）最小二乘解，我们方法实际上尤其简单：第一步：在方程组（1）两边分别左乘以系数矩阵X转置矩阵即可得到正规方程（4）；第二步：求解正规方程组（4）；第三步：检验解有效性，也就是把（4）解代入到（3）中计算出对应J(a)，从而可深入得到（标准差）sqrt(J(a)/n)。第14页9.3多项式拟合问题我们在第一节中已经给出了一个系统方法进行试验数据处理，在那里我们只是简单地把试验数据构想为满足一个线性关系。假如我们更普通地把试验数据构想为满足一个多项式函数关系，那么对应问题,也就是求多项式系数，就是多项式拟合。类似地，我们也能够利用最小二乘法来求解多项式拟合问题。第15页1.问题提法假如我们得到了某个未知函数Y=f(x)n+1对观察值现在要找一个m次多项式Y=a0+a1x+…+amxm（1）来近似地表示这个未知函数。多项式拟合问题关键就是怎样今可能准确地确定出这个多项式系数a0,,a1,…,am,,对此，我们可用前面学过最小二乘法来求解。XX0

X1

…

Xn

Yy0

y1

…

Yn

第16页2.写出矛盾方程把n+1对观察值{(xk,yk)|k=0,1,…,n}分别代入到（1）中，我们能够得到下面线性方程组（2）普通说来，当n>m时，这个线性方程组误解，所以我们也把它称为矛盾方程第17页3.矛盾方程矩阵形式记则（2）可简写为Xa=Y第18页4.正规方程矩阵形式现在我们能够绕过利用二次函数求极值来推导正规方程，而是利用最小二乘法结果直接写出正规方程，它矩阵形式是XTXa=XTy记Q=XTX,P=Xty则正规方程可表为Qa=P第19页5.正规方程详细形式为了求Qa=P数值解，我们需要更详细地写出Q和P,由Q和P定义，我们不难得到第20页Q和P计算格式为了得到问题最小二乘接，我们需要写出P和Q数值形式，然后再求

Qa=p数值解。对于多项式拟和来说，矩阵Q中只有2m+1个不一样元素，利用这一点，能够大大降低我们计算量。给出了P和Q数值形式后，再求解Qa=p就不难了。第21页9.4普通曲线拟合问题对于普通工程问题来说，我们用多项式进行曲线拟协议常能够满足要求，但情形并不总是这么，为此我们简明介绍普通曲线拟合原理和方法。所谓多项式拟合，实际就是利用基函数{1,x,x2,…,xm}一个线性组合来近似表示某个未知函数。假如需要话，我们也能够采取其它基函数，比如{sin(nx),cos(nx),n=0,1,…},线性组合，从而形成更普通曲线拟合方法。第22页1.问题提法假设我们搜集到了两个相关变量x,y或某个函数y=f(x)n对观察值列表我们希望用m+1个基函数

0(x),

1(x),…,

m(x)线性组合y=a0

0(x)+a1

1(x)+…+am

m(x)（1）来近似地表示x,y间函数关系，我们问题就是怎样依据给出观察值列表最正确确实定出组合系数a0,a1,…,am。注：假如我们取

k(x)=xk,k=0,1,…,m，则得到就是上一节多项式拟合问题。XX0

X1

…

Xn

Yy0

y1

…

Yn

第23页2.矛盾方程形式与多项式拟合类似，把n+1对观察值{（xk,yk）,k=0,1,…,n}分别代入到（1）中，就能够得到下面线性方程组

（2）普通说来，当n>m时，上面问题无解，但我们一样能够求它最小二乘解。第24页3.正规方程矩阵形式记

则（2）可简写为Xa=Y对应地正规方程为（XTX）a=XTY记Q=XTX,P=XTY则正规方程可表为Qa=Y第25页4.正规方程矩阵形式为了求（2）最小二乘解，我们能够直接写出它正规方程（XTX）a=XTY记Q=XTX,P=XTY则正规方程可表为Qa=Y第26页5.正规方程详细形式为了求Qa=P数值解，我们需要更详细地