版權(quán)說明:本文檔由用戶提供并上傳,收益歸屬內(nèi)容提供方,若內(nèi)容存在侵權(quán),請進行舉報或認領(lǐng)
文檔簡介
多項式插值擬合第1頁,共85頁,2023年,2月20日,星期一5.1關(guān)于多項式MATLAB命令一個多項式的冪級數(shù)形式可表示為:也可表為嵌套形式或因子形式
N階多項式n個根,其中包含重根和復根。若多項式所有系數(shù)均為實數(shù),則全部復根都將以共軛對的形式出現(xiàn)第2頁,共85頁,2023年,2月20日,星期一冪系數(shù):在MATLAB里,多項式用行向量表示,其元素為多項式的系數(shù),并從左至右按降冪排列。
例:
被表示為>>p=[2145]>>poly2sym(p)ans=2*x^3+x^2+4*x+5Roots:多項式的零點可用命令roots求的。
例:>>r=roots(p)得到
r=0.2500+1.5612i0.2500-1.5612i-1.0000所有零點由一個列向量給出。第3頁,共85頁,2023年,2月20日,星期一Poly:
由零點可得原始多項式的各系數(shù),但可能相差一個常數(shù)倍。例:>>poly(r)ans=1.00000.50002.00002.5000
注意:若存在重根,這種轉(zhuǎn)換可能會降低精度。例:
>>r=roots([1-615-2015-61])r=1.0042+0.0025i1.0042-0.0025i1.0000+0.0049i1.0000-0.0049i0.9958+0.0024i0.9958-0.0024i舍入誤差的影響,與計算精度有關(guān)。第4頁,共85頁,2023年,2月20日,星期一polyval:
可用命令polyval計算多項式的值。例:計算y(2.5)
>>c=[3,-7,2,1,1];xi=2.5;yi=polyval(c,xi)yi=23.8125如果xi是含有多個橫坐標值的數(shù)組,則yi也為與xi長度相同的向量。>>c=[3,-7,2,1,1];xi=[2.5,3];>>yi=polyval(c,xi)yi=23.812576.0000第5頁,共85頁,2023年,2月20日,星期一polyfit:給定n+1個點將可以唯一確定一個n階多項式。利用命令polyfit可容易確定多項式的系數(shù)。
例:>>x=[1.1,2.3,3.9,5.1];>>y=[3.887,4.276,4.651,2.117];>>a=polyfit(x,y,length(x)-1)a=-0.20151.4385-2.74775.4370>>poly2sym(a)ans=-403/2000*x^3+2877/2000*x^2-27477/10000*x+5437/1000多項式為Polyfit的第三個參數(shù)是多項式的階數(shù)。第6頁,共85頁,2023年,2月20日,星期一多項式積分:
功能:求多項式積分調(diào)用格式:py=poly_itg(p)p:被積多項式的系數(shù)
py:求積后多項式的系數(shù)
poly_itg.mfunctionpy=poly_itg(p)n=length(p);py=[p.*[n:-1:1].^(-1),0]不包括最后一項積分常數(shù)第7頁,共85頁,2023年,2月20日,星期一多項式微分:Polyder:求多項式一階導數(shù)的系數(shù)。調(diào)用格式為:b=polyder(c)c為多項式y(tǒng)的系數(shù),b是微分后的系數(shù),其值為:
第8頁,共85頁,2023年,2月20日,星期一兩個多項式的和與差:命令poly_add:求兩個多項式的和,其調(diào)用格式為:c=poly_add(a,b)
多項式a減去b,可表示為:c=poly_add(a,-b)第9頁,共85頁,2023年,2月20日,星期一功能:兩個多項式相加調(diào)用格式:b=poly_add(p1,p2)b:求和后的系數(shù)數(shù)組poly_add.mfunctionp3=poly_add(p1,p2)n1=length(p1);n2=length(p2);ifn1==n2p3=p1+p2;endifn1>n2p3=p1+[zeros(1,n1-n2),p2];endifn1<n2p3=[zeros(1,n2-n1),p1]+p2;end第10頁,共85頁,2023年,2月20日,星期一m階多項式與n階多項式的乘積是d=m+n階的多項式:計算系數(shù)的MATLAB命令是:c=conv(a,b)多項式除多項式的除法滿足:其中是商,是除法的余數(shù)。多項式和可由命令deconv算出。例:[q,r]=deconv(a,b)第11頁,共85頁,2023年,2月20日,星期一例
>>a=[2,-5,6,-1,9];b=[3,-90,-18];>>c=conv(a,b)c=6-195432-4539-792-162>>[q,r]=deconv(c,b)q=2-56-19r=0000000>>poly2sym(c)ans=6*x^6-195*x^5+432*x^4-453*x^3+9*x^2-792*x-162
第12頁,共85頁,2023年,2月20日,星期一5.2插值
5.2.1Lagrange插值方法介紹對給定的n個插值點及對應(yīng)的函數(shù)值,利用構(gòu)造的n-1次Lagrange插值多項式,則對插值區(qū)間內(nèi)任意x的函數(shù)值y可通過下式求的:MATLAB實現(xiàn)第13頁,共85頁,2023年,2月20日,星期一functiony=lagrange(x0,y0,x)ii=1:length(x0);y=zeros(size(x));fori=iiij=find(ii~=i);y1=1;forj=1:length(ij),y1=y1.*(x-x0(ij(j)));endy=y+y1*y0(i)/prod(x0(i)-x0(ij));end算例:給出f(x)=ln(x)的數(shù)值表,用Lagrange計算ln(0.54)的近似值。>>x=[0.4:0.1:0.8];>>y=[-0.916291,-0.693147,-0.510826,-0.356675,-0.223144];>>lagrange(x,y,[0.54,0.55,0.78])ans=-0.6161-0.5978-0.2484(精確解-0.616143)第14頁,共85頁,2023年,2月20日,星期一5.2.2Hermite插值方法介紹
不少實際問題不但要求在節(jié)點上函數(shù)值相等,而且要求導數(shù)值也相等,甚至要求高階導數(shù)值也相等,滿足這一要求的插值多項式就是Hermite插值多項式。下面只討論函數(shù)值與一階導數(shù)值個數(shù)相等且已知的情況。已知n個插值點及對應(yīng)的函數(shù)值和一階導數(shù)值。則對插值區(qū)間內(nèi)任意x的函數(shù)值y的Hermite插值公式:第15頁,共85頁,2023年,2月20日,星期一MATLAB實現(xiàn)%hermite.mfunctiony=hermite(x0,y0,y1,x)n=length(x0);m=length(x);fork=1:myy=0.0;fori=1:nh=1.0;a=0.0;forj=1:nifj~=ih=h*((x(k)-x0(j))/(x0(i)-x0(j)))^2;a=1/(x0(i)-x0(j))+a;endendyy=yy+h*((x0(i)-x(k))*(2*a*y0(i)-y1(i))+y0(i));endy(k)=yy;end第16頁,共85頁,2023年,2月20日,星期一算例:對給定數(shù)據(jù),試構(gòu)造Hermite多項式求出sin0.34的近似值。>>x0=[0.3,0.32,0.35];>>y0=[0.29552,0.31457,0.34290];>>y1=[0.95534,0.94924,0.93937];>>y=hermite(x0,y0,y1,0.34)y=0.3335>>sin(0.34)%與精確值比較ans=0.3335第17頁,共85頁,2023年,2月20日,星期一>>x=[0.3:0.005:0.35];y=hermite(x0,y0,y1,x);
>>plot(x,y)
>>y2=sin(x);holdon
>>plot(x,y2,'--r')第18頁,共85頁,2023年,2月20日,星期一5.2.3Runge現(xiàn)象問題的提出:根據(jù)區(qū)間[a,b]上給出的節(jié)點做插值多項式p(x)的近似值,一般總認為p(x)的次數(shù)越高則逼近f(x)的精度就越好,但事實并非如此。反例:在區(qū)間[-5,5]上的各階導數(shù)存在,但在此區(qū)間上取n個節(jié)點所構(gòu)成的Lagrange插值多項式在全區(qū)間內(nèi)并非都收斂。取n=10,用Lagrange插值法進行插值計算。第19頁,共85頁,2023年,2月20日,星期一>>x=[-5:1:5];y=1./(1+x.^2);x0=[-5:0.1:5];>>y0=lagrange(x,y,x0);>>y1=1./(1+x0.^2);%繪制圖形>>plot(x0,y0,'--r')%插值曲線>>holdon>>plot(x0,y1,‘-b')%原曲線為解決Rung問題,引入分段插值。第20頁,共85頁,2023年,2月20日,星期一算法分析:所謂分段插值就是通過插值點用折線或低次曲線連接起來逼近原曲線。MATLAB實現(xiàn)可調(diào)用內(nèi)部函數(shù)。命令1
interp1
功能:一維數(shù)據(jù)插值(表格查找)。該命令對數(shù)據(jù)點之間計算內(nèi)插值。它找出一元函數(shù)f(x)在中間點的數(shù)值。其中函數(shù)f(x)由所給數(shù)據(jù)決定。
格式1yi=interp1(x,Y,xi)
%返回插值向量yi,每一元素對應(yīng)于參量xi,同時由向量x與Y的內(nèi)插值決定。參量x指定數(shù)據(jù)Y的點。若Y為一矩陣,則按Y的每列計算。算例
對于t,beta、alpha分別有兩組數(shù)據(jù)與之對應(yīng),用分段線性插值法計算當t=321,440,571時beta、alpha的值。5.2.4分段插值第21頁,共85頁,2023年,2月20日,星期一>>temp=[300,400,500,600]';>>beta=1000*[3.33,2.50,2.00,1.67]';>>alpha=10000*[0.2128,0.3605,0.5324,0.7190]';>>ti=[321,400,571]';>>propty=interp1(temp,[beta,alpha],ti);%propty=interp1(temp,[beta,alpha],ti,’linear’);>>[ti,propty]ans=1.0e+003*0.32103.15572.43820.40002.50003.60500.57101.76576.6489第22頁,共85頁,2023年,2月20日,星期一格式2yi=interp1(Y,xi)
%假定x=1:N,其中N為向量Y的長度,或者為矩陣Y的行數(shù)。格式3yi=interp1(x,Y,xi,method)%用指定的算法計算插值:
’nearest’:最近鄰點插值,直接完成計算;’linear’:線性插值(缺省方式),直接完成計算;’spline’:三次樣條函數(shù)插值?!痗ubic’:分段三次Hermite插值。其它,如’v5cubic’。對于超出x范圍的xi的分量,使用方法’nearest’、’linear’、’v5cubic’的插值算法,相應(yīng)地將返回NaN。對其他的方法,interp1將對超出的分量執(zhí)行外插值算法。yi=interp1(x,Y,xi,method,'extrap')yi=interp1(x,Y,xi,method,extrapval)%確定超出x范圍的xi中的分量的外插值extrapval,其值通常取NaN或0。
第23頁,共85頁,2023年,2月20日,星期一算例>>year=1900:10:2010;>>product=[75.995,91.972,105.711,123.203,131.669,...150.697,179.323,203.212,226.505,249.633,256.344,267.893];>>p1995=interp1(year,product,1995)p1995=252.9885>>x=1900:1:2010;>>y=interp1(year,…product,x,'cubic');>>plot(year,…product,'o',x,y)第24頁,共85頁,2023年,2月20日,星期一例:已知的數(shù)據(jù)點來自函數(shù)
根據(jù)生成的數(shù)據(jù)進行插值處理,得出較平滑的曲線
直接生成數(shù)據(jù)。>>x=0:.12:1;>>y=(x.^2-3*x…+5).*exp(-5*x…).*sin(x);>>plot(x,y,x,y,'o')第25頁,共85頁,2023年,2月20日,星期一調(diào)用interp1()函數(shù):>>x1=0:.02:1;y0=(x1.^2-3*x1+5).*exp(-5*x1).*sin(x1);>>y1=interp1(x,y,x1);y2=interp1(x,y,x1,'cubic');>>y3=interp1(x,y,x1,'spline');y4=interp1(x,y,x1,'nearest');>>plot(x1,[y1',y2',y3',y4'],':',x,y,'o',x1,y0)誤差分析>>[max(abs(y0(1:49)…-y2(1:49))),max(abs(y0-y3)),max(abs(y0-y4))]ans=0.01770.00860.1598第26頁,共85頁,2023年,2月20日,星期一>>x0=-1+2*[0:10]/10;>>y0=1./(1+25*x0.^2);>>x=-1:.01:1;>>y=lagrange(x0,y0,x);%Lagrange插值>>ya=1./(1+25*x.^2);>>plot(x,ya,x,y,':')例第27頁,共85頁,2023年,2月20日,星期一>>y1=interp1(x0,y0,x,'cubic');y2=interp1(x0,y0,x,'spline');>>plot(x,ya,x,y1,':',x,y2,'--')第28頁,共85頁,2023年,2月20日,星期一命令2
interp2
功能
二維數(shù)據(jù)內(nèi)插值(表格查找)
格式1
ZI=interp2(X,Y,Z,XI,YI)
%返回矩陣ZI,其元素包含對應(yīng)于參量XI與YI(可以是向量、或同型矩陣)的元素。參量X與Y必須是單調(diào)的,且相同的劃分格式,就像由命令meshgrid生成的一樣。若Xi與Yi中有在X與Y范圍之外的點,則相應(yīng)地返回NaN。格式2
ZI=interp2(Z,XI,YI)%缺省地,X=1:n、Y=1:m,其中[m,n]=size(Z)。再按第一種情形進行計算。格式3ZI=interp2(X,Y,Z,XI,YI,method)%用指定的算法method計算二維插值:’linear’:雙線性插值算法(缺省算法);’nearest’:最臨近插值;’spline’:三次樣條插值;’cubic’:雙三次插值。第29頁,共85頁,2023年,2月20日,星期一算例:>>years=1950:10:1990;>>service=10:10:30;>>wage=[150.697199.592187.625179.323195.072250.287203.212179.092322.767226.505153.706426.730249.633120.281598.243];>>w=interp2(service,years,wage,15,1975)w=190.6288第30頁,共85頁,2023年,2月20日,星期一例>>[x,y]=meshgrid(-3:.6:3,-2:.4:2);>>z=(x.^2-2*x).*exp…(-x.^2-y.^2-x.*y);>>surf(x,y,z),axis([-3,3,-2,2,-0.7,1.5])第31頁,共85頁,2023年,2月20日,星期一選較密的插值點,用默認的線性插值算法進行插值>>[x1,y1]=meshgrid(-3:.2:3,-2:.2:2);>>z1=interp2(x,y,z,x1,y1);>>surf(x1,y1,z1),axis([-3,3,-2,2,-0.7,1.5])第32頁,共85頁,2023年,2月20日,星期一立方和樣條插值:>>z1=interp2(x,y,z,x1,y1,'cubic');>>z2=interp2(x,y,z,x1,y1,'spline');>>surf(x1,y1,z1),axis([-3,3,-2,2,-0.7,1.5])>>figure;surf(x1,y1,z2),axis([-3,3,-2,2,-0.7,1.5])第33頁,共85頁,2023年,2月20日,星期一算法誤差的比較>z=(x1.^2-2*x1).*exp(-x1.^2-y1.^2-x1.*y1);>>surf(x1,y1,abs(z-z1)),axis([-3,3,-2,2,0,0.08])>>figure;surf(x1,y1,abs(z-z2)),axis([-3,3,-2,2,0,0.025])第34頁,共85頁,2023年,2月20日,星期一二維一般分布數(shù)據(jù)的插值功能:可對非網(wǎng)格數(shù)據(jù)進行插值格式:z=griddata(x0,y0,z0,x,y,’method’)’v4’:MATLAB4.0提供的插值算法,公認效果較好;’linear’:雙線性插值算法(缺省算法);’nearest’:最臨近插值;’spline’:三次樣條插值;’cubic’:雙三次插值。例:在x為[-3,3],y為[-2,2]矩形區(qū)域隨機選擇一組坐標,用’v4’與’cubic’插值法進行處理,并對誤差進行比較。第35頁,共85頁,2023年,2月20日,星期一>>x=-3+6*rand(200,1);y=-2+4*rand(200,1);>>z=(x.^2-2*x).*exp(-x.^2-y.^2-x.*y);>>[x1,y1]=meshgrid(-3:.2:3,-2:.2:2);>>z1=griddata(x,y,z,x1,y1,'cubic');>>surf(x1,y1,z1),axis([-3,3,-2,2,-0.7,1.5])>>z2=griddata(x,y,z,x1,y1,'v4');>>figure;surf(x1,y1,z2),axis([-3,3,-2,2,-0.7,1.5])第36頁,共85頁,2023年,2月20日,星期一誤差分析>>z0=(x1.^2-2*x1).*exp(-x1.^2-y1.^2-x1.*y1);>>surf(x1,y1,abs(z0-z1)),axis([-3,3,-2,2,0,0.15])>>figure;surf(x1,y1,abs(z0-z2)),axis([-3,3,-2,2,0,0.15])第37頁,共85頁,2023年,2月20日,星期一例:
在x為[3,3],y為[-2,2]矩形區(qū)域隨機選擇一組坐標中,對分布不均勻數(shù)據(jù),進行插值分析。>>x=-3+6*rand(200,1);y=-2+4*rand(200,1);>>z=(x.^2-2*x).*exp(-x.^2-y.^2-x.*y);%生成已知數(shù)據(jù)>>plot(x,y,'x')%樣本點的二維分布>>figure,plot3(x,y,z,'x'),axis([-3,3,-2,2,-0.7,1.5]),grid第38頁,共85頁,2023年,2月20日,星期一去除在(-1,-1/2)點為圓心,以0.5為半徑的圓內(nèi)的點。>>x=-3+6*rand(200,1);y=-2+4*rand(200,1);%重新生成樣本點>>z=(x.^2-2*x).*exp(-x.^2-y.^2-x.*y);>>ii=find((x+1).^2+(y+0.5).^2>0.5^2);%找出滿足條件的點坐標>>x=x(ii);y=y(ii);z=z(ii);plot(x,y,'x')>>t=[0:.1:2*pi,2*pi];x0=-1+0.5*cos(t);y0=-0.5+0.5*sin(t);>>line(x0,y0)%在圖形上疊印該圓,可見,圓內(nèi)樣本點均已剔除第39頁,共85頁,2023年,2月20日,星期一用新樣本點擬合出曲面>>[x1,y1]=meshgrid(-3:.2:3,-2:.2:2);>>z1=griddata(x,y,z,x1,y1,'v4');>>surf(x1,y1,z1),axis([-3,3,-2,2,-0.7,1.5])第40頁,共85頁,2023年,2月20日,星期一誤差分析>>z0=(x1.^2-2*x1).*exp(-x1.^2-y1.^2-x1.*y1);>>surf(x1,y1,abs(z0-z1)),axis([-3,3,-2,2,0,0.1])>>contour(x1,y1,abs(z0-z1),30);holdon,plot(x,y,‘x’);line(x0,y0)%誤差的二維等高線圖第41頁,共85頁,2023年,2月20日,星期一命令3
interp3
三維網(wǎng)格生成用meshgrid()函數(shù),調(diào)用格式:[x,y,z]=meshgrid(x1,y1,z1)
其中x1,y1,z1為這三維所需要的分割形式,應(yīng)以向量形式給出,返回x,y,z為網(wǎng)格的數(shù)據(jù)生成,均為三維數(shù)組。
griddata3()三維非網(wǎng)格形式的插值擬合命令4
interpnn維網(wǎng)格生成用ndgrid()函數(shù),調(diào)用格式:
[x1,x2,…,xn]=ndgrid[v1,v2,…,vn]griddatan()n維非網(wǎng)格形式的插值擬合interp3()、interpn()調(diào)用格式同interp2()函數(shù)一致;griddata3()、griddatan()調(diào)用格式同griddata()函數(shù)一致。第42頁,共85頁,2023年,2月20日,星期一例:通過函數(shù)生成一些網(wǎng)格型樣本點,試根據(jù)樣本點進行擬合,并給出擬合誤差。>>[x,y,z]=meshgrid(-1:0.2:1);[x0,y0,z0]=meshgrid(-1:0.05:1);>>V=exp(x.^2.*z+y.^2.*x+z.^2.*y…).*cos(x.^2.*y.*z+z.^2.*y.*x);>>V0=exp(x0.^2.*z0+y0.^2.*x0…+z0.^2.*y0).*cos(x0.^2.*y0.*z0+z0.^2.*y0.*x0);>>V1=interp3(x,y,z,V,x0,y0,z0,'spline');err=V1-V0;max(err(:))ans=0.0419第43頁,共85頁,2023年,2月20日,星期一5.2.5樣條插值的MATLAB表示第44頁,共85頁,2023年,2月20日,星期一定義一個三次樣條函數(shù)類:
S=csapi(x,y)
其中x=[x1,x2,….,xn],y=[y1,y2,…,yn]為樣本點。S返回樣條函數(shù)對象的插值結(jié)果,包括子區(qū)間點、各區(qū)間點三次多項式系數(shù)等??捎胒nplt()繪制出插值結(jié)果,其調(diào)用格式:
fnplt(S)對給定的向量xp,可用fnval()函數(shù)計算,其調(diào)用格式:
yp=fnval(S,xp)其中得出的yp是xp上各點的插值結(jié)果。第45頁,共85頁,2023年,2月20日,星期一例:>>x0=[0,0.4,1,2,pi];y0=sin(x0);>>sp=csapi(x0,y0),fnplt(sp,':');holdon,sp=form:'pp'breaks:[00.4000123.1416]coefs:[4x4double]pieces:4order:4dim:1>>ezplot('sin(t)',[0,pi]);plot(x0,y0,'o')>>sp.coefsans=-0.16270.00760.99650-0.1627-0.18760.92450.38940.0244-0.48040.52380.84150.0244-0.4071-0.36370.9093第46頁,共85頁,2023年,2月20日,星期一在(0.4000,1)區(qū)間內(nèi),插值多項式可以表示為:第47頁,共85頁,2023年,2月20日,星期一例點,用三次樣條插值的方法對這些數(shù)據(jù)進行擬合>>x=0:.12:1;y=(x.^2-3*x+5).*exp(-5*x).*sin(x);>>sp=csapi(x,y);fnplt(sp)>>c=[sp.breaks(1:4)'sp.breaks(2:5)'sp.coefs(1:4,:),...sp.breaks(5:8)'sp.breaks(6:9)'sp.coefs(5:8,:)]c=Columns1through700.120024.7396-19.35884.515100.48000.12000.240024.7396-10.45260.93770.30580.60000.24000.36004.5071-1.5463-0.50220.31050.72000.36000.48001.91390.0762-0.67860.23580.8400
第48頁,共85頁,2023年,2月20日,星期一Columns8through120.6000-0.24040.7652-0.57760.15880.7200-0.47740.6787-0.40430.10010.8400-0.45590.5068-0.26210.06050.9600-0.45590.3427-0.16010.0356第49頁,共85頁,2023年,2月20日,星期一格式
S=csapi({x1,x2,…,xn},z)處理多個自變量的網(wǎng)格數(shù)據(jù)三次樣條插值類:第50頁,共85頁,2023年,2月20日,星期一>>x0=-3:.6:3;y0=-2:.4:2;[x,y]=ndgrid(x0,y0);%注意這里只能用ndgrid,否則生成的z矩陣順序有問題>>z=(x.^2-2*x).*…exp(-x.^2-y.^2-x.*y);>>sp=csapi({x0,y0},z);>>fnplt(sp);例第51頁,共85頁,2023年,2月20日,星期一函數(shù)spline功能三次樣條數(shù)據(jù)插值格式
yy=spline(x,y,xx)
例:對離散分布在y=exp(x)sin(x)函數(shù)曲線上的數(shù)據(jù)點進行樣條插值計算:
>>x=[024581212.817.219.920];>>y=exp(x).*sin(x);>>xx=0:.25:20;>>yy=spline(x,y,xx);>>plot(x,y,'o',xx,yy)第52頁,共85頁,2023年,2月20日,星期一第53頁,共85頁,2023年,2月20日,星期一5.3數(shù)據(jù)擬合用插值的方法對一函數(shù)進行近似,要求所得到的插值多項式經(jīng)過已知插值節(jié)點;在n比較大的情況下,插值多項式往往是高次多項式,這也就容易出現(xiàn)振蕩現(xiàn)象(龍格現(xiàn)象),即雖然在插值節(jié)點上沒有誤差,但在插值節(jié)點之外插值誤差變得很大,從“整體”上看,插值逼近效果將變得“很差”。所謂數(shù)據(jù)擬合是求一個簡單的函數(shù),例如是一個低次多項式,不要求通過已知的這些點,而是要求在整體上“盡量好”的逼近原函數(shù)。這時,在每個已知點上就會有誤差,數(shù)據(jù)擬合就是從整體上使誤差,盡量的小一些。第54頁,共85頁,2023年,2月20日,星期一5.3.1多項式擬合n次多項式:曲線與數(shù)據(jù)點的殘差為:殘差的平方和為:為使其最小化,可令R關(guān)于的偏導數(shù)為零,即:第55頁,共85頁,2023年,2月20日,星期一或或矩陣形式:第56頁,共85頁,2023年,2月20日,星期一多項式擬合MATLAB命令:polyfit
格式:p=polyfit(x,y,n)第57頁,共85頁,2023年,2月20日,星期一>>x0=0:.1:1;y0=(x0.^2-3*x0+5).*exp(-5*x0).*sin(x0);>>p3=polyfit(x0,y0,3);vpa(poly2sym(p3),10)%可以如下顯示多項式ans=2.839962923*x^3-4.789842696*x^2+1.943211631*x+.5975248921e-1例第58頁,共85頁,2023年,2月20日,星期一繪制擬合曲線:>>x=0:.01:1;ya=(x.^2-3*x+5).*exp(-5*x).*sin(x);>>y1=polyval(p3,x);plot(x,y1,x,ya,x0,y0,'o')第59頁,共85頁,2023年,2月20日,星期一就不同的次數(shù)進行擬合:>>p4=polyfit(x0,y0,4);y2=polyval(p4,x);>>p5=polyfit(x0,y0,5);y3=polyval(p5,x);>>p8=polyfit(x0,y0,8);y4=polyval(p8,x);>>plot(x,ya,x0,y0,'o',x,y2,x,y3,x,y4)第60頁,共85頁,2023年,2月20日,星期一擬合最高次數(shù)為8的多項式:>>vpa(poly2sym(p8),5)ans=-8.2586*x^8+43.566*x^7-101.98*x^6+140.22*x^5-125.29*x^4+74.450*x^3-27.672*x^2+4.9869*x+.42037e-6Taylor冪級數(shù)展開:>>symsx;y=(x^2-3*x+5)*exp(-5*x)*sin(x);>>vpa(taylor(y,9),5)ans=5.*x-28.*x^2+77.667*x^3-142.*x^4+192.17*x^5-204.96*x^6+179.13*x^7-131.67*x^8多項式表示數(shù)據(jù)模型是不唯一的,即是兩個多項式函數(shù)完全不同。在某一區(qū)域內(nèi)其曲線可能特別近似。第61頁,共85頁,2023年,2月20日,星期一多項式擬合的效果并不一定總是很精確的。>>x0=-1+2*[0:10]/10;y0=1./(1+25*x0.^2);>>x=-1:.01:1;ya=1./(1+25*x.^2);>>p3=polyfit(x0,y0,3);>>y1=polyval(p3,x);>>p5=polyfit(x0,y0,5);>>y2=polyval(p5,x);>>p8=polyfit(x0,y0,8);>>y3=polyval(p8,x);>>p10=polyfit(x0,y0,10);>>y4=polyval(p10,x);>>plot(x,ya,x,y1,x,y2,'-.',x,y3,'--',x,y4,':')例第62頁,共85頁,2023年,2月20日,星期一用Taylor冪級數(shù)展開效果將更差。>>symsx;y=1/(1+25*x^2);p=taylor(y,x,10)p=1-25*x^2+625*x^4-15625*x^6+390625*x^8多項式擬合效果>>x1=-1:0.01:1;>>ya=1./(1+25*x1.^2);>>y1=subs(p,x,x1);>>plot(x1,ya,'--‘,x1,y1)第63頁,共85頁,2023年,2月20日,星期一5.3.2函數(shù)線性組合的曲線擬合方法第64頁,共85頁,2023年,2月20日,星期一該方程的最小二乘解為:其中第65頁,共85頁,2023年,2月20日,星期一例第66頁,共85頁,2023年,2月20日,星期一>>x=[0,0.2,0.4,0.7,0.9,0.92,0.99,1.2,1.4,1.48,1.5]';>>y=[2.88;2.2576;1.9683;1.9258;2.0862;2.109;2.1979;2.5409;2.9627;3.155;3.2052];>>A=[ones(size(x)),exp(-3*x),cos(-2*x).*exp(-4*x),x.^2];>>c=A\y;c1=c'c1=1.22002.3397-0.67970.8700第67頁,共85頁,2023年,2月20日,星期一圖形顯示>>x0=[0:0.01:1.5]';>>A1=[ones(size(x0))exp(-3*x0),cos(-2*x0).*exp(-4*x0)x0.^2];>>y1=A1*c;>>plot(x0,y1,x,y,'x')第68頁,共85頁,2023年,2月20日,星期一數(shù)據(jù)分析>>x=[1.1052,1.2214,1.3499,1.4918,1.6487,1.8221,2.0138,...2.2255,2.4596,2.7183,3.6693];>>y=[0.6795,0.6006,0.5309,0.4693,0.4148,0.3666,0.3241,...0.2864,0.2532,0.2238,0.1546];>>plot(x,y,x,y,'*')例第69頁,共85頁,2023年,2月20日,星期一分別對x,y進行對數(shù)變換:>>x1=log(x);y1=log(y);plot(x1,y1)第70頁,共85頁,2023年,2月20日,星期一
>>A=[x1',ones(size(x1'))];c=[A\y1']‘c=-1.2339-0.2630>>exp(c(2))ans=0.7687第71頁,共85頁,2023年,2月20日,星期一>>x=[0:0.1:1]';y=(x.^2-3*x+5).*exp(-5*x).*sin(x);n=8;A=[];>>fori=1:n+1,A(:,i)=x.^(n+1-i);end>>c=A\y;vpa(poly2sym(c),5)ans=-8.2586*x^8+43.566*x^7-101.98*x^6+140.22*x^5-125.29*x^4+74.450*x^3-27.672*x^2+4.9869*x+.42037e-6例第72頁,共85頁,2023年,2月20日,星期一5.3.3最小二乘曲線擬合第73頁,共85頁,2023年,2月20日,星期一格式
溫馨提示
- 1. 本站所有資源如無特殊說明,都需要本地電腦安裝OFFICE2007和PDF閱讀器。圖紙軟件為CAD,CAXA,PROE,UG,SolidWorks等.壓縮文件請下載最新的WinRAR軟件解壓。
- 2. 本站的文檔不包含任何第三方提供的附件圖紙等,如果需要附件,請聯(lián)系上傳者。文件的所有權(quán)益歸上傳用戶所有。
- 3. 本站RAR壓縮包中若帶圖紙,網(wǎng)頁內(nèi)容里面會有圖紙預(yù)覽,若沒有圖紙預(yù)覽就沒有圖紙。
- 4. 未經(jīng)權(quán)益所有人同意不得將文件中的內(nèi)容挪作商業(yè)或盈利用途。
- 5. 人人文庫網(wǎng)僅提供信息存儲空間,僅對用戶上傳內(nèi)容的表現(xiàn)方式做保護處理,對用戶上傳分享的文檔內(nèi)容本身不做任何修改或編輯,并不能對任何下載內(nèi)容負責。
- 6. 下載文件中如有侵權(quán)或不適當內(nèi)容,請與我們聯(lián)系,我們立即糾正。
- 7. 本站不保證下載資源的準確性、安全性和完整性, 同時也不承擔用戶因使用這些下載資源對自己和他人造成任何形式的傷害或損失。
最新文檔
- 肺結(jié)核患者疼痛管理的觀察與護理策略
- 生活護理學習資料中心
- 跨境電商獨立站域名2025年爭議解決協(xié)議
- 初中政治考試內(nèi)容及答案
- 2025-2026人教版小學二年級語文上冊期末卷子
- 藥理麻醉藥試題及答案
- 2025-2026人教版五年級語文上學期模擬卷
- 腸道膽汁酸代謝與NASH進展
- 寢室衛(wèi)生獎罰制度
- 養(yǎng)老院清潔衛(wèi)生制度
- 2026年上半年眉山天府新區(qū)公開選調(diào)事業(yè)單位工作人員的參考題庫附答案
- 水產(chǎn)養(yǎng)殖技術(shù)手冊
- 英國汽車工業(yè)市場分析現(xiàn)狀供需格局投資前景未來規(guī)劃研究報告
- 2025年及未來5年市場數(shù)據(jù)中國吸塑、注塑行業(yè)發(fā)展前景預(yù)測及投資戰(zhàn)略數(shù)據(jù)分析研究報告
- 眼科醫(yī)療風險防范培訓
- 物流金融理論與實務(wù)課件
- 海內(nèi)外云廠商發(fā)展與現(xiàn)狀(三):資本開支壓力與海外云廠需求情況拆解-國信證券
- 2025年社區(qū)網(wǎng)格員招錄考試真題庫(含答案)
- GB/T 46510-2025玩具水基材料中游離甲醛的測定高效液相色譜法
- 溴化鋰清洗施工方案
- 第四方支付業(yè)務(wù)合規(guī)指引
評論
0/150
提交評論