昆明理工數(shù)值分析大作業(yè)最小二乘法

1、數(shù)值分析實(shí)驗(yàn)報(bào)告課題八 曲線擬合的最小二乘法一、問題提出從隨機(jī)的數(shù)據(jù)中找出其規(guī)律性， 給出其近似表達(dá)式的問題， 在生產(chǎn)實(shí)踐和科學(xué)實(shí)驗(yàn)中大量存在，通常利用數(shù)據(jù)的最小二乘法求得擬合曲線。在某冶煉過程中， 根據(jù)統(tǒng)計(jì)數(shù)據(jù)的含碳量與時(shí)間關(guān)系， 試求含碳量 y與時(shí)間 t 的擬合曲線。二、實(shí)驗(yàn)要求t(分）0510152025303540455055y(×10-4)01.272.162.863.443.874.154.374.514.584.024.641、用最小二乘法進(jìn)行曲線擬合；2、近似解析表達(dá)式為( t) = a1t + a2t 2 + a3t 3 ；3、打印出擬合函數(shù)(t)，并打印出 (tj

2、 )與 y(tj)的誤差， j = 1,2,",12 ；4、另外選取一個(gè)近似表達(dá)式，嘗試擬合效果的比較；5、 * 繪制出曲線擬合圖。三、實(shí)驗(yàn)?zāi)康?、掌握曲線擬合的最小二乘法；2、最小二乘法亦可用于解超定線代數(shù)方程組；3、探索擬合函數(shù)的選擇與擬合精度間的關(guān)系。四、實(shí)驗(yàn)原理最小二乘法擬合在函數(shù)的最佳平方逼近中f(x)a,b,對(duì)已知函數(shù)f(x)的一組離散數(shù)據(jù)(xi,yi),i=0,1,m，yi=f(xi)，求函數(shù)擬合S*(x),記誤差i=S*(xi)-yi 要求一個(gè)函數(shù)與所給數(shù)據(jù)的曲線擬合，這里，要求一個(gè)函數(shù)與所給數(shù)據(jù)擬合，若記誤差，設(shè)是上線性無關(guān)函數(shù)族，在中找一函數(shù),使誤差平方和， （4

3、.1）這里. （4.2） 這就是一般的最小二乘逼近，用幾何語(yǔ)言說，就稱為曲線擬合的最小二乘法。 用最小二乘法求擬合曲線時(shí)，首先要確定的形式。這不單純是數(shù)學(xué)問題，還與所研究問題的運(yùn)動(dòng)規(guī)律及所得觀測(cè)數(shù)據(jù)有關(guān)；通常要從問題的運(yùn)動(dòng)規(guī)律或給定數(shù)據(jù)描圖，確定的形式，并通過實(shí)際計(jì)算選出較好的結(jié)果這點(diǎn)將從下面的例題得到說明。的一般表達(dá)式為（4.2）式表示的線性形式。若是次多項(xiàng)式，就是次多項(xiàng)式。為了使問題的提法更有一般性，通常在最小二乘法中都考慮為加權(quán)平方和 （4.3）這里是上的權(quán)函數(shù)，它表示不同點(diǎn)處的數(shù)據(jù)比重不同，例如，可表示在點(diǎn)處重復(fù)觀測(cè)的次數(shù)，用最小二乘法求擬合曲線的問題時(shí)，就是在形如（4.2）式的中求一

4、函數(shù)，使（4.3）式取得最小。它轉(zhuǎn)化為求多元函數(shù) （4.4）的極小點(diǎn)的問題。由求多元函數(shù)極值的必要條件，有 五、實(shí)驗(yàn)設(shè)計(jì)和實(shí)驗(yàn)步驟本實(shí)驗(yàn)利用最小二乘法和題目所給函數(shù)形式對(duì)所給數(shù)據(jù)進(jìn)行擬合。而后選用不同函數(shù)形式對(duì)所給數(shù)據(jù)進(jìn)行擬合得到擬合數(shù)據(jù)和所給數(shù)據(jù)的相對(duì)誤差。而后在正交多項(xiàng)式函數(shù)組的基礎(chǔ)上建立相對(duì)最佳基函數(shù)選擇指標(biāo)P并且記錄了相對(duì)最佳擬合函數(shù)。現(xiàn)將實(shí)驗(yàn)步驟陳述如下:（1） 建立正交多項(xiàng)式函數(shù)族建立的程序思路和擬合函數(shù)系數(shù)求解。（2） 編輯程序進(jìn)行計(jì)算記錄實(shí)驗(yàn)結(jié)果。（3） 對(duì)所得結(jié)果進(jìn)行總結(jié)和分析，找出存在的問題，探索擬合函數(shù)的選擇與擬合精度間的關(guān)系。評(píng)價(jià)本次實(shí)驗(yàn)的結(jié)果。六、編程思路實(shí)驗(yàn)使用ma

5、tlab工具編寫了計(jì)算使用正交函數(shù)族=span(0(x),1(x)n(x)進(jìn)行最小二乘法擬合的擬合函數(shù)各項(xiàng)系數(shù)的函數(shù)程序。這是比較簡(jiǎn)單的。并計(jì)算了其誤差。現(xiàn)將編程思路陳述如下：一，按照不同和擬合多項(xiàng)式求擬合函數(shù)及其誤差：（1） 利用已知條件生成相應(yīng)正交多項(xiàng)式函數(shù)族=span(0(x),1(x)n(x)。（2） 求格拉姆矩陣（即求其對(duì)角線元素）G= （3） 求（k(x)，f(x)）（4） 對(duì)應(yīng)除以相應(yīng)格拉姆矩陣對(duì)角線元素得到擬合函數(shù)系數(shù)向量a.（5） 生成擬合函數(shù)繪圖并且計(jì)算原來擬合數(shù)據(jù)中對(duì)應(yīng)x點(diǎn)的函數(shù)值，求其誤差。（6） 改變所選擇擬合函數(shù)的類型進(jìn)行擬合如上過程求其擬合函數(shù)誤差及相應(yīng)圖像。（4

6、） 探索擬合函數(shù)的選擇與擬合精度間的關(guān)系。利用已經(jīng)建立的最小二乘法擬合函數(shù)來建立相應(yīng)的最佳擬合函數(shù)尋找函數(shù)，是一個(gè)比價(jià)復(fù)雜的過程，這里我們僅僅使用上面的正交多項(xiàng)式生成的函數(shù)族=span(0(x),1(x)n(x)來尋找適宜的基函數(shù)值個(gè)數(shù)，通過誤差是否達(dá)到某一個(gè)極值來反應(yīng)其適宜度的變化情況，現(xiàn)將該程序的建立思路簡(jiǎn)述如下：（1） 利用擬合數(shù)據(jù)和擬合函數(shù)值的相對(duì)誤差建立判斷指標(biāo)P,（2） 利用判斷指標(biāo)P隨著基函數(shù)個(gè)數(shù)的變化找出擬合數(shù)據(jù)范圍內(nèi)的最大值，最小值和峰值。（3） 記錄峰值處的擬合函數(shù)個(gè)數(shù)。（4） 在一定范圍內(nèi)畫出判斷指標(biāo)P隨基函數(shù)個(gè)數(shù)變化的函數(shù)圖。七、程序建立及實(shí)驗(yàn)結(jié)果利用正交多項(xiàng)式求最小二

7、乘法擬合函數(shù):function w=zjdxsnh(n,x,y);m=length(x);mf=zeros(1,n);mp=zeros(1,n);P=ones(n,m);aw=0;af=zeros(1,n-2);bf=zeros(1,n-2);for i=1:m; aw=x(i)+aw;endaf(1)=aw/m;for i=1:m P(2,i)=x(i)-af(1);endfor i=3:n; d=i-2; sump1=sum(P(i-1,:).2); sump2=sum(P(i-2,:).2); sump3=sum(P(i-1).2).*x); af(i-1)=sump3/sump1; b

8、f(i-2)=sump1/sump2; P(i,:)=x.*P(i-1,:)-af(i-1)*P(i-1,:)-bf(i-2)*P(i-2,:);endfor j=1:n; for k=1:m; mf(j)=mf(j)+P(j,k)2; mp(j)=mp(j)+P(j,k)*y(k); endend%正交多項(xiàng)式的系數(shù)計(jì)算A=zeros(n,n+2);A(:,3)=1;A(2,4)=af(1);A(:,1)=zeros();A(:,2)=zeros();for i=3:n; d=i-1; for j=4:(n+2); A(i,j)=A(i-1,j-1)-af(i-2)*A(i-1,j)-bf(i

9、-2)*A(i-2,j-2); endend%最終系數(shù)計(jì)算AY=zeros(1,n);a=mp./mf;AX=zeros(n,n+2);for z=1:n AX(z,:)=a(z)*A(z,:);endfor b=1:n; for t=1:b; AY(b)=AX(n-t+1,b+3-t)+AY(b); endendw=AY;由于使用的是課本上的正交多項(xiàng)式進(jìn)行擬合，其系數(shù)的計(jì)算會(huì)帶來巨大的誤差，故最終的擬合多項(xiàng)式和原來的擬合數(shù)據(jù)發(fā)生了很大誤差，這說明利用以上這一思路和實(shí)驗(yàn)設(shè)計(jì)來求解題目是不適合的。即：采用生成多項(xiàng)式的方法和思路去尋找最小二乘法擬合是不合理的。至于后面探索精度和擬合多項(xiàng)式類型之間關(guān)系的工作至此無法做起！由此本實(shí)驗(yàn)轉(zhuǎn)換思路利用matlab強(qiáng)大的擬合功能來進(jìn)行：所使用的代碼如下：x=05 10152025303540455055;y=01.272.162.863.443.874.154.374.514.584.024.64;plot(x,y,'*');hold onp=polyfit(x,y,n);xx=0:1:55;yy=polyval(p,xx);plot(xx,yy);分別采用5次多項(xiàng)式，4次多項(xiàng)式和3次多項(xiàng)式所得