【文章內(nèi)容簡介】 得到平滑的數(shù)據(jù)，往往需要采用樣條插值。 時(shí)至今日，隨著電子計(jì)算機(jī)的普及，插值法的應(yīng)用范圍已涉及到了生產(chǎn)、科研、的各個(gè)領(lǐng)域。特別是由于航空、造船、精密機(jī)械加工等實(shí)際問題的需要，更使得插值法在實(shí)踐與理論上顯得尤其重要并得到了進(jìn)一步發(fā)展，尤其是近幾十年發(fā)展起來的樣條(Spline)插值，更獲得了廣泛的應(yīng)用。 最小二乘法簡介在研究兩個(gè)變量之間的關(guān)系時(shí)，可以用回歸分析的方法進(jìn)行分析。當(dāng)確定了描述兩個(gè)變量之間的回歸模型后，就可以使用最小二乘法估計(jì)模型中的參數(shù)，進(jìn)而建立經(jīng)驗(yàn)方程。簡單地說，最小二乘的思想就是要使得觀測點(diǎn)和估計(jì)點(diǎn)的距離的平方和達(dá)到最小。里的“二乘”指的是用平方來度量觀測點(diǎn)與估計(jì)點(diǎn)的遠(yuǎn)近（在古漢語中“平方”稱為“二乘”），“最小”指的是參數(shù)的估計(jì)值要保證各個(gè)觀測點(diǎn)與估計(jì)點(diǎn)的距離的平方和達(dá)到最小。例如，對于回歸模型 ，若，…，為收集到的觀測數(shù)據(jù)，則應(yīng)該用來估計(jì)，這里是的估計(jì)值。這樣點(diǎn)的估計(jì)就是，它們之間距離的平方就是，進(jìn)而最小二乘估計(jì)量就是使得 (33)達(dá)到最小值的參數(shù)。特別當(dāng)各個(gè)和相應(yīng)的估計(jì)值相等，即時(shí)，最 (34)達(dá)到最小值的參數(shù)。 如果我們能夠在固定解釋變量值的前提下觀測預(yù)報(bào)變量，就認(rèn)為解釋變量的觀測值和估計(jì)值相等，從而可以通過(2)，人們常忽略“各個(gè)和相應(yīng)的估計(jì)值相等”的條件，而把(2)式的最小值點(diǎn)稱為參數(shù)的最小二乘估計(jì)量，其原因有二：其一是不知道最小二乘方法的原理；其二是找不到估計(jì)量的合理數(shù)學(xué)表達(dá)式，也就無法通過(1)式求最小二乘估計(jì)量，只好用(2)式的最小值點(diǎn)作為參數(shù)的估計(jì)。在教科書中，已知(x1，y1)，(x2，y2)，…，()是變量X和Y的一組觀測數(shù)據(jù)，要估計(jì)的是回歸直線方程y=b0＋b1x中參數(shù)b0，b1的值。所以這時(shí)目標(biāo)函數(shù)為。于是這時(shí)的最小二乘法就是尋求b0，b1的值，使在各點(diǎn)處的偏差－(b0＋b1xi)（=1，2，…，n），有意思的事情是：估計(jì)得到的直線=b0＋b1x一定經(jīng)過觀測數(shù)據(jù)點(diǎn)的中心(，)（，）。進(jìn)一步，若觀測數(shù)據(jù)全部落在某一直線上。關(guān)于最小二乘估計(jì)的計(jì)算，涉及更多的數(shù)學(xué)知識，這里不想詳述。其一般的過程是用目標(biāo)函數(shù)對各bi求偏導(dǎo)數(shù)，并令其等于0，并創(chuàng)設(shè)了解線性方程組的消元法——高斯消元法。 插值法和最小二乘法的比較在科學(xué)研究與工程技術(shù)中，常常需要從一組測量數(shù)據(jù)()(=0,1,……,n)出發(fā)，尋找變量x與y的函數(shù)關(guān)系的近似表達(dá)式，且是從給定的一組實(shí)驗(yàn)數(shù)據(jù)出發(fā)，尋求已知函數(shù)的一個(gè)逼近函數(shù)，使得逼近函數(shù)從總體上與已知函數(shù)的偏差按某種方法度量能達(dá)到最小，而又不需要通過全部的點(diǎn)(),這就是最小二乘曲線擬合。從計(jì)算的角度看，最小二乘法與插值法類似，都是處理數(shù)據(jù)的算法。但從創(chuàng)設(shè)的思想看，二者卻有本質(zhì)的不同。前者尋求一條曲線，使其與觀測數(shù)據(jù)“最接近”，其目的是代表觀測數(shù)據(jù)的趨勢；后者則是使曲線嚴(yán)格通過給定的觀測數(shù)據(jù)，其目的是通過來自函數(shù)模型的數(shù)據(jù)來近似刻畫該函數(shù)。在觀測數(shù)據(jù)帶有測量誤差的情況下，就會使得這些觀測數(shù)據(jù)偏離函數(shù)曲線，結(jié)果使得與觀測數(shù)據(jù)保持一致的插值法不如最小二乘法得到的曲線更符合客觀實(shí)際。因此用最小二乘法分析E和t的非線性關(guān)系是更符合實(shí)際的方法。 最小二乘法的應(yīng)用在《90國際溫標(biāo)通用熱電偶分度表手冊》中給出了原函數(shù)E=f(t)，在一定溫度范圍內(nèi)確定步長，均勻取點(diǎn)，得到一組數(shù)據(jù)(t0，E0)、(t1,E1)……(n=3)，這些就是最小二乘法中所提到的理論值。利用這組數(shù)據(jù)分段擬合，求出反函數(shù)t=g(E)，例如選擇(n=3)為擬合模型多項(xiàng)式(根據(jù)情況選擇多項(xiàng)式的階數(shù)，一般選擇的階數(shù)為35)，我們的目的是求多項(xiàng)式系數(shù)a0,a1,a2,a3。根據(jù)最小二乘定義，應(yīng)使各節(jié)點(diǎn)處的誤差平方和最小，即r=最小。所以,應(yīng)滿足條件(=0,1,2,3)，經(jīng)推導(dǎo)得到方程組：= （35）簡記為A a=b。由于一組數(shù)據(jù)(t0，E0 )(t1，E1)……(n=3)均已知，相當(dāng)于未知數(shù)為 (=0，1，2，3)，求出 (=0，1，2，3)即得出最小二乘擬合公式。我們用高斯約當(dāng)消去法求解線性方程組。高斯(Gauss)消去法是大家所熟悉的方法，雖然它是一種古老的方法，但用計(jì)算機(jī)求解線性方程組的實(shí)踐表明，它仍然是直接法中最常用和最有效的方法之一。其基本思想是逐步的消去未知數(shù)，把原方程組轉(zhuǎn)化為等價(jià)的三角方程組，這樣就極易求解，通過回代過程即可逐一求出各未知數(shù)。我們采用的是較高斯(Gauss)消去法更為簡單的高斯—約當(dāng)消去法。高斯—約當(dāng)消去法(Gauss—Jordan，簡稱約當(dāng)消去法)是一種無回代過程的消去法，其基本思想是對方程的增廣矩陣[A︱b]進(jìn)行初等變換，將A矩陣各列的非主元素全部化為零，主對角線上的元素化為1。Jordan法與Gauss法的區(qū)別在于Jordan法不僅將主元以下的未知數(shù)系數(shù)化為0，主元之上的系數(shù)同樣化為0，并且主元素自身歸1。原A矩陣位置變換為單位矩陣，而常數(shù)項(xiàng)位置上的值就是所求系數(shù)a (=0，1，2，3)，因此不需要回代過程。第四章 軟件編程 程序設(shè)計(jì)由前面介紹的數(shù)值分析方法方案知：最小二乘擬合生成系數(shù)矩陣方程和高斯約當(dāng)消去法解矩陣方程子程序是主程序中必須不斷調(diào)用高斯約當(dāng)消去法解矩陣方程的子程序，求解出多項(xiàng)式的系數(shù)。在本程序中對由E=f(t)精確公式取一組數(shù)據(jù)時(shí)的溫度范圍、步長、以及生成的擬合多項(xiàng)式的階數(shù)均使用了宏定義，也就是說此程序具有極大的通用性，適合各種類型的熱電偶、熱電阻。將用此程序計(jì)算所得的分段擬合多項(xiàng)式的系數(shù)存人微處理器的ROM區(qū)內(nèi)，可以使智能儀表完成對各類熱電偶、熱電阻溫度傳感器的精確測量。程序流程見圖41。在設(shè)計(jì)過程中應(yīng)注意以下幾點(diǎn)：是9次冪多項(xiàng)式（針對0—1372C，K偶），所以在式41,42中，n的取值為9，而不是擬合多項(xiàng)式的階數(shù)；(n)，其中n便是有效數(shù)字的位數(shù)；在解方程組時(shí)，只需n+1個(gè)點(diǎn)，所以在確定溫度范圍之后，應(yīng)該在這個(gè)范圍內(nèi)均勻取五個(gè)點(diǎn)。在求得系數(shù)之后，即得到t=g[E]，然后從分度表中以一定的間隔取E1,E2,……，將E的值代入t=g[E]中，得到溫度的估計(jì)值T，再用T與t比較，︱tT︱≤；以K偶為例，利用所編程序求得K偶在一般測量溫度范圍為50C —1300C的分段擬合多項(xiàng)式的系數(shù)。K偶E—t關(guān)系如下：當(dāng)溫度范圍為0—1372C時(shí)， (41)當(dāng)溫度范圍在270—0C時(shí)， (42)式中，E為電動勢，單位為mV；t是以ITS90為依據(jù)的溫度，與攝氏溫度t一樣；a,a,和C是9次冪多項(xiàng)式加上指數(shù)表達(dá)式的有關(guān)系數(shù) 。在用C語言表示時(shí)，應(yīng)利用一個(gè)for循環(huán)先表示出；如果沒有直接給數(shù)組賦值，那么，先將數(shù)組中所以元素設(shè)置為0；在使用pow()時(shí)，必須使用頭文件include ；A a = b中，A的每一個(gè)元素分別是由內(nèi)積得來，而可分別由一個(gè)矩陣的不同行表示；其中。b的每一個(gè)元素是由內(nèi)積而來，其中；得到A，b后，將它們組成增廣矩陣[A︱b]；在用高斯約當(dāng)消去法的過程中，可以先將每一列的元素變?yōu)?，前提是這些元素不為零，再將一些元素變?yōu)榱?。例如，將第二行到最后一行的第一個(gè)元素變?yōu)榱?，可先將第一列的所有元素變?yōu)?，然后用第二行到最后一行的元素減去第一行的元素，對應(yīng)相減。后面的變法以此類推。此程序是由幾個(gè)循環(huán)組成，較復(fù)雜，編程時(shí)應(yīng)注意每一個(gè)循環(huán)的適用范圍切記主對角線上的元素要變?yōu)?。且最后一列的元素就是所求的系數(shù)。開始宏定義溫度區(qū)間、步長length、階數(shù)n。利用E=f（t），計(jì)算出一組數(shù)據(jù)[t1，E1],[t2,E2]……[]利用數(shù)組表示矩陣A，b，c調(diào)用高斯約當(dāng)消去法解矩陣方程的子程序，求解出多項(xiàng)式的系數(shù)誤差△t≤輸出多項(xiàng)式系數(shù)a，并利用求出表52中的T結(jié)束NY圖41 流程圖 程序說明此程序適用于0C —1372C求多項(xiàng)式的系數(shù)a0,a1,a2,a3,a4；可將溫度劃分為90~250，250~650，650~1372三段。在編程的過程中發(fā)現(xiàn)用公式，計(jì)算t[n+1],E[n+1]與分度表中的值對比存在一定的誤差。為了避免誤差，我們可以直接給出t[n+1],E[n+1]中元素的值（可以直接在分度表中查找，若n=4，則分別在90~250，250~650，650~1372溫度范圍內(nèi)分別取5個(gè)點(diǎn)，且是均勻的）。例如，給出t[n+1]={650,830,1010,1190,1370}, E[n+1]= {,}，運(yùn)行程序就可以得到所求系數(shù)，然后用E0[12]={,,,}（從分度表中找出的）和公式計(jì)算出溫度T，再將T和t進(jìn)行比較。若不滿足要求，則縮小溫度的范圍，以便得到更高的精度。（后來補(bǔ)充的）分度表中給出的系數(shù)a0是錯(cuò)的，正確的是a0=，現(xiàn)在我們檢驗(yàn)誤差時(shí)不需要再從分度表中一個(gè)一個(gè)查值了，可以直接定義步長length=5，通過公式計(jì)算出一組(t0,E0)、(t1,E1)……(tn,En)，再將得到的E的值代入中，得到估計(jì)值T。對于90~250，由于溫度小于0時(shí)，其系數(shù)Ci和公式與溫度大于0是不同的，所以用2中的方法編程簡便。仿真的結(jié)果是按照2中的方法進(jìn)行的。附錄中的程序是更改后的程序，不用一個(gè)個(gè)取點(diǎn)。第五章 仿真結(jié)果與結(jié)論 仿真結(jié)果的說明 以下仿真結(jié)果只適用于K偶650~1372，其余溫度段250~650，90~250，需要更改t[n+1],E[n+1