最小二乘法曲線擬合原理和matlab實現(xiàn)

1、.wd曲線擬合curve-fitting：工程實踐中，用測量到的一些離散的數(shù)據(jù)求一個近似的函數(shù)來擬合這組數(shù)據(jù)，要求所得的擬合曲線能最好的反映數(shù)據(jù)的 基本趨勢即使最好地逼近，而不必滿足插值原那么。因此沒必要取=，只要使盡可能地小。原理：給定數(shù)據(jù)點。求近似曲線。并且使得近似曲線與的偏差最小。近似曲線在該點處的偏差，i=1,2,.,m。 常見的曲線擬合方法: 1.使偏差絕對值之和最小 2.使偏差絕對值最大的最小 3.使偏差平方和最小 最小二乘法： 按偏差平方和最小的原那么選取擬合曲線，并且采取二項式方程為擬合曲線的方法,稱為最小二乘法。推導過程： 1. 設擬合多項式為： 2. 各點到這條曲線的距離之

2、和，即偏差平方和如下： 3. 問題轉化為求待定系數(shù).對等式右邊求偏導數(shù)，因而我們得到了： . 4、 把這些等式化簡并表示成矩陣的形式，就可以得到下面的矩陣： 5. 將這個范德蒙得矩陣化簡后可得到: 6. 也就是說X*A=Y，那么A = (X'*X)-1*X'*Y，便得到了系數(shù)矩陣A，同時，我們也就得到了擬合曲線。MATLAB實現(xiàn)：MATLAB提供了polyfit函數(shù)命令進展最小二乘曲線擬合。調用格式：p=polyfit(x,y,n) p,s= polyfit(x,y,n) p,s,mu=polyfit(x,y,n)x,y為數(shù)據(jù)點，n為多項式階數(shù)，返回p為冪次從高到低的多項式系數(shù)

3、向量p。x必須是單調的。矩陣s包括R對x進展QR分解的三角元素、df(自由度、normr(殘差用于生成預測值的誤差估計。 p,s,mu=polyfit(x,y,n)在擬合過程中，首先對x進展數(shù)據(jù)標準化處理，以在擬合中消除量綱等影響，mu包含標準化處理過程中使用的x的均值和標準差。polyval( )為多項式曲線求值函數(shù)，調用格式： y=polyval(p,x) y,DELTA=polyval(p,x,s)y=polyval(p,x)為返回對應自變量x在給定系數(shù)P的多項式的值。y,DELTA=polyval(p,x,s) 使用polyfit函數(shù)的選項輸出s得出誤差估計Y DELTA。它假設pol

4、yfit函數(shù)數(shù)據(jù)輸入的誤差是獨立正態(tài)的，并且方差為常數(shù)。那么Y DELTA將至少包含50%的預測值。如下給定數(shù)據(jù)的擬合曲線：x=0.5,1.0,1.5,2.0,2.5,3.0，y=1.75,2.45,3.81,4.80,7.00,8.60。解：MATLAB程序如下：x=0.5,1.0,1.5,2.0,2.5,3.0;y=1.75,2.45,3.81,4.80,7.00,8.60;p=polyfit(x,y,2)x1=0.5:0.05:3.0;y1=polyval(p,x1);plot(x,y,'*r',x1,y1,'-b')運行結果如圖1計算結果為：p =0.5

5、614 0.8287 1.1560即所得多項式為y=0.5614x2+0.08287x+1.15560 圖1 最小二乘法曲線擬合例如比照檢驗擬合的有效性：例：在0,區(qū)間上對正弦函數(shù)進展擬合，然后在0,2區(qū)間畫出圖形，比較擬合區(qū)間和非擬合區(qū)間的圖形，考察擬合的有效性。在MATLAB中輸入如下代碼：clearx=0:0.1:pi;y=sin(x);p,mu=polyfit(x,y,9)x1=0:0.1:2*pi;y1=sin(x1);%實際曲線y2=polyval(p,x1);%根據(jù)由區(qū)間0到pi上進展擬合得到的多項式計算0到2pi上的函數(shù)值， %需要注意的是polyval返回的函數(shù)值在pi到2p

6、i上并沒有進展擬合plot(x1,y2,'k*',x1,y1,'k-')運行結果：p = 0.0000 0.0000 -0.0003 0.0002 0.0080 0.0002 -0.1668 0.0000 1.0000 0.0000mu = R: 10x10 double df: 22 normr: 1.6178e-07MATLAB的最優(yōu)化工具箱還提供了lsqcurvefit函數(shù)命令進展最小二乘曲線擬合(Solve nonlinear curve-fitting (data-fitting) problems in least-squares sense)。調用

7、格式：x = lsqcurvefit(fun,x0,xdata,ydata)x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub)x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub,options)x = lsqcurvefit(problem)x,resnorm = lsqcurvefit(.)x,resnorm,residual = lsqcurvefit(.)x,resnorm,residual,exitflag = lsqcurvefit(.)x,resnorm,residual,exitflag,output = lsqcur

8、vefit(.)x,resnorm,residual,exitflag,output,lambda = lsqcurvefit(.)x,resnorm,residual,exitflag,output,lambda,jacobian = x0為初始解向量;xdata，ydata為滿足關系ydata=F(x, xdata)的數(shù)據(jù);lb、ub為解向量的下界和上界 ，假設沒有指定界，那么lb= ，ub= ;options為指定的優(yōu)化參數(shù);fun為擬合函數(shù)，其定義方式為:x = lsqcurvefit(myfun,x0,xdata,ydata)，其中myfun已定義為 function F = myf

9、un(x,xdata)F = % 計算x處擬合函數(shù)值fun的用法與前面一樣;resnorm=sum (fun(x,xdata)-ydata).2)，即在x處殘差的平方和;residual=fun(x,xdata)-ydata，即在x處的殘差;exitflag為終止迭代的條件;output為輸出的優(yōu)化信息;lambda為解x處的Lagrange乘子;jacobian為解x處擬合函數(shù)fun的jacobian矩陣。例：lsqcurvefit()優(yōu)化程序Data = . 0.0000 5.8955 0.1000 3.5639 0.2000 2.5173 0.3000 1.9790 0.4000 1.8

10、990 0.5000 1.3938 0.6000 1.1359 0.7000 1.0096 0.8000 1.0343 0.9000 0.8435 1.0000 0.6856 1.1000 0.6100 1.2000 0.5392 1.3000 0.3946 1.4000 0.3903 1.5000 0.5474 1.6000 0.3459 1.7000 0.1370 1.8000 0.2211 1.9000 0.1704 2.0000 0.2636;t = Data(:,1);y = Data(:,2);% axis(0 2 -0.5 6) plot(t,y,'ro')tit

11、le('Data points')%We would like to fit the function y = c(1)*exp(-lam(1)*t) + c(2)*exp(-lam(2)*t) to the data%The lsqcurvefit function solves this type of problem easily.%To begin, define the parameters in terms of one variable x:%x(1) = c(1)%x(2) = lam(1)%x(3) = c(2)%x(4) = lam(2)%Then defi

12、ne the curve as a function of the parameters x and the data t:F = (x,xdata)x(1)*exp(-x(2)*xdata) + x(3)*exp(-x(4)*xdata);x0 = 1 1 1 0;x,resnorm,exitflag,output = lsqcurvefit(F,x0,t,y)hold onplot(t,F(x,t)hold offFsumsquares = (x)sum(F(x,t) - y).2);opts = optimset('LargeScale','off');x

13、unc,ressquared,eflag,outputu = .fminunc(Fsumsquares,x0,opts)fprintf('There were %d iterations using fminunc,' . ' and %d using lsqcurvefit.n', . outputu.iterations,output.iterations)fprintf('There were %d function evaluations using fminunc,' . ' and %d using lsqcurvefit.&

14、#39;, . outputu.funcCount,output.funcCount)type fitvectorx02 = 1 0;F2 = (x,t) fitvector(x,t,y);x2,resnorm2,exitflag2,output2 = lsqcurvefit(F2,x02,t,y)fprintf('There were %d function evaluations using the 2-d ' . 'formulation, and %d using the 4-d formulation.', . output2.funcCount,ou

15、tput.funcCount)x0bad = 5 1 1 0;xbad,resnormbad,exitflagbad,outputbad = . lsqcurvefit(F,x0bad,t,y)hold onplot(t,F(xbad,t),'g')legend('Data','Global fit','Bad local fit','Location','NE')hold offfprintf('The residual norm at the good ending point is %

16、f,' . ' and the residual norm at the bad ending point is %f.', . resnorm,resnormbad)displayEndOfDemoMessage(mfilename)擬合效果如下：直線的最小二乘擬合：ya+bx式中有兩個待定參數(shù)，a代表截距，b代表斜率。對于等精度測量所得到的N組數(shù)據(jù)xi，yi，i1，2，N，xi值被認為是準確的，所有的誤差只聯(lián)系著yi。下面利用最小二乘法把觀測數(shù)據(jù)擬合為直線。 用最小二乘法估計參數(shù)時，要求觀測值yi的偏差的加權平方和為最小。對于等精度觀測值的直線擬合來說，可使下式的值最?。荷鲜椒謩e對a、b求偏導得：整理后得到方程組：解上述方程組便可求得直線參數(shù)a和b的最正確估計值。1、 可看成是一階多項式擬合，跟前面曲線擬合方法一樣。2、 利用linefit()函數(shù)進展最小二乘的直線擬合 使用： clear x=0.5 1 1.5 2 2.5 3; y= 1.75 2.45 3.81 4.8 8 8.6; k,b=linefit(x,y) %得到斜率k和常數(shù)b y1=polyval(k,b,x); plot(x,y1,k-,x,y,k