數(shù)學(xué)實(shí)驗(yàn)6-1多項(xiàng)式、插值及數(shù)據(jù)擬合

1、第六講 多項(xiàng)式、插值與數(shù)據(jù)擬合（上） 多項(xiàng)式MATLAB命令 插值 Lagrange插值 Hermite插值 Runge現(xiàn)象和分段插值 分段插值 樣條插值的MATLAB表示 數(shù)據(jù)擬合多項(xiàng)式擬合函數(shù)線性組合的曲線擬合方法最小二乘曲線擬合B樣條函數(shù)及其MATLAB表示5.1 關(guān)于多項(xiàng)式MATLAB命令 一個(gè)多項(xiàng)式的冪級(jí)數(shù)形式可表示為： 也可表為嵌套形式 或因子形式 N階多項(xiàng)式n個(gè)根，其中包含重根和復(fù)根。若多項(xiàng)式所有系數(shù)均為實(shí)數(shù)，則全部復(fù)根都將以共軛對(duì)的形式出現(xiàn) 1121nnnnyc xc xc xc1231()nnyc xc xc xc xc112()()()nyc xrxrxr 冪系數(shù)：在MAT

2、LAB里，多項(xiàng)式用行向量表示，其元素為多項(xiàng)式的系數(shù)，并從左至右按降冪排列。 例： 被表示為 p=2 1 4 5 poly2sym(p) ans = 2*x3+x2+4*x+5 Roots: 多項(xiàng)式的零點(diǎn)可用命令roots求的。 例： r=roots(p) 得到 r = 0.2500 + 1.5612i 0.2500 - 1.5612i -1.0000 所有零點(diǎn)由一個(gè)列向量給出。32245yxxx Poly: 由零點(diǎn)可得原始多項(xiàng)式的各系數(shù)，但可能相差一個(gè)常數(shù)倍。 例： poly(r)ans = 1.0000 0.5000 2.0000 2.5000 注意：若存在重根，這種轉(zhuǎn)換可能會(huì)降低精度。 例

3、： r=roots(1 -6 15 -20 15 -6 1) r = 1.0042 + 0.0025i 1.0042 - 0.0025i 1.0000 + 0.0049i 1.0000 - 0.0049i 0.9958 + 0.0024i 0.9958 - 0.0024i舍入誤差的影響，與計(jì)算精度有關(guān)。665432(1)615201561yxxxxxxx polyval: 可用命令polyval計(jì)算多項(xiàng)式的值。 例： 計(jì)算y(2.5) c=3,-7,2,1,1； xi=2.5; yi=polyval(c,xi) yi = 23.8125如果xi是含有多個(gè)橫坐標(biāo)值的數(shù)組，則yi也為與xi長(zhǎng)度相同

4、的向量。 c=3,-7,2,1,1; xi=2.5,3; yi=polyval(c,xi)yi = 23.8125 76.00004323721yxxxx polyfit:給定n+1個(gè)點(diǎn)將可以唯一確定一個(gè)n階多項(xiàng)式。利用命令polyfit可容易確定多項(xiàng)式的系數(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.2015 1.4385 -2.7477 5.4370 poly2sym(a) ans = -403/2000*x3+2877/2000*x2-27477/10000*x+5

5、437/1000 多項(xiàng)式為Polyfit的第三個(gè)參數(shù)是多項(xiàng)式的階數(shù)。320.20151.43852.74775.4370yxxx多項(xiàng)式積分： 功能：求多項(xiàng)式積分 調(diào)用格式：py=poly_itg(p) p:被積多項(xiàng)式的系數(shù) py：求積后多項(xiàng)式的系數(shù) poly_itg.m function py=poly_itg(p) n=length(p); py=p.*n:-1:1.(-1),0不包括最后一項(xiàng)積分常數(shù)12121212nnnnncccYydxxxxcxcnn1121nnnnyc xc xc xc多項(xiàng)式微分： Polyder: 求多項(xiàng)式一階導(dǎo)數(shù)的系數(shù)。 調(diào)用格式為: b=polyder(c )

6、c為多項(xiàng)式y(tǒng)的系數(shù)，b是微分后的系數(shù)，其值為： 1121nnnnyc xc xc xc1212(1)nnnync xnc xc12,(1),nncncc兩個(gè)多項(xiàng)式的和與差： 命令poly_add:求兩個(gè)多項(xiàng)式的和,其調(diào)用格式為: c= poly_add(a,b) 多項(xiàng)式a減去b，可表示為: c= poly_add(a,-b)1121nnbnnyb xb xb xb1121mmammya xa xa xa 功能：兩個(gè)多項(xiàng)式相加 調(diào)用格式：b=poly_add(p1,p2) b:求和后的系數(shù)數(shù)組 poly_add.mfunction p3=poly_add(p1,p2)n1=length(p1);

7、n2=length(p2);if n1=n2 p3=p1+p2;endif n1n2 p3=p1+zeros(1,n1-n2),p2;endif n1 a=2,-5,6,-1,9; b=3,-90,-18; c=conv(a,b)c = 6 -195 432 -453 9 -792 -162 q,r=deconv(c,b)q = 2 -5 6 -1 9r = 0 0 0 0 0 0 0 poly2sym(c) ans = 6*x6-195*x5+432*x4-453*x3+9*x2-792*x-162 5.2 插值5.2.1 Lagrange插值 方法介紹 對(duì)給定的n個(gè)插值點(diǎn) 及對(duì)應(yīng)的函數(shù)值

8、，利用構(gòu)造的n-1次Lagrange插值多項(xiàng)式，則對(duì)插值區(qū)間內(nèi)任意x的函數(shù)值y可通過(guò)下式求的： MATLAB實(shí)現(xiàn)12,nx xx12,ny yy11( )()nnjkkjkjjkxxy xyxxfunction y=lagrange(x0,y0,x)ii=1:length(x0); y=zeros(size(x);for i=ii ij=find(ii=i); y1=1; for j=1:length(ij), y1=y1.*(x-x0(ij(j); end y=y+y1*y0(i)/prod(x0(i)-x0(ij);end 算例：給出f(x)=ln(x)的數(shù)值表，用Lagrange計(jì)算ln

9、(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)11( )()nnjkkjkjj kxxy xyxx5.2.2 Hermite插值 方法介紹 不少實(shí)際問(wèn)題不但要求在節(jié)點(diǎn)上函數(shù)值相等，而且要求導(dǎo)數(shù)值也相等，甚至要求高階導(dǎo)數(shù)值也相等，滿足這一要求的插值多項(xiàng)式就是Hermite插值多項(xiàng)式。下面只討論函數(shù)值與一階導(dǎo)數(shù)值個(gè)數(shù)相等且已知的情況。 已知n

10、個(gè)插值點(diǎn) 及對(duì)應(yīng)的函數(shù)值 和一階導(dǎo)數(shù)值 。則對(duì)插值區(qū)間內(nèi)任意x的函數(shù)值y的Hermite插值公式：12,nx xx12,ny yy12,ny yy1211( )()(2)1 () ;niiiiiiinnjiijjijijjijiy xhxxa yyyxxhaxxxx其 中 MATLAB實(shí)現(xiàn)% hermite.mfunction y=hermite(x0,y0,y1,x)n=length(x0); m=length(x);for k=1:m yy=0.0; for i=1:n h=1.0; a=0.0; for j=1:n if j=i h=h*(x(k)-x0(j)/(x0(i)-x0(j)2

11、; a=1/(x0(i)-x0(j)+a; end end yy=yy+h*(x0(i)-x(k)*(2*a*y0(i)-y1(i)+y0(i); end y(k)=yy;end1211( )()(2)1 () ;niiiiiiinnjiijjijijj ij iy xh xxayyyx xhaxxxx其中 算例:對(duì)給定數(shù)據(jù),試構(gòu)造Hermite多項(xiàng)式求出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

12、 sin(0.34) 與精確值比較ans = 0.3335 x=0.3:0.005:0.35;y=hermite(x0,y0,y1,x); plot(x,y) y2=sin(x); hold on plot(x,y2,-r)5.2.3 Runge現(xiàn)象 問(wèn)題的提出：根據(jù)區(qū)間a,b上給出的節(jié)點(diǎn)做插值多項(xiàng)式p(x)的近似值，一般總認(rèn)為p(x)的次數(shù)越高則逼近f(x)的精度就越好，但事實(shí)并非如此。 反例： 在區(qū)間-5,5上的各階導(dǎo)數(shù)存在，但在此區(qū)間上取n個(gè)節(jié)點(diǎn)所構(gòu)成的Lagrange插值多項(xiàng)式在全區(qū)間內(nèi)并非都收斂。 取n=10,用Lagrange插值法進(jìn)行插值計(jì)算。21()1fxx x=-5:1:5;

13、 y=1./(1+x.2); x0=-5:0.1:5; y0=lagrange(x,y,x0); y1=1./(1+x0.2);繪制圖形 plot(x0,y0,-r)插值曲線 hold on plot(x0,y1,-b)原曲線 為解決Rung問(wèn)題，引入分段插值。 算法分析：所謂分段插值就是通過(guò)插值點(diǎn)用折線或低次曲線連接起來(lái)逼近原曲線。 MATLAB實(shí)現(xiàn) 可調(diào)用內(nèi)部函數(shù)。 命令1 interp1 功能 ： 一維數(shù)據(jù)插值（表格查找）。該命令對(duì)數(shù)據(jù)點(diǎn)之間計(jì)算內(nèi)插值。它找出一元函數(shù)f(x)在中間點(diǎn)的數(shù)值。其中函數(shù)f(x)由所給數(shù)據(jù)決定。 格式1 yi = interp1(x,Y,xi) %返回插值向量

14、yi，每一元素對(duì)應(yīng)于參量xi，同時(shí)由向量x與Y的內(nèi)插值決定。參量x指定數(shù)據(jù)Y的點(diǎn)。若Y為一矩陣，則按Y的每列計(jì)算。 算例 對(duì)于t，beta 、alpha分別有兩組數(shù)據(jù)與之對(duì)應(yīng)，用分段線性插值法計(jì)算當(dāng)t=321, 440, 571時(shí)beta 、alpha的值。5.2.4 分段插值 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(

15、temp,beta,alpha,ti ,linear); ti,proptyans = 1.0e+003 * 0.3210 3.1557 2.4382 0.4000 2.5000 3.6050 0.5710 1.7657 6.6489 格式2 yi = interp1(Y,xi) %假定x=1:N，其中N為向量Y的長(zhǎng)度，或者為矩陣Y的行數(shù)。 格式3 yi = interp1(x,Y,xi,method) %用指定的算法計(jì)算插值： nearest：最近鄰點(diǎn)插值，直接完成計(jì)算； linear：線性插值（缺省方式），直接完成計(jì)算； spline：三次樣條函數(shù)插值。 cubic： 分段三次Hermit

16、e插值。 其它，如v5cubic 。 對(duì) 于 超 出 x 范 圍 的 x i 的 分 量 ， 使 用 方法nearest、linear、v5cubic的插值算法，相應(yīng)地將返回NaN。對(duì)其他的方法，interp1將對(duì)超出的分量執(zhí)行外插值算法。 yi = interp1(x,Y,xi,method,extrap) yi = interp1(x,Y,xi,method,extrapval) %確定超出x范圍的xi中的分量的外插值extrapval，其值通常取NaN或0。 算例 year=1900:10:2010; product=75.995,91.972,105.711,123.203,131.6

17、69,. 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) 例:已知的數(shù)據(jù)點(diǎn)來(lái)自函數(shù)已知的數(shù)據(jù)點(diǎn)來(lái)自函數(shù)根據(jù)生成的數(shù)據(jù)進(jìn)行插值處理，得出較平滑的曲線根據(jù)生成的數(shù)據(jù)進(jìn)行插值處理，得出較平滑的曲線直接生成數(shù)據(jù)。直接生成數(shù)據(jù)。 x=0:.12:1; y=(x.2-3*x+5).*exp

18、(-5*x).*sin(x); plot(x,y,x,y,o) 調(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

19、= 0.0177 0.0086 0.1598 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,:) 例 y1=interp1(x0,y0,x,cubic); y2=interp1(x0,y0,x,spline); plot(x,ya,x,y1,:,x,y2,-) 命令2 interp2 功能 二維數(shù)據(jù)內(nèi)插值（表格查找） 格式1 ZI = interp2(X,Y,Z,XI,YI) %返回矩陣ZI，其元素包含對(duì)應(yīng)于參量

20、XI與YI（可以是向量、或同型矩陣）的元素。參量X與Y必須是單調(diào)的，且相同的劃分格式，就像由命令meshgrid生成的一樣。若Xi與Yi中有在X與Y范圍之外的點(diǎn)，則相應(yīng)地返回NaN。 格式2 ZI = interp2(Z,XI,YI) %缺省地，X=1:n、Y=1:m，其中m,n=size(Z)。再按第一種情形進(jìn)行計(jì)算。 格式3 ZI = interp2(X,Y,Z,XI,YI,method) %用指定的算法method計(jì)算二維插值： linear：雙線性插值算法（缺省算法）； nearest：最臨近插值； spline：三次樣條插值； cubic：雙三次插值。 算例： years=1950:

21、10:1990; service=10:10:30; wage = 150.697 199.592 187.625 179.323 195.072 250.287 203.212 179.092 322.767 226.505 153.706 426.730 249.633 120.281 598.243; w = interp2(service,years,wage,15,1975)w = 190.6288 例 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

22、,-0.7,1.5) 選較密的插值點(diǎn)，用默認(rèn)的線性插值算法進(jìn)行插值 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) 立方和樣條插值： 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) 算法誤差的比較 z

23、=(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)二維一般分布數(shù)據(jù)的插值 功能:可對(duì)非網(wǎng)格數(shù)據(jù)進(jìn)行插值 格式:z=griddata(x0,y0,z0,x,y,method) v4 :MATLAB4.0提供的插值算法，公認(rèn)效果較好； linear：雙線性插值算法（缺省算法）； nearest：最臨近插值； spline：三次樣條插值； cubic：雙三次插值。 例： 在x為-3，3，

24、y為2，2矩形區(qū)域隨機(jī)選擇一組坐標(biāo)，用 v4 與cubic插值法進(jìn)行處理，并對(duì)誤差進(jìn)行比較。222( , )(2 )xyxyzf x yxx e 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)

25、,axis(-3,3,-2,2,-0.7,1.5)誤差分析 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) 例： 在x為3，3，y為2，2矩形區(qū)域隨機(jī)選擇一組坐標(biāo)中，對(duì)分布不均勻數(shù)據(jù)，進(jìn)行插值分析。 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(

26、x,y,x) % 樣本點(diǎn)的二維分布 figure, plot3(x,y,z,x), axis(-3,3,-2,2,-0.7,1.5)，grid222( , )(2 )xyxyzf x yxx e 去除在(-1,-1/2)點(diǎn)為圓心，以0.5為半徑的圓內(nèi)的點(diǎn)。 x=-3+6*rand(200,1); y=-2+4*rand(200,1); % 重新生成樣本點(diǎn) z=(x.2-2*x).*exp(-x.2-y.2-x.*y); ii=find(x+1).2+(y+0.5).20.52); % 找出滿足條件的點(diǎn)坐標(biāo) x=x(ii); y=y(ii); z=z(ii); plot(x,y,x) t=0:.

27、1:2*pi,2*pi; x0=-1+0.5*cos(t); y0=-0.5+0.5*sin(t); line(x0,y0) % 在圖形上疊印該圓，可見(jiàn)，圓內(nèi)樣本點(diǎn)均已剔除 用新樣本點(diǎn)擬合出曲面 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) 誤差分析 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(

28、x1,y1,abs(z0-z1),30); hold on, plot(x,y,x); line(x0,y0) 誤差的二維等高線圖 命令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 interpn n維網(wǎng)格生成用ndgrid( )函數(shù),調(diào)用格式: x1,x2,xn=ndgridv1,v2,vn griddatan( ) n維非網(wǎng)格形式的插值擬合in

29、terp3 ( )、 interpn( )調(diào)用格式同interp2( )函數(shù)一致；函數(shù)一致；griddata3( )、 griddatan( )調(diào)用格式同griddata( )函數(shù)函數(shù)一致。一致。 例： 通過(guò)函數(shù)生成一些網(wǎng)格型樣本點(diǎn)，試根據(jù)樣本點(diǎn)進(jìn)行擬合，并給出擬合誤差。 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

30、+z0.2.*y0.*x0); V1=interp3(x,y,z,V,x0,y0,z0,spline); err=V1-V0; max(err(:)ans = 0.041922222( , , )cos()x z y x z yV x y zex yzz yx5.2.5樣條插值的MATLAB表示定義一個(gè)三次樣條函數(shù)類(lèi)： S=csapi(x,y) 其中x=x1,x2,.,xn, y=y1,y2,yn為樣本點(diǎn)。S返回樣條函數(shù)對(duì)象的插值結(jié)果，包括子區(qū)間點(diǎn)、各區(qū)間點(diǎn)三次多項(xiàng)式系數(shù)等。 可用 fnplt()繪制出插值結(jié)果，其調(diào)用格式： fnplt(S) 對(duì)給定的向量xp，可用fnval()函數(shù)計(jì)算，其調(diào)

31、用格式： yp=fnval(S,xp)其中得出的yp是xp上各點(diǎn)的插值結(jié)果。 例： x0=0,0.4,1,2,pi; y0=sin(x0); sp=csapi(x0,y0), fnplt(sp,:); hold on,sp = form: pp breaks: 0 0.4000 1 2 3.1416 coefs: 4x4 double pieces: 4 order: 4 dim: 1 ezplot(sin(t),0,pi); plot(x0,y0,o) sp.coefsans = -0.1627 0.0076 0.9965 0 -0.1627 -0.1876 0.9245 0.3894 0.

32、0244 -0.4804 0.5238 0.8415 0.0244 -0.4071 -0.3637 0.9093在(0.4000, 1)區(qū)間內(nèi)，插值多項(xiàng)式可以表示為：322( )0.1627(0.4)0.1876(0.4)0.9245(0.4)0.3894S xxxx 例點(diǎn)，用三次樣條插值的方法對(duì)這些數(shù)據(jù)進(jìn)行擬合 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.brea

33、ks(6:9) sp.coefs(5:8,:) c = Columns 1 through 7 0 0.1200 24.7396 -19.3588 4.5151 0 0.4800 0.1200 0.2400 24.7396 -10.4526 0.9377 0.3058 0.6000 0.2400 0.3600 4.5071 -1.5463 -0.5022 0.3105 0.7200 0.3600 0.4800 1.9139 0.0762 -0.6786 0.2358 0.8400 Columns 8 through 12 0.6000 -0.2404 0.7652 -0.5776 0.1588

34、 0.7200 -0.4774 0.6787 -0.4043 0.1001 0.8400 -0.4559 0.5068 -0.2621 0.0605 0.9600 -0.4559 0.3427 -0.1601 0.0356 格式 S=csapi(x1,x2,xn,z)處理多個(gè)自變量的網(wǎng)格數(shù)據(jù)三次樣條插值類(lèi)： x0=-3:.6:3; y0=-2:.4:2; x,y=ndgrid(x0,y0); % 注意這里只能用 ndgrid,否則生成的 z 矩陣順序有問(wèn)題 z=(x.2-2*x).*exp(-x.2-y.2-x.*y); sp=csapi(x0,y0,z); fnplt(sp);例函數(shù)函數(shù)sp

35、line 功能功能 三次樣條數(shù)據(jù)插值三次樣條數(shù)據(jù)插值 格式格式 yy = spline(x,y,xx) 例:對(duì)離散分布在y=exp(x)sin(x)函數(shù)曲線上的數(shù)據(jù)點(diǎn)進(jìn)行樣條插值計(jì)算： x = 0 2 4 5 8 12 12.8 17.2 19.9 20; y = exp(x).*sin(x); xx = 0:.25:20; yy = spline(x,y,xx); plot(x,y,o,xx,yy)5.3 數(shù)據(jù)擬合 用插值的方法對(duì)一函數(shù)進(jìn)行近似用插值的方法對(duì)一函數(shù)進(jìn)行近似,要求所得到的要求所得到的插值多項(xiàng)式經(jīng)過(guò)已知插值節(jié)點(diǎn)插值多項(xiàng)式經(jīng)過(guò)已知插值節(jié)點(diǎn);在在n比較大的情比較大的情況下況下,插值多

36、項(xiàng)式往往是高次多項(xiàng)式插值多項(xiàng)式往往是高次多項(xiàng)式,這也就容這也就容易出現(xiàn)振蕩現(xiàn)象（龍格現(xiàn)象），即雖然在插值易出現(xiàn)振蕩現(xiàn)象（龍格現(xiàn)象），即雖然在插值節(jié)點(diǎn)上沒(méi)有誤差節(jié)點(diǎn)上沒(méi)有誤差,但在插值節(jié)點(diǎn)之外插值誤差變但在插值節(jié)點(diǎn)之外插值誤差變得很大得很大,從從“整體整體”上看上看,插值逼近效果將變得插值逼近效果將變得“很差很差”。 所謂數(shù)據(jù)擬合是求一個(gè)簡(jiǎn)單的函數(shù)所謂數(shù)據(jù)擬合是求一個(gè)簡(jiǎn)單的函數(shù),例如是一個(gè)例如是一個(gè)低次多項(xiàng)式低次多項(xiàng)式,不要求通過(guò)已知的這些點(diǎn)不要求通過(guò)已知的這些點(diǎn),而是要而是要求在整體上求在整體上“盡量好盡量好”的逼近原函數(shù)。這時(shí)的逼近原函數(shù)。這時(shí),在在每個(gè)已知點(diǎn)上就會(huì)有誤差每個(gè)已知點(diǎn)上就會(huì)有誤

37、差,數(shù)據(jù)擬合就是從整體數(shù)據(jù)擬合就是從整體上使誤差上使誤差,盡量的小一些。盡量的小一些。5.3.1 多項(xiàng)式擬合 n次多項(xiàng)式： 曲線與數(shù)據(jù)點(diǎn)的殘差為： 殘差的平方和為： 為使其最小化，可令R關(guān)于 的偏導(dǎo)數(shù)為零，即：1121( )nnng xc xc xc( ),1,2,iiiryg xiL21LiiRrjc0,1,2,1jRjnc 或 或矩陣形式：1221111(),1,2,1nLLnj knkijiijiixcxy kn 221111112111211110111.LLLLnnnniiiiiiiiiLLLnnniiiiiiinLLLniiiiiixxxx ycxxxyccxxy多項(xiàng)式擬合MATL

38、AB命令：polyfit格式：p=polyfit(x,y,n) x0=0:.1:1; y0=(x0.2-3*x0+5).*exp(-5*x0).*sin(x0); p3=polyfit(x0,y0,3); vpa(poly2sym(p3),10) % 可以如下顯示多項(xiàng)式ans =2.839962923*x3-4.789842696*x2+1.943211631*x+.5975248921e-1例例 繪制擬合曲線： 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) 就不同

39、的次數(shù)進(jìn)行擬合： 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) 擬合最高次數(shù)為8的多項(xiàng)式： vpa(poly2sym(p8),5)ans =-8.2586*x8+43.566*x7-101.98*x6+140.22*x5-125.29*x4+74.450*x3-27.672*x2+4.9869*x+.42037e-6 Taylor冪級(jí)數(shù)展開(kāi)：

40、 syms x; y=(x2-3*x+5)*exp(-5*x)*sin(x); vpa(taylor(y,9),5)ans =5.*x-28.*x2+77.667*x3-142.*x4+192.17*x5-204.96*x6+179.13*x7-131.67*x8 多項(xiàng)式表示數(shù)據(jù)模型是不唯一的，即是兩個(gè)多項(xiàng)式函數(shù)完全不同。在某一區(qū)域內(nèi)其曲線可能特別近似。多項(xiàng)式擬合的效果并不一定總是很精確的。 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

41、); 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,:)例例 用Taylor冪級(jí)數(shù)展開(kāi)效果將更差。 syms x; y=1/(1+25*x2); p=taylor(y,x,10)p =1-25*x2+625*x4-15625*x6+390625*x8 多項(xiàng)式擬合效果 x1=-1:0.01:1; ya=1./(1+25*x1.2); y

42、1=subs(p,x,x1); plot(x1,ya,-,x1,y1)5.3.2 函數(shù)線性組合的曲線擬合方法該方程的最小二乘解為:其中例例 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=Ay; c1=cc1 = 1.2200 2.3397 -0.6797 0.8700 圖形顯示 x0=

43、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) 數(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,*)例例 分別對(duì)x,y進(jìn)行對(duì)數(shù)變換： x1=log(x); y1=log(y); plot(x1,y1) A=x1, ones(size(x1); c=Ay1c = -1.2339 -0.2630 exp(c(2)ans = 0.7687 x=0:0.1:1; y=(x.2-3*x+5).*exp(-5*x).*sin(x); n=8; A=; for i=1:n+1, A(:,i)=x.(n+1-i); end c=Ay; vpa(poly2sym(c),5)ans =-8.2586*x