數(shù)學建模插值法與曲線擬合講課

關于數(shù)學建模插值法與曲線擬合講課第1頁，講稿共65頁，2023年5月2日，星期三

一、問題的提出

在生產(chǎn)和實驗中，關于函數(shù)f(x)，經(jīng)常存在兩種情況：（1）其表達式不便于計算；（2）無表達式.

而只有函數(shù)在給定點的函數(shù)值，怎樣預測其它點的函數(shù)值？xx0x1x2…xnyy0y1y2…yn第2頁，講稿共65頁，2023年5月2日，星期三飛機機翼制造

下表給出的x、y數(shù)據(jù)位于機翼端面的輪廓線上，Y1和Y2分別對應輪廓的上下線。假設需要得到x坐標每改變0.1時的y坐標，試完成加工所需數(shù)據(jù)，畫出曲線.x035791112131415Y101.82.22.73.03.12.92.52.01.6Y201.21.72.02.02.01.81.21.01.6第3頁，講稿共65頁，2023年5月2日，星期三山體地貌要在某山區(qū)方圓大約27平方公里范圍內(nèi)修建一條公路，從山腳出發(fā)經(jīng)過一個居民區(qū)，再到達一個礦區(qū)。橫向縱向分別每隔400米測量一次，得到一些地點的高程：試做出該山區(qū)的地貌圖.第4頁，講稿共65頁，2023年5月2日，星期三船在該海域會擱淺嗎？---作業(yè)

在某海域測得一些點(x,y)處的水深z由下表給出，船的吃水深度為5英尺，在矩形區(qū)域（75，200）*（-50，150）里的哪些地方船要避免進入.第5頁，講稿共65頁，2023年5月2日，星期三

水深和流速的問題

在水文數(shù)據(jù)測量中，不同水深的流速是不同的.水文數(shù)據(jù)的測量時天天進行的，為了減少測量的工作，希望得到確定的水深和水流之間的關系.為此測量了一系列不同水深和流速值.下表給出了對某河流的測量數(shù)據(jù)，其中水深和流速根據(jù)適當?shù)膯挝贿M行了規(guī)范化，共10個值.第6頁，講稿共65頁，2023年5月2日，星期三美國人口問題據(jù)美國人口普查局數(shù)據(jù)：從1790每隔10年至2000年的總人口（單位：百萬）如下示

t=1790:10:2000;p=[3.9,5.3,7.2,9.6,12.9,17.1,23.1,31.4,38.6,50.2,62.9,76,92,105.7,122.8,131.7,150.7,179,205,226.5,251.4,281.422];

預測2001，2002年的美國人口數(shù)？并與調(diào)查數(shù)據(jù)285.318，288.369比較，選擇擬合較好的模型。第7頁，講稿共65頁，2023年5月2日，星期三農(nóng)作物施肥效果分析1992年A題

在農(nóng)業(yè)生產(chǎn)試驗研究中，對某地區(qū)土豆的產(chǎn)量與化肥的關系做了一實驗，得到了氮肥、磷肥的施肥量與土豆產(chǎn)量的對應關系如下表：

1.根據(jù)上表數(shù)據(jù)分別給出土豆產(chǎn)量與氮、磷肥的關系式。

2.施肥問題優(yōu)化策略氮肥量（公斤/公頃）03467101135202259336404471土豆產(chǎn)量（公斤）15.1821.3625.7232.293439.4543.1543.4640.8330.75磷肥量（公斤/公頃）024497398147196245294342土豆產(chǎn)量（公斤）33.4632.4736.0637.964140.141。342.240.442.7第8頁，講稿共65頁，2023年5月2日，星期三

配藥方案---作業(yè)

一種新藥用于臨床之前,必須設計給藥方案.在快速靜脈注射的給藥方式下,所謂給藥方案是指,每次注射劑量多大,間隔時間多長.

藥物進入機體后隨血液輸送到全身,在這個過程中不斷地被吸收,分布,代謝,最終排出體外.藥物在血液中的濃度,即單位體積血液中的藥物含量,稱血藥濃度.在最簡單的一室模型中,將整個機體看作一個房室,稱中心室,室內(nèi)的血藥濃度是均勻的.快速靜脈注射后,濃度立即上升;然后逐漸下降.當濃度太低時,達不到預期的治療效果;血藥濃度太高,又可能導致藥物中毒或副作用太強.臨床上,每種藥物有一個最小有效濃度c1和一個最大治療濃度c2.設計給藥方案時,要使血藥濃度保持在c1-c2之間.設本題所研究藥物的最小有效濃度c1=10,最大治療濃度c2=25().第9頁，講稿共65頁，2023年5月2日，星期三

顯然,要設計給藥方案,必須知道給藥后血藥濃度隨時間變化的規(guī)律.為此,從實驗和理論兩方面著手.在實驗方面,對某人用快速靜脈注射方式一次注入該藥物300mg后,在一定時刻t(小時)采集血樣,測得血藥濃度c.如表:血藥濃度c(t)的測試數(shù)據(jù)

t0.250.511.523468c19.2118.1515.3614.1012.899.327.455.243.01問題：1.在快速靜脈注射的給藥方式下，研究血藥濃度（單位體積血液中的藥物含量）的變化規(guī)律；2.給定藥物的最小有效濃度和最大治療濃度，設計給藥方案：每次注射劑量多大；間隔時間多長？

配藥方案第10頁，講稿共65頁，2023年5月2日，星期三二、問題的解決（1）插值法；（2）曲線擬合法．

1、問題的抽象xx1x2…xmyy1y2…ym構造一個簡單易于計算的近似函數(shù)

p(x)f(x)

（精確函數(shù)）。2、構造近似函數(shù)，p(x)的方法有兩種：在實驗中經(jīng)常給出一組離散點，第11頁，講稿共65頁，2023年5月2日，星期三插值法定義：當精確函數(shù)y=f(x)非常復雜或未知時，在一系列節(jié)點

x0…xn

處測得函數(shù)值y0

=f(x0),…，yn

=f(xn)，由此構造一個簡單易算的近似函數(shù)p(x)

f(x)，滿足條件p(xi)=f(xi)(i=0,…n)，（插值條件）這里的p(x)

稱為f(x)的插值函數(shù)；構造插值函數(shù)的方法為插值法。第12頁，講稿共65頁，2023年5月2日，星期三曲線擬合

但是不要求使p(xi)=yi,而只要p(xi)yi

總體上盡可能小。這種構造近似函數(shù)p(x)

的方法稱為曲線擬合法，p(x)稱為擬合函數(shù)。定義：當精確函數(shù)y=f(x)非常復雜或未知時，在一系列節(jié)點x0…xn

處,測得函數(shù)值y0

,…，yn

，由此構造一個簡單易算的近似函數(shù)p(x)

f(x)，第13頁，講稿共65頁，2023年5月2日，星期三插值與擬合的相同點都需要根據(jù)已知數(shù)據(jù)構造函數(shù)?？墒褂玫玫胶瘮?shù)計算未知點的函數(shù)值。xx1x2…xmyy1y2…ym求一個簡單易算的近似函數(shù)

p(x)f(x)

。第14頁，講稿共65頁，2023年5月2日，星期三插值與擬合的不同點插值:過節(jié)點;；擬合:不過點,整體近似；第15頁，講稿共65頁，2023年5月2日，星期三插值法拉格朗日插值牛頓插值三次埃爾米特插值法分段線性插值分段三次埃爾米特插值法三次樣條插值第16頁，講稿共65頁，2023年5月2日，星期三1、拉格朗日插值公式（１）定義對給定的n+1個節(jié)點x0,x1,x2,…,xn及對應的函數(shù)值y0,

y1,y2,…,yn，構造一個n次插值多項式：即為拉格朗日插值公式，其中插值基函數(shù)第17頁，講稿共65頁，2023年5月2日，星期三拉格朗日插值的matlab實現(xiàn)functiony=lagrange(x0,y0,x)%x0插值節(jié)點，y0插值節(jié)點處的函數(shù)值，x要計算函數(shù)值的點；n=length(x0);%計算x0的長度m=length(x);%計算x的長度fori=1:ms=0;z=x(i);

fork=1:np=1.0;forj=1:n ifj~=kp=p*(z-x0(j))/(x0(k)-x0(j));%計算插值基函數(shù)

endends=p*y0(k)+s;endy(i)=s;%計算在x(i)處的函數(shù)值（拉格朗日）end第18頁，講稿共65頁，2023年5月2日，星期三2、牛頓插值法牛頓插值公式：Nn(x)=f(x0)+f[x0,x1](x-x0)+…+f[x0,x1,x2,…,xn](x-x0)(x-x1)…(x-xn)其中：

f[x0,x1]一階差商

f[x0,x1,x2,…,xn]n階差商注：牛頓插值法與拉格朗日插值法，同一個多項式，不同的表達方式，但是計算量不一樣，牛頓插值法的計算量小。第19頁，講稿共65頁，2023年5月2日，星期三龍格現(xiàn)象Runge在上個世紀初發(fā)現(xiàn)：在[-5,5]上用n+1個等距節(jié)點作n次插值多項式Pn(x)，當在n→∞時，插值多項式Pn(x)在區(qū)間中部趨于f(x)=1/(1+x2),但對于3.63≤∣x∣≤1的x，Pn(x)嚴重發(fā)散。用圖形分析問題。

第20頁，講稿共65頁，2023年5月2日，星期三forn=10:2:20%從10等份到20等份x0=[-5:10/n:5];%插值節(jié)點y0=1./(1+x0.^2);%插值節(jié)點處的精確函數(shù)值x=[-5:0.1:5];%要進行計算函數(shù)值的點y=lagrange(x0,y0,x);%調(diào)用函數(shù)計算x點的函數(shù)值plot(x0,y0,‘*’,x,1./(1+x.^2),‘r’,x,y)%繪制圖形pause%等待，按任意鍵end第21頁，講稿共65頁，2023年5月2日，星期三3、分段低次插值法（1）分段線性插值

定義：已知n+1個不同節(jié)點x0,x1,…,xn,構造分段多項式I(x)，使之滿足l

I(x)在[a，b]上連續(xù)；l

I(xk)=yk；l

I(x)在[xi,xi+1]上是一次多項式；

I(x)=第22頁，講稿共65頁，2023年5月2日，星期三（2）分段三次埃爾米特插值法定義：已知n+1個不同節(jié)點x0,x1,…,xn,構造分段多項式I(x)，使之滿足：l

I(x)在[a，b]上二階連續(xù)導數(shù)；l

I(xk)=yk，I’(xk)=y’k，；l

I(x)在[xi,xi+1]上是三次次多項式。第23頁，講稿共65頁，2023年5月2日，星期三4、三次樣條插值法

對于給定n+1個不同節(jié)點x0,x1,…,xn及函數(shù)值y0,y1,…,yn，其中a=x0<x1<…<xn=b，構造三次樣條插值函數(shù)S(x)。

S(x)稱為三次樣條函數(shù)時需滿足：l

S(x)在[a,b]上二階導數(shù)連續(xù)；l

S(xk)=yk(k=0,1,…,n)；l

每個子區(qū)間[xk,xk+1]上S(x)是三次多項式(k=0,1,…,n)。第24頁，講稿共65頁，2023年5月2日，星期三插值法的matlab實現(xiàn)—一維插值

命令：interp1(x0,y0,x，’method’)

其中：x0：插值節(jié)點；

y0：插值節(jié)點處的函數(shù)值；

x：要計算函數(shù)值的點；

method：

linear

：分段線性插值；

cubic：分段三次埃爾米特插值；

spline：三次樣條插值。第25頁，講稿共65頁，2023年5月2日，星期三插值法的應用一水庫上游河段降暴雨，根據(jù)預報測算上游流入水庫的流量為Q(t)(102立方米/秒)

：

t(時)81216243044485660Q（t）3654789210135251613

通過這個預報值，分別用不同的數(shù)值方法插值法來估計14和20時上游流入水庫的流量。第26頁，講稿共65頁，2023年5月2日，星期三二維插值的MATLAB實現(xiàn)

在MATLAB中，二維插值命令常用的有兩個，

1、一個是網(wǎng)格節(jié)點插值：

z=interp2(x0,y0,z0,x,y,’method’)

其中，

z：被插值點處的函數(shù)值；

x0,y0,z0：插值節(jié)點，x0,y0為向量，z0是矩陣，其列數(shù)等于x0的長度，行數(shù)等于y0的長度；

x,y：要計算函數(shù)值的點；

interp1(x0,y0,x，’method’)第27頁，講稿共65頁，2023年5月2日，星期三山體地貌要在某山區(qū)方圓大約27平方公里范圍內(nèi)修建一條公路，從山腳出發(fā)經(jīng)過一個居民區(qū)，再到達一個礦區(qū)。橫向縱向分別每隔400米測量一次，得到一些地點的高程：試做出該山區(qū)的地貌圖，并對幾種插值法進行比較.第28頁，講稿共65頁，2023年5月2日，星期三程序設計：clearx0=[1200:400:4000];y0=[1200:400:3600];z0=[11301250128012301040900500700;13201450142014001300700900850;139015001500140090011001060950;15001200110013501450120011501010;15001200110015501600155011801070;15001550160015501600160016001550;1480150015501510143013001200980];xi=1200:10:4000;%加密數(shù)據(jù)點yi=1200:10:3600;zil=interp2(x0,y0,z0,xi',yi,'linear');%線性插值zic=interp2(x0,y0,z0,xi',yi,'cubic');%三次插值zis=interp2(x0,y0,z0,xi',yi,'spline');%樣條插值subplot(2,2,1)mesh(x0,y0,z0)subplot(2,2,2)mesh(xi,yi,zil)subplot(2,2,3)mesh(xi,yi,zic)subplot(2,2,4)mesh(xi,yi,zis)第29頁，講稿共65頁，2023年5月2日，星期三第30頁，講稿共65頁，2023年5月2日，星期三二維插值的MATLAB實現(xiàn)2、另一個是離散數(shù)據(jù)節(jié)點的插值命令：

z=griddata(x0,y0,z0,x,y,’method’)

其中，

z：被插值點處的函數(shù)值；

x0,y0,z0：插值節(jié)點，x0,y0,z0均為向量；

x,y：被插值點；

method：插值方法，包括：

'linear'——線性插值；

'cubic'——三次插值；第31頁，講稿共65頁，2023年5月2日，星期三船在該海域會擱淺嗎？

在某海域測得一些點(x,y)處的水深z由下表給出，船的吃水深度為5英尺，在矩形區(qū)域（75，200）*（-50，150）里的哪些地方船要避免進入.第32頁，講稿共65頁，2023年5月2日，星期三解決問題的步驟：1.作出測量點的分布圖；2.求出矩形區(qū)域（75，200）*（-50，150）的細分網(wǎng)格節(jié)點之橫、縱坐標向量；3.利用MATLAB中的散點插值函數(shù)求網(wǎng)格節(jié)點的水深；4.作出海底曲面圖形和等高線圖；5.作出水深小于5的海域范圍.第33頁，講稿共65頁，2023年5月2日，星期三程序clearx0=[129140103.588185.5195105157.5107.57781162162117.5];y0=[7.5141.52314722.5137.585.5-6.5-81356.5-66.584-33.5];z0=[-4-8-6-8-6-8-8-9-9-8-8-9-4-9];subplot(2,2,1)plot(x0,y0,'+');%作出測量點的分布圖；x=75:1:200;%加密y=-50:1:150;[x,y]=meshgrid(x,y);z=griddata(x0,y0,z0,x,y,'cubic');subplot(2,2,2)mesh(x,y,z),%用插值方法求出網(wǎng)格節(jié)點處的z坐標矩陣,繪制出三維圖形subplot(2,2,3)meshc(x,y,z),%繪制等高線subplot(2,2,4)contour(x,y,z,[-5-5]);%水深5英尺處海底曲面的等高線gridon第34頁，講稿共65頁，2023年5月2日，星期三第35頁，講稿共65頁，2023年5月2日，星期三擬合的標準（1）用各點誤差絕對值的和表示（2）用各點誤差按絕對值的最大值表示（3）用各點誤差的平方和表示第36頁，講稿共65頁，2023年5月2日，星期三最小二乘擬合式中R2稱為均方誤差。由于計算均方誤差的最小值的原則容易實現(xiàn)而被廣泛采用。按均方誤差達到極小構造擬合曲線的方法稱為最小二乘法。第37頁，講稿共65頁，2023年5月2日，星期三++++++++++++++++++++++++++++++p=a1+a2xp=a1+a2x+a3x2p=a1+a2x+a3x2p=a1+a2/xp=aebxp=ae-bx將數(shù)據(jù)(xi,yi)i=1,…,n

作圖，通過直觀判斷確定p(x)：第38頁，講稿共65頁，2023年5月2日，星期三MATLAB---曲線擬合工具箱Matlab有一個功能強大的曲線擬合工具箱（CurveFittingToolbox

）cftool，使用方便，能實現(xiàn)多種類型的線性、非線性曲線擬合。調(diào)用：cftool界面如下所示第39頁，講稿共65頁，2023年5月2日，星期三第40頁，講稿共65頁，2023年5月2日，星期三“Data”按鈕數(shù)據(jù)的選取點擊“Data”按鈕，彈出“Data”窗口；利用Xdata和Ydata的下拉菜單讀入數(shù)據(jù)x,y，可修改數(shù)據(jù)集名“Datasetname”，然后點擊“Createdataset”按鈕，退出“Data”窗口，返回工具箱界面，這時會自動畫出數(shù)據(jù)集的曲線圖；

第41頁，講稿共65頁，2023年5月2日，星期三第42頁，講稿共65頁，2023年5月2日，星期三“Fitting”按鈕曲線擬合點擊“Fitting”按鈕，彈出“Fitting”窗口；點擊“Newfit”按鈕，可修改擬合項目名稱“Fitname”，通過“Dataset”下拉菜單選擇數(shù)據(jù)集，然后通過下拉菜單“Typeoffit”選擇擬合曲線的類型。第43頁，講稿共65頁，2023年5月2日，星期三SSEThesumofsquaresduetoerror.Thisstatisticmeasuresthedeviationoftheresponsesfromthefittedvaluesoftheresponses.Avaluecloserto0indicatesabetterfit.偏差平方和，越接近0越好第44頁，講稿共65頁，2023年5月2日，星期三R-square

Thecoefficientofmultipledetermination.Thisstatisticmeasureshowsuccessfulthefitisinexplainingthevariationofthedata.Avaluecloserto1indicatesabetterfit.復相關系數(shù)平方(決定系數(shù))，越接近1越好第45頁，講稿共65頁，2023年5月2日，星期三AdjustedR-squareThedegreeoffreedomadjustedR-square.Avaluecloserto1indicatesabetterfit.Itisgenerallythebestindicatorofthefitqualitywhenyouaddadditionalcoefficientstoyourmodel.修正的復相關系數(shù)平方，越接近1越好第46頁，講稿共65頁，2023年5月2日，星期三AdjustedR-square下列公式中的m為擬合函數(shù)中待估參數(shù)個數(shù)，如：對一元一次多項式擬合，f(x)=a+bx，此時m=2，n為數(shù)據(jù)點個數(shù)。該修正類似修正的樣本方差使其為總體方差的無偏估計。第47頁，講稿共65頁，2023年5月2日，星期三RMSETherootmeansquarederror.Avaluecloserto0indicatesabetterfit.偏差平方的均值的算術平方根，越接近0越好第48頁，講稿共65頁，2023年5月2日，星期三曲線擬合好壞如何評價首要指標是目標函數(shù)誤差最小（擬合度最大）；其次是應考慮關鍵點的吻合，這些關鍵點包括：初始點（有時是原點）、拐點、峰值點、極值點、中間點、漸近點、終值點等，在這些關鍵點上，數(shù)據(jù)觀察值點與函數(shù)值點應盡可能一致；再次是擬合的模型應盡可能簡單（模型的形式簡單，參數(shù)數(shù)少）。第49頁，講稿共65頁，2023年5月2日，星期三實踐中如何選擇模型？在數(shù)據(jù)擬合實踐中，理性模型畢竟是少數(shù)，大多數(shù)的情形是根據(jù)數(shù)據(jù)的趨勢尋找合適的模型，有時好幾個模型對數(shù)據(jù)都有較好的擬合，但通過對關鍵點的比較總會找到一種最合適的模型。在選擇不同的模型時，合理性和可解釋性是首要考慮的因素。第50頁，講稿共65頁，2023年5月2日，星期三美國人口問題據(jù)美國人口普查局數(shù)據(jù)：從1790每隔10年至2000年的總人口（單位：百萬）如下示

t=1790:10:2000;p=[3.9,5.3,7.2,9.6,12.9,17.1, 23.1,31.4,38.6,50.2,62.9,76,92, 105.7,122.8,131.7,150.7,179,205,226.5,251.4,281.422];第51頁，講稿共65頁，2023年5月2日，星期三美國人口問題t=1790:10:2000;p=[3.9,5.3,7.2,9.6,12.9,17.1, 23.1,31.4,38.6,50.2,62.9,76,92, 105.7,122.8,131.7,150.7,179,205,226.5,251.4,281.422];預測2001，2002年的美國人口數(shù)？并與調(diào)查數(shù)據(jù)285.318，288.369比較，選擇擬合較好的模型。第52頁，講稿共65頁，2023年5月2日，星期三Matlab求解在命令窗口輸入命令cftool回車，得擬合的圖形用戶界面第53頁，講稿共65頁，2023年5月2日，星期三第54頁，講稿共65頁，2023年5月2日，星期三第55頁，講稿共65頁，2023年5月2日，星期三結果分析LinearmodelPoly2:f(x)=p1*x^2+p2*x+p3Coefficients(with95%confidencebounds):p1=0.006757(0.006369,0.007144)p2=-24.32(-25.78,-22.85)p3=2.188e+004(2.049e+004,2.327e+004)Goodnessoffit:SSE:184.5