1、Mathematics Laboratory阮小娥博士Experiments in Mathematics赵小艳数学试验数学试验办公地址:理科楼办公地址:理科楼214第1页试验试验13 13 人口预测与数据拟合人口预测与数据拟合2 2、了解利用最小二乘法进行、了解利用最小二乘法进行数据拟合基本思想,掌握用数数据拟合基本思想,掌握用数据拟正当寻找最正确拟合曲线据拟正当寻找最正确拟合曲线方法。方法。3、了解多元函数极值在数、了解多元函数极值在数据拟正当中应用。据拟正当中应用。实验目1 1、学会用、学会用MATLABMATLAB软件进行数软件进行数据拟合。据拟合。4、经过对实际问题进行分、经过对实际
2、问题进行分析研究,初步掌握建立数据析研究,初步掌握建立数据拟合数学模型方法。拟合数学模型方法。第2页据人口统计年鉴,知我国从据人口统计年鉴,知我国从19491949年至年至19941994年人口数据资料以下:年人口数据资料以下:(人口数单位为:百万人口数单位为:百万)(1 1)在直角坐标系上作出人口数图象。)在直角坐标系上作出人口数图象。(2 2)建立人口数与年份函数关系,并估算)建立人口数与年份函数关系,并估算19991999年人年人口数。口数。试验问题试验问题年份年份19491954 1959 1964 1969人口数人口数 541.67602.66 672.09 704.99 806.7
3、1 年份年份 1974 1979 1984 1989 1994人口数人口数 908.59 975.42 1034.751106.761176.74 第3页怎样确定怎样确定a,b?线性模型线性模型第4页1 曲线拟合问题提法曲线拟合问题提法:已知一组(二维)数据,即平面上已知一组(二维)数据,即平面上n个点个点),(iiyx,ixni,2,1L=互不相同,寻求一个函数(曲线)互不相同,寻求一个函数(曲线))(xfy=,使使)(xf在观察点在观察点x x处所取得值处所取得值f(x)f(x)分别与观察值分别与观察值y y在某种在某种 xy0+一、曲线拟合一、曲线拟合准则下最为靠近,即曲线拟合得最好,如
4、图准则下最为靠近,即曲线拟合得最好,如图第5页从几何上讲,并不要求曲线严格经过已知点,但从几何上讲,并不要求曲线严格经过已知点,但要求曲线在各数据点和已知数据点之间总体误差要求曲线在各数据点和已知数据点之间总体误差最小,通常称为最小,通常称为数据拟合。数据拟合。到达最小。到达最小。最小二乘准则最小二乘准则 而我们经常是确定而我们经常是确定f(x)使得偏差平方和使得偏差平方和第6页数据插值数据插值已知一组(二维)数据,即平面上已知一组(二维)数据,即平面上n个点个点),(iiyx,ixni,2,1L=互不相同,寻求一个函数(曲线)互不相同,寻求一个函数(曲线))(xfy=数据插值数据插值第7页.
5、用什么样曲线拟合已知数据用什么样曲线拟合已知数据?惯用曲线函数系惯用曲线函数系ri(x)类型:类型:)画图观察)画图观察)理论分析)理论分析指数曲线:指数曲线:双曲线(一支双曲线(一支):):多项式:多项式:直线:直线:第8页比如比如指数函数拟合指数函数拟合三角函数拟合三角函数拟合多项式拟合多项式拟合第9页 拟合函数组中系数确定第10页4 4 用用matlabmatlab软件进行数据拟合软件进行数据拟合(1)lsqcurvefit命令命令-最小二乘拟合最小二乘拟合a=lsqcurvefit(fun,x0,xdata,ydata)a,resnorm=lsqcurvefit(fun,x0,xdat
6、a,ydata)是依据给定数据是依据给定数据xdata,ydata,按照函数文按照函数文件件funfun给定函数,以给定函数,以x0 x0为初值做最小二乘拟为初值做最小二乘拟合,返回函数中系数向量合,返回函数中系数向量a a和残差平方和和残差平方和resnorm。第11页例首先编写函数文件首先编写函数文件function y=f(a,x)f=a(1)*exp(x)+a(2)*x.2+a(3)*x.3保留为保留为f.mf.m,其次调用该函数,其次调用该函数x=0:0.1:1;y=3.1,3.27,3.81,4.5,5.18,6,7.05,8.56,9.69,11.25,13.17;a0=0 0
7、0;x,resnorm=lsqcurvefit(f,a0,x,y)第12页也能够用也能够用inlineinline命令定义函数命令定义函数x=0:0.1:1;y=3.1,3.27,3.81,4.5,5.18,6,7.05,8.56,9.69,11.25,13.17;f=inline(a(1)*exp(x)+a(2)*x.2+a(3)*x.3,a,x);a0=0 0 0;a,resnorm=lsqcurvefit(f,a0,x,y)plot(x,y,*)hold ong=a(1)*exp(x)+a(2)*x.2+a(3)*x.3;plot(x,g,r-)第13页 a=polyfit(xdata,
8、ydata,n)其中其中n n表示多项式最高阶数表示多项式最高阶数 xdata,ydata 为要拟合数据,它是用向量为要拟合数据,它是用向量方式输入。方式输入。输出参数输出参数a为拟合多项式为拟合多项式 y=anxn+a1x+a0系数系数a=an,a1,a0。多项式在多项式在x x处值处值y y可用下面程序计算。可用下面程序计算。y=polyval(a,x)因为高次多项式曲线改变不稳定,所以多项式次数选取不因为高次多项式曲线改变不稳定,所以多项式次数选取不宜过高宜过高。(2)(2)polyfit命令命令-多项式曲线拟合多项式曲线拟合第14页比如clear;clc;x=0:0.1:1;y=-0.
9、447 1.978 3.28 6.16 7.08 7.34 7.66,9.56,9.48,9.3,11.2;plot(x,y,k.,markersize,25);axis(0 1.3-2 16);p3=polyfit(x,y,3)p6=polyfit(x,y,6)t=0:0.01:1.2;s=polyval(p3,t);s1=polyval(p6,t);hold onplot(t,s,r-,linewidth,2);plot(t,s1,b-,linewidth,2);grid第15页第16页编写程序调用编写程序调用matlab命令命令x=1949:5:1994;y=541.67,602.66,
10、672.09,704.99,806.71,908.59,975.42,1034.75,1106.76,1176.74;plot(x,y,r*,linewidth,2)gridf=inline(a(1)+a(2)*x,a,x);a0=0 5;a,resnorm=lsqcurvefit(f,a0,x,y)hold ong=a(1)+a(2)*x;plot(x,g,b-,linewidth,2)二、人口预测线性模型二、人口预测线性模型第17页或者调用或者调用M函数函数function f=nihe(a,x)f=a(1)+a(2)*x;保留成保留成nihe.m,在新窗口编写程序在新窗口编写程序x=19
11、49:5:1994;y=541.67,602.66,672.09,704.99,806.71,908.59,975.42,1034.75,1106.76,1176.74;a0=10 10;a,resnorm=lsqcurvefit(nihe,a0,x,y)注意:该命令与初值相关系。注意:该命令与初值相关系。第18页也能够直接编写程序以下也能够直接编写程序以下:clc;clf;x=1949:5:1994;y=541.67,602.66,672.09,704.99,806.71,908.59,975.42,1034.75,1106.76,1176.74;plot(x,y,r*,linewidth,
12、2)grida11=10;a12=sum(x);a21=a12;a22=sum(x.2);d1=sum(y);d2=sum(x.*y);A=a11,a12;a21,a22;D=d1;d2;ab=inv(A)*Dplot(x,g,b-,linewidth,2)t=1949:5:;g=ab(1)+ab(2)*t;hold onplot(t,g,b-,linewidth,2)y=ab(1)+ab(2)*y=ab(1)+ab(2)*y=ab(1)+ab(2)*axis(1945 500 1450)plot(,1295.3,g*,linewidth,2)plot(,1306.28,g*,linewidt
13、h,2)plot(,1370.5,g*,linewidth,2)第19页计算得计算得从而得到人口数与年份函数关系为从而得到人口数与年份函数关系为线性预测模型线性预测模型年份年份预测(百万)1266.61339.11411.7真实值(百万)1295.331306.281370.5并预测2000,2005,年人口第20页仿真结果表明:仿真结果表明:线性模型在短线性模型在短期内基本上能期内基本上能比较准确地反比较准确地反应人口自然增应人口自然增加规律,但长加规律,但长久预测误差较久预测误差较大。大。第21页 英国统计学家英国统计学家MalthusMalthus于于17981798年提出了一个关于年提
14、出了一个关于生物种群繁殖生物种群繁殖指数增加模型指数增加模型:假设种群数量增加率:假设种群数量增加率与该时刻种群个体数量成正比。与该时刻种群个体数量成正比。三、人口预测三、人口预测Malthus模型模型基本假设基本假设 :人口人口(相对相对)增加率增加率 r 是常数是常数x(t)时刻时刻t 人口人口,t=0时人口数为时人口数为x0指数增加模型指数增加模型实际中,惯用实际中,惯用第22页解:解:取得最小值取得最小值.其中其中,表示人口数量表示人口数量。表示年份表示年份,解方程组解方程组:即得参数即得参数值值.使得使得问题转化为求参数问题转化为求参数第23页计算得计算得从而得到人口数与年份函数关系
15、为从而得到人口数与年份函数关系为指数预测模型指数预测模型年份年份预测(百万)1363.6 1488.8 1625.4真实值(百万)1295.331306.281370.5并预测2000,2005,年人口第24页仿真结果表明:仿真结果表明:人口增加指数人口增加指数模型在较短期模型在较短期内基本上能比内基本上能比较准确地反应较准确地反应人口自然增加人口自然增加规律,但长久规律,但长久预测误差很大。预测误差很大。第25页四、人口预测四、人口预测Logistic模型模型假如人口增加符合假如人口增加符合Malthus模型,则当模型,则当人口增加到一定数量后,增加率下降原因:人口增加到一定数量后,增加率下
16、降原因:资源、环境等原因对人口增加阻滞作用资源、环境等原因对人口增加阻滞作用 1838 1838年,荷兰生物学家年,荷兰生物学家VerhulstVerhulst对对MalthusMalthus模型作模型作了深入分析后指出:了深入分析后指出:造成上述不符合实际情况主要造成上述不符合实际情况主要原因是未能考虑原因是未能考虑“密度制约密度制约”原因。原因。即最终造成地球上人口爆炸,这与实际是不相符。即最终造成地球上人口爆炸,这与实际是不相符。且阻滞作用随人口数量增加而变大且阻滞作用随人口数量增加而变大r是是x减函数减函数第26页假设假设r固有增加率固有增加率(x很小时很小时)k人口容量(资源、环境能
17、容纳最大数量)人口容量(资源、环境能容纳最大数量)第27页结果显示:Logistic模型在长久预测时基本上能比较准确地反应人口自然增加规律,在2005,年较为吻合。第28页五、人口预测多项式模型五、人口预测多项式模型仿真结果表明仿真结果表明,人口模型用低阶多项式拟合能比人口模型用低阶多项式拟合能比较准确地反应人口自然增加规律,而高阶多项式较准确地反应人口自然增加规律,而高阶多项式拟合预测效果很差。拟合预测效果很差。青色青色-二次多项式二次多项式蓝色蓝色三次多项式三次多项式紫红色紫红色-四次多项四次多项式式第29页例例2:2:海底光缆线长度预测模型海底光缆线长度预测模型某某一一通通信信企企业业在
18、在一一次次施施工工 中中,需需要要在在水水面面宽宽为为2 20 0m m河河沟沟底底沿沿直直线线走走向向铺铺设设一一条条沟沟底底光光缆缆.在在铺铺设设光光缆缆之之前前需需要要对对沟沟底底地地形形做做初初B2468101214161820986420ADC探测到一组等分点位置深度数据以下表所表示探测到一组等分点位置深度数据以下表所表示.25步探测步探测,从而预计所需光缆长度从而预计所需光缆长度,为工程预算提为工程预算提供依据供依据.基本情况如图所表示基本情况如图所表示.第30页10.9310.809.818.867.957.959.1510.2211.2912.6113.322019181716
19、15141312111013.2812.2611.1810.139.058.027.967.968.969.01深度深度(m)9876543210分点分点2121个等分点处深度个等分点处深度(1)(1)预测经过这条河沟所需光缆长度近似值预测经过这条河沟所需光缆长度近似值.(2)(2)作出铺设沟底光缆曲线图作出铺设沟底光缆曲线图.第31页解:解:用用12次多项式函数拟合光缆走势曲线图以下次多项式函数拟合光缆走势曲线图以下仿真结果表仿真结果表明明,拟合曲拟合曲线能较准确线能较准确地反应光缆地反应光缆走势图走势图.The length of the label is L=26.3809(m)假设所铺
20、设光缆足够柔软假设所铺设光缆足够柔软,在铺设过程中光缆触地走在铺设过程中光缆触地走势光滑势光滑,紧贴地面紧贴地面,而且忽略水流对光缆冲击而且忽略水流对光缆冲击.第32页%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,100);P
21、=9.01,8.96,7.96,7.97,8.02,9.05,10.13,11.18,12.26,13.2P=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.88,13.32,12.61,11.29,10.22,9.15,7.90,7.95,8.86,9.81,10.80,10.93;0,10.93;a,s=polyfit(t,P,12);a,s=polyfit(t,P,12);yy=polyval(a,x);yy=polyval(
22、a,x);disp(yy=);disp(yy);disp(yy=);disp(yy);plot(x,yy,r*-,t,P,b+-);plot(x,yy,r*-,t,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(The length of the label is L=);disp(L);disp(The length of the label is L=);disp(L);第33页上机任务上机任务 李继成(书)李继成(书)Page 194 练习练习1 第第1、2题题 Page 198 练习练习2 第第2、3题题 学号末尾是单号者做第学号末尾是单号者做第1个题,双号者做第个题,双号者做第2个题个题第34页第35页