matlab插值計算

1、插值方法晚上做一個曲線擬合，結(jié)果才開始用最小二乘法擬合時，擬合出來的東西太難看了！于是嘗試用其他方法。經(jīng)過一番按圖索驥，終于發(fā)現(xiàn)做曲線擬合的話，采用插值法是比較理想的方法。尤其是樣條插值，插完后線條十分光滑。方法付后，最關(guān)鍵的問題是求解時要積分，放這里想要的時候就可以直接過來拿，不用死去搜索啦。呵呵插值方法的Matlab實現(xiàn)一維數(shù)據(jù)插值MATLAB中用函數(shù)interp1來擬合一維數(shù)據(jù)，語法是YI = INTERP1(X,Y,XI,方法)其中（X， Y） 是已給的數(shù)據(jù)點， XI 是插值點， 其中方法主要有      'linea

2、r'   -線性插值，默認(rèn)      'pchip'   -逐段三次Hermite插值      'spline'   -逐段三次樣條函數(shù)插值其中最后一種插值的曲線比較平滑例：x=0:.12:1; x1=0:.02:1;%（其中x=0:.12:1表示顯示的插值點，x1=0:.02:1表示插值的步長）y=(x.2-3*x+5).*exp(-5*x).*sin(x);plot(x,y,'o&#

3、39;); hold on;y1=interp1(x,y,x1,'spline');plot(x1,y1,':')如果要根據(jù)樣本點求函數(shù)的定積分，而函數(shù)又是比較光滑的，則可以用樣條函數(shù)進行插值后再積分，在MATLAB中可以編寫如下程序：function y=quadspln(x0,y0,a,b)f=inline('interp1(x0,y0,x,''spline'')','x','x0','y0');y=quadl(f,a,b,1e-8,x0,y0); 現(xiàn)求

4、sin(x)在區(qū)間0，pi上的定積分，只取5點x0=0,0.4,1,2,pi;y0=sin(x0);I=quadspln(x0,y0,0,pi)結(jié)果得到的值為 2.01905， 精確值為2求一段matlab插值程序懸賞分：20 - 解決時間：2009-12-26 19:57 已知5個數(shù)據(jù)點：x=0.25 0.5 0.75 1 y=0 0.3104 0.6177 0.7886 1 ，求一段matlab插值程序，求過這5個數(shù)據(jù)點的插值多項式，并在x-y坐標(biāo)中畫出y=f（x）圖形，并且求出f（x）與x軸圍成圖形的面積（積分），不勝感激！ 使用Lagrange 插值多項式的方法：首先把下面的代碼復(fù)制到

5、M文件中，保存成lagranfunction C,L=lagran(X,Y)% input - X is a vector that contains a list of abscissas% - Y is a vector that contains a list of ordinates% output - C is a matrix that contains the coefficients of the lagrange interpolatory polynomial%- L is a matrix that contains the lagrange coefficients p

6、olynomialw=length(X);n=w-1;L=zeros(w,w);for k=1:n+1V=1;for j=1:n+1if k=jV=conv(V,poly(X(j)/(X(k)-X(j);endendL(k,:)=V;endC=Y*L;然后在命令窗口中輸入以下內(nèi)容： x=0 0.25 0.5 0.75 1;y=0 0.3104 0.6177 0.7886 1;lagran(x,y)ans = 3.3088 -6.3851 3.3164 0.7599 0得到的數(shù)據(jù)就是多項式各項的系數(shù)，注意最后一個是常數(shù)項，即x0，所以表達式為:f=3.3088*x.4-6.3851*x.3+3.

7、3164*x.2 +0.7599*x求面積就是積分求解>> f=(x)3.3088*x.4-6.3851*x.3+3.3164*x.2 +0.7599*x;>> quad(f,0,1)ans = 0.5509 這些點肯定是通過這個多項式的！ MATLAB插值與擬合 §1曲線擬合實例：溫度曲線問題氣象部門觀測到一天某些時刻的溫度變化數(shù)據(jù)為：t012345678910T1315171416192624262729試描繪出溫度變化曲線。曲線擬合就是計算出兩組數(shù)據(jù)之間的一種函數(shù)關(guān)系，由此可描繪其變化曲線及估計非采集數(shù)據(jù)對應(yīng)的變量信息。曲線擬合有多種方式，下面

8、是一元函數(shù)采用最小二乘法對給定數(shù)據(jù)進行多項式曲線擬合，最后給出擬合的多項式系數(shù)。1.線性擬合函數(shù)：regress()調(diào)用格式：  b=regress(y,X)                     b,bint,r,rint,stats= regress(y,X)          &#

9、160;          b,bint,r,rint,stats= regress(y,X,alpha)說明：b=regress(y,X)返回X與y的最小二乘擬合值，及線性模型的參數(shù)值、。該函數(shù)求解線性模型：y=X+是p´1的參數(shù)向量；是服從標(biāo)準(zhǔn)正態(tài)分布的隨機干擾的n´1的向量；y為n´1的向量；X為n´p矩陣。bint返回的95%的置信區(qū)間。r中為形狀殘差，rint中返回每一個殘差的95%置信區(qū)間。Stats向量包含R2統(tǒng)計量、回歸的F值和p值。例1：設(shè)y的值為

10、給定的x的線性函數(shù)加服從標(biāo)準(zhǔn)正態(tài)分布的隨機干擾值得到。即y=10+x+ ；求線性擬合方程系數(shù)。程序： x=ones(10,1) (1:10)'      y=x*10;1+normrnd(0,0.1,10,1);      b,bint=regress(y,x,0.05)結(jié)果：  x =     1     1     1 

11、60;   2     1     3     1     4     1     5     1     6     1     7    

12、; 1     8     1     9     1    10y =   10.9567   11.8334   13.0125   14.0288   14.8854   16.1191   17.1189   17.9962 

13、60; 19.0327   20.0175b =              9.9213              1.0143bint =            9.7889   10.0537  

14、          0.9930    1.0357即回歸方程為：y=9.9213+1.0143x2.多項式曲線擬合函數(shù)：polyfit( )調(diào)用格式：  p=polyfit(x,y,n)                     p,s= polyfit(x,y,n)說明：x,y

15、為數(shù)據(jù)點，n為多項式階數(shù)，返回p為冪次從高到低的多項式系數(shù)向量p。矩陣s用于生成預(yù)測值的誤差估計。(見下一函數(shù)polyval)例2：由離散數(shù)據(jù)x0.1.2.3.4.5.6.7.8.91y.3.511.41.61.9.6.4.81.52擬合出多項式。程序：              x=0:.1:1;            y=.3 .5 1 1.4

16、 1.6 1.9 .6 .4 .8 1.5 2;            n=3;            p=polyfit(x,y,n)            xi=linspace(0,1,100);%linspace用于創(chuàng)建向量，如：x

17、=linspace(a1,a2,a3);a1為第一個元素，a2為最末一個元素，a3表示x共有a3個元素，每個元素間距相等。 z=polyval(p,xi); %多項式求值            plot(x,y,'o',xi,z,'k:',x,y,'b')            legend('

18、;原始數(shù)據(jù)','3階曲線')結(jié)果：p =   16.7832  -25.7459   10.9802   -0.0035多項式為：16.7832x3-25.7459x2曲線擬合圖形：如果是n=6，則如下圖： 也可由函數(shù)給出數(shù)據(jù)。例3：x=1:20,y=x+3*sin(x)程序：       x=1:20;       y=x+3*sin(x);  

19、60;    p=polyfit(x,y,6)       xi=linspace(1,20,100);       z=polyval(p,xi);     %多項式求值函數(shù)       plot(x,y,'o',xi,z,'k:',x,y,'b')     

20、  legend('原始數(shù)據(jù)','6階曲線')結(jié)果：p =0.0000   -0.0021    0.0505   -0.5971    3.6472   -9.7295   11.3304 再用10階多項式擬合      程序：x=1:20;y=x+3*sin(x);p=polyfit(x,y,10)xi=linspace(1,20,100);z=pol

21、yval(p,xi);plot(x,y,'o',xi,z,'k:',x,y,'b')legend('原始數(shù)據(jù)','10階多項式')結(jié)果：p =  Columns 1 through 7    0.0000   -0.0000    0.0004   -0.0114    0.1814   -1.8065   11.2360  C

22、olumns 8 through 11  -42.0861   88.5907  -92.8155   40.2671 可用不同階的多項式來擬合數(shù)據(jù)，但也不是階數(shù)越高擬合的越好。3.         多項式曲線求值函數(shù)：polyval( )調(diào)用格式：  y=polyval(p,x)             &

23、#160;       y,DELTA=polyval(p,x,s)說明：y=polyval(p,x)為返回對應(yīng)自變量x在給定系數(shù)P的多項式的值。y,DELTA=polyval(p,x,s) 使用polyfit函數(shù)的選項輸出s得出誤差估計Y DELTA。它假設(shè)polyfit函數(shù)數(shù)據(jù)輸入的誤差是獨立正態(tài)的，并且方差為常數(shù)。則Y DELTA將至少包含50%的預(yù)測值。(未完)4.         多項式曲線擬合的評價和置信區(qū)間函數(shù)：polyconf( )調(diào)用格式

24、：  Y,DELTA=polyconf(p,x,s)                     Y,DELTA=polyconf(p,x,s,alpha)說明：Y,DELTA=polyconf(p,x,s)使用polyfit函數(shù)的選項輸出s給出Y的95%置信區(qū)間Y DELTA。它假設(shè)polyfit函數(shù)數(shù)據(jù)輸入的誤差是獨立正態(tài)的，并且方差為常數(shù)。1-alpha為置信度。例4：給出上面例1的預(yù)

25、測值及置信度為90%的置信區(qū)間。程序：   x=0:.1:1;        y=.3 .5 1 1.4 1.6 1.9 .6 .4 .8 1.5 2        n=3;        p,s=polyfit(x,y,n)        alpha=0.0

26、5;       Y,DELTA=polyconf(p,x,s,alpha)       結(jié)果：   p =   16.7832  -25.7459   10.9802   -0.0035s =   R: 4x4 double  df: 7normr: 1.1406Y =  Columns 1 through 9

27、60;  -0.0035    0.8538    1.2970    1.4266    1.3434    1.1480    0.9413    0.8238    0.8963  Columns 10 through 11    1.2594    2.01405.

28、60;        穩(wěn)健回歸函數(shù)：robust( )穩(wěn)健回歸是指此回歸方法相對于其他回歸方法而言，受異常值的影響較小。調(diào)用格式：  b=robustfit(x,y)                     b,stats=robustfit(x,y)     

29、0;               b,stats=robustfit(x,y,wfun,tune,const)說明：b返回系數(shù)估計向量；stats返回各種參數(shù)估計；wfun指定一個加權(quán)函數(shù)；tune為調(diào)協(xié)常數(shù)；const的值為on(默認(rèn)值)時添加一個常數(shù)項；為off 時忽略常數(shù)項。例5：演示一個異常數(shù)據(jù)點如何影響最小二乘擬合值與穩(wěn)健擬合。首先利用函數(shù)y=10-2x加上一些隨機干擾的項生成數(shù)據(jù)集，然后改變一個y的值形成異常值。調(diào)用不同的擬合函數(shù)，通過圖形觀查

30、影響程度。程序：x=(1:10);y=10-2*x+randn(10,1);y(10)=0;bls=regress(y,ones(10,1) x) %線性擬合brob=robustfit(x,y) %穩(wěn)健擬合scatter(x,y)hold onplot(x,bls(1)+bls(2)*x,:)plot(x,brob(1)+brob(2)*x,r)結(jié)果 ： bls =                    8.445

31、2                   -1.4784brob =                   10.2934            

32、0;      -2.0006 分析：穩(wěn)健擬合(實線)對數(shù)據(jù)的擬合程度好些，忽略了異常值。最小二乘擬合(點線)則受到異常值的影響，向異常值偏移。6.         向自定義函數(shù)擬合對于給定的數(shù)據(jù)，根據(jù)經(jīng)驗擬合為帶有待定常數(shù)的自定義函數(shù)。所用函數(shù)：nlinfit( )調(diào)用格式：  beta,r,J=nlinfit(X,y,fun,beta0)說明：beta返回函數(shù)fun中的待定常數(shù)；r表示殘差；J表示雅可比矩陣。X,y為數(shù)據(jù)；fun自定義函數(shù)；b

33、eta0待定常數(shù)初值。例6：在化工生產(chǎn)中獲得的氯氣的級分y隨生產(chǎn)時間x下降，假定在x8時，y與x之間有如下形式的非線性模型：      現(xiàn)收集了44組數(shù)據(jù)，利用該數(shù)據(jù)通過擬合確定非線性模型中的待定常數(shù)。x            y                

34、0;  x            y                   x            y8       &

35、#160;    0.49               16           0.43               28     

36、60;     0.418            0.49               18           0.46       &#

37、160;       28           0.4010           0.48               18      &#

38、160;    0.45               30           0.4010           0.47         

39、      20           0.42               30           0.4010        

40、   0.48               20           0.42               30        &

41、#160;  0.3810           0.47               20           0.43           

42、;    32           0.4112           0.46               20          

43、; 0.41               32           0.4012           0.46            

44、60;  22           0.41               34           0.4012           0.45

45、0;              22           0.40               36           0.41

46、12           0.43               24           0.42             

47、60; 36           0.3614           0.45               24           0.40 

48、0;             38           0.4014           0.43               24

49、0;          0.40               38           0.4014           0.43   &#

50、160;           26           0.41               40           0.3616  &#

51、160;        0.44               26           0.40               42 

52、0;         0.3916           0.43               26           0.41    &#

53、160;  首先定義非線性函數(shù)的m文件：fff6.mfunction yy=model(beta0,x)  a=beta0(1);  b=beta0(2);  yy=a+(0.49-a)*exp(-b*(x-8);       程序：x=8.00 8.00 10.00 10.00 10.00 10.00 12.00 12.00 12.00 14.00 14.00 14.00.      16.00 16.00 16.00 18.00 18.00

54、20.00 20.00 20.00 20.00 22.00 22.00 24.00.       24.00 24.00 26.00 26.00 26.00 28.00 28.00 30.00 30.00 30.00 32.00 32.00.     34.00 36.00 36.00 38.00 38.00 40.00 42.00'   y=0.49 0.49 0.48 0.47 0.48 0.47 0.46 0.46 0.45 0.43 0.45 0.43 0.4

55、3 0.44 0.43.     0.43 0.46 0.42 0.42 0.43 0.41 0.41 0.40 0.42 0.40 0.40 0.41 0.40 0.41 0.41.     0.40 0.40 0.40 0.38 0.41 0.40 0.40 0.41 0.38 0.40 0.40 0.39 0.39'     beta0=0.30 0.02;betafit = nlinfit(x,y,model ,beta0)結(jié)果：betafit = 

56、               0.38960.1011       即：a=0.3896 ，b=0.1011 擬合函數(shù)為：yy=0.3896+0.1004*exp(-0.1011*(x-8)excell數(shù)據(jù)調(diào)入matlab中的方法：import datadatax=data(:,所調(diào)出數(shù)據(jù)的列數(shù))如果想調(diào)出某列數(shù)據(jù)中間隔相等的數(shù)據(jù)則：x=data(1:2:10,colume2)則調(diào)用原始數(shù)據(jù)的第二列中前十個

57、數(shù)據(jù)中1、3、5、7、9行數(shù)據(jù)。§2 插值問題在應(yīng)用領(lǐng)域中，由有限個已知數(shù)據(jù)點，構(gòu)造一個解析表達式，由此計算數(shù)據(jù)點之間的函數(shù)值，稱之為插值。實例：海底探測問題某公司用聲納對海底進行測試，在5×5海里的坐標(biāo)點上測得海底深度的值，希望通過這些有限的數(shù)據(jù)了解更多處的海底情況。并繪出較細致的海底曲面圖。一、一元插值一元插值是對一元數(shù)據(jù)點(xi,yi)進行插值。1  線性插值：由已知數(shù)據(jù)點連成一條折線，認(rèn)為相臨兩個數(shù)據(jù)點之間的函數(shù)值就在這兩點之間的連線上。一般來說，數(shù)據(jù)點數(shù)越多，線性插值就越精確。調(diào)用格式：yi=interp1(x,y,xi,linear)  %線

58、性插值zi=interp1(x,y,xi,spline)  %三次樣條插值wi=interp1(x,y,xi,cubic)  %三次多項式插值說明：yi、zi、wi為對應(yīng)xi的不同類型的插值。x、y為已知數(shù)據(jù)點。例1：已知數(shù)據(jù)：x0.1.2.3.4.5.6.7.8.91y.3.511.41.61.9.6.4.81.52求當(dāng)xi=0.25時的yi的值。程序：x=0:.1:1;y=.3 .5 1 1.4 1.6 1 .6 .4 .8 1.5 2;yi0=interp1(x,y,0.025,'linear')xi=0:.02:1;yi=interp1(x,y,xi

59、,'linear');zi=interp1(x,y,xi,'spline');wi=interp1(x,y,xi,'cubic');plot(x,y,'o',xi,yi,'r+',xi,zi,'g*',xi,wi,'k.-')legend('原始點','線性點','三次樣條','三次多項式')結(jié)果：yi0 =  0.3500 要得到給定的幾個點的對應(yīng)函數(shù)值，可用：xi = 0.2500  0

60、.3500  0.4500yi=interp1(x,y,xi,'spline')結(jié)果：yi =1.2088  1.5802  1.3454(二) 二元插值二元插值與一元插值的基本思想一致，對原始數(shù)據(jù)點(x,y,z)構(gòu)造見世面函數(shù)求出插值點數(shù)據(jù)(xi,yi,zi)。一、單調(diào)節(jié)點插值函數(shù)，即x,y向量是單調(diào)的。調(diào)用格式1：zi=interp2(x,y,z,xi,yi,linear)liner 是雙線性插值 (缺省)調(diào)用格式2：zi=interp2(x,y,z,xi,yi,nearest)nearest 是最近鄰域插值調(diào)用格式3：zi=interp2(x

61、,y,z,xi,yi,spline)spline是三次樣條插值說明：這里x和y是兩個獨立的向量，它們必須是單調(diào)的。z是矩陣，是由x和y確定的點上的值。z和x,y之間的關(guān)系是z(i,:)=f(x,y(i) z(:,j)=f(x(j),y) 即：當(dāng)x變化時，z的第i行與y的第i個元素相關(guān)，當(dāng)y變化時z的第j列與x的第j個元素相關(guān)。如果沒有對x,y賦值，則默認(rèn)x=1:n, y=1:m。n和m分別是矩陣z的行數(shù)和列數(shù)。例2：已知某處山區(qū)地形選點測量坐標(biāo)數(shù)據(jù)為：x=0  0.5  1  1.5  2  2.5  3  3.5 

62、; 4  4.5  5y=0  0.5  1  1.5  2  2.5  3  3.5  4  4.5  5  5.5  6海拔高度數(shù)據(jù)為：z=89 90 87 85 92 91 96 93 90 87 82   92 96 98 99 95 91 89 86 84 82 84   96 98 95 92 90 88 85 84 83 81 85   80 81 82 89 95 96 93 92 89 86 86   82 85 87 98 99 96 97 88 85 82 83   82 85 89 94 95 93 92 91 86 84 8