1、第五章 数据拟合,1、拟合的概念,2、拟合的方法和MATLAB编程实现,3、调用MATLAB命令实现拟合,4、范例,在一项工程实践中,通过观测,得到了一个离散的函数关系(xi,yi) i=1,2,n。由于工程的需要,我们希望揭示出反映这组离散数据的一个解析的函数关系。用几何术语来表达:根据平面上的观测点,要求确定一个函数曲线y=f(x), 使曲线尽量接近这些点。实现这个愿望的方法简称为曲线拟合(fitting a curve).,5.1 导 言,在生产实践和科学实验中,经常会遇到大量的不同类型的数据(data).这些数据提供了有用的信息,可以帮助我们认识事物的内在规律、研究事物之间的关系等.,
2、曲线拟合是根据实验获得的数据,建立自变量与因变量之间的函数关系,为进一步的深入研究提供线索。,5.2 引例:浓度变化规律 在化学反应中,为研究某化合物的浓度随时间的变化规律,测得一组数据如表5.1,曲线拟合与函数插值都是要根据一组二维数据(xi,yi)构造一个函数作为近似。由于近似的要求不同(是否通过观测点),二者在数学方法上有所不同。希望通过本章的学习,读者能掌握一些曲线拟合的基本方法,弄清楚曲线拟合与插值方法之间的区别,学会使用MATLAB软件进行曲线拟合。,表 5.1,t 时间 1 2 3 4 5 6 7 8 y浓度 4 6.4 8.0 8.4 9.28 9.5 9.7 9.86t时间
3、9 10 11 12 13 14 15 16 y浓度 10 10.2 10.32 10.42 10.5 10.55 10.58 10.6,表5.1中的数据反映了浓度随时间变化的函数关系,它是一种离散关系.若需要推断第20、40分钟的浓度值,就要用一个解析的函数y=f(t)来拟合表5.1中的离散数据,然后再算浓度f(20),f(40)。首先将这些离散数据描绘在直角坐标系下,得到散点图。然后观察浓度与时间之间呈现什么规律。,图5.1,浓度 y 随时间 t 呈“抛物线(二次函数)状”变化.,根据散点图,可以认为y与t的函数为y=a+bt+ct2,其中a,b,c为待定,称为参数。参数的选择需要科学的方
4、法和实验修正。,提示,函数形式确定以后,关键是要确定函数中含有的待定参数a,b,c.常用的方法是最小二乘法(method of least squares),下面介绍该方法的基本原理。,5.3 最小二乘法,平面上的点 (xi,yi) i=1,2,n。揭示出一个离散的函数关系; 设连续可微的函数y=f(x)很接近上述离散的函数关系。一般来说,最小,因此,我们的愿望是:如何选取 f(x) 的系数使得,,yi f(xi) i=1,2,n。,这就是函数拟合的最小二乘法问题。这一愿望的几何意义如图,设对 数据点,我们选定了二次多项式,代入 xi 于是有,寻求最合适的参数,使得,最小,即求,记,数学分析的
5、解法:,求偏导数,确定极小值,得到关于a,b,c的函数。,令3个偏导数为零,有,得到关于a,b,c 的线性方程组,MATLAB实现: x=1:16; y=4 6.4 8.0 8.4 9.28 9.5 9.7 9.86 10 10.2 10.32 10.42 10.5 10.55 10.58 10.6; a11=16; a12=sum(x); a21=a12; a13=sum(x.2);,a31=a13; a22=a13; a23=sum(x.3); a32=a23; a33=sum(x.4); A=a11 a12 a13;a21 a22 a23;a31 a32 a33; b1=sum(y);
6、b2=sum(y.*x); b3=sum(y.*x2); b=b1;b2;b3; X=Ab 得常数项4.3252 , 一次项系数1.0711, 二次项系数 0.0445。,5.4 曲线拟合的MATLAB实现,MATLAB软件提供基本的曲线拟合函数的命令。,多项式函数拟合:a=polyfit(xdata,ydata,n) 其中(xdata,ydata)为观测数据,n表示多项式的最高阶数。 输出参数a =a1,an,an+1 与多项式f(x)=a1xn+anx+an+1对应,例如:引例中的问题 t=1:16; y=4 6.4 8.0 8.4 9.28 9.5 9.7 9.86 10 10.2 10
7、.32 10.42 10.5 10.55 10.58 10.6; a=polyfit(t,y,2) a = -0.0445 1.0711 4.3252 即拟合函数为f(t)=-0.0445t2+1.0711t+4.3252 这一结果与我们采用数学分析的方法得到的结果一致。,对函数的拟合效果如何检测呢?仍然以图形来检测,我们将散点与拟合曲线画在一个画面上即可看出。,xi=linspace(0,16,160);yi=polyval(a,xi); %plot(x,y,o,xi,yi) % 图略,右图是以8次多项式拟合的效果 a=polyfit(t,y,8); xi=linspace(0,16,160
8、); yi=polyval(a,xi); plot(t,y,o,xi,yi, g),给出多项式的系数a,计算多项式函数的值,一般的曲线拟合:p=lsqcurvefit(Fun,p0,xdata,ydata)(xdata,ydata)是观测数据。对于这组观测数据我们选择了自认 是拟合效果比较好的函数形式f(x),其中参数以字母表示,取值待定 .我们把这个函数形式写入名为Fun的M文件.,例如:拟合函数为 y=ae-bx+ce-dx,编写M文件ex.mfunction y=ex(p, x) y= p(1)*exp( -p(2)*x)+ p(3)*exp( -p(4)*x);,自变量写为x,在lsq
9、curvefit命令中由xdata传递。待定参数写为 p(1) , p(2) , 我们对待定参数应有一个大致的估计,体现在拟合命令lsqcurvefit中的初始向量p0。,调用后返回的p就是按照最小二乘原则求得的待定参数。这时 再把p的分量对位代入函数形式的相应位置,就得到了完整的拟合函数。,minpsum( fun(p,xdata) ydata ).2,lsqcurvefit()命令的求解原理是,在所有 可能的 参数p 中挑选 使sum 最小的,函数法 则、待 定参数 p,观测到 的函数 值序列,观测到 的自变 量序列,例如: x=0:0.1:1; y= 4.0000 2.8297 2.01
10、83 1.4524 1.0550 0.7739 0.5733 0.4290 0.3242 0.2473 0.1903; 绘图认识观测数据体现的函数关系: plot(x,y),若要求解点x处的函数值可用程序f=ex(p,x)计算。 如 x=0:0.1:1; y=ex(p,x); plot(x,y),选择函数形式 y=ae-bx+ce-dx,编写M文件ex.m function y=ex(p, x) y=p(1)*exp(-p(2)*x)+p(3)*exp(-p(4)*x);,调用 p=lsqcurvefit(ex,1 2 1 4,x,y); 拟合函数为 y= p(1)*exp(-p(2)*x)+
11、p(3)*exp(-p(4)*x),x=0 0.3142 0.6283 0.9425 1.2566 1.57081.8850 2.1991 2.5133 2.8274 3.1416 ;,y=0 1.6180 1.9021 0.6180 -1.1756 -2.0000 -1.1756 0.6180 1.9021 1.6180 0.0000;,课堂练习:对下述观测数据,给出拟合函数,提示: 1.绘图估计解析关系 2.建立解析关系的M文件 3.调用lsqcurvefit命令 4.针对解析关系绘图,与观测点对比,评价拟合函数优劣。,5.6 实 验,5.6.1 实验一:拟合Malthus人口指数增长模型
12、中的参数,1790-1980年间美国每隔10年的人口记录如表5.3:,表 5.3,年份 1790 1800 1810 1820 1830 1840 1850 人口(*106) 3.9 5.3 7.2 9.6 12.9 17.1 23.2,年份 1860 1870 1880 1890 1900 1910 1920 人口(*106) 31.4 38.6 50.2 62.9 76.0 92.0 106.5 年份 1930 1940 1950 1960 1970 1980 人口(*106) 123.2 131.7 150.7 179.3 204.0 226.5,用以上数据检验马尔萨斯(Malthus)
13、人口指数增长模型(求解下面的微分方程),根据检验结果进一步讨论马尔萨斯人口模型的改进.,Malthus模型的基本假设是单位时间上的人口的增长率为常数,记为r.记时刻t的人口为x(t)(即x(t)为模型的状态变量),且初始时刻的人口为x0,于是得到如下微分方程,提示,需要先求微分方程的解,再用拟合模型的参数。,5.6.2 实验二:旧车价格预测,某年美国旧车价格的调查资料如表5.4,其中xi表示轿车的使用年数,yi表示相应的平均价格。试分析用什么形式的曲线来拟合上述的数据,并预测使用4.5年后轿车的平均价格大致为多少。,表5.4,xi 1 2 3 4 5 6 7 8 9 10 yi 2615 19
14、43 1494 1087 765 538 484 290 226 204,5.6.3 实验三:经济增长模型,增加生产、发展经济所依靠的主要动力有,在研究国民经济产值与这些因素的数量关系时,由于技术水平不像资金、劳动力那样容易定量化,作为初步的模型,可认为技术水平不变,只讨论产值和资金、劳动力之间的关系,即给出,增加投资 增加劳动力 以及技术革新等因素,产值= f (资金,劳动力),用 Q,K,L 分别表示产值、资金、劳动力. 即寻求,Q= Q (K,L),Cobb-Douglas生产函数,式中 要由经济统计数据确定.现有美国马萨诸塞州1900-1926年上述三个经济指数的统计数据,如表5.5,
15、试用数据拟合的方法,求出(*)式中的参数,(*),两边同乘以L得,首先分析一个劳动力分摊的资金K/L与他创造的产值Q/L。,一般化,就得到,通常Q/L是随着K/L的增加而递增的,但增速递减。因此 有,表5.5,由于(*)式对参数 是非线性的,因此,可以有两种方式进行拟合,一是直接使用MATLAB软件中的曲线或曲面拟合命令,另一个是将非线性 函数转化成线性函数的形式,使用线性函数拟合.,提示,参 考 文 献,1 郭锡伯,徐安农,高等数学实验课讲义,中国标准出版社,1998年. 2 李尚志等,数学建模竞赛教程,江苏教育出版社,1996年. 3 谭永基等,数学模型,复旦大学出版社,1998年.,关于
16、cobb-douglus生产模型的参数拟合,function f=stu12(p,KL) f=p(1)*(KL(1,:).p(2).*(KL(2,:).p(3);以stu12为名存盘退出,第一步 建立生产函数的M文件,(*),参照,k=1.04 1.06 1.16 1.22 1.27 1.37 1.44 1.53 1.57 2.05 2.51 2.63 2.74 2.82 3.24 3.24 3.61 4.10 4.36 4.77 4.75 4.54 4.54 4.58 4.58 4.58 4.54; l=1.05 1.08 1.18 1.22 1.17 1.30 1.39 1.47 1.31
17、 1.43 1.58 1.59 1.66 1.68 1.65 1.62 1.86 1.93 1.96 1.95 1.90 1.58 1.67 1.82 1.60 1.61 1.64; q=1.05 1.18 1.29 1.30 1.30 1.42 1.50 1.52 1.46 1.60 1.69 1.81 1.93 1.95 2.01 2.00 2.09 1.96 2.20 2.12 2.16 2.08 2.24 2.56 2.34 2.45 2.58;,第二步:建立主程序 Main.m,plot3(k,l,q,r+) hold on KL=k;l; p0=1 0.1 0.5; p=lsqcurvefit(stu12,p0,KL,q) q1=p(1)*k.p(2).*l.p(3); plot3(k,l,q1),注意:拟合效果较差。这是因为从,(*),引入独立参数 到,是一个纯数学的推广,没有直观的经济学意义。,% 空间散点图,请同学们以,及马萨诸塞州的数据K./L , Q./L去拟合参数,并观察 上述生产关系式的模型效果。,单位 劳动 力创 造的 产值,单位 劳动 力投 入的 资金,马萨诸塞州是美国发展历史最久的老工业基地,曾在二战后长期陷于严重的衰退之中。20世纪七八十年代,该州成功实现了经济转型,被称为“马萨诸塞奇迹”。其中勃兴的高技术产业的作用可谓居功至伟,