MATLAB的差值与曲线拟合.ppt

MATLAB的差值与曲线拟合,MATLAB数据处理与应用2011-2012学年选修课,第八讲,王文健,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,2,主要内容,插值运算曲线拟合分段函数拉格拉日插值牛顿插值Hermite插值,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,3,MATLAB插值与曲线拟合,简介把试验数据绘制成图是点状分布，将这些点练成曲线，然后转化为有意义的数学函数，才能对数据做对比和分析上述过程涉及插值(interpolation)和曲线拟合(curve-fitting)插值过程中认为数据是准确的，求取描述点之间的数据曲线拟合中，假定已知曲线的规律，做曲线的最佳逼近，不需要经过所有数据点,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,4,MATLAB插值与曲线拟合,插值计算假设试验得到一组数据形态为：f(xk)；其中k=1,2,n；且x1=axn=c如果某些点xi不属于上述的xi，但是a=xi=c，要估计这些点的函数值f(xi)就需要做插值运算根据原始数据所描述的函数的复杂程度，存在有：一维插值二维插值Spline插值,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,5,MATLAB插值与曲线拟合,一维插值一维插值可以是线性的，也可以是三次多项式或spline插值一维线性插值是假设两个数据中的变化为线性关系，因此可由已知点的坐标(f(a),a)和(f(c),c)计算b点的函数值f(b)一维线性插值是最简单的插值，适用范围很小可用较为复杂的三次多项式或spline来近似找到原函数f(x),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,6,MATLAB插值与曲线拟合,一维插值MATLAB的一维插值函数为interp1(a,c,b)和interp1(a,c,b,method)其中a,c为已知数据，而b为要插值的数据点，method为预先设定的插值方法，分别为线性(linear)、三次多项式(cubic)和spline如果数据的变化较大，以spline插值所形成的曲线最为平滑，效果最好三次多项式所得的曲线平滑度介于linear和spline之间,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,7,MATLAB插值与曲线拟合,一维插值举例假设汽车引擎在额定转速下温度与时间的测量值关系如下00；120；260；368；477；5110；用一维插值法估计时间为2.6s和4.9s时的温度值。分别用三种方法估计其值x=0 1 2 3 4 5;y=0 20 60 68 77 110;y1=interp1(x,y,2.6)y2=interp1(x,y,2.6 4.9)y3=interp1(x,y,2.6 4.9,cubic)y3=interp1(x,y,2.6 4.9,spline),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,8,MATLAB插值与曲线拟合,一维插值举例画出插值后的时间与温度变化曲线图。t=0 1 2 2.6 3 4 4.9 5;y1=0 20 60 64.8 68 77 106.7 110;y2=0 20 60 65.2 68 77 105.6 110;y3=0 20 60 67.3 68 77 105.2 110;plot(t,y1,t,y2,t,y3),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,9,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,10,MATLAB插值与曲线拟合,二维插值假设试验所得到的一组数据形态为z=f(xk,yk)，期中k=1，2，n；如果某些点(xi,yi)不属于上述点，要估计这些点的函数值f(xi,yi)就需要进行二维擦绘制运算，二维插值相当于二元函数运算MATLAB中二维插值函数是interp2(x,y,z,xi,yi)，期中xi,yi为要插值的数据点，x，y，z为已知数据通过设定interp2(x,y,z,xi,yi,method)也可以设定不同的插值方法，有linear、cubic和spline三种,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,11,MATLAB插值与曲线拟合,二维插值举例假设汽车的引擎转速、温度与时间的测量值关系如下表，估计时间为2.6s和速度为2500r时的温度。,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,12,MATLAB插值与曲线拟合,二维插值举例应用函数interp2函数分别用linear、cubic、spline来估计时间2.6s和转速2500r的温度值。t=1 2 3 4 5;speed=2000 3000 4000;temp=20 110 176;60 180 220;68 240 349;77 310 450;110 405 503;temp1=interp2(speed,t,temp,2500,2.6)temp2=interp2(speed,t,temp,2500,2.6,cubic)temp3=interp2(speed,t,temp,2500,2.6,spline),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,13,MATLAB插值与曲线拟合,spline插值函数spline有两种应用形式，spline(x,y,xi)和spline(x,y)spline(x,y,xi)与interp1(x,y,xi,spline)效果一样举例x=0 1 2 3 4 5;y=0 20 60 68 77 110;y1=spline(x,y,2.6)y1=spline(x,y,2.6,4.9),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,14,MATLAB插值与曲线拟合,spline插值举例：在0,10区间上按正弦规律y=sin(x)取10个点yi，再在区间0,10上取41个点，以yi为已知数据，对这41个点做spline插值运算，得到函数值并绘图，比较其值与正弦波偏离程度。x=0:10;y=sin(x);xx=0:0.25:10;yy=spline(x,y,xx);plot(x,y,o,xx,yy,xx,sin(xx);,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,15,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,16,MATLAB插值与曲线拟合,spline插值采用一维线性插值运算x=0:10;y=sin(x);xx=0:0.25:10;yy=interp1(x,y,xx);plot(x,y,o,xx,yy,xx,sin(xx);,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,17,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,18,MATLAB插值与曲线拟合,spline插值采用一维三次多项式插值运算x=0:10;y=sin(x);xx=0:0.25:10;yy=interp1(x,y,xx,cubic);plot(x,y,o,xx,yy,xx,sin(xx);,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,19,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,20,MATLAB插值与曲线拟合,spline插值向量t=1900:10:1990，p=75.995 91.972 105.711 123.203 131.669 150.697 179.323 203.212 226.505 249.633表示从1900年到1990年人口普查的数字，利用spline插值法预测2000年的人口数，并绘制人口变化趋势图。t=1900:10:1990;p=75.995 91.972 105.711 123.203 131.669 150.697 179.323 203.212 226.505 249.633;y=spline(t,p,2000)plot(t 2000,p y,o),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,21,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,22,MATLAB插值与曲线拟合,spline插值画一个圆并在数据点y(:,2)，y(:,6)处标上字母o。x=pi*0:0.5:2;y=0 1 0-1 0 1 0;1 0 1 0-1 0 1;pp=spline(x,y);yy=ppval(pp,linspace(0,2*pi,101);%计算分算函数的值plot(yy(1,:),yy(2,:),-b,y(1,2:5),y(2,2:5),or),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,23,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,24,MATLAB插值与曲线拟合,曲线拟合曲线拟合(curve-fitting)顾名思义就是用近似的曲线方程来代表一组离散的数据曲线拟合与插值有许多相似之处，二者最大区别再于曲线拟合要找出一个曲线函数式，而插值仅求出与数据点对应的函数值而已曲线拟合运用最小二乘法原理如果拟合的曲线限定为多项式就称为多项式最小二乘法曲线拟合MATLAB提供的曲线拟合有多项式最小二乘拟合、普通最小二乘拟合等,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,25,MATLAB插值与曲线拟合,线性回归假设有一组数据形态为y=y(x)，其中x=0 1 2 3 4 5，y=0 20 60 68 77 110，用一个最简单的方程式来表示这组数据，一阶线性方程式最简单随意假设其方程式为y=20 xx=0 1 2 3 4 5;y=0 20 60 68 77 110;y1=20*xsum_sq=sum(y-y1).2)plot(x,y1,x,y,o),title(Linear estimate),grid,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,26,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,27,MATLAB插值与曲线拟合,线性回归从图中看出，所假设的方程与实际数据的差距很大，从理论上讲，误差平方和最小才是代表数据的最佳方程式，以此标准寻找数据方程式的方法称为最小平方误差法或线性回归MATLAB中的polyfit()函数提供了从一阶到高阶多项式的回归运算，语法为coef=polyfit(x,y,n)，期中x、y为输入数据，n为多项式的阶数，coef为输出数组，数值为多项式的系数项,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,28,MATLAB插值与曲线拟合,线性回归假设一组数据x=0 1 2 3 4 5和y=0 20 60 68 77 110，期中y是x的函数，使用最小二乘法建立一阶线性回归方程用polyfit函数求出系数x=0 1 2 3 4 5;y=0 20 60 68 77 110;coef=polyfit(x,y,1),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,29,MATLAB插值与曲线拟合,线性回归画出方程的直线及x和y的数据点a1=coef(1);a0=coef(2);ybest=a1*x+a0;sum_sq=sum(y-ybest).2)plot(x,ybest,x,y,o);title(Linear regression estimate),grid通过比较发现，最小平方误差法获得的误差平方和要小于y=20 x的值,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,30,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,31,MATLAB插值与曲线拟合,多项式回归函数polyfit(x,y,n)中当n大于等于2时，称为多项式回归高阶多项式回归有时并不如低阶的更能代表原始数据对上述数据的x和y分别求二、三、四、五阶多项式回归方程，画出其拟合曲线二阶：x=0 1 2 3 4 5;y=0 20 60 68 77 110;coef=polyfit(x,y,2);ybest=coef(1)*x.2+coef(2)*x+coef(3);plot(x,ybest,x,y,o),title(Linear regression estimate),grid,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,32,MATLAB插值与曲线拟合,多项式回归三阶：x=0 1 2 3 4 5;y=0 20 60 68 77 110;coef=polyfit(x,y,3);ybest=coef(1)*x.3+coef(2)*x.2+coef(3)*x+coef(4);plot(x,ybest,x,y,o),title(Linear regression estimate),grid三阶的拟合曲线最为光滑,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,33,MATLAB插值与曲线拟合,多项式回归四阶：x=0 1 2 3 4 5;y=0 20 60 68 77 110;coef=polyfit(x,y,4);ybest=coef(1)*x.4+coef(2)*x.3+coef(3)*x.2+coef(4)*x+coef(5);plot(x,ybest,x,y,o),title(Linear regression estimate),grid,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,34,MATLAB插值与曲线拟合,多项式回归五阶：x=0 1 2 3 4 5;y=0 20 60 68 77 110;coef=polyfit(x,y,5);ybest=coef(1)*x.5+coef(2)*x.4+coef(3)*x.3+coef(4)*x.2+coef(5)*x+coef(6);plot(x,ybest,x,y,o),title(Linear regression estimate),grid五阶以上的多项式回归方程会通过所有的原始数据点，但曲振荡剧烈,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,35,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,36,MATLAB插值与曲线拟合,多次式拟合及函数的计算polyval也可用来做多项式函数计算，由polyfit计算出多项式的系数后，再由polyval计算出多项式的函数值语法polyval(coef,x)，期中coef为多项式的各个系数所构成的数组x=0:0.1:1;y=-.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2;coef=polyfit(x,y,2)xi=0:0.01:1;yi=polyval(coef,xi);plot(x,y,o,xi,yi);,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,37,MATLAB插值与曲线拟合,多次式拟合及函数的计算举例：对数据x=0 1 2 3 4 5和y=0 20 60 68 77 110做二阶到九阶的多项式回归分析，然后计算多项式的函数值，绘制其曲线x=0 1 2 3 4 5;y=0 20 60 68 77 110;newx=0:0.05:5;for n=2:9f(:,n)=polyval(polyfit(x,y,n),newx);plot(newx,f(:,n),x,y,o);hold on;end,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,38,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,39,MATLAB插值与曲线拟合,分段多项式函数函数pp=mkpp(breaks,coefs)是区间breaks上根据多项式的系数coefs构建一个分段多项式ppbreaks是一个长度为L+1的向量，元素是每一个区间的起点与终点，共有L个区间coefs是一个L(k+1)的矩阵函数ppval()可估算分段多项式pp在特殊点xx上的值v=ppval(pp,xx),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,40,MATLAB插值与曲线拟合,分段多项式函数举例在区间-8-4 0 4 8上分别画出首尾相接的y=-0.5x+1；y=0.5x-1；y=-0.5x+1；y=0.5x-1四条线段breaks=-8-4 0 4 8;coefs=-0.5 1;0.5-1;-0.5 1;0.5-1;pp=mkpp(breaks,coefs)x=-8:0.1:8;y=ppval(pp,x);plot(x,y,-b);,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,41,MATLAB插值与曲线拟合,分段多项式函数二次多项式y=-x2/4+x，在区间-8-4 0 4 8上分别画出y;-y;y;-y的曲线breaks=-8-4 0 4 8;coefs=-1/4 1 0;pp=mkpp(breaks,coefs;-coefs;coefs;-coefs)x=-8:0.1:8;y=ppval(pp,x);plot(x,y,K-);,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,42,MATLAB插值与曲线拟合,分段多项式函数在区间-8-4上应用函数pp=mkpp(breaks,coefs)构建二次多项式函数y=-0.25x2+xcoefs=-0.25 1 0;pp=mkpp(-8-4,coefs)x=-8:0.1:-4;y=ppval(pp,x);plot(x,y);,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,43,MATLAB插值与曲线拟合,分段多项式函数画出分段函数在连续三个周期内的图像x=0 0.075 0.1 0.175 0.2 0.275 0.3;y=40/3 0;-40 1;40/3 0;-40 1;40/3 0;-40 1;pp=mkpp(x,y);xx=0:0.001:0.3;yy=ppval(pp,xx);plot(xx,yy,b-);,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,44,MATLAB插值与曲线拟合,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,45,MATLAB插值与曲线拟合,获得分段多项式的详细资料使用函数breaks,coefs,l,k,d=unmkpp(pp)可以获得分段多项式的详细资料Breaks是断点数，coefs是系数，l是段数，k是项数，d是维数，pp是由于spine或mkpp所创建为分段二项多项式-0.25x4+x和-(-0.25x4+x)创建一个描述y=-0.25 1 0;pp=mkpp(-8-4 0 4 8,y;-y;y;-y);breaks,coefs,l,k,d=unmkpp(pp),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,46,MATLAB插值与曲线拟合,拉格拉日插值法通过平面上不同两点可以去定一条直线经过这两点，这就是拉格拉日插值问题，对于不在同一条直线的三个点得到哦的插值多项式为抛物线举例：已知cos30=sqrt(3)/2、cos45=sqrt(2)/2、cos60=0.5、cos90=0，利用拉格拉日插值法计算cos(-40)、cos47、cos53、cos79、cos174的值，并给出插值多项式步骤：1.编写拉格拉日插值方法函数文件lag.m2.编写主函数文件lag_main.m3.命令窗口运行lag_main,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,47,MATLAB插值与曲线拟合,拉格拉日插值法x=pi/4,pi/6,pi/3,pi/2;y=cos(pi/4),cos(pi/6),cos(pi/3),cos(pi/2);t=-40*pi/180,47*pi/180,53*pi/180,79*pi/180,174*pi/180;du=-40 47 53 79 174yt=lag(x,y,t)yreal=cos(-40*pi/180)cos(47*pi/180)cos(53*pi/180)cos(79*pi/180)cos(174*pi/180)dy=yt-yrealyt=lag(x,y)ezplot(yt,-pi/4,pi)hold onezplot(cos(t),-pi/4,pi);,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,48,MATLAB插值与曲线拟合,牛顿插值法牛顿插值法主要是基于插商的概念进行的，当增加插值基点时原来的计算结果对后面的计算仍然有用举例：已知零阶Bvessel函数f(x)在若干点处的函数值如下，X 1.0 1.3 1.6 1.9 2.2f(x)0.7651977 0.6200860 0.4554022 0.2818186 0.1103623x=1.0 1.3 1.6 1.9 2.2;y=0.7651977 0.6200860 0.4554022 0.2818186 0.1103623;yt=niudun(x,y,1.5),Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,49,MATLAB插值与曲线拟合,牛顿插值法5阶勒让德多项式为P5=(63X5-70X3+15X)/8，在-1,1上任意取9个点，如x=-1-0.3 0.4-0.7 0-0.4 0.8 0.7 1，不管取点顺序1.根据5阶勒让德多项式，计算这9点的函数值；2.依据这9点作牛顿插值，画出插值多项式图像，与勒让德多项式图像相比较；3.计算x=0.24-0.46 0.83时的近似结果，与真实结果相比。步骤：1.编写牛顿插值方法函数文件niudun.m2.编写主函数文件niudun_main.m3.命令窗口运行niudun_main,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,50,MATLAB插值与曲线拟合,牛顿插值法clearsyms xfx=(63*x5-70*x3+15*x)/8;v=-1,1,-1,1;subplot(1,2,1);ezplot(fx),grid onaxis(v)x0=-1-0.3 0.4-0.7 0-0.4 0.8 0.7 1;y0=subs(fx,x0);yt=niudun(x0,y0);subplot(1,2,2);ezplot(yt),axis(v),grid t=0.24-0.46 0.83legen5=subs(fx,t)chazhi=niudun(x0,y0,t)wucha=legen5-chazhi,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,51,MATLAB插值与曲线拟合,Hermite插值法前面的插值要求多项式在插值点处的取值和函数值相等，当不仅要求插值多项式在插值点与函数值相同，而且还要求导数相同时，就是Hermite插值基本思想：满足多项式的点的一阶导数与函数的一阶导数值相等举例：人为的制造一个函数g=x2sinx，计算出这个函数在一下点列中的函数值与导数值，X=-3-2-1 0 1 2 3，利用上面的三组数据插值计算以下点的函数近似值，P=-5.1-3.2-1.4 1.6 5.3，算出这些点的真实函数值，并比较不同点的误差,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,52,MATLAB插值与曲线拟合,Hermite插值法步骤：1.编制Hermite插值函数M文件2.编写主函数文件herm_main1文件3.命令窗口中运行命令herm_main1,Tribology Research Institute SOUTHWEST JIAOTONG UNIVERSITY,53,MATLAB插值与曲线拟合,插值中的龙格现象当利用上述插值方法进行计算时，是否N越大误差就会越趋于0，即是否阶数越高插值效果越好，对与sinx这样的函数由于所有导数有相同的常数界，答案是肯定的，但对于一般函数，不一定阶数越高插值效果越好，这就是插值过程中的龙格现象注意：对于实际工程与应用的插值问题实际上很少采用很高阶插值方法对于插值中的龙格现象，一般可采用分段插值的方法进行消除,谢 谢,