數(shù)學建模插值與擬合及MATLAB實現(xiàn)學習教案

1、會計學1數(shù)學數(shù)學(shxu)建模插值與擬合及建模插值與擬合及MATLAB實實現(xiàn)現(xiàn)第一頁，共50頁。 我們經(jīng)常會遇到大量的數(shù)據(jù)需要(xyo)處理，而處理數(shù)據(jù)的關鍵就在于這些算法，例如數(shù)據(jù)擬合、參數(shù)估計、插值等數(shù)據(jù)處理算法。此類問題在MATLAB中有很多現(xiàn)成的函數(shù)可以調用，熟悉MATLAB，這些方法都能游刃有余的用好。 一、概述(i sh)第1頁/共50頁第二頁，共50頁。 數(shù)據(jù)擬合(n h)在很多賽題中有應用，與圖形處理有關的問題很多與插值和擬合(n h)有關系，例如98年美國賽A題，生物組織切片的三維插值處理，94年A題逢山開路，山體海拔高度的插值計算，2003年吵的沸沸揚揚的“非典”問題也要

2、用到數(shù)據(jù)擬合(n h)算法，觀察數(shù)據(jù)的走向進行處理， 2005年的雨量預報的評價的插值計算。2001年的公交車調度擬合(n h)問題，2003年的飲酒駕車擬合(n h)問題。第2頁/共50頁第三頁，共50頁。預測點和實測(sh c)點的圖形插值后的圖形(txng)第3頁/共50頁第四頁，共50頁。喝兩瓶酒的擬合(n h)曲線喝1-5瓶酒的擬合(n h)曲線2004年全國大學生數(shù)學建模競賽題（酒后駕車）中給出某人在短時間內喝下啤酒后，間隔一定的時間測量他的血液中酒精含量（毫克百毫升）第4頁/共50頁第五頁，共50頁。 在實際中，常常要處理由實驗或測量所得到(d do)的一些離散數(shù)據(jù)。插值與擬合方

3、法就是要通過這些數(shù)據(jù)去確定某一類已知函數(shù)的參數(shù)或尋求某個近似函數(shù)，使所得到(d do)的近似函數(shù)與已知數(shù)據(jù)有較高的擬合精度。 如果要求這個近似函數(shù)（曲線或曲面）經(jīng)過所已知的所有數(shù)據(jù)點，則稱此類問題為插值問題。 （不需要函數(shù)表達式）二、第5頁/共50頁第六頁，共50頁。 如果不要求近似函數(shù)通過所有數(shù)據(jù)點，而是要求它能較好地反映數(shù)據(jù)變化規(guī)律的近似函數(shù)的方法稱為數(shù)據(jù)擬合(n h)。（必須有函數(shù)表達式） 近似函數(shù)不一定（曲線或曲面）通過所有的數(shù)據(jù)點。 第6頁/共50頁第七頁，共50頁。1、聯(lián)系都是根據(jù)實際中一組已知數(shù)據(jù)來構造一個能夠反映數(shù)據(jù)變化規(guī)律的近似函數(shù)的方法。2、區(qū)別插值問題(wnt)不一定得到

4、近似函數(shù)的表達形式，僅通過插值方法找到未知點對應的值。數(shù)據(jù)擬合要求得到一個具體的近似函數(shù)的表達式。 三、插值與擬合的區(qū)別(qbi)和聯(lián)系第7頁/共50頁第八頁，共50頁。四、插值的使用(shyng)及求解 當數(shù)據(jù)量不夠，需要補充，且認定(rndng)已有數(shù)據(jù)可信時, 通常利用函數(shù)插值方法。 實際問題當中碰到的函數(shù) f (x) 是各種各樣的，有的表達式很復雜，有的甚至給不出數(shù)學的式子，只提供了一些離散數(shù)據(jù)，警如，某些點上的函數(shù)值和導數(shù)值。 4.1 引言(ynyn)第8頁/共50頁第九頁，共50頁。 選用(xunyng)不同類型的插值函數(shù)，逼近的效果就不同，一般有：（1）拉格朗日插值（lagran

5、ge插值）（2）分段線性插值（3）Hermite（4）三次樣條插值。 4.2 插值方法(fngf)第9頁/共50頁第十頁，共50頁。 Matlab 實現(xiàn)：實現(xiàn)分段線性插值不需要編制函數(shù)(hnsh)程序，它自身提供了內部的功能函數(shù)(hnsh)interp1(一維插值) interp2(二維) interp3(三維) intern(n維) 4.3 MATLAB實現(xiàn)(shxin)插值第10頁/共50頁第十一頁，共50頁。用MATLAB作插值計算(j sun)一維插值函數(shù)(hnsh)：yi=interp1(x，y，xi，method)插值方法被插值點插值節(jié)點xi處的插值結果nearest 最鄰近(ln

6、 jn)插值；linear 線性插值；spline 三次樣條插值；cubic 立方插值； 缺省時 分段線性插值 注意：所有的插值方法都要求x是單調的，并且xi不能夠超過x的范圍第11頁/共50頁第十二頁，共50頁。例：從1點12點的11小時內，每隔1小時測量一次溫度，測得的溫度的數(shù)值依次(yc)為：5，8，9，15，25，29，31，30，22，25，27，24試估計每隔1/10小時的溫度值hours=1:12;temps=5 8 9 15 25 29 31 30 22 25 27 24;h=1:0.1:12;t=interp1(hours,temps,h,spline); plot(hour

7、s,temps,k+,h,t,b,hours,temps,r:) %作圖xlabel(Hour),ylabel(Degrees Celsius)第12頁/共50頁第十三頁，共50頁。xy機翼下輪廓線例 已知飛機(fij)下輪廓線上數(shù)據(jù)如下，求x每改變時的y值第13頁/共50頁第十四頁，共50頁。返回(fnhu)第14頁/共50頁第十五頁，共50頁。 要求x0,y0單調；x，y可取為矩陣，或x取行向量，y取為列向量，x,y的值分別(fnbi)不能超出x0,y0的范圍z=interp2(x0,y0,z0,x,y,method)被插值點插值方法用MATLAB作網(wǎng)格(wn )節(jié)點數(shù)據(jù)的插值插值節(jié)點被插

8、值點的函數(shù)值nearest 最鄰近(ln jn)插值；linear 雙線性插值；cubic 雙三次插值；缺省時 雙線性插值.第15頁/共50頁第十六頁，共50頁。例：測得平板表面35網(wǎng)格點處的溫度分別為： 82 81 80 82 84 79 63 61 65 81 84 84 82 85 86 試作出平板表面的溫度分布(fnb)曲面z=f(x,y)的圖形輸入以下輸入以下(yxi)命令：命令：x=1:5;y=1:3;temps=82 81 80 82 84;79 63 61 65 81;84 84 82 85 86;mesh(x,y,temps)1.先在三維坐標畫出原始數(shù)據(jù)，畫出粗糙先在三維坐標

9、畫出原始數(shù)據(jù)，畫出粗糙(cco)的溫度的溫度分布曲線圖分布曲線圖.2以平滑數(shù)據(jù),在 x、y方向上每隔個單位的地方進行插值.第16頁/共50頁第十七頁，共50頁。再輸入以下命令再輸入以下命令:xi=1:0.2:5;yi=1:0.2:3;zi=interp2(x,y,temps,xi,yi,cubic);mesh(xi,yi,zi)畫出插值后的溫度畫出插值后的溫度(wnd)分布曲面圖分布曲面圖. 第17頁/共50頁第十八頁，共50頁。 通過此例對最近鄰點插值、雙線性插值方法和雙三次插值方法的插值效果進行(jnxng)比較第18頁/共50頁第十九頁，共50頁。輸入(shr)原始數(shù)據(jù):第19頁/共50

10、頁第二十頁，共50頁。figure(1);meshz(x,y,z)xlabel(X),ylabel(Y),zlabel(Z) xi=1200:50:4000;yi=1200:50:3600; figure(2)z1i=interp2(x,y,z,xi,yi,nearest);surfc(xi,yi,z1i)xlabel(X),ylabel(Y),zlabel(Z) figure(3)z2i=interp2(x,y,z,xi,yi);surfc(xi,yi,z2i)xlabel(X),ylabel(Y),zlabel(Z) figure(4)z3i=interp2(x,y,z,xi,yi,cub

11、ic);surfc(xi,yi,z3i)xlabel(X),ylabel(Y),zlabel(Z) figure(5)subplot(1,3,1),contour(xi,yi,z1i,10,r);subplot(1,3,2),contour(xi,yi,z2i,10,r);subplot(1,3,3),contour(xi,yi,z3i,10,r);注:surfc 函數(shù)功能 在矩形區(qū)域(qy)內顯示三維帶陰影曲面圖,且在曲面下面畫出等高線。meshc函數(shù) 生成帶等高線網(wǎng)線圖 meshz函數(shù) 生成帶垂簾的網(wǎng)線圖數(shù)據(jù)(shj)插值與圖形繪制:第20頁/共50頁第二十一頁，共50頁。 插值函數(shù)(hn

12、sh)griddata格式為: cz =griddata（x，y，z，cx，cy，method）用MATLAB作散點數(shù)據(jù)(shj)的插值計算 要求(yoqi)cx取行向量，cy取為列向量被插值點插值方法插值節(jié)點被插值點的函數(shù)值nearest最鄰近插值linear 雙線性插值cubic 雙三次插值v4- MATLAB提供的插值方法缺省時, 雙線性插值第21頁/共50頁第二十二頁，共50頁。 例 在某海域測得一些點(x,y)處的水深z由下表給出，船的吃水深度為5英尺(yngch)，在矩形區(qū)域（75，200）（-50，150）里的哪些地方船要避免進入第22頁/共50頁第二十三頁，共50頁。4.作出水

13、深小于5的海域(hiy)范圍,即z=5的等高線.2.在矩形區(qū)域(qy)(75,200)(-50,150)進行插值。 1.輸入(shr)插值基點數(shù)據(jù) 3. 作海底曲面圖 第23頁/共50頁第二十四頁，共50頁。%程序(chngx)一：插值并作海底曲面圖 x =129.0 140.0 103.5 88.0 185.5 195.0 105.5 157.5 107.5 77.0 81.0 162.0 117.5 ;y = 7.5 141.5 23.0 147.0 22.5 137.5 85.5 -6.5 -81 3.0 56.5 -66.5 84.0 -33.5 ;z = 4 8 6 8 6 8 8

14、9 9 8 8 9 4 9 ;x1=75:1:200;y1=-50:1:150;x1,y1=meshgrid(x1,y1);z1=griddata(x,y,z,x1,y1,v4);meshc(x1,y1,z1) 第24頁/共50頁第二十五頁，共50頁。海底(hi d)曲面圖第25頁/共50頁第二十六頁，共50頁。%程序(chngx)二：插值并作出水深小于5的海域范圍。x1=75:1:200;y1=-50:1:150;x1,y1=meshgrid(x1,y1);z1=griddata(x,y,z,x1,y1,v4); %插值z1(z1=5)=nan; %將水深大于5的置為nan，這樣繪圖就不會顯

15、示出來meshc(x1,y1,z1) 第26頁/共50頁第二十七頁，共50頁。水深小于5的海域(hiy)范圍第27頁/共50頁第二十八頁，共50頁。5.1 引言 對于情況較復雜的實際問題（因素不易化簡，作用機理不詳）可直接使用(shyng)數(shù)據(jù)組建模，尋找簡單的因果變量之間的數(shù)量關系， 從而對未知的情形作預報。這樣組建的模型為擬合模型。 擬合模型的組建主要是處理好觀測數(shù)據(jù)的誤差，使用(shyng)數(shù)學表達式從數(shù)量上近似因果變量之間的關系。擬合模型的組建是通過對有關變量的觀測數(shù)據(jù)的觀察、分析和選擇恰當?shù)臄?shù)學表達方式得到的。 五、擬合(n h)的使用及求解第28頁/共50頁第二十九頁，共50頁。5

16、.2 擬合(n h)模型的分類 5.2.1 直線擬合5.2.2 曲線擬合 觀察數(shù)據(jù)修勻 對于已給一批實測數(shù)據(jù)，由于實測方法、實驗環(huán)境等一些外界因素的影響，不可避免地會產(chǎn)生隨機干擾和誤差。我們自然希望(xwng)根據(jù)數(shù)據(jù)分布的總趨勢去剔除觀察數(shù)據(jù)中的偶然誤差，這就是所謂的數(shù)據(jù)修勻（或稱數(shù)據(jù)平滑）問題。 第29頁/共50頁第三十頁，共50頁。直 線 擬 合 問 題 引 例 1溫度t(C) 20.5 32.7 51.0 73.0 95.7電阻R() 765 826 873 942 1032已知熱敏電阻數(shù)據(jù)：求60C時的電阻(dinz)R 設 R=at+ba,b為待定系數(shù)(xsh)第30頁/共50頁第

17、三十一頁，共50頁。曲 線 擬 合 問 題 引 例 2 t (h) 0.25 0.5 1 1.5 2 3 4 6 8c (g/ml) 19.21 18.15 15.36 14.10 12.89 9.32 7.45 5.24 3.01已知一室模型快速靜脈注射下的血藥濃度數(shù)據(jù)(t=0注射300mg)求血藥濃度隨時間(shjin)的變化規(guī)律c(t).在直角坐標(zh jio zu bio)系下作圖如下(plot)0( )e,ktc tcc k為待定系數(shù)第31頁/共50頁第三十二頁，共50頁。曲 線 擬 合 問 題 的 提 法已知一組（二維）數(shù)據(jù)，即平面上 n個點（xi,yi) i=1,n, 尋求(x

18、nqi)一個函數(shù)（曲線）y=f(x), 使 f(x) 在某種準則下與所有數(shù)據(jù)點最為接近，即曲線擬合得最好 +xyy=f(x)(xi,yi)ii 為點（xi,yi) 與曲線(qxin) y=f(x) 的距離第32頁/共50頁第三十三頁，共50頁。擬合(n h)與插值的關系 函數(shù)插值與曲線擬合都是要根據(jù)一組數(shù)據(jù)構造一個函數(shù)作為近似，由于(yuy)近似的要求不同，二者在數(shù)學方法上是完全不同的 實例：下面數(shù)據(jù)是某次實驗所得，希望得到(d do)X和 f之間的關系？x124791 21 31 51 7f1 .53 .96 .611 .71 5 .61 8 .81 9 .62 0 .62 1 .1MATL

19、AB(cn)問題：給定一批數(shù)據(jù)點，需確定滿足特定要求的曲線或曲面解決方案：若不要求曲線（面）通過所有數(shù)據(jù)點，而是要求它反映對象整體的變化趨勢，這就是數(shù)據(jù)擬合，又稱曲線擬合或曲面擬合若要求所求曲線（面）通過所給所有數(shù)據(jù)點，就是插值問題；第33頁/共50頁第三十四頁，共50頁。最臨近(ln jn)插值、線性插值、樣條插值與曲線擬合結果：第34頁/共50頁第三十五頁，共50頁。曲線擬合問題(wnt)最常用的解法線性最小二乘法的基本思路第一步:先選定(xun dn)一組函數(shù) r1(x), r2(x), ,rm(x), mn, 令 f(x)=a1r1(x)+a2r2(x)+ +amrm(x) （1）其中

20、 a1,a2, ,am 為待定系數(shù) 第二步: 確定a1,a2, ,am 的準則(zhnz)（最小二乘準則(zhnz)）：使n個點（xi,yi) 與曲線 y=f(x) 的距離i 的平方和最小 記 )2()()(),(211211221iiknimkkininiiimyxrayxfaaaJ 問題歸結為，求 a1,a2, ,am 使 J (a1,a2, ,am) 最小第35頁/共50頁第三十六頁，共50頁。用MATLAB作線性最小二乘擬合(n h)1. 作多項式f(x)=a1xm+ +amx+am+1擬合(n h),可利用已有程序:a=polyfit(x,y,m)輸入同長度的數(shù)組x，y擬合多項式次數(shù)

21、2.多項式在x處的值y可用以下命令(mng lng)計算： y=polyval（a，x）第36頁/共50頁第三十七頁，共50頁。1）輸入以下命令）輸入以下命令(mng lng)： x=0:0.1:1; y=-0.447 1.978 3.28 6.16 7.08 7.34 . 7.66 9.56 9.48 9.30 11.2; A=polyfit(x,y,2) z=polyval(A,x); plot(x,y,k+,x,z,r) %作出數(shù)據(jù)點和擬合曲線的圖形作出數(shù)據(jù)點和擬合曲線的圖形2）計算結果： 用多項式擬合(n h)的命令0317.01293.208108.9)(2xxxf第37頁/共50頁

22、第三十八頁，共50頁。如何預報人口的增長 人口的增長是當前世界上引起普遍關注的問題，并且我們會發(fā)現(xiàn)在不同的刊物預報同一時間的人口數(shù)字不相同，這顯然是由于用了不同的人口模型計算的結果。 我國是世界第一人口大國，基本上地球每九個人中就有一個中國人。有效地控制我國人口的增長是使我過全面進入小康社會、到21世紀中葉建成富強民主文明的社會主義國家的需要。而有效控制人口增長的前提是要認識人口數(shù)量的變化規(guī)律，建立(jinl)人口模型，作出較準確的預報。 例:如何(rh)預報人口的增長第38頁/共50頁第三十九頁，共50頁。例如：1949年1994年我國人口數(shù)據(jù)資料如下： 年 份xi 1949 1954 19

23、59 1964 1969 1974 1979 1984 1989 1994 人口數(shù)yi 5.4 6.0 6.7 7.0 8.1 9.1 9.8 10.3 11.3 11.8 建模分析我國人口增長(zngzhng)的規(guī)律,預報1999年我國人口數(shù)。 第39頁/共50頁第四十頁，共50頁。模型(mxng)一：假設人口隨時間線性地增加 模型：參數(shù)估計觀測值的模型： 擬合的精度(jn d): 誤差平方和。 可以算出：a = -283.2320 b模型：y = + x bxaynibxayii,2, 1,)(22iiiibxayeQ第40頁/共50頁第四十一頁，共50頁。則可看成是線性方程,用 poly

24、fit命令(mng lng)計算得：模型(mxng)二：指數(shù)增長模型(mxng) bxaey bxay lnln可變?yōu)閤ey0179. 033. 2YA=+BXa=2.33, b則所求模型(mxng)為:第41頁/共50頁第四十二頁，共50頁。程序如下(rxi)：x=1949 1954 1959 1964 1969 1974 1979 1984 1989 1994; y=5.4 6.0 6.7 7.0 8.1 9.1 9.8 10.3 11.3 11.8 ; a=polyfit(x,y,1); x1=1949:10:1994; y1=a(2)+a(1)*x1; b=polyfit(x,log(

25、y),1); y2=exp(b(2)*exp(b(1)*x1); plot(x,y,*) hold on plot(x1,y1,-r) hold on plot(x1,y2,-k) legend(原曲線,模型一曲線,模型二曲線) 第42頁/共50頁第四十三頁，共50頁。結論的比較結論的比較(bjio)(bjio)如下表：如下表： 年年 份份 xi 1949 1954 1959 1964 1969 1974 1979 1984 1989 1994 xi 1949 1954 1959 1964 1969 1974 1979 1984 1989 1994 人口數(shù)人口數(shù) yi 5.4 6.0 6.7

26、7.0 8.1 9.1 9.8 10.3 11.3 11.8 yi 5.4 6.0 6.7 7.0 8.1 9.1 9.8 10.3 11.3 11.8 模型一值模型一值 5.24 5.97 6.70 7.43 8.16 8.90 9.62 10.36 11.09 11.82 5.24 5.97 6.70 7.43 8.16 8.90 9.62 10.36 11.09 11.82 誤誤 差差 0.16 0.03 0.00 -0.43 -0.06 0.20 0.18 -0.06 0.01 -0.02 0.16 0.03 0.00 -0.43 -0.06 0.20 0.18 -0.06 0.01

27、-0.02 模型二值模型二值 5.55 6.06 6.62 7.23 7.90 8.64 9.44 10.31 11.26 12.31 5.55 6.06 6.62 7.23 7.90 8.64 9.44 10.31 11.26 12.31 誤差誤差 第43頁/共50頁第四十四頁，共50頁。結果分析： （1） Q1 = 0.2915 0.7437 = Q2 . 線性模型更適合中國人口(rnku)的增長。 （2） 預報：1999 年 12.55 億，13.43 億 （3） 統(tǒng)計年鑒(ninjin)： 2005 年 13.3 億， 2010 年 14 億 模型 I 2005 年13.43 億， 2010 年14.16 億 模型 II 2005 年14.94 億， 2010 年 16.33 億第44頁/共50頁第四十五頁，共50頁。1. lsqcurvefit1. lsqcurvefit已知數(shù)據(jù)已知數(shù)據(jù)(shj)(shj)點：點： xdata= xdata=（xdata1,xdata2,xdatanxdata1,xdata2,xdatan）, , ydata= ydata=（ydata1,ydata2,ydatanydata1,ydata2,ydatan） 用MATLAB作非線性最小二乘擬合(n h) MATLAB提供了求非線性最小二乘擬合的函數(shù)：ls