最小二乘法曲线拟合原理及matlab实现.docx

最小二乘法曲线拟合原理及matlab实现.docx

《最小二乘法曲线拟合原理及matlab实现.docx》由会员分享，可在线阅读，更多相关《最小二乘法曲线拟合原理及matlab实现.docx（11页珍藏版）》请在冰豆网上搜索。

最小二乘法曲线拟合原理及matlab实现

曲线拟合（curve-fitting）：

工程实践中，用测量到的一些离散的数据

求一个近似的函数

来拟合这组数据，要求所得的拟合曲线能最好的反映数据的基本趋势（即使

最好地逼近

，而不必满足插值原则。

因此没必要取

=

，只要使

尽可能地小）。

原理：

给定数据点

。

求近似曲线

。

并且使得近似曲线与

的偏差最小。

近似曲线在该点处的偏差

，i=1,2,...,m。

常见的曲线拟合方法:

1.使偏差绝对值之和最小

2.使偏差绝对值最大的最小

3.使偏差平方和最小

最小二乘法：

按偏差平方和最小的原则选取拟合曲线，并且采取二项式方程为拟合曲线的方法,称为最小二乘法。

推导过程：

1.设拟合多项式为：

2.各点到这条曲线的距离之和，即偏差平方和如下：

3.问题转化为求待定系数

...

对等式右边求

偏导数，因而我们得到了：

.......

4、把这些等式化简并表示成矩阵的形式，就可以得到下面的矩阵：

5.将这个范德蒙得矩阵化简后可得到:

6.也就是说X*A=Y，那么A=（X'*X）-1*X'*Y，便得到了系数矩阵A，同时，我们也就得到了拟合曲线。

MATLAB实现：

MATLAB提供了polyfit（）函数命令进行最小二乘曲线拟合。

调用格式：

p=polyfit（x,y,n）

[p,s]=polyfit（x,y,n）

[p,s,mu]=polyfit（x,y,n）

x,y为数据点，n为多项式阶数，返回p为幂次从高到低的多项式系数向量p。

x必须是单调的。

矩阵s包括R（对x进行QR分解的三角元素）、df（自由度）、normr（残差）用于生成预测值的误差估计。

[p,s,mu]=polyfit（x,y,n）在拟合过程中，首先对x进行数据标准化处理，以在拟合中消除量纲等影响，mu包含标准化处理过程中使用的x的均值和标准差。

polyval（）为多项式曲线求值函数，调用格式：

y=polyval（p,x）

[y,DELTA]=polyval（p,x,s）

y=polyval（p,x）为返回对应自变量x在给定系数P的多项式的值。

[y,DELTA]=polyval（p,x,s）使用polyfit函数的选项输出s得出误差估计YDELTA。

它假设polyfit函数数据输入的误差是独立正态的，并且方差为常数。

则YDELTA将至少包含50%的预测值。

如下给定数据的拟合曲线：

x=[,,,,,]，

y=[,,,,,]。

解：

MATLAB程序如下：

x=[,,,,,];

y=[,,,,,];

p=polyfit（x,y,2）

x1=:

:

;

y1=polyval（p,x1）;

plot（x,y,'*r',x1,y1,'-b'）

运行结果如图1

计算结果为：

p=

即所得多项式为y=^2++

图1最小二乘法曲线拟合示例

对比检验拟合的有效性：

例：

在[0,π]区间上对正弦函数进行拟合，然后在[0,2π]区间画出图形，比较拟合区间和非拟合区间的图形，考察拟合的有效性。

在MATLAB中输入如下代码：

clear

x=0:

:

pi;

y=sin（x）;

[p,mu]=polyfit（x,y,9）

x1=0:

:

2*pi;

y1=sin（x1）;%实际曲线

y2=polyval（p,x1）;%根据由区间0到pi上进行拟合得到的多项式计算0到2pi上的函数值，

%需要注意的是polyval（）返回的函数值在pi到2pi上并没有进行拟合

plot（x1,y2,'k*',x1,y1,'k-'）

运行结果：

p=

mu=

R:

[10x10double]

df:

22

normr:

MATLAB的最优化工具箱还提供了lsqcurvefit（）函数命令进行最小二乘曲线拟合（Solvenonlinearcurve-fitting（data-fitting）problemsinleast-squaressense）。

调用格式：

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]=lsqcurvefit（...）

[x,resnorm,residual,exitflag,output,lambda]=lsqcurvefit（...）

[x,resnorm,residual,exitflag,output,lambda,jacobian]=

x0为初始解向量;xdata，ydata为满足关系ydata=F（x,xdata）的数据;

lb、ub为解向量的下界和上界，若没有指定界，则lb=[]，ub=[];

options为指定的优化参数;

fun为拟合函数，其定义方式为:

x=lsqcurvefit（@myfun,x0,xdata,ydata），

其中myfun已定义为functionF=myfun（x,xdata）

F=…%计算x处拟合函数值fun的用法与前面相同;

resnorm=sum（（fun（x,xdata）-ydata）.^2），即在x处残差的平方和;

residual=fun（x,xdata）-ydata，即在x处的残差;

exitflag为终止迭代的条件;

output为输出的优化信息;

lambda为解x处的Lagrange乘子;

jacobian为解x处拟合函数fun的jacobian矩阵。

例：

lsqcurvefit（）优化程序

Data=...

[

];

t=Data（:

1）;

y=Data（:

2）;

%axis（[026]）

plot（t,y,'ro'）

title（'Datapoints'）

%Wewouldliketofitthefunctiony=c

（1）*exp（-lam

（1）*t）+c

（2）*exp（-lam

（2）*t）tothedata

%Thelsqcurvefitfunctionsolvesthistypeofproblemeasily.

%Tobegin,definetheparametersintermsofonevariablex:

%x

（1）=c

（1）

%x

（2）=lam

（1）

%x（3）=c

（2）

%x（4）=lam

（2）

%Thendefinethecurveasafunctionoftheparametersxandthedatat:

F=@（x,xdata）x

（1）*exp（-x

（2）*xdata）+x（3）*exp（-x（4）*xdata）;

x0=[1110];

[x,resnorm,~,exitflag,output]=lsqcurvefit（F,x0,t,y）

holdon

plot（t,F（x,t））

holdoff

Fsumsquares=@（x）sum（（F（x,t）-y）.^2）;

opts=optimset（'LargeScale','off'）;

[xunc,ressquared,eflag,outputu]=...

fminunc（Fsumsquares,x0,opts）

fprintf（['Therewere%diterationsusingfminunc,'...

'and%dusinglsqcurvefit.\n'],...

fprintf（['Therewere%dfunctionevaluationsusingfminunc,'...

'and%dusinglsqcurvefit.'],...

typefitvector

x02=[10];

F2=@（x,t）fitvector（x,t,y）;

[x2,resnorm2,~,exitflag2,output2]=lsqcurvefit（F2,x02,t,y）

fprintf（['Therewere%dfunctionevaluationsusingthe2-d'...

'formulation,and%dusingthe4-dformulation.'],...

x0bad=[5110];

[xbad,resnormbad,~,exitflagbad,outputbad]=...

lsqcurvefit（F,x0bad,t,y）

holdon

plot（t,F（xbad,t）,'g'）

legend（'Data','Globalfit','Badlocalfit','Location','NE'）

holdoff

fprintf（['Theresidualnormatthegoodendingpointis%f,'...

'andtheresidualnormatthebadendingpointis%f.'],...

resnorm,resnormbad）

displayEndOfDemoMessage（mfilename）

拟合效果如下：

直线的最小二乘拟合：

y＝a+bx

式中有两个待定参数，a代表截距，b代表斜率。

对于等精度测量所得到的N组数据（xi，yi），i＝1，2……，N，xi值被认为是准确的，所有的误差只联系着yi。

下面利用最小二乘法把观测数据拟合为直线。

用最小二乘法估计参数时，要求观测值yi的偏差的加权平方和为最小。

对于等精度观测值的直线拟合来说，可使下式的值最小：

上式分别对a、b求偏导得：

整理后得到方程组：

解上述方程组便可求得直线参数a和b的最佳估计值。

1、可看成是一阶多项式拟合，跟前面曲线拟合方法一样。

2、利用linefit（）函数进行最小二乘的直线拟合

使用：

clear

x=[123];

y=[8];

[k,b]=linefit（x,y）%得到斜率k和常数b

y1=polyval（[k,b],x）;

plot（x,y1,’k-’,x,y,’k*’）

MATLAB一元到多元线性回归方程的计算和检验_XX文库

研究生数理统计）多元线性回归及显著性检验Matlab程序（完美版）_XX文库