matlab數(shù)值分析-插值與曲線擬合

1、1?MATLAB數(shù)據(jù)處理與應(yīng)用數(shù)據(jù)處理與應(yīng)用?2021-2021學(xué)年選修課學(xué)年選修課西南交通大學(xué)摩擦學(xué)研究所2主要內(nèi)容插值運(yùn)算曲線擬合分段函數(shù)拉格拉日插值牛頓插值Hermite插值西南交通大學(xué)摩擦學(xué)研究所3MATLAB插值與曲線擬合插值與曲線擬合 簡介把試驗(yàn)數(shù)據(jù)繪制成圖是點(diǎn)狀分布，將這些點(diǎn)練成曲線，然后轉(zhuǎn)化為有意義的數(shù)學(xué)函數(shù)，才能對數(shù)據(jù)做比照和分析上述過程涉及插值(interpolation)和曲線擬合(curve-fitting)插值過程中認(rèn)為數(shù)據(jù)是準(zhǔn)確的，求取描述點(diǎn)之間的數(shù)據(jù)曲線擬合中，假定曲線的規(guī)律，做曲線的最正確逼近，不需要經(jīng)過所有數(shù)據(jù)點(diǎn)西南交通大學(xué)摩擦學(xué)研究所4MATLAB插值與曲線

2、擬合插值與曲線擬合插值計(jì)算假設(shè)試驗(yàn)得到一組數(shù)據(jù)形態(tài)為：f(xk)；其中k=1,2,n；且x1=axn=c如果某些點(diǎn)xi不屬于上述的xi ，但是a=xi=c，要估計(jì)這些點(diǎn)的函數(shù)值f(xi)就需要做插值運(yùn)算根據(jù)原始數(shù)據(jù)所描述的函數(shù)的復(fù)雜程度，存在有：一維插值二維插值Spline插值西南交通大學(xué)摩擦學(xué)研究所5MATLAB插值與曲線擬合插值與曲線擬合一維插值一維插值可以是線性的，也可以是三次多項(xiàng)式或spline插值一維線性插值是假設(shè)兩個數(shù)據(jù)中的變化為線性關(guān)系，因此可由點(diǎn)的坐標(biāo)(f(a),a)和(f(c),c)計(jì)算b點(diǎn)的函數(shù)值f(b)一維線性插值是最簡單的插值，適用范圍很小可用較為復(fù)雜的三次多項(xiàng)式或sp

3、line來近似找到原函數(shù)f(x)()()()(afcfacabafbf西南交通大學(xué)摩擦學(xué)研究所6MATLAB插值與曲線擬合插值與曲線擬合一維插值MATLAB的一維插值函數(shù)為interp1(a,c,b)和interp1(a,c,b,method)其中a,c為數(shù)據(jù)，而b為要插值的數(shù)據(jù)點(diǎn)，method為預(yù)先設(shè)定的插值方法，分別為線性(linear)、三次多項(xiàng)式(cubic)和spline如果數(shù)據(jù)的變化較大，以spline插值所形成的曲線最為平滑，效果最好三次多項(xiàng)式所得的曲線平滑度介于linear和spline之間西南交通大學(xué)摩擦學(xué)研究所7MATLAB插值與曲線擬合插值與曲線擬合一維插值舉例假設(shè)汽車引

4、擎在額定轉(zhuǎn)速下溫度與時間的測量值關(guān)系如下00； 120； 260； 368； 477； 5110；用一維插值法估計(jì)時間為和時的溫度值。分別用三種方法估計(jì)其值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)西南交通大學(xué)摩擦學(xué)研究所8MATLAB插值與曲線擬合插值與曲線擬合一維插值舉例畫出插值后的時間與溫度變化曲線圖。t=0 1 2 2.6 3 4 4.9 5;y1=0 20 6

5、0 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)西南交通大學(xué)摩擦學(xué)研究所9MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所10MATLAB插值與曲線擬合插值與曲線擬合二維插值假設(shè)試驗(yàn)所得到的一組數(shù)據(jù)形態(tài)為z=f(xk,yk)，期中k=1，2，n；如果某些點(diǎn)(xi,yi)不屬于上述點(diǎn)，要估計(jì)這些點(diǎn)的函數(shù)值f(xi,yi)就需要進(jìn)行二維擦繪制運(yùn)算，二維插值相當(dāng)于二元函數(shù)運(yùn)算MATLAB中二維插值函數(shù)是interp2(x

6、,y,z,xi,yi)，期中xi,yi為要插值的數(shù)據(jù)點(diǎn)，x，y，z為數(shù)據(jù)通過設(shè)定interp2(x,y,z,xi,yi,method)也可以設(shè)定不同的插值方法，有l(wèi)inear、cubic和spline三種西南交通大學(xué)摩擦學(xué)研究所11MATLAB插值與曲線擬合插值與曲線擬合二維插值舉例假設(shè)汽車的引擎轉(zhuǎn)速、溫度與時間的測量值關(guān)系如下表，估計(jì)時間為和速度為2500r時的溫度。T/SpeedTime/s2000r3000r4000r00001201101762601802203682403494773104505110405503西南交通大學(xué)摩擦學(xué)研究所12MATLAB插值與曲線擬合插值與曲線擬合二維

7、插值舉例應(yīng)用函數(shù)interp2函數(shù)分別用linear、cubic、spline來估計(jì)時間和轉(zhuǎn)速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)西南交通大學(xué)摩擦學(xué)研究所13MATLAB插值與曲線擬合插

8、值與曲線擬合spline插值函數(shù)spline有兩種應(yīng)用形式，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)西南交通大學(xué)摩擦學(xué)研究所14MATLAB插值與曲線擬合插值與曲線擬合spline插值舉例：在0,10區(qū)間上按正弦規(guī)律y=sin(x)取10個點(diǎn)yi，再在區(qū)間0,10上取41個點(diǎn)，以yi為數(shù)據(jù)，對這41個點(diǎn)做spline插值運(yùn)算，得到函數(shù)值并繪圖，比較其值

9、與正弦波偏離程度。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);西南交通大學(xué)摩擦學(xué)研究所15MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所16MATLAB插值與曲線擬合插值與曲線擬合spline插值采用一維線性插值運(yùn)算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);西南交通大學(xué)摩擦學(xué)研究所17MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所18MATLAB

10、插值與曲線擬合插值與曲線擬合spline插值采用一維三次多項(xiàng)式插值運(yùn)算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);西南交通大學(xué)摩擦學(xué)研究所19MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所20MATLAB插值與曲線擬合插值與曲線擬合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年到19

11、90年人口普查的數(shù)字，利用spline插值法預(yù)測2000年的人口數(shù)，并繪制人口變化趨勢圖。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)西南交通大學(xué)摩擦學(xué)研究所21MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所22MATLAB插值與曲線擬合插值與曲線擬合spline插值畫一個圓并在數(shù)據(jù)點(diǎn)y(:,2)，y(:,6)處標(biāo)上字母o。x=pi*0:0.5:2;y=0 1

12、0 -1 0 1 0;1 0 1 0 -1 0 1;pp=spline(x,y);yy=ppval(pp,linspace(0,2*pi,101);%計(jì)算分算函數(shù)的值plot(yy(1,:),yy(2,:),-b,y(1,2:5),y(2,2:5),or)西南交通大學(xué)摩擦學(xué)研究所23MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所24MATLAB插值與曲線擬合插值與曲線擬合曲線擬合曲線擬合(curve-fitting)顧名思義就是用近似的曲線方程來代表一組離散的數(shù)據(jù)曲線擬合與插值有許多相似之處，二者最大區(qū)別再于曲線擬合要找出一個曲線函數(shù)式，而插值僅求出與數(shù)據(jù)點(diǎn)對應(yīng)的函數(shù)值而已曲

13、線擬合運(yùn)用最小二乘法原理如果擬合的曲線限定為多項(xiàng)式就稱為多項(xiàng)式最小二乘法曲線擬合MATLAB提供的曲線擬合有多項(xiàng)式最小二乘擬合、普通最小二乘擬合等西南交通大學(xué)摩擦學(xué)研究所25MATLAB插值與曲線擬合插值與曲線擬合線性回歸假設(shè)有一組數(shù)據(jù)形態(tài)為y=y(x)，其中x=0 1 2 3 4 5，y=0 20 60 68 77 110，用一個最簡單的方程式來表示這組數(shù)據(jù)，一階線性方程式最簡單隨意假設(shè)其方程式為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 est

14、imate),grid西南交通大學(xué)摩擦學(xué)研究所26MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所27MATLAB插值與曲線擬合插值與曲線擬合線性回歸從圖中看出，所假設(shè)的方程與實(shí)際數(shù)據(jù)的差距很大，從理論上講，誤差平方和最小才是代表數(shù)據(jù)的最正確方程式，以此標(biāo)準(zhǔn)尋找數(shù)據(jù)方程式的方法稱為最小平方誤差法或線性回歸MATLAB中的polyfit()函數(shù)提供了從一階到高階多項(xiàng)式的回歸運(yùn)算，語法為coef=polyfit(x,y,n)，期中x、y為輸入數(shù)據(jù)，n為多項(xiàng)式的階數(shù)，coef為輸出數(shù)組，數(shù)值為多項(xiàng)式的系數(shù)項(xiàng)西南交通大學(xué)摩擦學(xué)研究所28MATLAB插值與曲線擬合插值與曲線擬合線性回歸假

15、設(shè)一組數(shù)據(jù)x=0 1 2 3 4 5和y=0 20 60 68 77 110，期中y是x的函數(shù)，使用最小二乘法建立一階線性回歸方程用polyfit函數(shù)求出系數(shù)x=0 1 2 3 4 5;y=0 20 60 68 77 110;coef=polyfit(x,y,1)西南交通大學(xué)摩擦學(xué)研究所29MATLAB插值與曲線擬合插值與曲線擬合線性回歸畫出方程的直線及x和y的數(shù)據(jù)點(diǎn)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),gri

16、d通過比較發(fā)現(xiàn)，最小平方誤差法獲得的誤差平方和要小于y=20 x的值西南交通大學(xué)摩擦學(xué)研究所30MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所31MATLAB插值與曲線擬合插值與曲線擬合多項(xiàng)式回歸函數(shù)polyfit(x,y,n)中當(dāng)n大于等于2時，稱為多項(xiàng)式回歸高階多項(xiàng)式回歸有時并不如低階的更能代表原始數(shù)據(jù)對上述數(shù)據(jù)的x和y分別求二、三、四、五階多項(xiàng)式回歸方程，畫出其擬合曲線二階：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,yb

17、est,x,y,o),title(Linear regression estimate),grid西南交通大學(xué)摩擦學(xué)研究所32MATLAB插值與曲線擬合插值與曲線擬合多項(xiàng)式回歸三階：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三階的擬合曲線最為光滑西南交通大學(xué)摩擦學(xué)研究所33MATLAB插值與曲線擬合插值與曲線擬合多項(xiàng)式回

18、歸四階：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西南交通大學(xué)摩擦學(xué)研究所34MATLAB插值與曲線擬合插值與曲線擬合多項(xiàng)式回歸五階：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

19、.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五階以上的多項(xiàng)式回歸方程會通過所有的原始數(shù)據(jù)點(diǎn)，但曲振蕩劇烈西南交通大學(xué)摩擦學(xué)研究所35MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所36MATLAB插值與曲線擬合插值與曲線擬合屢次式擬合及函數(shù)的計(jì)算polyval也可用來做多項(xiàng)式函數(shù)計(jì)算，由polyfit計(jì)算出多項(xiàng)式的系數(shù)后，再由polyval計(jì)算出多項(xiàng)式的函數(shù)值語法polyval(coef,x)，期中coef為多項(xiàng)式的各

20、個系數(shù)所構(gòu)成的數(shù)組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);西南交通大學(xué)摩擦學(xué)研究所37MATLAB插值與曲線擬合插值與曲線擬合屢次式擬合及函數(shù)的計(jì)算舉例：對數(shù)據(jù)x=0 1 2 3 4 5和y=0 20 60 68 77 110做二階到九階的多項(xiàng)式回歸分析，然后計(jì)算多項(xiàng)式的函數(shù)值，繪制其曲線x=0 1 2 3 4 5;y=0 20 60 68 77 110;ne

21、wx=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西南交通大學(xué)摩擦學(xué)研究所38MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所39MATLAB插值與曲線擬合插值與曲線擬合分段多項(xiàng)式函數(shù)函數(shù)pp=mkpp(breaks,coefs)是區(qū)間breaks上根據(jù)多項(xiàng)式的系數(shù)coefs構(gòu)建一個分段多項(xiàng)式ppbreaks是一個長度為L+1的向量，元素是每一個區(qū)間的起點(diǎn)與終點(diǎn)，共有L個區(qū)間coefs是一個L(k+1)的矩陣函數(shù)ppval()可估算分段多項(xiàng)式pp

22、在特殊點(diǎn)xx上的值v=ppval(pp,xx)西南交通大學(xué)摩擦學(xué)研究所40MATLAB插值與曲線擬合插值與曲線擬合分段多項(xiàng)式函數(shù)舉例在區(qū)間-8 -4 0 4 8上分別畫出首尾相接的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);西南交通大學(xué)摩擦學(xué)研究所41MATLAB插值與曲線擬合插值與曲線擬合分段多項(xiàng)式函數(shù)二次多項(xiàng)式y(tǒng)=-x2/4+x，在區(qū)間-8 -4 0 4 8上分別畫

23、出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-);西南交通大學(xué)摩擦學(xué)研究所42MATLAB插值與曲線擬合插值與曲線擬合分段多項(xiàng)式函數(shù)在區(qū)間-8 -4上應(yīng)用函數(shù)pp=mkpp(breaks,coefs)構(gòu)建二次多項(xiàng)式函數(shù)2+xcoefs=-0.25 1 0;pp=mkpp(-8 -4,coefs)x=-8:0.1:-4;y=ppval(pp,x);plot(x,y);西南交通大學(xué)摩擦學(xué)研究所43M

24、ATLAB插值與曲線擬合插值與曲線擬合分段多項(xiàng)式函數(shù)畫出分段函數(shù)在連續(xù)三個周期內(nèi)的圖像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-); 1 . 0 ,075. 0 , 140075. 0 , 0 ,3/40 xxy西南交通大學(xué)摩擦學(xué)研究所44MATLAB插值與曲線擬合插值與曲線擬合西南交通大學(xué)摩擦學(xué)研究所45MATLAB插值與曲線擬合插值與曲線擬合獲得分段多項(xiàng)式的詳細(xì)資料使用函

25、數(shù)breaks,coefs,l,k,d=unmkpp(pp)可以獲得分段多項(xiàng)式的詳細(xì)資料Breaks是斷點(diǎn)數(shù)，coefs是系數(shù)，l是段數(shù)，k是項(xiàng)數(shù)，d是維數(shù)，pp是由于spine或mkpp所創(chuàng)立為分段二項(xiàng)多項(xiàng)式4+x和4+x)創(chuàng)立一個描述y=-0.25 1 0;pp=mkpp(-8 -4 0 4 8,y;-y;y;-y);breaks,coefs,l,k,d=unmkpp(pp)西南交通大學(xué)摩擦學(xué)研究所46MATLAB插值與曲線擬合插值與曲線擬合拉格拉日插值法通過平面上不同兩點(diǎn)可以去定一條直線經(jīng)過這兩點(diǎn)，這就是拉格拉日插值問題，對于不在同一條直線的三個點(diǎn)得到哦的插值多項(xiàng)式為拋物線舉例：cos3

26、0=sqrt(3)/2、cos45=sqrt(2)/2、 、cos90=0，利用拉格拉日插值法計(jì)算cos(-40)、cos47、cos53、cos79、cos174的值，并給出插值多項(xiàng)式步驟：1.編寫拉格拉日插值方法函數(shù)文件2.編寫主函數(shù)文件3.命令窗口運(yùn)行l(wèi)ag_main西南交通大學(xué)摩擦學(xué)研究所47MATLAB插值與曲線擬合插值與曲線擬合拉格拉日插值法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=-

27、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);西南交通大學(xué)摩擦學(xué)研究所48MATLAB插值與曲線擬合插值與曲線擬合牛頓插值法牛頓插值法主要是基于插商的概念進(jìn)行的，當(dāng)增加插值基點(diǎn)時原來的計(jì)算結(jié)果對后面的計(jì)算仍然有用舉例：零階Bvessel函數(shù)f(x)在假設(shè)干點(diǎn)處的函數(shù)值如下，f(xx=1.0

28、 1.3 1.6 1.9 2.2;y=0.7651977 0.6202160 0.4554022 0.2818186 0.1103623;yt=niudun(x,y,1.5)西南交通大學(xué)摩擦學(xué)研究所49MATLAB插值與曲線擬合插值與曲線擬合牛頓插值法5階勒讓德多項(xiàng)式為P5=(63X5-70X3+15X)/8，在-1,1上任意取9個點(diǎn)，如x=-1 -0.3 0.4 -0.7 0 -0.4 0.8 0.7 1，不管取點(diǎn)順序1.根據(jù)5階勒讓德多項(xiàng)式，計(jì)算這9點(diǎn)的函數(shù)值；2.依據(jù)這9點(diǎn)作牛頓插值，畫出插值多項(xiàng)式圖像，與勒讓德多項(xiàng)式圖像相比較；3.計(jì)算x=0.24 -0.46 0.83時的近似結(jié)果，與真實(shí)結(jié)果相比。步驟：1.編寫牛頓插值方法函數(shù)文件2.編寫主函數(shù)文件3.命令窗口運(yùn)行niudun_main西南交通大學(xué)摩擦學(xué)研究所50MATLAB插值與曲線擬合插值與曲線擬合牛頓插值法cle