第六章MATLAB解曲線擬合問題

1、黑龍江八一農(nóng)墾大學數(shù)學建模選修課-1 -黑龍江八一農(nóng)墾大學數(shù)學建模選修課-# -黑龍江八一農(nóng)墾大學數(shù)學建模選修課第六章 MATLA解曲線擬合問題-# -黑龍江八一農(nóng)墾大學數(shù)學建模選修課-# -黑龍江八一農(nóng)墾大學數(shù)學建模選修課曲線擬合就是計算出兩組數(shù)據(jù)之間的一種函數(shù)關系，由此可描繪其變化曲線及估計非采集數(shù)據(jù)對應的變量信息。曲線擬合有多種方式，下面是一元函數(shù)采用最小二乘法對給定數(shù)據(jù)進行多項式曲線擬合， 最后給出擬合的多項式系數(shù)。 1線性擬合一、數(shù)學模型y=X B + B是p 1的參數(shù)向量；是服從標準正態(tài)分布的隨機干擾的n 1的向量；y為n 1的向量；X為n p矩陣。二、 求解線性擬合函數(shù)regre

2、ss調(diào)用格式：b=regress(y,X)b,bi nt,r,ri nt,stats= regress(y,X)b,b in t,r,ri nt,stats= regress(y,X,alpha)說明：b=regress(y,X) 返回X處y的最小二乘擬合值。 bint返回3的95%的置信區(qū)間。r 中為形狀殘差，rint中返回每一個殘差的 95%置信區(qū)間。Stats向量包含R2統(tǒng)計量、回歸的F 值和p值。三、舉例例1:設y的值為給定的x的線性函數(shù)加服從標準正態(tài)分布的隨機干擾值得到。即y=10x+ ；求線性擬合方程系數(shù)。程序： 編寫M文件mainG.m如下： x=ones(10,1) (1:10

3、) %ones(m,n)產(chǎn)生一個 m行n列的元素全為1的矩陣%(1： 10)產(chǎn)生一個10行1列的元素值從110的矩陣%A B將矩陣A和B拼接成新矩陣或者寫成A ；By=x*10;1+normrnd(0,0.1,10,1)%R= normrnd(MU,SIGMA,m,n)依據(jù)參數(shù) MU SIGMA生成一個隨機數(shù)，m和n是R的行數(shù)和列數(shù).b,b in t=regress(y,x,0.05)結(jié)果：x =111 21314151 6171 8191 10-# -黑龍江八一農(nóng)墾大學數(shù)學建模選修課y = 10.956711.833413.012514.028814.885416.119117.118917

4、.996219.032720.0175b = 9.92131.0143bint =9.7889 10.05370.99301.0357即回歸方程為：y=9.9213+1.0143x 2多項式曲線擬合二、求解多項式曲線擬合函數(shù)ployfit調(diào)用格式：p=polyfit(x,y, n)p,s= polyfit(x,y, n)說明：x,y為數(shù)據(jù)點，n為多項式階數(shù)，p為幕次從高到低的多項式系數(shù)向量。矩陣 于生成預測值的誤差估計。(見下一函數(shù)polyval)三、舉例例2:由離散數(shù)據(jù)擬合出多項式。x0.1.2.3.4.5.6.7.8.91y.3.511.41.61.9.6.4.81.52程序：x=0:.1

5、:1;y=.3 .5 1 1.4 1.6 1.9 .6 .4 .8 1.5 2；n=3;p=polyfit(x,y, n)xi=lin space(0,1,100);z=polyval(p,xi);%多項式求值plot(x,y,o,xi,z,k:,x,y,b)legend(原始數(shù)據(jù),3 階曲線)結(jié)果：p = 16.7832 -25.7459 10.9802 -0.003532多項式為：16.7832x -25.7459x +10.9802x-0.0035曲線擬合圖形：例 3: x=1:20,y=x+3*sin(x)程序：-3 -黑龍江八一農(nóng)墾大學數(shù)學建模選修課y=x+3*s in (x);p=

6、polyfit(x,y,6)xi=1i nspace(1,20,100);z=poyval(p,xi);%多項式求值函數(shù)plot(x,y, o ,xi,z, k: ,x,y, b) legend( 原始數(shù)據(jù)，6階曲線)結(jié)果：-9.7295 11.3304圖10階曲線p =0.0000 -0.00210.0505 -0.59713.6472圖6階曲線 例4:再用10階多項式擬合 程序：x=1:20;y=x+3*s in (x);p=polyfit(x,y,10)xi=lin space(1,20,100);z=polyval(p,xi);plot(x,y,o,xi,z,k:,x,y,b) leg

7、end(原始數(shù)據(jù),10階多項式)11.2360結(jié)果：p = Columns 1 through 70.0000 -0.00000.0004 -0.01140.1814 -Colu mns 8 through 11-42.0861 88.5907 -92.8155 40.2671說明：可用不同階的多項式來擬合數(shù)據(jù)，但也不是階數(shù)越高擬合的越好。 3多項式曲線求值函數(shù)polyval()調(diào)用格式：y=polyval(p,x)y,DELTA=polyval(p,x,s)說明：y=polyval(p,x) 為返回對應自變量x在給定系數(shù)P的多項式的值。y,DELTA=polyval(p,x,s) 使用pol

8、yfit 函數(shù)的選項輸出 s得出誤差估計 丫士 DELTA它 假設polyfit 函數(shù)數(shù)據(jù)輸入的誤差是獨立正態(tài)的，并且方差為常數(shù)。則丫-DELTA將至少包含50%勺預測值。 4多項式曲線擬合的評價和置信區(qū)間函數(shù)polyconf()調(diào)用格式：Y,DELTA=polyco nf(p,x,s)Y,DELTA=polyco nf(p,x,s,alpha)說明：Y,DELTA=polyconf(p,x,s)使用polyfit 函數(shù)的選項輸出 s給出Y的95%置信區(qū)間Y - DELTA它假設polyfit 函數(shù)數(shù)據(jù)輸入的誤差是獨立正態(tài)的，并且方差為常數(shù)。1-alpha為置信度。例5:給出上面例1的預測值及

9、置信度為90%勺置信區(qū)間。程序： x=0:.1:1;y=.3 .5 1 1.4 1.6 1.9 .6 .4 .8 1.5 2n=3;p,s=polyfit(x,y, n)alpha=0.05;Y,DELTA=polyco nf(p,x,s,alpha)結(jié)果： p = 16.7832 -25.7459 10.9802 -0.0035s = R: 4x4 doubledf: 7n ormr: 1.1406Y = Colu mns 1 through 7-0.00350.85381.29701.42661.34341.14800.9413Colu mns 8 through 110.82380.89

10、631.25942.0140DELTA = Colu mns 1 through 71.36391.15631.15631.15891.13521.12021.1352Colu mns 8 through 111.15891.15631.15631.3639 5穩(wěn)健回歸函數(shù)：robust()穩(wěn)健回歸是指此回歸方法相對于其他回歸方法而言，受異常值的影響較小。調(diào)用格式：b=robustfit(x,y)b,stats=robustfit(x,y) b,stats=robustfit(x,y, wfun ,tune, con st )說明：b返回系數(shù)估計向量；stats返回各種參數(shù)估計；wfun指定一

11、個加權函數(shù)；tune 為調(diào)協(xié)常數(shù)；con st的值為on(默認值)時添加一個常數(shù)項；為off時忽略常數(shù)項。例6:演示一個異常數(shù)據(jù)點如何影響最小二乘擬合值與穩(wěn)健擬合。首先利用函數(shù)y=10-2x加上一些隨機干擾的項生成數(shù)據(jù)集，然后改變一個y的值形成異常值。調(diào)用不同的擬合函數(shù)，通過圖形觀查影響程度。程序：x=(1:10);y=10-2*x+ra ndn( 10,1);y(10)=0； bls=regress(y,o nes(10,1) x) %線性擬合brob=robustfit(x,y) % 穩(wěn)健擬合 scatter(x,y) hold onplot(x,bls(1)+bls(2)*x, :) p

12、lot(x,brob(1)+brob(2)*x, r ) 結(jié)果:bls = 8.4452-1.4784brob =10.2934-2.0006分析：穩(wěn)健擬合(實線)對數(shù)據(jù)的擬合程度好些,忽略了異常值。最小二乘擬合(點線)則受到異常值的影響，向異常值偏移。 6向自定義函數(shù)擬合對于給定的數(shù)據(jù)，根據(jù)經(jīng)驗擬合為帶有待定常數(shù)的自定義函數(shù)。二、向自定義函數(shù)擬合所用函數(shù)：nli nfit()調(diào)用格式：beta,r,J=nlinfit(X,y, fun ,betao)說明：beta返回函數(shù)fun中的待定常數(shù)；r表示殘差；J表示雅可比矩陣。X,y為數(shù)據(jù)； fun 自定義函數(shù)；betaO待定常數(shù)初值。例7:在化工

13、生產(chǎn)中獲得的氯氣的級分y隨生產(chǎn)時間x下降，假定在x 8時，y與x之間有如下形式的非線性模型：y = a (0.49_a)e43現(xiàn)收集了44組數(shù)據(jù)，利用該數(shù)據(jù)通過擬合確定非線性模型中的待定常數(shù)。xyxyxy80.49160.43280.4180.49180.46280.40100.48180.45300.40100.47200.42300.40100.48200.42300.38100.47200.43320.41120.46200.41320.40120.46220.41340.40120.45220.40360.41120.43240.42360.36140.45240.40380.4014

14、0.43240.40380.40140.43260.41400.36160.44260.40420.39160.43260.41首先定義非線性函數(shù)的m 文件：fff6.mfun cti on yy=model(betaO,x)a=betaO(1);b=beta0(2);yy=a+(0.49-a)*exp(-b*(x-8);程序：x=8.00 8.00 10.00 10.00 10.00 10.00 12.00 12.00 12.00 14.00 14.00 14.00. 16.00 16.00 16.00 18.00 18.00 20.00 20.00 20.00 20.00 22.00 22.00 24.00. 24.00 24.00 26.00 26.00 26.00 28.00 28.00 30.00 30.00 30.00 32.00 32.00. 34.00 36.00 36.00 38.00 38.00 40.00 42.00;y=0.49 0.49 0.48 0.47 0.48 0.47 0.46 0.46 0.45 0.43 0.45 0.43 0.43 0.44 0.43. 0.43 0.46 0.42 0.42 0.43 0.41 0.41 0.40 0.42 0.40