數(shù)值計(jì)算_第6章曲線擬合的最小二乘法

1、第6章 曲線擬合的最小二乘法6.1擬合曲線通過(guò)觀察或測(cè)量得到一組離散數(shù)據(jù)序列-,當(dāng)所得數(shù)據(jù)比較準(zhǔn)確時(shí),可構(gòu)造插值函數(shù)廠:門逼近客觀存在的函數(shù) 廠Fl罔，構(gòu)造的原則是要求插值函數(shù)通過(guò)這些 數(shù)據(jù)點(diǎn)，即一；.一_二。此時(shí)，序列- _： r.J.與兒是相等的。如果數(shù)據(jù)序列（心必）廠L曲，含有不可避免的誤差（或稱 噪音”，如圖6.1 所示；如果數(shù)據(jù)序列無(wú)法同時(shí)滿足某特定函數(shù)，如圖6.2所示，那么，只能要求所做逼近函數(shù)廠丫最優(yōu)地靠近樣點(diǎn)，即向量 - -7： 匸匸與，.珂、的誤差或距離最小。按：一與？之間誤差最小原則作為最優(yōu)”標(biāo)準(zhǔn)構(gòu)造的逼近函數(shù)，稱為擬合函數(shù)。圖6.1含有噪聲”的數(shù)據(jù)圖6.2 一條直線公路與

2、多個(gè)景點(diǎn)插值和擬合是構(gòu)造逼近函數(shù)的兩種方法。插值的目標(biāo)是要插值函數(shù)盡量靠近離散點(diǎn)；擬合的目標(biāo)是要離散點(diǎn)盡量靠近擬合函數(shù)。向量：與匸之間的誤差或距離有各種不同的定義方法。例如：用各點(diǎn)誤差絕對(duì)值的和表示：用各點(diǎn)誤差按模的最大值表示:用各點(diǎn)誤差的平方和表示：小音或nf(6.1)其中二稱為均方誤差，由于計(jì)算均方誤差的最小值的方法容易實(shí)現(xiàn)而被廣泛采用。按均方誤差達(dá)到極小構(gòu)造擬合曲線的方法稱為最小二乘法。本章主要講述用最小二乘法構(gòu)造擬合曲線的方法。在運(yùn)籌學(xué)、統(tǒng)計(jì)學(xué)、逼近論和控制論中，最小二乘法都是很重要的求解方法。例如，它 是統(tǒng)計(jì)學(xué)中估計(jì)回歸參數(shù)的最基本方法。關(guān)于最小二乘法的發(fā)明權(quán)，在數(shù)學(xué)史的研究中尚未定

3、論。有材料表明高斯和勒讓德分別獨(dú)立地提出這種方法。勒讓德是在1805年第一次公開(kāi)發(fā)表關(guān)于最小二乘法的論文，這時(shí)高斯指出，他早在1795年之前就使用了這種方法。但數(shù)學(xué)史研究者只找到了高斯約在1803年之前使用了這種方法的證據(jù)。在實(shí)際問(wèn)題中，怎樣由測(cè)量的數(shù)據(jù)設(shè)計(jì)和確定最貼近”的擬合曲線？關(guān)鍵在選擇適當(dāng)?shù)臄M合曲線類型，有時(shí)根據(jù)專業(yè)知識(shí)和工作經(jīng)驗(yàn)即可確定擬合曲線類型；在對(duì)擬合曲線一無(wú)所知的情況下，不妨先繪制數(shù)據(jù)的粗略圖形，或許從中觀測(cè)出擬合曲線的類型；更一般地，對(duì) 數(shù)據(jù)進(jìn)行多種曲線類型的擬合，并計(jì)算均方誤差，用數(shù)學(xué)實(shí)驗(yàn)的方法找出在最小二乘法意義 下的誤差最小的擬合函數(shù)。例如，某風(fēng)景區(qū)要在已有的景點(diǎn)之間

4、修一條規(guī)格較高的主干路，景點(diǎn)與主干路之間由各具特色的支路聯(lián)接。設(shè)景點(diǎn)的坐標(biāo)為點(diǎn)列- - ;設(shè)主干路為一條直線廠門一:二，即擬合函數(shù)是一條直線。通過(guò)計(jì)算均方誤差二：最小值而確定直線方程（見(jiàn)圖6.2）。6.2線性擬合和二次擬合函數(shù)線性擬合給定一組數(shù)據(jù)，做擬合直線 匸二廠-二，均方誤差為（6.2）是二元函數(shù)，Q（afb）的極小值要滿足整理得到擬合曲線滿足的方程:2-1S-1*1 3-17U-i丿(6.3)稱式（6.3）為擬合曲線的法方程。用消元法或克萊姆法則解出方程: 12 1*kt-1 J例6.1下表為P. Sale及R. Dybdall在某處作的魚類抽樣調(diào)查，表中為魚的數(shù)量，：為魚的種類。請(qǐng)用線

5、性函數(shù)擬合魚的數(shù)量和種類的函數(shù)關(guān)系。13151621222325293031361110111212131312141617404255606264707210013013142214212124172334解：設(shè)擬合直線-：；1；，并計(jì)算得下表:編號(hào)xyxy2 x113111431692151015022531611176256421122524415221226448421130344420169009563441891361640將數(shù)據(jù)代入法方程組（6.3）中，得到:解方程得：=8.2084 , :. = 0.1795擬合直線為： 處） =8.2084 + 0.1795 X二次擬合函數(shù)給

6、定數(shù)據(jù)序列一，用二次多項(xiàng)式函數(shù)擬合這組數(shù)據(jù)。- 2設(shè)1./，作出擬合函數(shù)與數(shù)據(jù)序列的均方誤差:Q（砌厲）=（應(yīng)）-耳）二丫（嗎+昭+的彳（6.4）由多元函數(shù)的極值原理,-的極小值滿足整理得二次多項(xiàng)式函數(shù)擬合的法方程:mU1iL2-1i 4al二i-1-J1如丿W 丿i-1*(6.5)解此方程得到在均方誤差最小意義下的擬合函數(shù)二I。方程組（6.5）稱為多項(xiàng)式擬合的法方程，法方程的系數(shù)矩陣是對(duì)稱的。當(dāng)擬保多項(xiàng)式階-時(shí)，法方程的系數(shù)矩陣是病態(tài)的，在計(jì)算中要用雙精度或一些特殊算法以保護(hù)解的準(zhǔn)確性。例6.2給定一組數(shù)據(jù)，如下表。用二次多項(xiàng)式函數(shù)擬合的這組數(shù)據(jù)。32101234230125解：設(shè)心；昭由計(jì)

7、算得下表:34129362781224488161331311000000011111112244881635159452781139287096將數(shù)據(jù)代入式（6.5），相應(yīng)的法方程為：解方程得：=0.66667, 1 = 1.39286 , ； = 0.13095=0.66667 1.39286 ： 0.13095 丄擬合曲線的均方誤差:工卄工（嗆）-訂=3.09524結(jié)果見(jiàn)圖6.3。圖6.3擬合曲線與數(shù)據(jù)序6.3解矛盾方程組在6.2節(jié)中用最小二乘法構(gòu)造擬合函數(shù)，本節(jié)中用最小二乘法求解線性矛盾方程組的方法構(gòu)造擬合函數(shù)。給定數(shù)據(jù)序列，作擬合多項(xiàng)式 卩二一一；，如果 要求一4 J過(guò)這些點(diǎn)，那么有

8、矛盾方程組:即:入(6.6)我們作一輔助函數(shù)這是自變量為（-I的多元函數(shù)，要使# define MAX_N 25 typedef struct tagPOINT/ (x_i, y_i)的最大維數(shù)點(diǎn)的結(jié)構(gòu) double x ；double y; POINT ；int main ( )int m ；int i ；POINT points MAX_N ；static double u11 ， u12， u21， u22， c1，c2；double a， b，tmp；printf ( n“Input n value ： ”；) / 輸入點(diǎn)數(shù) mscanf ( “% ，d”&m) ；if (m MAX_

9、N )printf (“The input n is larger then MAX_，N please redefine the MAX_N. n”；)return 1 ；if (m v = 0) printf (“Please input a number between 1 and % dn.”， MAX_N ) ；return 1 ； printf (“ Now in put t(e_i, y_i) , i=0 ,，% d: n”，m 1);/ 輸入點(diǎn)for (i=0 ； iv = m； i+)scanf ( “% ，lf ”&tmp) ；printf i.x=tmp ；scanf (

10、“ l, ”tmp);printf i.y=tmp ;/ scanf (“%lf %lf, & printf i.x , & printf i.y);for (i=0 ; i v = n; i+)/ 列出方程 U(a , b)T = cu21+ = points i.x ;u22+ = points i.x I points i.x ;c1+ = points i.y ;c2+ = points i.y 1 points i.y ;u12 = u21;u11 = m;/求解a = ( cl 1 u22 - c2 1 u12) / (u111 u22- u12 u21);b = ( cl *u2

11、1 c2u11)/ (u211 u12 u22 u11);printf(“Solve)(x) = %f +%f x n” a, b);/ 輸出return 0; / End of File計(jì)算實(shí)例用線性函數(shù)擬合下列數(shù)據(jù)（例 6.1 ）。13151621222325293031361110111212131312141617404255606264707210013013142214212124172334程序輸入輸出in put m value : 21Now in put the (x_i , y_i) , i=0 ,，20:13 11 15 10 16 11 21 12 22 12 23

12、13 25 13 29 1230 14 31 16 36 17 40 13 42 14 55 22 60 14 62 2164 21 70 24 72 17 100 23 130 34Solve： p(x) = 8.208408 + 0.179522 x程序6.2用形如 P=/的函數(shù)擬合定數(shù)據(jù)二-1宀一。算法描述輸入：值，及-L：# include v math.h# define MAX_N 25/ (x_i , y_i)的最大維數(shù)double x ；double y ； POINT ；int main ( ) int m ；int i ；POINT points MAX_N ；static

13、 double u11， u12，u21，u22，c1，c2；double A ， B ，tmp ；printf ( n“Input n value ： ”；)/ 輸入點(diǎn)數(shù) mscanf ( “% ，d”&m) ；if (m MAX_N )printf (“The input n is larger then MAX_，N please redefine the MAX_N. n”；)return 1 ；if (m v = 0)printf (“Please input a number between 1 and % dn.”， MAX_N ) ；return 1 ；printf (“ No

14、w in put t(e_i, y_i) , i=0 ,，% d: n”，m 1);/ 輸入點(diǎn)for (i=0 ; i v = m ; i+) seanf (“ % l, &tmp);printf i.x = tmp ;scanf (“ l, ”tmp);printf i.y = tmp ;/ scanf (“%lf %lf, & printf i.x , & printf i.y);for (i=0 ; i v = n; i+)/ 列出方程 X(a , b)T = Y其中Y為對(duì)原方程取對(duì)數(shù)后的數(shù)據(jù)u21+ = points i.x ;u22+ = points i.x * points i.x ;c1+ = log (points i.)y ;c2+ = points i.y 丨 log ( points i.y);u12 = u21;u11 = m;/求解a = ( cl 1 u22 - c2 u12) / (u111 u22- u12 u21);b = ( cl u21 - c2*u11) / (u21 * u12- u22u11);printf(“ Solve)(x) = % lfexp (%lf x) n, exp(A) , B);/輸出 /return 0;End of Fi