数据建模常规方法的Matlab实现课件

,*,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,Matlab在数学建模中的应用II,数据常规方法的Matlab实现,Matlab在数学建模中的应用II数据常规方法的Matl,1,1.1 数据的读入与写出,1.2 数据拟合方法,1.3 数据拟合应用实例,1.4 数据的可视化,数据建模常规方法的MATLAB实现,1.1 数据的读入与写出数据建模常规方法的MATLAB实,2,1.1,数据的读入与写出,数学建模不可避免地要用到大量的数据，,最简单的方法,是复制、粘贴，但是不方便。,另一种方法,是与Excel和记事本（*.dat或者*.txt的文件）进行交互。,（算 法）,程 序,数据读入,输出,图形,数据,1.1 数据的读入与写出 数学建模不可避免地,3,1.1.1 Excel与MATLAB的数据交互,首先要安装Excel和MATLAB。,第一，打开Excel的工具宏安全性安全级（中）,第二，打开Excel的工具加载宏浏览,安装MATLAB的目录toolbox exlink excllink.xla,确定，得到如下的工具条（即可使用）：,1.1.1 Excel与MATLAB的数据交互,4,1.1.2 记事本与MATLAB的数据交互,%从记事本t.txt中读取数据,load(filename.*) %记事本中记录的全是数据,name,type, x,y,answer,=textread(,t.txt,%s Type %n %f %n %s,2,),%将Matlab数据写入记事本,fid= fopen(tp.txt,wt);,%文件扩展名可以为：*.dat或*.xls,fprintf(fid,This is the database of class 1.n);,name=Sally;types=1;x=3.1;y=45;answer=Yes;,fprintf(fid,% s Type %u %f %u %s n,name,types,x,y,answer);,name=Tom;types=2;x=2.5;y=20;answer=No;,fprintf(fid,% s Type %u %f %u %s n,name,types,x,y,answer);,fclose(fid);,每一列要保存的变量名,保存的文件名,读取的格式,读取2行数据,%从记事本t.txt中读取数据,load(filename.*) %记事本中记录的全是数据,name,type, x,y,answer,=textread(,t.txt,%s Type %n %f %n %s,2,),%从记事本t.txt中读取数据,load(filename.*) %记事本中记录的全是数据,name,type, x,y,answer,=textread(,t.txt,%s Type %n %f %n %s,2,),%从记事本t.txt中读取数据,load(filename.*) %记事本中记录的全是数据,name,type, x,y,answer,=textread(,t.txt,%s Type %n %f %n %s,2,),%从记事本t.txt中读取数据,load(filename.*) %记事本中记录的全是数据,name,type, x,y,answer,=textread(,t.txt,%s Type %n %f %n %s,2,),1.1.2 记事本与MATLAB的数据交互%从记事本t.tx,5,1.2,数据拟合方法,1.2.1 多项式拟合,一般情况下，数据点较少的用,插值法,。数据点较多的，只需要考察数据的总体变化趋势的，用,拟合法,。最常用的拟合方法是,最小二乘拟合法,。,1.多项式拟合指令,polyfit(X,Y,N)：多项式拟合，返回降幂排列的多项式系数。,polyval(P,xi)：计算多项式的值。,说明：X,Y是数据点的值，N是拟合的最高次幂，P是返回的多项式系数，xi是要求的点的横坐标。,1.2 数据拟合方法1.2.1 多项式拟合,6,例1：对所给数据进行多项式拟合,x,1,2,3,4,5,6,7,8,9,y,9,7,6,3,-1,2,5,7,20,解:Matlab运行程序:,x=1 2 3 4 5 6 7 8 9;,y=9 7 6 3 -1 2 5 7 20;,P=polyfit(x,y,3),xi=0:0.2:10;,yi=polyval(P,xi);,plot(xi,yi,x,y,r*);,运行结果:P = 0.1481 -1.4030 1.8537 8.2698(降幂排),例1：对所给数据进行多项式拟合x123456789y9763,7,画出数据点：, x=1 2 3 4 5 6 7 8 9;y=9 7 6 3 -1 2 5 7 20;plot(x,y,r*);,在图形窗口中单击Tools-Basic Fitting，打开对话框。,2.图形窗口的多项式拟合,画出数据点： 2.图形窗口的多项式拟合,8,数据建模常规方法的Matlab实现课件,9,例2：在某次阻尼振荡实验中测得18组数据点。,1.2.2 指定函数拟合,x,0,0.4,1.2,2,2.8,3.6,4.4,5.2,6,y,1,0.85,0.29,-0.27,-0.53,-0.4,-0.12,0.17,0.28,x,7.2,8,9.2,10.4,11.6,12.4,13.6,14.4,15,y,0.15,-0.03,-0.15,-0.071,0.059,0.08,0.032,-0.015,-0.02,解,：(1),画出散点图：,x=0;0.4;1.2;2;2.8;3.6;4.4;5.2;6;7.2;8;9.2;10.4;11.6;12.4;13.6;14.4;15;,y=1;0.85;0.29;-0.27;-0.53;-0.4;-0.12;0.17;0.28;0.15;-0.03;-0.15;-0.071;0.059;0.08;0.032;-0.015;-0.02;,plot(x,y,r*),例2：在某次阻尼振荡实验中测得18组数据点。,x,0,0.4,1.2,2,2.8,3.6,4.4,5.2,6,y,1,0.85,0.29,-0.27,-0.53,-0.4,-0.12,0.17,0.28,x,7.2,8,9.2,10.4,11.6,12.4,13.6,14.4,15,y,0.15,-0.03,-0.15,-0.071,0.059,0.08,0.032,-0.015,-0.02,例2：在某次阻尼振荡实验中测得18组数据点。1.2.2 指定,10,(2)对应函数形式为 ,进行拟合.,程序如下: syms t,x=0;0.4;1.2;2;2.8;3.6;4.4;5.2;6;7.2;8;9.2;10.4;11.6;12.4;13.6;14.4;15;,y=1;0.85;0.29;-0.27;-0.53;-0.4;-0.12;0.17;0.28;0.15;-0.03;-0.15;-0.071;0.059;0.08;0.032;-0.015;-0.02;,f=fittype(a*cos(k*t)*exp(w*t),independent,t,coefficients,a,k,w);,cfun=fit(x,y,f) %显示拟合函数,xi=0:.1:20;,yi=cfun(xi);,plot(x,y,r*,xi,yi,b-),说明:fittype是自定义拟合函数,cfun=fit(x,y,f)是根据自定义的拟合函数f来拟合数据x,y.此处数据必须是列向量的形式.,(2)对应函数形式为,11,MATLAB主窗口左下角：,starttoolboxesCurveFittingCurve Fitting Tool (cftool),1.2.3 曲线拟合工具箱,MATLAB主窗口左下角：1.2.3 曲线,12,练习:,1、有一只对温度敏感的电阻，已经测得一组温度t和电阻R数据：,20.5,32.7,51.0,73.0,95.7,765,826,873,942,1032,拟合引例中电阻R与温度t之间的关系R=at+b。,2、x为0到1之间的数，间隔0.1，y为-0.447，1.978，3.28，6.16，7.08，7.34，7.66，9.56，9.48，9.30，11.2，,分别用二次、五次和十次拟合曲线来拟合这组数据，并进行比较。,练习:1、有一只对温度敏感的电阻，已经测得一组温度t和电阻R,13,3、已知在某实验中测的某质点的位移和速度随时间的变化如下，求质点的速度与位移随时间的变化曲线以及位移随速度的变化曲线。,t,0,0.5,1.0,1.5,2.0,2.5,3.0,v,0,0.4794,0.8415,0.9975,0.9093,0.5985,0.1411,s,1,1.5,2,2.5,3,3.5,4,3、已知在某实验中测的某质点的位移和速度随时间的变化如下，求,14,R与t大致呈直线关系，即,1、解答：,R与t大致呈直线关系，即1、解答：,15,数据建模常规方法的Matlab实现课件,16,2、Matlab运行程序：,2、Matlab运行程序：,17,从图形上看，次数越高拟合程度越好。,从图形上看，次数越高拟合程度越好。,18,3、Matlab运行程序：,3、Matlab运行程序：,19,1.3,数据拟合应用实例,1.3.1 人口预测模型,某地区1971-2000年人口数据，给出该地区人口增长的数学模型。,根据所给数据作出散点图。人口随时间的变化是非线性的，存在一条与x轴平行的渐近线，因此用Logistic曲线模型进行拟合。,1.3 数据拟合应用实例1.3.1 人口预测模型,20,1.3.2 薄膜渗透率的测定,1.3.2 薄膜渗透率的测定,21,数据建模常规方法的Matlab实现课件,22,数据建模常规方法的Matlab实现课件,23,数据建模常规方法的Matlab实现课件,24,下面用Matlab软件进行计算。,下面用Matlab软件进行计算。,25,1.,lsqcurvefit,已知,数据点,：,xdata=（xdata,1,xdata,2,xdata,n,）,ydata=（ydata,1,ydata,2,ydata,n,）,补充：用,MATLAB,作非线性最小二乘拟合,MATLAB,提供了两个求非线性最小二乘拟合的函数：,lsqcurvefit,和,lsqnonlin,两个命令都要先建立M文件,fun.m,，在其中定义函数,f,(,x,),，但两者定义,f,(,x,)的方式是不同的。,lsqcurvefit,用以求含参量,x,（向量）的向量值函数,F(x,xdata)=（F（x,xdata,1,）,F（x，xdata,n,）,T,中的参变量x(向量),使得,1. lsqcurvefit补充：用MATLAB作非线性最小,26,fun是一个事先建立的定义函数,F(x,xdata) 的M文件,自变量为,x和xdata,说明:,x=lsqcurvefit(fun,x0,xdata,ydata,options);,迭代初值,已知数据点,所有优化函数,fun是一个事先建立的定义函数F(x,xdata) 的M文件,27,lsqnonlin,用以求含参量,x,（向量）的向量值函数,f,(,x,)=(,f,1,(,x,),f,2,(,x,),f,n,(,x,),T,中的参量,x,，使得,最小,其中,f,i,（,x,）=,f,（,x,，,xdata,i,，,ydata,i,）,=,F,(,x,xdata,i,)-,ydata,i,2.,lsqnonlin,已知数据点：,xdata=（xdata,1,，xdata,2,，xdata,n,）,ydata=（ydata,1,，ydata,2,，ydata,n,）,lsqnonlin用以求含参量x（向量）的向,28,输入格式为：,1),x=,lsqnonlin,（fun，x0）,；,2),x=lsqnonlin（fun，x0，options）,；,3),x= lsqnonlin（fun，x0，optionsgrad）,；,4),x，options=lsqnonlin （fun，x0，）,；,5),x，options，funval=lsqnonlin（funx0，）,；,说明：,x=,lsqnonlin,（fun，x0，options）,；,fun,是一个事先建立的定义函数,f,(,x,)的M文件，,自变量为,x,迭代初值,所有优化函数,输入格式为：说明：x= lsqnonlin （fun，x,29,MATLAB(liti21),1）编写M文件,curvefun1.m,function,f=curvefun1(,x,tdata,),f=x(1)+x(2)*exp(-0.02*x(3)*tdata),;,%,其中 x(1)=a; x(2)=b；x(3)=k;,2）输入命令,tdata=100:100:1000,;,cdata=1e-03*4.54,4.99,5.35,5.65,5.90,6.10,6.26,6.39,6.50,6.59;,x0=0.2,0.05,0.05;,x=lsqcurvefit (curvefun1,x0,tdata,cdata),f= curvefun1(x,tdata),plot(tdata,cdata,o,tdata,f,r-),F,(,x,，,tdata,)= ，,x,=(,a,，,b,，,k,),解法1,.,用命令,lsqcurvefit,MATLAB(liti21) 1）编写M文件 curvef,30,3）运算结果为,：,f =0.0043 0.0051 0.0056 0.0059 0.0061,0.0062 0.0062 0.0063 0.0063 0.0063,x = 0.0063 -0.0034 0.2542,4）结论,：,a,=0.0063,b,=-0.0034,k,=0.2542,3）运算结果为：4）结论：a=0.0063, b=-0.00,31,MATLAB(liti22),解法,2,：,用命令,lsqnonlin,f,(,x,)=,F,(x,tdata,ctada,)=,x,=（,a,，,b,，,k,）,1）,编写M文件,curvefun2.m,function,f=curvefun2(,x,),tdata=100:100:1000;,cdata=1e-03*4.54,4.99,5.35,5.65,5.90,6.10,6.26,6.39,6.50,6.59;,f=x(1)+x(2)*exp(-0.02*x(3)*tdata)- cdata;,2）输入命令:,x0=0.2,0.05,0.05;,x=lsqnonlin(,curvefun2,x0),f= curvefun2(x),函数,curvefun2,的自变量是,x,，,cdata,和,tdata,是已知参数，故应将,cdata tdata,的值写在,curvefun2.m,中,MATLAB(liti22) 解法2 ： 用命令lsq,32,3）运算结果为,f,=1.0e-003 *(0.2322 -0.1243 -0.2495 -0.2413 -0.1668 -0.0724 0.0241 0.1159 0.2030 0.2792,x =0.0063 -0.0034 0.2542,可以看出,两个命令的计算结果是相同的,.,4）结论,：,即拟合得,a,=,0.0063,b,=-0.0034,k,=0.2542,3）运算结果为可以看出,两个命令的计算结果是相同的.4）,33,背景,年 1625 1830 1930 1960 1974 1987 1999,人口(亿) 5 10 20 30 40 50 60,世界人口增长概况,中国人口增长概况,年 1908 1933 1953 1964 1982 1990 1995 2000,人口(亿) 3.0 4.7 6.0 7.2 10.3 11.3 12.0 13.0,研究人口变化规律,控制人口过快增长,1.3.3,如何预报人口的增长,背景 年 1625 1830,34,指数增长模型马尔萨斯提出 (,1798,),常用的计算公式,x,(,t,) ,时刻,t,的,人口,基本假设,: 人口(相对)增长率,r,是常数,今年人口,x,0, 年增长率,r,k,年后人口,随着时间增加，人口按指数规律无限增长.,与常用公式的一致,rt,e,x,t,x,0,),(,=,?,指数增长模型马尔萨斯提出 (1798)常用的计算公式x(,35,指数增长模型的应用及局限性,与19世纪以前欧洲一些地区人口统计数据吻合.,适用于19世纪后迁往加拿大的欧洲移民后代.,可用于短期人口增长预测.,不符合19世纪后多数地区人口增长规律.,不能预测较长期的人口增长过程.,19世纪后人口数据,人口增长率,r,不是常数(逐渐下降),指数增长模型的应用及局限性 与19世纪以前欧洲一些地区人口统,36,阻滞增长模型(,Logistic,模型),人口增长到一定数量后，增长率下降的原因：,资源、环境等因素对人口增长的阻滞作用,且阻滞作用随人口数量增加而变大,假设,r,固有增长率(,x,很小时),x,m,人口容量（资源、环境能容纳的最大数量）,r,是,x,的减函数,阻滞增长模型(Logistic模型)人口增长到一定数量后，增,37,dx,/,dt,x,0,x,m,x,m,/2,x,m,t,x,0,x,(,t,),S,形曲线,x,增加先快后慢,x,0,x,m,/2,阻滞增长模型(,Logistic,模型),指数增长模型,dx/dtx0xmxm/2xmtx0x(t)S形曲线, x,38,参数估计,用指数增长模型或阻滞增长模型作人口,预报，必须先估计模型参数,r,或,r, x,m,.,根据统计数据利用,线性最小二乘法,作拟合,阻滞增长模型(,Logistic,模型),例：美国人口数据(百万),1860 1870 1880 1960 1970 1980 1990 2000,31.4 38.6 50.2 179.3 204.0 226.5 251.4 281.4,t,x,数据(,t,x,),数据(,x,y,),用最小二乘法估计,r,s,r,x,m,参数估计用指数增长模型或阻滞增长模型作人口根据统计数据利用线,39,模型检验,用模型计算2000年美国人口,误差不到3%,阻滞增长模型(,Logistic,模型),r,=0.2557,x,m,=392.1,用美国18601990年数据(去掉个别异常数据),与实际数据(2000年为,281.4,)比较,1790年为零点,=,274.5,模型检验用模型计算2000年美国人口误差不到3%阻滞增长模型,40,Logistic,模型的应用,模型应用,加入2000年人口数据后重新估计模型参数,r,=0.2490,x,m,=434.0,x,(2010)=306.0,预报美国,2010,年的人口,经济领域中的增长规律(耐用消费品的售量).,种群数量模型 (鱼塘中的鱼群, 森林中的树木).,预报人口的增长,指数增长模型,阻滞增长模型,参数估计和模型检验是建模的重要步骤.,线性最小二乘法是参数估计的基本方法.,修改假设,Logistic 模型的应用模型应用加入2000年人口数据后,41,1.4.1 地形地貌图形的绘制,1.4,数据的可视化,对某地地貌测量所得结果（相对某高度），x，y方向均从110，用这些数据尽量准确地绘制出该地区的地形图。,关键是要将未测量的高度用数据插值的方法求出来，然后用MATLAB绘制出来。,1.4.1 地形地貌图形的绘制1.4 数据的可视化,42,补充：用,MATLAB,作插值计算,一维插值函数：,yi=interp1(x，y，xi，method),插值方法,被插值点,插值节点,xi,处的插值结果,nearest,最邻近插值；,linear,线性插值；,spline,三次样条插值；,cubic,立方插值；,缺省时 分段线性插值,注意：所有的插值方法,都要求,x,是单调的，并且,xi,不,能够超过,x,的范围,补充：用MATLAB作插值计算一维插值函数：yi=inter,43,要求,x0,y0,单调；,x,，,y,可取,为矩阵，或,x,取行向量，,y,取为列向量，,x,y,的值分别不能超出,x0,y,0,的范围,z=interp2(x0,y0,z0,x,y,method),被插值点,插值方法,用,MATLAB,作网格节点数据的插值,插值节点,被插值点的函数值,nearest,最邻近插值；,linear,双线性插值；,cubic,双三次插值；,缺省时 双线性插值.,要求x0,y0单调；x，y可取为矩阵，或x取行向量，,44,程序如下：,x,y=meshgrid(1:10); %构造测量网格,h=0,0.02,-0.12,0,-2.09,0,-0.58,-0.08,0,0;,0.02,0,0,-2.38,0,-4.96,0,0,0,-0.1;,0,0.1,1,0,-3.04,0,-0.53,0,0.1,0;,0,0,0,3.52,0,0,0,0,0,0;,-0.43,-1.98,0,0,0,0.77,0,2.17,0,0;,0,0,-2.29,0,0.69,0,2.59,0,0.3,0;,-0.09,-0.31,0,0,0,4.27,0,0,0,-0.01;,0,0,0,5.13,7.4,0,1.89,0,0.04,0;,0.1,0,0.58,0,0,1.75,0,-0.11,0,0;,0,-0.01,0,0,0.3,0,0,0,0,0.01; %测量数据点,xi,yi=meshgrid(1:0.1:10); %构造插值网格,hi=interp2(x,y,h,xi,yi,spline); %二维插值命令,surf(hi); %绘制地形图,xlabel(x);ylabel(y),zlabel(h);,程序如下：x,y=meshgrid(1:10); %构造,45,1.4.2 车灯光源投影区域的绘制CUMCM2002A,2002年CUMCM的A题，其中的一个要求时绘制投影区域。先建立车灯投影的数学模型，再根据模型绘出投影效果图。,根据得到的线光源长度，用投点法可以画出测试屏上的反射光亮区。MATLAB程序如下：,1.4.2 车灯光源投影区域的绘制CUMCM2002A,46,MATLAB程序如下：,p=0.03;x=25.0216;,for y1=-0.002:0.0004:0.002,y0=(-0.036:0.001:0.036)*ones(1,73);,z0=ones(73,1)*(-0.036:0.001:0.036);,x0=(y0.2+z0.2)/(2*p);,xn=(p3+4*x0*2*p.*x0+p*(-4*y1*y0+3*2*p*x0)./(2*(p2+2*p*x0);,yn=(2*p*x0.*y0+p2*(-y1+y0)+y1*(y0.2-z0.2)./(p2+2*p*x0);,zn=(p2+2*p*x0+2*y1*y0).*z0./(p2+2*p*x0);,y=y0+(yn-y0).*(x-x0)./(xn-x0);,z=z0+(zn-z0).*(x-x0)./(xn-x0);,plot(y,z,b.),xlabel(y);ylabel(z);,hold on,end,MATLAB程序如下：p=0.03;x=25.0216;,47,