Thursday, February 5, 2009

非线性最小二乘曲线拟和问题

非线性最小二乘曲线拟和问题

MATLAB优化工具箱中,使用lsqnonlin函数可以进行非线性最小二乘问题以及非线性曲线拟和问题。其调用格式如:
x = lsqnonlin(fun,x0)
x = lsqnonlin(fun,x0,lb,ub)
x = lsqnonlin(fun,x0,lb,ub,options)
[x,resnorm] = lsqnonlin(...)
[x,resnorm,residual] = lsqnonlin(...)
[x,resnorm,residual,exitflag] = lsqnonlin(...)
[x,resnorm,residual,exitflag,output] = lsqnonlin(...)
[x,resnorm,residual,exitflag,output,lambda] = lsqnonlin(...)
[x,resnorm,residual,exitflag,output,lambda,jacobian] = lsqnonlin(...)
使用非线性最小二乘求解函数lsqnonlin或者lsqcurvefit同样可以进行最小二乘的曲线拟合。Lsqcurvefit函数调用格式如下:
x = lsqcurvefit(fun,x0,xdata,ydata)
x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub)
x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub,options)
[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] = lsqcurvefit(...)

以下通过实例比较二者的区别:
首先需要绘制数据点图形,确定拟和数据点的大致曲线形式。
xdata =[0,4,6,10,12,16,22,26,30];
ydata =[-33.8980,-13.3859,-6.2068,1.9980,3.0236,-1.0788,-22.6163,-47.2308,-80.0499];
plot(xdata,ydata,'bo');
title('original data points');

Lsqnonlin函数
使用最小二乘法进行非线性曲线的拟和,首先需要建立目标函数的M文件:
function [feval]=curvefit_fun(x,xdata,ydata)
yval=x(1)*(xdata-x(2)).^2+x(3);
feval=yval-ydata;

在命令窗口中运行以下程序段:
x0=[0;0;0];
[x,resnorm,residual,exitflag,output] = lsqnonlin(@(x)curvefit_fun(x,xdata,ydata),x0)
x1=linspace(min(xdata),max(xdata),1000);
y1=x(1)*(x1-x(2)).^2+x(3);
hold on
plot(x1,y1,'r-');
legend('数据点','lsqnonlin函数拟和曲线')

lsqnonlin函数拟和曲线结果如图11.10所示。多项式曲线的系数向量及拟和曲线在散点处的留数及留数平方和数据如下:

x =
-0.2564
12.0000
3.0232
resnorm =
9.5716e-007
residual =
1.0e-003 *
-0.1805 -0.4245 -0.3716 -0.4164 -0.4139 -0.3594 -0.2535 0.0665 0.2526

运行结果如图:



Lsqcurvefit函数曲线拟和
首先,同样需要建立拟和曲线的目标函数,注意不同函数,目标函数的定义不完全相同:

function [feval]=curvefit_fun(x,xdata)
feval=x(1)*(xdata-x(2)).^2+x(3);

MATLAB的命令行窗口中输入以下程序段:

x0=[0;0;0];
[x,resnorm,residual] = lsqcurvefit(@(x,xdata)curvefit_fun(x,xdata),x0,xdata,ydata)
x1=linspace(min(xdata),max(xdata),1000);
y1=x(1)*(x1-x(2)).^2+x(3);
hold on
plot(x1,y1,'r-')
legend('数据点','lsqcurvefit拟和曲线')

运行结果如下:
x =
-0.2564
12.0000
3.0232
resnorm =
9.5716e-007
residual =
1.0e-003 *
-0.1805 -0.4245 -0.3716 -0.4164 -0.4139 -0.3594 -0.2535 0.0665 0.2526

运行结果如图:

More details could be found in my published book:
MATLAB编程基础与典型应用
北京:人民邮电出版社,2008
ISBN:978-7-115-17932-6/TP

Pls contact me with Email:lhd06@mails.tsinghua.edu.cn

No comments: