官术网_书友最值得收藏!

任務(wù)二 查表

知識目標(biāo)

掌握查表的兩種方法插值和曲線擬合,掌握線性內(nèi)插法原理,了解其他插值原理和方法;了解曲線擬合的原理和方法。

能力目標(biāo)

能用Excel軟件采用單元格驅(qū)動、ExcelVBA自定義插值函數(shù)進(jìn)行插值,能用Excel軟件進(jìn)行曲線擬合。

在生產(chǎn)和科學(xué)實驗中,經(jīng)常采用表格的形式來反映實際生產(chǎn)數(shù)據(jù)和實驗觀測數(shù)據(jù),如表1-7為碳酸氫鈉溶解度實驗數(shù)據(jù)表,這種用數(shù)據(jù)表格形式給出的因變量溶解度(Ct)與自變量溫度(t)之間的關(guān)系,即函數(shù)Ct=ft),通常稱作列表函數(shù)。

表1-7 NaHCO3溶解度Ct(g/100g)數(shù)據(jù)表

列表函數(shù)通常具有以下特點:①自變量與函數(shù)值一一對應(yīng)(不允許多值);②函數(shù)值具有相當(dāng)可靠的精確度;③自變量與函數(shù)間的解析表達(dá)式可能不清楚,或者函數(shù)關(guān)系的解析表達(dá)式非常復(fù)雜不便于計算。

思考:1.何謂插值?常用的插值方法有哪幾種?

  2.何謂線性插值、Lagrange插值、三次樣條插值?

  3.何謂曲線擬合?

由表1-7可知,表格中只提供了幾個溫度為整數(shù)值時所對應(yīng)的溶解度數(shù)據(jù),當(dāng)我們需要某一表格中未給出溫度(如26.8℃)的溶解度數(shù)據(jù)時,該如何處理呢?

目前,處理表格函數(shù)的方法可分為兩大類:①插值;②曲線擬合。插值就是通過列表函數(shù)中若干點數(shù)據(jù)構(gòu)造一個比較簡單的函數(shù)來近似表達(dá)原表格函數(shù)進(jìn)行數(shù)據(jù)處理的方法,顯然選用不同的插值函數(shù),就有不同的計算方法和計算結(jié)果。曲線擬合也稱經(jīng)驗建模,通常是將表格中的全部數(shù)據(jù)作為處理對象,將其描述成數(shù)學(xué)表達(dá)式的方法。

本任務(wù)將著重介紹插值中的線性插值法和Lagrange多項式插值法及借助于Excel進(jìn)行曲線擬合的方法。

一、采用常用插值方法處理表格函數(shù)

為便于寫出計算通式,將表1-7中的實驗數(shù)據(jù)抽象為表1-8中的表格函數(shù)。

表1-8 表格函數(shù)

(一)線性插值法

線性插值就是將表格函數(shù)中的相鄰兩點之間的函數(shù)關(guān)系視為直線關(guān)系,即通過兩點的數(shù)據(jù)構(gòu)造一直線方程來處理位于兩點之間數(shù)據(jù)關(guān)系的方法,以點1、點2為例,兩點之間的因變量(y)的計算式為:

線性插值又稱兩點插值,其幾何意義見圖1-9。可見,線性插值將原函數(shù)(曲線)fx)近似為直線px)處理,原函數(shù)fx)在點Q的值yQ被現(xiàn)函數(shù)px)在點Q’的值所替代,這勢必導(dǎo)致誤差的存在。因此,線性插值是一種比較近似的插值處理方法,只有當(dāng)原函數(shù)fx)近似為直線或插值區(qū)間(x1x2)比較小時才可應(yīng)用,否則誤差較大。

圖1-9 線性插值法示意圖

【例1-3】 利用表1-7采用線性插值方法查溫度為26.8℃時NaHCO3的溶解度?

分析:該題在利用公式(1-2)時,公式中自變量x對應(yīng)于本題中的溫度t,因變量y對應(yīng)于系統(tǒng)中的溶解度Ct

解:由題意可得如下數(shù)據(jù)表

由式(1-2)得:

(二)拉格朗日(Lagrange)分段插值法

下面以一元三點拉格朗日插值法為例說明分段插值法的要點。

已知n+1個數(shù)據(jù)點的數(shù)值,其順序依次為x0x1,…,xn。插值點的數(shù)值為x,則

xxi時,使用結(jié)點xi-1xixi+1i=1,2,…,n-1);

x>xn-1時,使用結(jié)點xn-2xn-1xn

Lagrange插值計算式如下:

【例1-4】 利用Lagrange插值方法重新計算例1-3中溫度為26.8℃時溶解度值。

解:由式(1-3)得:

此結(jié)果與線性內(nèi)插值不完全一致,相對而言,查表時采用拉格朗日分段插值法的精度要高于線性插值,尤其是當(dāng)表格函數(shù)偏離直線關(guān)系時此法更適用。

*(三)其他插值法

前面所述的線性插值法和分段拉格朗日插值法只是眾多插值法中的兩種,工程上還存在著多種其他類型的插值法,限于篇幅以下只簡單介紹一下其他幾種常用的插值方法的名稱和定義,詳細(xì)內(nèi)容可參閱有關(guān)書籍。

(1)差商與牛頓插值法

拉格朗日插值法定義具有直觀性,計算機(jī)程序也很簡明,這是它的優(yōu)點。但如精度不滿足要求,需要增加插值節(jié)點時原來計算出的數(shù)據(jù)不能利用,必須重新計算。牛頓插值公式能克服這一缺點。

(2)差分與等距節(jié)點插值法

上面所介紹的插值法是采用節(jié)點任意分布的插值公式,但實際應(yīng)用時經(jīng)常遇到等距節(jié)點的情形,這時插值公式可進(jìn)一步簡化,計算也簡單得多。這就是差分與等距節(jié)點插值法的優(yōu)點。

(3)分段插值法

由于有些函數(shù)會在某些取值范圍發(fā)生突變,因此,對于處于此種類型下的插值,并不意味插值點越多精度越高,在實際進(jìn)行插值計算時,通常將插值范圍分為若干段,然后在每個分段上使用低階插值(如線性插值或二次插值),這就是所謂的分段插值法。

(4)三次樣條插值

在工程上經(jīng)常要求通過平面上n+1個已知點作一條連接光滑的曲線。譬如船體放樣與機(jī)翼設(shè)計均要求曲線不僅連續(xù)而且處處光滑。就高速飛機(jī)機(jī)翼的設(shè)計來說,要求盡可能采用流線型,使氣流沿機(jī)冀表面能形成平滑的流線,以減少空氣的阻力。解決此類問題,當(dāng)節(jié)點很多時,構(gòu)造一個高次插值多項式是不理想的,可能出現(xiàn)龍格(C.Runge)現(xiàn)象(所謂龍格現(xiàn)象即插值多項式在插值區(qū)間發(fā)生劇烈振蕩,出現(xiàn)函數(shù)不收斂的現(xiàn)象)。所以,在工程上進(jìn)行放樣時,描圖員常用富有彈性的細(xì)長木條稱樣條,把它用壓鐵固定在樣點上,其他地方讓它自由彎曲,然后畫下長條的曲線,稱為樣條曲線。該曲線可以看成由一段一段的三次多項式曲線拼合而成,在拼接處,不僅函數(shù)自身是連續(xù)的,而且它的一階和二階導(dǎo)數(shù)也是連續(xù)的,這種對描圖員描出的樣條曲線進(jìn)行數(shù)學(xué)模擬得出的函數(shù)叫做樣條插值函數(shù)。此插值函數(shù)需通過求解一個三對角矩陣(采用追趕法)才能求得。

思考:如何采用Excel進(jìn)行線性插值、Lagrange插值?

二、基于Excel單元格的常用插值法

Excel具有強(qiáng)大的計算功能,同樣可用來處理上述插值問題。下面采用單元格驅(qū)動法以例1-3的線性插值和例1-4的拉格朗日(Lagrange)分段插值為例加以說明。

步驟1:打開Excel表,在單元格B19、C19、D19、E19、B21、D21、E21依次輸入已知各數(shù)據(jù)點和待插數(shù)據(jù)點,如圖1-10所示。

步驟2:在待插值單元格C21中輸入線性插值計算公式(1-2),即:

=B21+(D21-B21)/(D19-B19)*(C19-B19),按“Enter”鍵,得計算值10.62;

若采用Lagrange插值計算法,則在單元格C23中輸入:

=B21*(C19-D19)*(C19-E19)/(B19-D19)/(B19-E19)+D21*(C19-B19)*(C19-E19)/(D19-B19)/(D19-E19)+E21*(C19-B19)*(C19-D19)/(E19-B19)/(E19-D19)

按“Enter”鍵,得計算值10.61。

由圖1-10可知,插值計算結(jié)果與前面例題計算結(jié)果完全一致。

圖1-10 采用Excel單元格驅(qū)動線性內(nèi)插與Lagrange內(nèi)插示意圖

思考:如何利用ExcelVBA開發(fā)線性插值和Lagrange插值自定義函數(shù)?

三、采用Excel曲線擬合法處理表格函數(shù)

線性插值處理表格函數(shù)時需要相鄰兩點數(shù)據(jù)且要求數(shù)據(jù)密度較大或相互之間有直線關(guān)系或近似直線關(guān)系,Lagrange分段插值需相鄰三點數(shù)據(jù),改變插值要求,常常需要重新選擇插值區(qū)間。當(dāng)需同時處理若干插值點時,計算工作量相對相大。當(dāng)然,通過編寫VBA自定義函數(shù)(參閱本任務(wù)技能拓展部分)也能利用插值方法實現(xiàn)對離散數(shù)據(jù)表連續(xù)查詢數(shù)據(jù)的功能。

若能將表格中離散的數(shù)據(jù)描述成數(shù)學(xué)表達(dá)式,則可很方便地實現(xiàn)連續(xù)插值。在數(shù)學(xué)上利用回歸分析可實現(xiàn)將離散型數(shù)據(jù)描述成函數(shù)式,這種方法也稱曲線擬合或經(jīng)驗建模(常用的方法有線性與非線性最小二乘法,具體可參閱有關(guān)數(shù)學(xué)資料,本任務(wù)的知識拓展部分簡單介紹了線性二乘法的相關(guān)知識)。Excel軟件已具備了部分曲線擬合功能,擬合的表達(dá)式類型有:①線性;②指數(shù);③對數(shù);④多項式;⑤冪;⑥移動平均六種類型,以下通過具體例子詳細(xì)說明Excel曲線擬合過程。

【例1-5】 將表1-7的溶解度數(shù)據(jù),用Excel擬合成函數(shù)關(guān)系Ct=ft),并利用此函數(shù)計算溫度為26.8℃時的溶解度值。

步驟1:將表1-7中的數(shù)據(jù)按列依次輸入如圖1-11所示A27~I(xiàn)28的區(qū)域。

圖1-11 例1-5數(shù)據(jù)表及擬合數(shù)據(jù)散點圖

步驟2:選中表中數(shù)據(jù)區(qū),即單元格A27~I(xiàn)28范圍,選擇“插入”?“散點圖”?如圖1-11。單擊“第一張散點圖”,出現(xiàn)如圖1-12所示的散點圖。

圖1-12 例1-5散點圖

步驟3:選中數(shù)據(jù)點?右擊數(shù)據(jù)點,出現(xiàn)如圖1-13所示的右鍵菜單。

圖1-13 添加趨勢線右鍵菜單

步驟4:單擊如圖1-13“添加趨勢線”,出現(xiàn)如圖1-14所示的“設(shè)置趨勢線格式”對話框。在“趨勢預(yù)測/回歸分析類型”中選擇“多項式”,在其后的“順序”選擇“3”,在“顯示公式”、“顯示R平方值”之前打勾,?如圖1-15所示的擬合公式,即:

圖1-14 “設(shè)置趨勢線格式”對話框

圖1-15 例1-5數(shù)據(jù)表擬合的初次表達(dá)式

y=2E-06x3+0.0004x2+0.126x+6.8595

式中“2E-06”是Excel科學(xué)記數(shù)表達(dá)形式,即2×10-6

注意:

①在圖1-14中,Excel默認(rèn)的趨勢線是“線性”,在選擇時可根據(jù)曲線形狀和“R平方值”來選擇回歸分析類型,R2值越趨近于1,所擬合的方程與數(shù)據(jù)點之間的關(guān)系越接近,利用此方程計算的結(jié)果與原數(shù)據(jù)點的誤差就越小。

②從圖1-15中可看出,Excel初次擬合的表達(dá)式中自變量、應(yīng)變量分別用xy表示,可通過修改成為tCt,使之與實際情況相符。

③初次表達(dá)式中,當(dāng)系數(shù)較小時,Excel默認(rèn)采用科學(xué)記數(shù)格式且只保留1位整數(shù),直接應(yīng)用此表達(dá)式時可能會帶來一定的誤差,此時可通過改變數(shù)據(jù)表達(dá)方式進(jìn)行改善。

步驟5:選中圖1-15的表達(dá)式,右擊?選中右鍵菜單“設(shè)置趨勢線標(biāo)簽格式”,如圖1-16?“設(shè)置趨勢線標(biāo)簽格式”對話框,如圖1-17,選擇其中“科學(xué)記數(shù)”,在“小數(shù)位數(shù)”中輸入你想要的位數(shù),如“5”?保留了5位小數(shù)的表達(dá)式。

圖1-16 例1-5數(shù)據(jù)表擬合表達(dá)式“設(shè)置趨勢線標(biāo)簽格式”

圖1-17 例1-5數(shù)據(jù)表最終擬合表達(dá)式

由圖1-17得擬合的表達(dá)式為:

Ct=2.31664×10-6t3+3.96321×10-4t2+0.126676t+6.85955

t=26.8℃代入上式得:

Ct=2.31664×10-6×26.83+3.96321×10-4×26.82+0.126676×26.8+6.85955

 =10.58(g/100g)

由此可見,計算結(jié)果與例1-3、例1-4十分接近。同時通過例1-3~例1-5的計算表明,即使同一數(shù)據(jù)來源處理相同的數(shù)據(jù),當(dāng)采用不同的處理方法時,其最終的結(jié)果不完全一致,具體采用何種方法,需認(rèn)真考慮,正確選擇。

用擬合表達(dá)式重新計算表1-7中的實驗點的溶解度數(shù)據(jù)并列于同一表格中,如表1-9所示。

表1-9 例1-5擬合值與原值的對比表

由表1-9可看出,每一點的擬合值與實驗值相對誤差均很小,最大誤差值不超過0.71%,工程上完全可用此擬合的表達(dá)式代替表1-7。

【例1-6】 如表1-10所示為一彈簧荷重與彈簧伸長之間的關(guān)系實驗數(shù)據(jù),試用線性內(nèi)插法、Lagrange分段內(nèi)插法、Excel曲線擬合法求當(dāng)荷重為11.24kg時彈簧長度?

表1-10 彈簧荷重與彈簧伸長的關(guān)系

解:(1)線性內(nèi)插法

由式(1-2)得:

(2)Lagrange分段內(nèi)插法

由式(1-3)得:

(3)Excel曲線擬合法

處理步驟同例1-5相同,具體步驟省略,其最終擬合結(jié)果如圖1-18。

圖1-18 例1-6數(shù)據(jù)表及最終擬合直線圖

根據(jù)力學(xué)知識,在彈簧的彈性限度內(nèi),符合虎克定律,彈簧伸長y與荷重x成正比,即yx的線性函數(shù),因此擬合時選擇“線性”。

由圖1-18得擬合的直線方程為:

y=0.62275x+30.00578

將荷重x=11.24kg代入上式得:

y=0.62275×11.24+30.00578=37.01(cm)

將擬合結(jié)果與實驗值比較,如表1-11。

表1-11 例1-6擬合計算與實驗值對比表

可見,擬合結(jié)果與實驗值相當(dāng)符合。因此,針對此題而言,方法(3)的計算結(jié)果優(yōu)于方法(1)和方法(2)。

【例1-7】 某化學(xué)反應(yīng)速率常數(shù)k與熱力學(xué)溫度T的實驗數(shù)據(jù)如表1-12所示,試用線性內(nèi)插法、Lagrange分段內(nèi)插法、Excel曲線擬合法求當(dāng)T=388.5K時的反應(yīng)速率?

表1-12 例1-7反應(yīng)速率與熱力學(xué)溫度的實驗數(shù)據(jù)表

解:(1)線性內(nèi)插法

由式(1-2)得:

(2)Lagrange分段內(nèi)插法

由式(1-3)得:

(3)Excel曲線擬合法

處理步驟同例1-5相同,具體步驟省略,其最終擬合結(jié)果如圖1-19。

圖1-19 例1-7數(shù)據(jù)表最終擬合曲線圖

由圖1-19得擬合的方程為:

k=2.05974×10-13e0.0667705T

T=388.5K代入上式得:k=2.05974×10-13e0.0667705×388.5=3.7980×10-2(min-1

由反應(yīng)動力學(xué)可知,反應(yīng)動力學(xué)方程通常表示為,其中反應(yīng)速率常數(shù)k與溫度T的關(guān)系符合阿累尼烏斯(S.A.Arrhenius)方程,即k=k0exp[-E/(RT)]的形式,顯然T-k之間的關(guān)系不符合直線關(guān)系。

本例若要使擬合方程符合阿累尼烏斯方程,可將原數(shù)據(jù)表1-12中T-k之間的關(guān)系轉(zhuǎn)換為1/T-k之間的關(guān)系,再利用Excel進(jìn)行曲線擬合,具體結(jié)果如圖1-20所示。

圖1-20 例1-7數(shù)據(jù)表采用1/T-k關(guān)系擬合曲線圖

由圖1-20得擬合的方程為:

k=3.27786×109e-9771.68/T

T=388.5K代入上式得:k=3.27786×109e-9771.68/388.5=3.9091×10-2(min-1

將兩種曲線擬合結(jié)果與原實驗數(shù)據(jù)進(jìn)行比較,計算結(jié)果如表1-13所示。

表1-13 例1-7兩種擬合計算結(jié)果與實驗數(shù)據(jù)比較表

例1-7小結(jié):

①由表1-13可看出,若已知自變量與因變量之間的關(guān)系,采用正確的擬合方程(1/T-k關(guān)系),所得結(jié)果與實驗值之間的誤差極小;若以1/T-k關(guān)系擬合計算k值為基準(zhǔn),則將其余三種方法的計算結(jié)果與之比較的誤差值列于表1-14。

表1-14 例1-7四種計算結(jié)果之間的比較表

②由表1-14看出,線性插值法計算的誤差較大,Lagrange分段插值的計算結(jié)果優(yōu)于直接采用T-k關(guān)系擬合曲線的計算結(jié)果。

因此,若能已知原表格函數(shù)的函數(shù)類型,再與Excel曲線擬合相結(jié)合,可以得到更好的擬合效果。

四、技能拓展——ExcelVBA自定義插值函數(shù)插值

從Office97開始,微軟為所有的Office組件引入了統(tǒng)一的應(yīng)用程序自動化語言——Visual Basic For Application(VBA),并提供了VBA的IDE(Integrated Development Environment)環(huán)境。VBA集成開發(fā)環(huán)境是進(jìn)行VBA程序和代碼編寫的地方,同一版本的Office共享同一IDE。VBA代碼和Excel文件是保存在一起的,可以通過打開VBA的IDE環(huán)境進(jìn)行程序設(shè)計和代碼編寫,以下以Excel2007為例介紹線性插值和Lagrange插值自定義函數(shù)LineIn、LagrangeIn的具體開發(fā)步驟。

步驟1:打開Excel?開發(fā)工具?Visual Basic,如圖1-21所示。

圖1-21 啟動Excel中VisualBasic編輯器

步驟2:選擇菜單“插入”?“模塊”,如圖1-22所示。在未插入模塊之前,“過程”是灰色的,不可用。

圖1-22 插入模塊

步驟3:選擇菜單“插入”?“過程”?“添加過程”對話框,如圖1-23。在名稱輸入框中輸入過程名:LineIn(意為線性插值,此過程名可根據(jù)函數(shù)的用途由用戶自己取名),在“類型”中選擇“函數(shù)”,“范圍”中兩者皆可選擇,此處選擇“公共的”。按“確定”按鈕,出現(xiàn)圖1-24的模塊1(代碼)編輯窗口。

圖1-23 “添加過程”對話框

圖1-24 模塊1(代碼)編輯窗口

步驟4:在函數(shù)LineIn的括號內(nèi)依次輸入計算用的變量:Y1、Y2、X1、X2、X,在函數(shù)體部分輸入具體的計算公式,完整的代碼如下:

Public Function LineIn(Y2,Y1,X2,X1,X)As Double

LineIn=Y1+(Y2-Y1)/(X2-X1)*(X-X1)     '線性內(nèi)插計算公式

End Function

步驟5:利用已開發(fā)的線性插值自定義函數(shù)即可在Excel中進(jìn)行插值計算,如圖1-25在單元格C22中輸入:=LINEIN(D21,B21,D19,B19,C19),按“Enter”鍵,可得與圖1-10單元格C21相同的計算結(jié)果。

圖1-25 自定義線性插值函數(shù)LineIn的計算示意

Lagrange自定義插值函數(shù)LagrangeIn的開發(fā)步驟與LineIn自定義函數(shù)的完全相同,其完整代碼如下:

Public Function Lagrange In(Y2,Y1,Y0,X2,X1,X0,X)AsDouble   'Lagrange內(nèi)插函數(shù)

Dim A1,A2,A3AsDouble

A1=Y0*(X-X1)*(X-X2)/(X0-X1)/(X0-X2)

A2=Y1*(X-X0)*(X-X2)/(X1-X0)/(X1-X2)

A3=Y2*(X-X0)*(X-X1)/(X2-X0)/(X2-X1)

LagrangeIn=A1+A2+A3

End Function

上述Lagrange自定義插值函數(shù)開發(fā)完成后,在圖1-25的單元格C24中輸入:

=LagrangeIn(E21,D21,B21,E19,D19,B19,C19),即可得結(jié)果10.61。顯然,采用該法與在單元格中直接輸入計算公式相比,既簡單又可避免輸入錯誤。

注意:在函數(shù)的括號內(nèi)輸入的變量不分大小寫,輸入順序也無規(guī)定,只要在使用此函數(shù)時變量的順序與此函數(shù)的變量順序一致即可。

五、知識拓展——線性最小二乘法

(一)關(guān)聯(lián)函數(shù)的選擇和線性化

實測數(shù)據(jù)關(guān)聯(lián)成數(shù)學(xué)模型的方法一般有以下幾種:

①具有一定的理論依據(jù),可直接根據(jù)機(jī)理選擇關(guān)聯(lián)函數(shù)的形式。

如反應(yīng)動力學(xué)方程通常表示為,其中反應(yīng)速率常數(shù)k與溫度T的關(guān)系符合阿累尼烏斯(S.A.Arrhenius)方程,k=k0exp[-E/(RT)]的形式。此法的關(guān)鍵在于確定上述公式中knk0E等未知系數(shù),以使模型密切逼近實測數(shù)據(jù)。這種模型稱為半經(jīng)驗?zāi)P停ぷ饕c在于參數(shù)估計。

②尚無任何理論依據(jù),但已有一些經(jīng)驗公式可選擇。

很多物性數(shù)據(jù)如熱容、密度、飽和蒸氣壓等與溫度的關(guān)系常表示為:

фT)=b0+b1T+b2T2+b3T3+b4lnT+b5/T

當(dāng)然不一定上述公式中六個系數(shù)都很重要,有的物性也只取前三、四項即可滿足精度要求,這樣可使模型簡單化。

③沒有任何經(jīng)驗可循的情況。

對于此類情況,通常只能將實驗數(shù)據(jù)畫出圖形與已知函數(shù)圖形進(jìn)行比較,選擇圖形接近的函數(shù)形式作擬合模型。

不論上述哪種情況,在選定關(guān)聯(lián)函數(shù)的形式之后,就是如何根據(jù)實驗數(shù)據(jù)去確定所選關(guān)聯(lián)函數(shù)中的待定系數(shù),最常用的方法是最小二乘法,具體可參閱有關(guān)數(shù)學(xué)專著。對于一些相對比較簡單的函數(shù)類型,通常采用線性最小二乘法,以下介紹線性二乘法原理及其應(yīng)用。

思考:何謂最小二乘法?

一元線性模型:

Y=A+BX  (1-4)

多元線性模型:

Y=B0+B1X1+B2X2+……  (1-5)

對于一些非線性模型,應(yīng)事先將其變換成線性形式,即線性化處理,然后再用線性最小二乘法進(jìn)行關(guān)聯(lián)。

表1-15列出了化工中常用的幾種函數(shù)類型及線性化的方法。此表中所列均為單變量問題,經(jīng)線性化處理后的線性模型均可統(tǒng)一用式(1-4)表示。

表1-15 常用函數(shù)線性化方法

對于多變量函數(shù)關(guān)系

y=fx1x2,…)

若采用冪函數(shù)的乘積作為關(guān)聯(lián)函數(shù),即將上面的函數(shù)關(guān)系寫成如下形式

可作如下線性化處理,令

Y=lnyX1=lnx1X2=lnx2,…

B0=lnaB1=bB2=c,…

經(jīng)線性化處理后的模型即式(1-5)。

對于一元非線性化方程,如:

y=a+bx+cx2+dx3…  (1-7)

Y=yX1=xX2=x2,…

B0=aB1=bB2=c,…

經(jīng)線性化處理后的模型也為式(1-5)。

(二)線性最小二乘法

關(guān)聯(lián)函數(shù)的形式確定之后,如何由實驗數(shù)據(jù)比較精確地去確定關(guān)聯(lián)函數(shù)中的待定系數(shù)仍是一個重要問題,最常用的方法就是線性最小二乘法,以下先以例1-6的表1-10數(shù)據(jù)加以說明。

將表1-10數(shù)據(jù)點畫在圖1-26上,可以看出盡管荷重與伸長兩者呈直線關(guān)系,但不同的人可能所畫的直線并不嚴(yán)格在一條直線上,說明由于讀數(shù)或其他影響因素造成數(shù)據(jù)包含有隨機(jī)誤差。

圖1-26 荷重與彈簧伸長長度關(guān)系

根據(jù)力學(xué)上的虎克定律,彈簧伸長y應(yīng)該與荷重x成正比,即yx的線性函數(shù),通過實驗確定比例系數(shù)(彈簧的彈性系數(shù))。一般地直線方程模型表示為

y'=a+bx  (1-8)

如果用直尺將圖1-26上的點連成直線,由于9個點不在一直線上,所以可以畫出多條直線。也即式(1-8)線性模型中參數(shù)ab可以有多種取值,于是產(chǎn)生這樣一個問題,圖1-26眾多的連線中哪一條直線最能體現(xiàn)物理現(xiàn)象的本質(zhì)呢?換句話說線性模型式(1-8)中截距a和斜率b取什么值為最佳選擇?為說明這個問題,這里引入“殘差”的概念。

設(shè)有n對實驗數(shù)據(jù)(xiyi)(i=1,2,…,n),需要尋找一個近似函數(shù)模型y'=fx)來擬合這一組數(shù)據(jù)。令第i點實測函數(shù)值yi與模型計算值之差為殘差,即

顯然,δi刻劃了yi與回歸模型計算值的偏離程度。如果每一個點的殘差δi=0,說明實驗數(shù)據(jù)(xiyi)完全可用直線擬合,但出于存在實驗誤差,δi=0是不可能的。也就是說最佳的ab應(yīng)使δi的和最小。但用δi的和最小原則估計參數(shù)ab,在應(yīng)用上不很方便,所以,一般采用最小二乘法,其原理可以這樣描述:所謂最小二乘原理就是使殘差的平方和最小,即

用最小二乘原理選擇最佳擬合模型的物理意義是顯見的,即在上例中找一條直線,使它與各實測點的距離(即δi)平方加和最小,這樣得到的擬合方程就可完全替代原數(shù)據(jù)表中數(shù)據(jù)之間的關(guān)系。這種采用使擬合函數(shù)值與原函數(shù)值之間的差值即殘差平方和最小的方法來確定擬合函數(shù)的方法稱最小二乘法。

將式(1-9)代入式(1-10),可得

顯然,Qab的函數(shù)。由數(shù)學(xué)分析多元函數(shù)求極值的必要條件,使Q最小的ab必須滿足以下方程組

由式(1)得

其中

分別表示yixi的平均值。由式(2)可推得

由式(3)、式(4)經(jīng)整理可得回歸系數(shù)b的計算公式

具體計算時,先由式(1-11)或式(1-12)求得b,再代入式(3)得a

由于計算ab的公式中所有的量都可以從觀測數(shù)據(jù)得出,因此回歸直線方程

y'=a+bx

便可確定。

為了簡化公式,將上述公式中的(下同)均用∑代替。令任一數(shù)據(jù)點xi與其平均值之差稱為離差xi的離差的平方和記為lxx,即

同樣,將yi的離差的平方和記為lyy,即

xi的離差與yi的離差的乘積之和記為lxy,即

則式(1-12)可表示為

利用式(1-17)、式(1-13)計算表1-10中的數(shù)據(jù),列于表1-16中。

表1-16 彈簧荷重與彈簧伸長的關(guān)系采用最小二乘法計算數(shù)據(jù)

lxx=816-(72)2/9=240

lxy=2668.58-72×314.89/9=149.46

b=lxy/lxx=149.46/240=0.62275

a=314.89/9-0.62275×72/9=30.006

故回歸方程為

y'=30.006+0.62275x

可見,與圖1-18Excel擬合的曲線表達(dá)式完全一致。

對于非線性模型作了線性化處理后的回歸方程

Y'=A+BX

需注意,此模型并非最終模型,最后需恢復(fù)線性化處理前的模型原樣。

以例1-7表1-12數(shù)據(jù)為例,由反應(yīng)動力學(xué)知,反應(yīng)速率常數(shù)與熱力學(xué)溫度的關(guān)系一般服從阿累尼烏斯(S.A.Arrhenius)方程,即:

k=k0exp[-E/(RT)]

由于此式為非線性函數(shù),需進(jìn)行線性化處理,兩邊取對數(shù)后得

Y=lnkX=1/TA=lnk0B=-E/R

Y=A+BX

將原k-T數(shù)據(jù)組換算成X-Y數(shù)據(jù)組列于表1-17。

表1-17 例1-7反應(yīng)速率與熱力學(xué)溫度的實驗數(shù)據(jù)采用最小二乘法計算數(shù)據(jù)

將相關(guān)數(shù)據(jù)代入式(1-11)得:

代入式(1-13)得

A=-18.190/5-(-9771.68)×13.073×10-3/5=21.910

因為A=lnk0,所以k0=3.27786×109

所以原關(guān)聯(lián)方程為:k=3.27786×109exp(-9771.68/T

同樣與圖1-20由Excel擬合的完全一致。

(三)曲線擬合效果分析

1.線性相關(guān)系數(shù)與顯著性檢驗

需要指出,曲線擬合處理的是隨機(jī)變量問題,觀測值xy不存在確定性函數(shù)關(guān)系,而只是一種相關(guān)關(guān)系。線性最小二乘法只適宜處理變量xy具有相關(guān)的問題,但在線性最小二乘法應(yīng)用過程中,并不需要限制兩個變量之間一定具有線性相關(guān)關(guān)系,就是說即使平面圖上一堆完全雜亂無章的散點,也可用此方法給它們配一條直線方程模型。顯然這樣做是毫無意義的。只有當(dāng)兩個變量大致呈線性關(guān)系時才適宜用直線模型去擬合數(shù)據(jù),因此必須給出一個數(shù)量性指標(biāo)描述兩個變量線性關(guān)系的密切程度,該指標(biāo)稱相關(guān)系數(shù),通常記作r,其表達(dá)式為

圖1-27說明了r取各種不同數(shù)值時散點的分布情況。

圖1-27 r不同時散點分布情況

r=0,此時lxy=0,因此b=0,即根據(jù)最小二乘法確定的回歸直線平行于x軸,這說明y的變化與x無關(guān),此時xy毫無線性關(guān)系,通常這時散點分布是完全不規(guī)則的。如圖1-27中的(a)。

②0<|r|<1,這是絕大多數(shù)情形,xy之間存在一定的線性關(guān)系。當(dāng)|r|越接近于1,說明線性相關(guān)越大,也就是散點與回歸直線越靠近。當(dāng)r>0時,b>0,yx增加而增加,稱為正相關(guān)。當(dāng)r<0時,b<0,yx增加而減小,稱為負(fù)相關(guān)。如圖1-27中的(b)、(c)。

③|r|=1,所有數(shù)據(jù)都在回歸直線上,此時,xy完全相關(guān),實際上此時xy間存在確定的線性函數(shù)關(guān)系。如圖1-27中的(d)。

利用式(1-18)對表1-10中彈簧荷重與彈簧伸長長度之間關(guān)系的回歸方程進(jìn)行線性相關(guān)系數(shù)的計算。

由式(1-16)得

lyy=11110.42-(314.89)2/9=93.124

由此可見,此例變量間線性相關(guān)程度很好。

必須指出,相關(guān)系數(shù)只表示xy的線性關(guān)系的密切程度,當(dāng)r很小或為零時,并不表示xy不存在其他關(guān)系。如圖1-27中的(e),xy呈某種曲線關(guān)系。可以對它進(jìn)行線性化處理,變換成為Y=A+BX直線方程模型,此時可以用相關(guān)系數(shù)來討論新變量XY之間線性相關(guān)程度,但是,新變量XY線性相關(guān)程度并不能直接說明原始數(shù)據(jù)xy與非線性模型擬合效果的優(yōu)劣。因此,對于非線性模型擬合的效果常用另一指標(biāo)——相關(guān)指數(shù)來衡量,記作R2,Excel曲線擬合中提供的R2即相關(guān)指數(shù),其計算式如下:

式中 yi——未經(jīng)線性變換的原始數(shù)據(jù);

y'i——非線性模型的計算值;

y——原始數(shù)據(jù)的平均值。

顯然R2<1,R2值越接近于1,擬合曲線效果越好,當(dāng)R2=1時,說明yi趨于一致,實測點完全落在擬合曲線上。

利用表1-16的數(shù)據(jù)采用式(1-19)計算相關(guān)指數(shù)R2,其計算數(shù)據(jù)列于表1-18。

表1-18 例1-6曲線擬合相關(guān)指數(shù)計算數(shù)據(jù)表

其中=314.89/9=34.988

將表1-18中的數(shù)據(jù)代入式(1-19)得:

與圖1-18中由Excel擬合的結(jié)果R2=0.99949完全一致。

同樣,利用表1-17的數(shù)據(jù)采用式(1-19)也可進(jìn)行R2計算,其計算數(shù)據(jù)列于表1-19。

表1-19 例1-7曲線擬合相關(guān)指數(shù)計算數(shù)據(jù)表

其中=0.19645/5=3.9290×10-2

將表1-19中的數(shù)據(jù)代入式(1-19)得:

與圖1-20中由Excel擬合的結(jié)果R2=0.999998完全一致。

2.相關(guān)系數(shù)r與顯著性水平α

對于一個具體問題,只有當(dāng)相關(guān)系數(shù)r的絕對值大到一定程度時方可用回歸直線來近似表示xy之間的關(guān)系。表1-20給出了r的起碼值,它與觀測次數(shù)n及顯著性水平α有關(guān),當(dāng)|r|大于表中相應(yīng)的值時,所回歸的直線才有意義。舉例來說,當(dāng)n-2=3時,即用5個數(shù)據(jù)來回歸直線時,相關(guān)系數(shù)r至少為0.878,所得直線方程的置信度為95%。

表1-20 相關(guān)系數(shù)r與顯著性水平α的關(guān)系

主站蜘蛛池模板: 库尔勒市| 徐州市| 瑞金市| 巴彦淖尔市| 高平市| 贵港市| 额济纳旗| 泸溪县| 玉龙| 枞阳县| 湘乡市| 固原市| 肇东市| 田东县| 玉屏| 马边| 南投县| 昌江| 凤山市| 镇沅| 丘北县| 梓潼县| 景泰县| 贵德县| 乌苏市| 竹北市| 维西| 屯昌县| 云安县| 永泰县| 岳普湖县| 三台县| 彭水| 兴宁市| 苏尼特左旗| 洛隆县| 定陶县| 洛浦县| 凌云县| 油尖旺区| 南澳县|