插值擬合方法

第4講數(shù)據(jù)的插值擬合法4.1數(shù)據(jù)的擬合方法4.2數(shù)據(jù)的插值方法4.1數(shù)據(jù)的擬合方法已知n個數(shù)據(jù)點（xi,yi)，i=1,2,…,n，xi互不相同，如何尋求函數(shù)y=f(x)，使f(x)在某種準(zhǔn)則下與這n個點最接近？

擬合模型通過尋找簡單的因果變量之間的數(shù)量關(guān)系，對未知的情形作出預(yù)測與預(yù)報.Remark：不必要求近似函數(shù)的曲線或曲面通過所有的數(shù)據(jù)點.一、擬合原理問題二、擬合模型的分類與方法1.直線擬合：用一次函數(shù)或稱線性函數(shù)擬合數(shù)據(jù).2.曲線擬合：若直線擬合效果不佳，可提高擬合精度，用曲線擬合數(shù)據(jù).常用的是二次函數(shù)、三次函數(shù)等高次多項式，有時也會用到指數(shù)函數(shù)、對數(shù)函數(shù)、三角函數(shù)等.3.觀察數(shù)據(jù)修勻：根據(jù)數(shù)據(jù)分布的總趨勢去剔除觀察數(shù)據(jù)中的偶然誤差，即數(shù)據(jù)修勻（或稱數(shù)據(jù)光滑）問題.4.分段擬合：在不同段上用不同的低次多項式進行擬合.選擇何種曲線擬合最好呢？首先可在坐標(biāo)軸上畫出數(shù)據(jù)的散點圖，通過觀察選擇幾種合適的曲線分別擬合，通過比較，哪條曲線的最小二乘指標(biāo)J最小即為最好的擬合曲線.表示第i個點的擬合值與實際值的絕對誤差第i個點的擬合值第i個點的實際值一、模型假設(shè)假設(shè)所測數(shù)據(jù)均為抽樣數(shù)據(jù).二、模型的分析與建立先畫出散點圖，見下圖.問題一

【溫度與電阻的關(guān)系模型】有一個對溫度敏感的電阻，現(xiàn)測得一組溫度t與電阻R的數(shù)據(jù)。見表9.1t20.532.751.07395.7R7658268739421032試給出溫度與電阻間的函數(shù)關(guān)系，并計算溫度為60度時的電阻值表9.1t=[20.532.7517395.7];R=[7658268739421032];plot(t,R,'*')xlabel('t')ylabel('R')觀察上圖不難發(fā)現(xiàn)：散點基本上在一條直線上，因此可假設(shè)電阻與溫度滿足一次函數(shù)（或稱線性函數(shù)）.

設(shè)擬合函數(shù)為y=a1+a2xt=[20.532.7517395.7];R=[7658268739421032];A=[ones(5,1),t'];b=R';a1a2=A\ba1a2=702.09683.3987三、模型求解可用Matlab計算.解法一：使用regress()函數(shù)t=[20.532.7517395.7];R=[7658268739421032];x=[ones(5,1),t'];y=R';b=regress(y,x,0.05)b=702.09683.3987解法二：使用ployfit()函數(shù)

cleart=[20.532.7517395.7];R=[7658268739421032];x=t';y=R';p=polyfit(x,y,1)p=3.3987702.0968即a1=702.0968,a2=3.3987，因此擬合的函數(shù)為y=702.0968+3.3987x將x=60代入計算得y=906cleart=[20.532.7517395.7];R=[7658268739421032];plot(t,R,'*')xlabel('t');ylabel('R')holdont=20:0.1:100;plot(t,702.0968+3.3987.*t)實際值與擬合圖regress()

函數(shù)主要用于線性擬合，在擬合時進行顯著性檢驗，故稱為回歸函數(shù).Polyfit()函數(shù)主要是利用多項式擬合.它可以是線性或非線性.Remark:polyfit(x,y,m)表示用m次多項式擬合數(shù)據(jù)x,y.

regress()

函數(shù)與polyfit()

函數(shù)的區(qū)別問題2

【農(nóng)業(yè)生產(chǎn)實驗?zāi)Ｐ汀吭谘芯哭r(nóng)業(yè)生產(chǎn)的試驗中，為分析某地區(qū)土豆產(chǎn)量與化肥的關(guān)系，得到了每公頃地的氮肥的施肥量與土豆產(chǎn)量的對應(yīng)關(guān)系，見圖9.2氮肥量(kg)03467101135202259336404471土豆產(chǎn)量(kg)15.1821.3625.7232.2934.0339.4543.1543.4640.8330.75請根據(jù)表9.2的數(shù)據(jù)，給出土豆產(chǎn)量與氮肥施肥量之間的關(guān)系表9.2一、模型假設(shè)

假設(shè)實驗數(shù)據(jù)為抽樣數(shù)據(jù).

假設(shè)其他化肥用量不變.二、模型分析與建立先用Matlab畫出散點圖.clearx=[03467101135202259336404471];y=[15.1821.3625.7232.2934.0339.4543.1543.4640.8330.75];plot(x,y,'+')xlabel('x(danfeiliang)')ylabel('y(tudouchanliang)')可以看出散點圖呈二次曲線圖形，可取如下擬合函數(shù)其中x表示氮肥量，y表示土豆產(chǎn)量。a,b,c為待定系數(shù)。三、模型求解abc=-0.00030.197114.7416clearx=[03467101135202259336404471];y=[15.1821.3625.7232.2934.0339.4543.1543.4640.8330.75];abc=polyfit(x,y,2)x1=0:471;y1=polyval(abc,x1);plot(x,y,'*',x1,y1)xlabel('x(danfeiliang)')ylabel('y(tudouchanliang)')問題3

【血藥濃度模型】通過實驗測得一次性快速靜脈注射300mg藥物后的血藥濃度數(shù)據(jù)，見表9.3t(h)0.250.511.523468y(ug/ml)19.2118.1515.3614.1012.899.327.455.243.01求血藥深度隨時間的變化規(guī)律y(t)表9.3一、模型假設(shè)1.假設(shè)t=0時，y=0.2.假設(shè)實驗數(shù)據(jù)為抽樣數(shù)據(jù)，能反映血藥濃度與時間的關(guān)系二、模型分析與建立

利用Matlab畫出散點圖.cleart=[0.250.511.523468];y=[19.2118.1515.3614.1012.89 9.327.455.243.01];plot(t,y,'o')xlabel('t(shijian)')ylabel('y(nongdu)')由圖可見，散點圖大致呈負(fù)指數(shù)函數(shù)形態(tài)，可令其中a,b>0為待定系數(shù).三、模型求解解法一：將其線性化處理.將等式

兩邊同時取對數(shù)，得令則方程變?yōu)橛肕atlab求解如下:cleart=[0.250.511.523468];y=[19.2118.1515.3614.1012.89 9.327.455.243.01];Y=log(y);b=polyfit(t,Y,1)b=-0.23472.9943即b1=2.9943b2=-0.2347則血藥濃度與時間的關(guān)系為cleart=[0.250.511.523468];y=[19.2118.1515.3614.1012.89 9.327.455.243.01];Y=log(y);b=polyfit(t,Y,1);a=exp(b(2));b=-b(1);plot(t,y,'*');holdont=0:0.01:9;plot(t,a.*exp(-b.*t),'r')xlabel('t')ylabel('y')解法2

使用非線性擬合函數(shù)nlinfit().首先建立Fun7_1.m文件如下：functiony=Fun7_1(a,t)y=a(1)*exp(-a(2)*t);不妨取a=20,由t=1時，y=15.36算出b=0.264取a0=[200.264];用Matlab求解如下cleart=[0.250.511.523468];y=[19.2118.1515.3614.1012.89 9.327.455.243.01];a0=[20,0.264];a=nlinfit(t,y,'fun7_1',a0)a=20.24130.2420t=[0.250.511.523468];y=[19.2118.1515.3614.1012.89 9.327.455.243.01];plot(t,y,'*');holdona0=[20,0.264];a=nlinfit(t,y,'fun7_1',a0);t=0:0.01:9;plot(t,a(1).*exp(-a(2).*t))holdonplot(t,19.9714.*exp(-0.2347.*t),'r')xlabel('t');ylabel('y')nlinfit函數(shù)的用法

nlinfit函數(shù)采用的迭代法，其中a0為迭代初值.若所給初值離最優(yōu)解比較近，則迭代求出該最優(yōu)解的概率就很高.

如何估計初值，暫無確切方法.可以在得到解后畫出函數(shù)圖形，看看實驗點是否都在曲線附近.若相差太大，可以考慮重新給出初值重新計算.問題4一、模型假設(shè)假設(shè)氯氣與生產(chǎn)時間之間滿足其中a,b為待定系數(shù).二、模型分析與建立此問題實質(zhì)上是確定待定系數(shù)a,b的值.三、模型求解首先定義非線性函數(shù)fun7_2.m文件：functiony=fun7_2(beta0,x)a=beta0(1);b=beta0(2);y=a+(0.49-a)*exp(-b*(x-8));然后在命令窗口中輸入clearx=[88101010101212121214141416161618182020202022222424242626262828303030323234363638384042];y=[0.490.490.480.470.480.470.460.460.450.430.450.430.430.440.430.430.460.450.420.420.430.410.410.40.420.40.40.410.40.410.410.40.40.40.380.410.40.40.410.380.40.40.390.39];beta0=[0.300.20];ab=nlinfit(x,y,'fun7_2',beta0)ab=0.39040.1028即a=0.3904,b=0.1028.所以模型函數(shù)為

作業(yè)：用Matlab畫出相關(guān)圖形.問題5一、模型假設(shè)假設(shè)中國人口的變化滿足一定規(guī)律.二、模型分析、建立與求解首先畫出散點圖.clearx=1949:5:1994;y=[5.466.778.19.19.810.311.311.8];plot(x,y,'*')xlabel('x(nianfen)')ylabel('y(renkoushu)')三、模型求解模型一：用線性函數(shù)擬合.設(shè)我國人口數(shù)量滿足以下模型y=a+bx其中a,b為待定系數(shù).通過計算得y=-283.232+0.148x模型二：用指數(shù)函數(shù)擬合.設(shè)我國人口數(shù)量滿足以下模型y=aebx其中a,b為待定參數(shù).模型兩邊取對數(shù)得lny=lna+bxy=4.1444*10-15e0.0179xx=1949:5:1994;y=[5.466.778.19.19.810.311.311.8];a=polyfit(x,y,1);x1=1949:0.1:1994;y1=a(2)+a(1)*x1;b=polyfit(x,log(y),1);y2=exp(b(2))*exp(b(1)*x1);plot(x,y,'*')holdonplot(x1,y1)holdonplot(x1,y2,'r')xlabel('x(nianfen)')ylabel('y(renkoushu)')用兩個模型分別計算相應(yīng)年度的人口數(shù)及其誤差見下表.分別計算兩個模型的最小二乘指標(biāo)得從而線性模型更適合中國人口的增長.4.2數(shù)據(jù)的插值方法若知道函數(shù)y=f(x)在n個互異的點x1,x2,…,xn的函數(shù)值y1,y2,…,yn

.如何估計此函數(shù)在另一點a的函數(shù)值？一、插值原理問題考慮構(gòu)造一個過x1,x2,…,xn的次數(shù)不超過n的多項式y(tǒng)=Ln(x)，使其滿足Ln(xk)=yk,k=1,2,…,n

然后用Ln(a)作為準(zhǔn)確值L(a)的近似值.這種方法叫作插值.

Remark:插值方法要求近似函數(shù)經(jīng)過所有的已知點.插值方法一般有：拉格朗日（Lagrange）插值、分段線性插值、埃爾米特（Hermite）插值、樣條插值以及分形插值等.Matlab的插值函數(shù)分為內(nèi)部插值與外部插值.內(nèi)部插值：要求已知點x是單調(diào)的，且被插值點xi不能超過x的范圍.如interp1()、interp2()、interpn()等.而griddata()既可以計算內(nèi)部插值，也可以計算外部插值.二、插值方法數(shù)據(jù)擬合要求得到一個具體的近似函數(shù)表達(dá)式，而數(shù)據(jù)插值不一定得到近似函數(shù)的表達(dá)式，僅通過插值方法找到未知點對應(yīng)的近似函數(shù)值.數(shù)據(jù)擬合得到的近似函數(shù)表達(dá)式不一定要經(jīng)過所有已知點，而數(shù)據(jù)插值要求近似函數(shù)經(jīng)過已知的所有數(shù)據(jù)點.數(shù)據(jù)擬合與數(shù)據(jù)插值的區(qū)別一、模型假設(shè)假設(shè)飛機機翼截面下輪廓的變化是連續(xù)的.二、模型分析與建立畫出散點圖問題7從散點圖很難觀察出函數(shù)形式,用插值的方法.由于只有一個變量且是單調(diào)的，可用一元插值函數(shù)interp1()進行計算.三、模型求解用Matlab求解如下：x1=[035791112131415];y1=[01.21.722.121.81.211.6];x=0:0.1:15;y=interp1(x1,y1,x,'spline')plot(x1,y1,'.',x,y,'r')gridtitle('spline')xlabel('x')ylabel('y')插值后的結(jié)果如下圖顯示.一、模型假設(shè)假設(shè)電壓的變化是連續(xù)的.二、模型分析與建立畫出散點圖.問題8由散點圖可觀察出電容器兩端的電壓隨充電時間的關(guān)系大致呈對數(shù)關(guān)系，本題不要求給出函數(shù)形式，因此同樣使用插值的方法.三、模型求解用Matlab求解如下：插值后結(jié)果如圖所示.functionwenti8_2x1=[12346.5912];y1=[6.27.38.299.610.110.4];x=0:0.1:12;y=interp1(x1,y1,x,'spline')plot(x1,y1,'.',x,y,'r')gridtitle('spline')xlabel('x')ylabel('y')Remark:一元插值函數(shù)interp1的基本調(diào)用格式為interp1(x,y,cx,‘method’).

其中x,y分別表示已知數(shù)據(jù)點的橫、縱坐標(biāo).x必須單調(diào)，cx為需要插值的橫坐標(biāo)，method為可選參數(shù)，可為以下四個值之一：

1）nearest——最近鄰點插值

2）linear——線性插值（可缺省）

3）spline——三次樣條插值

4）cubic——三次插值一、模型假設(shè)假設(shè)溫度的變化是連續(xù)的.二、模型分析與建立首先畫出散點圖.問題9從散點圖難以觀察出函數(shù)形式,下面用插值的方法.三、模型求解用Matlab求解如下.x1=1:1:12;y1=[589152529313022252724];x=[3.26.57.111.7];y=interp1(x1,y1,x)plot(x1,y1,'*',x,y,'O')xlabel('x')ylabel('y')y=10.200030.000030.900024.9000運行結(jié)果如圖所示.問題10一、模型假設(shè)1.假設(shè)河床的深度是連續(xù)變化的.2.假設(shè)緊靠河床鋪設(shè)電纜，即電纜長度等于河床長度.二、模型分析與建立畫出河床觀測的散點圖x=5:5:100;y=[2.412.962.152.653.124.235.126.215.68...4.223.913.262.852.353.023.634.123.462.080];y1=10-y;plot(x,y1,'*')axis([0100010])xlabel('x')ylabel('y')gridon

三、模型求解利用分段線性插值，并可以在此基礎(chǔ)上利用梯形法求積分命令trapz計算河床面積.同時利用每段連續(xù)線長度之和來近似河床曲線長度.x=5:5:100;y=[2.412.962.152.653.124.235.126.215.68...4.223.913.262.852.353.023.634.123.