37数值计算方法教案曲线拟合与函数逼近.docx
《37数值计算方法教案曲线拟合与函数逼近.docx》由会员分享，可在线阅读，更多相关《37数值计算方法教案曲线拟合与函数逼近.docx（17页珍藏版）》请在冰豆网上搜索。
37数值计算方法教案曲线拟合与函数逼近
第三章插值法与最小二乘法
3.7最小二乘法
一、教学目标及基本要求
通过对本节课的学习，使学生掌握数值逼近的拟合方法。
二、教学内容及学时分配
本章主要介绍数值分析的最小二乘法。
具体内容如下：
曲线拟合原理，最小二乘法。
三、教学重点难点
1．教学重点：
曲线拟合。
2.教学难点：
最小二乘法。
四、教学中应注意的问题
多媒体课堂教学为主。
适当提问，加深学生对概念的理解。
一．曲线拟合
1.问题提出：
已知多组数据，由此预测函数的表达式。
数据特点：
（1）点数较多。
（2）所给数据存在误差。
解决方法：
构造一条曲线反映所给数据点的变化总趋势，即所谓的“曲线拟合”。
2.直线拟合的概念
设直线方程为y=a+bx。
则残差为：
，其中。
残差是衡量拟合好坏的重要标志。
可以用MATLAB软件绘制残差的概念。
x=1:
6;
y=[3,4.5,8,10,16,20];
p=polyfit（x,y,1）;
xi=0:
0.01:
7;
yi=polyval（p,xi）;
plot（xi,yi,x,y,'o'）;
y1=polyval（p,x）;
holdon
fori=1:
6
plot（[i,i],[y（i）,y1（i）],'r'）;
end
可以绘制出如下图形：
三个准则：
（1）最小
（2）最小
（3）最小
3.最小二乘法的直线拟合
问题：
对于给定的数据点，求一次多项式y=a+bx，使得总误差Q最小。
其中。
根据
故有以下方程组（正则方程）：
例1.给定数据表，求最小二乘拟合一次多项式
xi
165
123
150
123
141
yi
187
126
172
125
148
解：
N=5，=702，=758，=99864，=108396。
则有方程组
解得a=-60.9392，b=1.5138，则一次多项式为y=-60.9392+1.5138b
用MATLAB计算并画图如下：
x=[165,123,150,123,141];
y=[187,126,172,125,148];
A（1,1）=5;A（1,2）=sum（x）;A（1,3）=sum（y）;
A（2,1）=sum（x）;A（2,2）=sum（x.^2）;A（2,3）=x*y';
B=rref（A）;
a=B（1,3）;b=B（2,3）;
p=[b,a];
%以上四行，可以用一行命令p=polyfit（x,y,1）;替代。
xi=min（x）-1:
0.01:
max（x）+1;
yi=polyval（p,xi）;
plot（xi,yi,x,y,'o'）;
绘制如下图形
4.最小二乘法的多项式拟合
问题：
对于给定的数据点，求m次多项式（m<其中。
根据
故有正则方程：
当m=2时，有
例2.求数据表的最小二乘法拟合的二次多项式函数
xi
-1
-0.75
-0.5
-0.25
0
0.25
0.5
0.75
1
yi
50
40
25
20
18
21
35
56
66
在MATLAB命令窗口输入：
x=-1:
0.25:
1;
y=[50,40,25,20,18,21,35,56,66];
A（1,1）=length（x）;A（1,2）=sum（x）;A（1,3）=sum（x.^2）;A（1,4）=sum（y）;
A（2,1）=sum（x）;A（2,2）=sum（x.^2）;A（2,3）=sum（x.^3）;A（2,4）=y*x';
A（3,1）=sum（x.^2）;A（3,2）=sum（x.^3）;A（3,3）=sum（x.^4）;A（3,4）=y*（x.^2）';
B=rref（A）;
p=[B（3,4）,B（2,4）,B（1,4）];
%以上五行可以用p=polyfit（x,y,2）;替代
xi=min（x）-0.1:
0.01:
max（x）+0.1;
yi=polyval（p,xi）;
plot（xi,yi,x,y,'o'）;
可以绘制出如下图形：
例3.从三次多项式上找出21个点，然后对这21个点进行“差错处理”，得到新的21个点，根据新的21个点拟合一个新的3次多项式函数，然后和原函数进行比较。
解：
在MATLAB命令窗口输入：
p3=inline（'2.*x.^3-3.*x.^2+4.*x-5'）;
x=-10:
10;
y=p3（x）;
e=randn（1,length（x））*80;
y=y+e;
p=polyfit（x,y,3）;
xi=-10:
0.01:
10;
yi=polyval（p,xi）;
plot（xi,yi,x,y,'o'）;
holdon
fplot（p3,[-10,10],'r'）;
5.利用MATLAB的多项式拟合命令polyfit来实现多项式的插值
例1.过随机6个数据点，构造5次多项式函数。
解：
在MATLAB命令窗口输入：
x=1:
6;
y=round（10*randn（1,6））;
p=polyfit（x,y,length（x）-1）;
xi=1:
0.01:
6;
yi=polyval（p,xi）;
plot（xi,yi,x,y,'o'）;
可以得到以下图形：
6.利用最小二乘法解超定方程组
例1．解下列超定方程组
解：
设超定方程的解为。
方法一：
点到4条直线的距离平方分别为：
，，，
设，根据，有：
=0
=0
化简有：
解得
方法二：
最小二乘法：
点关于4条直线的残差平方和为：
+
根据，有：
=0
=0
化简有：
解得
用MATLAB命令有：
symsx0y0
f1=4*（2*x0+4*y0-11）+2*3*（3*x0-5*y0-3）+2*（x0+2*y0-6）+2*2*（2*x0+y0-7）
f2=8*（2*x0+4*y0-11）-2*5*（3*x0-5*y0-3）+4*（x0+2*y0-6）+2*（2*x0+y0-7）
解得：
f1=
36*x0-6*y0-102
f2
f2=
-6*x0+92*y0-96
继续在MATLAB命令窗口输入：
A=[36,-6,102;-6,92,96];
B=rref（A）
x0=B（1,3）
y0=B（2,3）
解得：
x0=
3.04029304029304
y0=
1.24175824175824
方法三：
最小二乘法（矩阵运算）
针对方程组的最小二乘近似解即为方程组的解
于是，在MATLAB命令窗口输入：
A=[2,4;3,-5;1,2;2,1];
b=[11;3;6;7];
x=inv（A'*A）*A'*b
计算结果为：
x=
3.04029304029304
1.24175824175824
方法四：
用MATLAB左除命令“\”
在MATLAB命令窗口输入：
A=[2,4;3,-5;1,2;2,1];
b=[11;3;6;7];
x=A\b
即可以得到答案
x=
3.04029304029304
1.24175824175824
可以看出用MATLAB的左除“\”命令计算得到的答案与最小二乘法得到的答案是一致的。
其实，MATLAB的左除“\”命令就是按照最小二乘法的原来来编写的。
另外，可以用MATLAB的ezplot命令绘制四条直线的图形
ezplot（'2*x+4*y=11'）;
holdon
ezplot（'3*x-5*y=3'）;
ezplot（'x+2*y=6'）;
ezplot（'2*x+y=7'）;
plot（2.99,1.30,'o'）;
A=[2,4;3,-5;1,2;2,1];
b=[11;3;6;7];
x=A\b
plot（x
（1）,x
（2）,'*'）;
绘制图形如下：
二．函数逼近
问题，已知函数f（x），求一个多项式函数在区间[a,b]上逼近f（x）。
解决方法：
函数的最佳平方逼近。
令
，使Q最小，则有
例1.求一次多项式在上逼近函数。
解：
构造直线为：
，，则有
，
解得：
a=0.6644389，b=0.1147707
在MATLAB命令窗口输入：
xi=0:
0.01:
pi/2;
yi=sin（xi）;
p=polyfit（xi,yi,1）;
pi=polyval（p,xi）;
plot（xi,yi,xi,pi）;
可以绘制以下图形：
作业：
（1）用最小二乘法求一个形如的经验公式，使它与下列数据表拟合。
xi
19
25
31
38
41
yi
19.0
32.3
49.0
73.3
97.8
解：
方法一，最小二乘法；方法二，用解超定方程组的思路来解题。
，根据，有：
在MATLAB命令窗口输入：
x=[19,25,31,38,44];
y=[19,32.3,49,73.3,97.8];
A=[5,sum（x.^2）,sum（y）;sum（x.^2）,sum（x.^4）,x.^2*y'];
B=rref（A）;
p=[B（2,3）,0,B（1,3）];
xi=min（x）:
0.01:
max（x）;
yi=polyval（p,xi）;
plot（xi,yi,x,y,'o'）;
绘制图形如下：
（2）已知数据表如下，试用二次多项式拟合。
xi
0
1
2
3
4
5
6
yi
15
14
14
14
14
15
16
（3）求一个形如（a，b为常数，a>0）的经验公式，使它能和下表数据拟合。
xi
1
1.25
1.5
1.75
2
yi
5.1
5.79
6.53
7.45
8.46
解：
公式可以变为：
lny=lna+bx，进一步可以写为Y=A+bx。
其中Y=lny，A=lna，对应表格为：
xi
1
1.25
1.5
1.75
2
yi
5.1
5.79
6.53
7.45
8.46
Yi
1.629
1.756
1.867
2.008
2.135
（4）求函数在区间[1/4,1]上的最小一次式。
小结
这节课我们主要学习了最小二乘法的原理及其上机实现过程，大家要求掌握最小二乘法的原理，并会设计算法，设计Matlab程序。
作业：
课后相应习题。
升级会员 