数学实验8-曲线拟合

1、Mathematics Laboratory阮小娥博士Experiments in Mathematics李换琴数学实验数学实验办公地址：理科楼办公地址：理科楼225电话：电话：82663174实验实验13 13 人口数量预测模型实验人口数量预测模型实验2 2、掌握在最小二乘意、掌握在最小二乘意义下数据拟合的理论和义下数据拟合的理论和方法方法. .1 1、学会用、学会用MATLABMATLAB软件软件进行数据拟合进行数据拟合3、通过对实际问题的、通过对实际问题的分析和研究，初步掌分析和研究，初步掌握建立数据拟合数学握建立数据拟合数学模型的方法模型的方法实验目的实验目的据人口统计年鉴，知我国从据

2、人口统计年鉴，知我国从19491949年至年至19941994年人口数据资料如下：年人口数据资料如下： ( (人口数单位为：百万人口数单位为：百万) )（1 1）在直角坐标系上作出人口数的图象。）在直角坐标系上作出人口数的图象。（2 2）建立人口数与年份的函数关系，并估算）建立人口数与年份的函数关系，并估算19991999年年的人口数。的人口数。实验问题实验问题年份年份19491954 1959 1964 1969人口数人口数 541.67602.66 672.09 704.99 806.71 年份年份 1974 1979 1984 1989 1994人口数人口数 908.59 975.42

3、1034.751106.761176.74 如何确定如何确定a,b?baxy 线性模型线性模型1 曲线拟合问题的提法曲线拟合问题的提法: 已知一组（二维）数据，即平面上的已知一组（二维）数据，即平面上的n 个点个点),(iiyx， ixni, 2 , 1L互不相同，寻求一个函数（曲线）互不相同，寻求一个函数（曲线）)(xfy ，使使)(xf在在某种准则下某种准则下与所有数据点最为接近，即曲线拟合与所有数据点最为接近，即曲线拟合得最好，如图得最好，如图: xy0+i ),(iiyx)(xfy 一、曲线拟合一、曲线拟合 niiiniiyxf1212)( 确定确定f(x)使得使得 达到最小达到最小

4、最小二乘准则最小二乘准则 . 用什么样的曲线拟合已知数据用什么样的曲线拟合已知数据?常用的曲线函数系类型：常用的曲线函数系类型：画图观察；画图观察；理论分析理论分析xaeay21 指数曲线：指数曲线： 21axay 双曲线（一支双曲线（一支):): 11 mmmaxaxayL多项式：多项式： 21axay 直线：直线： 拟合函数组中系数的确定拟合函数组中系数的确定达到最小。达到最小。 niiiniiyaxaaaJ12211221),( 21,aa,)(21为例为例以以axaxf 使得使得，即确定即确定21aa 0) 20) 2121121niiiiniiiyaxaxyaxa二、人口预测线性模型

5、二、人口预测线性模型对于开始提出的实验问题对于开始提出的实验问题, 代如数据，计算得代如数据，计算得27754,1521 aa从而得到人口数与年份的函数关系为从而得到人口数与年份的函数关系为2775415 xy把把x=1999代如，估算出代如，估算出1999年的人口数为年的人口数为 y=1252.1（百万）（百万）12.52亿亿1999年实际人口数量为年实际人口数量为.亿。亿。线性预测模型线性预测模型 英国人口学家英国人口学家MalthusMalthus根据百余年的人口统计资根据百余年的人口统计资料，于料，于17981798年提出了著名的人口自然增长的年提出了著名的人口自然增长的指数增指数增长

6、模型长模型。 btaetx 三、人口预测的三、人口预测的Malthus模型模型 0)0(,xxrxdtdxtrextx0)( 基本假设基本假设 : : 人口人口( (相对相对) )增长率增长率 r 是常数是常数x(t) 时刻时刻t 的的人口人口,t=0时人口数为时人口数为x0 btatx ln指数增长模型指数增长模型实际中，常用实际中，常用1. 1. 由前由前100100年的数据求出美国的人口增长年的数据求出美国的人口增长Malthus模型模型。2. 2. 预测后预测后100100年（每隔年（每隔1010年）的人口状况。年）的人口状况。3. 3. 根据预测的人口状况和实际的人口数量根据预测的人

7、口状况和实际的人口数量, ,讨论人讨论人口模型的改进情况。口模型的改进情况。美国美国17901790年年19801980年每隔年每隔1010年的人口记录年的人口记录226.5204.0179.3150.7131.7123.2106.592.076.062.9人口人口( (百万百万) )1980197019601950194019301920191019001890年份年份50.238.631.423.217.112.99.67.25.33.9人口人口( (百万百万) )1880187018601850184018301820181018001790年份年份例例解：解： btaetx 设设bta

8、x )ln(2101)ln(),(iiixbtabaJ 取得最小值取得最小值. .其中其中, , iitxx 表示人口数量表示人口数量。it表示年份表示年份, ,解方程组解方程组: 0ln210ln10211011012101101101iiiiiiiiiiitxbtatbJxbtaaJ即得参数即得参数ba,的值的值.ba,使得使得问题转化为求参数问题转化为求参数 prog41.m %prog41.m % % This program is to predict the number of This program is to predict the number of population

9、%population %format longformat longt1=1790;1800;1810;1820;1830;1840;1850;1860;187t1=1790;1800;1810;1820;1830;1840;1850;1860;1870;1880;0;1880;t2=1890;1900;1910;1920;1930;1940;1950;1960;197t2=1890;1900;1910;1920;1930;1940;1950;1960;1970;1980;0;1980;x1=3.9;5.3;7.2;9.6;12.9;17.1;23.2;31.4;38.6;50.x1=3.9

10、;5.3;7.2;9.6;12.9;17.1;23.2;31.4;38.6;50.2;2;x2=62.9;76.0;92.0;106.5;123.2;131.7;150.7;179.x2=62.9;76.0;92.0;106.5;123.2;131.7;150.7;179.3;204.0;226.5;3;204.0;226.5;lnx1=log(x1); lnx2=log(x2);lnx1=log(x1); lnx2=log(x2);a12=sum(t1);a11=10;a21=a12;a22=sum(t1.2); d1=sum(lnx1);d2=sum(lnx1.*t1);A=a11,a12

11、;a21,a22;D=d1;d2;ab=inv(A)*D;disp(a=);disp(ab(1);disp(b=);disp(ab(2);forfor i=1:10 i=1:10 xx1(i)=exp(ab(1)+ab(2) xx1(i)=exp(ab(1)+ab(2)* *t1(i);t1(i);endendforfor i=1:10 i=1:10 xx2(i)=exp(ab(1)+ab(2) xx2(i)=exp(ab(1)+ab(2)* *t2(i);t2(i);endendplot(t1,x1,rplot(t1,x1,r* *-,t1,xx1,b+-, -,t1,xx1,b+-, t2

12、,x2,gt2,x2,g* *-,t2,xx2,m+-);-,t2,xx2,m+-);a= -49.79535457790735b=0.02859807120038仿真结果表明：仿真结果表明： 人口增加的指人口增加的指数模型在短期数模型在短期内基本上能比内基本上能比较准确地反映较准确地反映人口自然增长人口自然增长的规律，但长的规律，但长期预测误差很期预测误差很大，需要修正大，需要修正预测模型。预测模型。拟合曲线拟合曲线原始数据曲线原始数据曲线四、人口预测的四、人口预测的Logistic模型模型人口增长到一定数量后，增长率下降的原因：人口增长到一定数量后，增长率下降的原因：资源、环境等因素对人口

13、增长的阻滞作用资源、环境等因素对人口增长的阻滞作用且阻滞作用随人口数量增加而变大且阻滞作用随人口数量增加而变大假设假设r固有增长率固有增长率(x很小时很小时)k人口容量（资源、环境能容纳的最大数量）人口容量（资源、环境能容纳的最大数量）)1()(kxrxr r是是x的减函数的减函数krs 0)( kr)1()(mxxrxxxrdtdx 000)()(xexkkxtxtr 例的例的LogisticLogistic模型留给同学们练习模型留给同学们练习五、多项式拟合的五、多项式拟合的Matlab指令指令a = polyfit（xdata，ydata，n）其中其中n n表示多项式的最高阶数表示多项式的

14、最高阶数 xdata，ydata 为要拟合的数据，它是为要拟合的数据，它是用向量的方式输入。用向量的方式输入。输出参数输出参数a为拟合多项式为拟合多项式 y = a1xn + + anx + an+1的的系数系数a = a1, , an, an+1。多项式在多项式在x x处的值处的值y y可用下面程序计算。可用下面程序计算。 y = polyval (a, x) 用多项式拟合人口模型用多项式拟合人口模型% % This program is to predict the model of population by This program is to predict the model of

15、 population by 4-degree polynomial% 4-degree polynomial% % %prog42.m%prog42.m%format longformat longt1=1790;1800;1810;1820;1830;1840;1850;1860;1870;1880;t1=1790;1800;1810;1820;1830;1840;1850;1860;1870;1880;t2=1890;1900;1910;1920;1930;1940;1950;1960;1970;1980;t2=1890;1900;1910;1920;1930;1940;1950;196

16、0;1970;1980;t=t1;t2;t=t1;t2;P1=3.9;5.3;7.2;9.6;12.9;17.1;23.2;31.4;38.6;50.2;P1=3.9;5.3;7.2;9.6;12.9;17.1;23.2;31.4;38.6;50.2;P2=62.9;76.0;92.0;106.5;123.2;131.7;150.7;179.3;204.0;2P2=62.9;76.0;92.0;106.5;123.2;131.7;150.7;179.3;204.0;226.5;26.5;P=P1;P2;P=P1;P2;n=4; n=4; % The degree of the fitting

17、polynomial% The degree of the fitting polynomial% a,sa,s=polyfit(t1,P1,n);=polyfit(t1,P1,n);y=y=polyval(a,tpolyval(a,t););% % a a is the coefficients vector from n-degree to 0- is the coefficients vector from n-degree to 0-degree%degree%plot(t,P,rplot(t,P,r* *-,-,t,y,bt,y,b+-);+-);23a =1.0e+006 * -0

18、.00000000000014 0.00000000107892 -0.00000304878595 0.00381927346813 -1.79012132225427仿真结果表明仿真结果表明, , 人口增加的模型用多项式拟合能人口增加的模型用多项式拟合能比较准确地反映人口自然增长的规律，对长期预比较准确地反映人口自然增长的规律，对长期预测具有指导意义。测具有指导意义。例例2: 2: 海底光缆线长度预测模型海底光缆线长度预测模型某一通信公司在一次某一通信公司在一次施工中施工中, ,需要在水面宽需要在水面宽为为2020m m的河沟底沿直线的河沟底沿直线走向铺设一条沟底光走向铺设一条沟底光缆缆.

19、 .在铺设光缆之前需在铺设光缆之前需要对沟底的地形做初要对沟底的地形做初Bix2468101214161820986420ADCih0 x20 x探测到一组等分点位置的深度数据如下表所示探测到一组等分点位置的深度数据如下表所示. .25步探测步探测, ,从而估计所需光缆的长度从而估计所需光缆的长度, ,为工程预算为工程预算提供依据提供依据. .基本情况如图所示基本情况如图所示. .10.9310.809.818.867.957.959.1510.2211.2912.6113.32201918171615141312111013.2812.2611.1810.139.058.027.967.96

20、8.969.01深度深度(m)9876543210分点分点2121个等分点处的深度个等分点处的深度(1) (1) 预测通过这条河沟所需光缆长度的近似值预测通过这条河沟所需光缆长度的近似值. .(2) (2) 作出铺设沟底光缆的曲线图作出铺设沟底光缆的曲线图. .024681012141618207891011121314解：解： 用用12次多项式函数拟合光缆走势的曲线图如下次多项式函数拟合光缆走势的曲线图如下仿真结果表仿真结果表明明,拟合曲线拟合曲线能较准确地能较准确地反映光缆的反映光缆的走势图走势图.The length of the label is L= 26.3809 (m)假设所铺设

21、的光缆足够柔软假设所铺设的光缆足够柔软,在铺设过程中光缆触地在铺设过程中光缆触地走势光滑走势光滑,紧贴地面紧贴地面,并且忽略水流对光缆的冲击并且忽略水流对光缆的冲击.% % prog45.m This program is to fit the data by polynomial %prog45.m This program is to fit the data by polynomial %format longformat longt=linspace(0,20,21);t=linspace(0,20,21);x=linspace(0,20,100);x=linspace(0,20,10

22、0);P=9.01,8.96,7.96,7.97,8.02,9.05,10.13,11.18,12.26,13.P=9.01,8.96,7.96,7.97,8.02,9.05,10.13,11.18,12.26,13.28,13.32,12.61,11.29,10.22,9.15,7.90,7.95,8.86,9.81,10.28,13.32,12.61,11.29,10.22,9.15,7.90,7.95,8.86,9.81,10.80,10.93;80,10.93; a,sa,s=polyfit(t,P,12);=polyfit(t,P,12);yyyy= =polyval(a,xpoly

23、val(a,x););plot(x,yy,rplot(x,yy,r* *-,-,t,P,bt,P,b+-);+-);L=0;L=0;forfor i=2:100 i=2:100 L=L+sqrt(x(i)-x(i-1)2+(yy(i)-yy(i-1)2); L=L+sqrt(x(i)-x(i-1)2+(yy(i)-yy(i-1)2);endenddisp(Thedisp(The length of the label is L=); length of the label is L=);disp(Ldisp(L););format longt=linspace(0,20,21);x=linsp

24、ace(0,20,100);P=9.01,8.96,7.96,7.97,8.02,9.05,10.13,11.18,12.26,13.28,13.32,12.61,11.29,10.22,9.15,7.90,7.95,8.86,9.81,10.80,10.93;n=input(n=) %通过键盘输入拟合次数a,s=polyfit(t,P,n);yy=polyval(a,x);p1=polyval(a,t);d=norm(P-p1) %计算拟合误差plot(x,yy,r*-,t,P,b+-);L=0;for i=2:100 L=L+sqrt(x(i)-x(i-1)2+(yy(i)-yy(i-1)

25、2);enddisp(The length of the label is L=);disp(L);x0 x1x2x3x4xg(x)(xfy 函函数数：:在一些点处测得函数值在一些点处测得函数值nixfyiiL2 , 1 , 0 )( 使使得得求求一一个个简简单单易易算算的的函函数数),(xgniyxgiiL2 , 1 , 0 )( ),()(,0 xfxgxxn 上上在在.)()(的的插插值值函函数数称称为为xfxg实验实验14 插值问题插值问题,表表达达式式未未知知或或非非常常复复杂杂 f(x)w插值函数有各种类型，如代数多项式，插值函数有各种类型，如代数多项式，三角函数，有理函数等。三角

26、函数，有理函数等。w当当p(x)为多项式时，称为（代数）插为多项式时，称为（代数）插值多项式。值多项式。一阶一阶 二阶二阶 三阶三阶nnyyyxfyxxxxLL1010)( ，设设nnnxaxaaxp L10)()1 . 1(2 , 1 , 0 )(niyxpiinL 满满足足,1,)(10处处的的函函数数值值个个互互异异点点上上在在已已知知nxxxnbaxfL 由克莱姆法则知方程组有唯一解，即满足由克莱姆法则知方程组有唯一解，即满足(1.1)的插值多项式存在且唯一。的插值多项式存在且唯一。Matlab指令：指令：yb=interp1(x,y,xb,method)详见课本详见课本204页页实验任务实验任务:观测观测序号序号12345678910X46495152545657585960Y40505563727077739093观测观测序号序号11121314151617181920X61626364666768717271Y9688991101131201271371321371、下表中，、下表中，X是华氏温度，是华氏温度，Y是一分钟内一只蟋是一分钟内一只蟋蟀的鸣叫次数，试用多项式模型拟合这些数据，蟀的鸣叫次数，试用多项式模型拟合这些数据，画出拟合曲线，分析你的拟合模型是否很好？画出拟合曲线，分析你的拟合模型是否很好？2、(1)