matlab非线性参数拟合估计_很好的参考材料

使⽤nlinfit、fminsearch在matlab中实现基于最⼩⼆乘法的⾮线性参数拟合（整理⾃⽹上资源）

最⼩⼆乘法在曲线拟合中⽐较普遍。拟合的模型主要有1.直线型2.多项式型3.分数函数型4.指数函数型5.对数线性型6.⾼斯函数型......

⼀般对于LS问题，通常利⽤反斜杠运算“\\”、fminsearch或优化⼯具箱提供的极⼩化函数求解。在Matlab中，曲线拟合⼯具箱也提供了曲线拟合的图形界⾯操作。在命令提⽰符后键⼊：cftool，即可根据数据，选择适当的拟合模型。“\\”命令

1.假设要拟合的多项式是：y=a+b*x+c*x^2.⾸先建⽴设计矩阵X：X=[ones(size(x)) x x^2];执⾏：para=X\\y

para中包含了三个参数：para(1)=a;para(2)=b;para(3)=c;这种⽅法对于系数是线性的模型也适应。2.假设要拟合：y=a+b*exp(x)+cx*exp(x^2)设计矩阵X为

X=[ones(size(x)) exp(x) x.*exp(x.^2)];para=X\\y

3.多重回归(乘积回归)

设要拟合：y=a+b*x+c*t，其中x和t是预测变量，y是响应变量。设计矩阵为X=[ones(size(x)) x t] %注意x,t⼤⼩相等！para=X\\ypolyfit函数

polyfit函数不需要输⼊设计矩阵，在参数估计中，polyfit会根据输⼊的数据⽣成设计矩阵。1.假设要拟合的多项式是：y=a+b*x+c*x^2p=polyfit(x,y,2)

然后可以使⽤polyval在t处预测：y_hat=polyval(p,t)

polyfit函数可以给出置信区间。[p S]=polyfit(x,y,2) %S中包含了标准差

[y_fit,delta] = polyval(p,t,S) %按照拟合模型在t处预测在每个t处的95%CI为：(y_fit-1.96*delta, y_fit+1.96*delta)2.指数模型也适应

假设要拟合：y = a+b*exp(x)+c*exp(x.?2)p=polyfit(x,log(y),2)fminsearch函数

fminsearch是优化⼯具箱的极⼩化函数。LS问题的基本思想就是残差的平⽅和(⼀种范数，由此，LS产⽣了许多应⽤)最⼩，因此可以利⽤fminsearch函数进⾏曲线拟合。假设要拟合：y = a+b*exp(x)+c*exp(x.?2)⾸先建⽴函数，可以通过m⽂件或函数句柄建⽴：x=[......]';y=[......]';

f=@(p,x) p(1)+p(2)*exp(x)+p(3)*exp(x.?2) %注意向量化:p(1)=a;p(2)=b;p(3)=c;

%可以根据需要选择是否优化参数%opt=options()p0=ones(3,1);%初值

para=fminsearch(@(p) (y-f(p,x)).^2,p0) %可以输出Hessian矩阵res=y-f(para,x)%拟合残差曲线拟合⼯具箱

提供了很多拟合函数，对⼤样本场合⽐较有效！⾮线性拟合nlinfit函数clear all;

x1=[0.4292 0.4269 0.381 0.4015 0.4117 0.3017]';x2=[0.00014 0.00059 0.0126 0.0061 0.00425 0.0443]';x=[x1 x2];

y=[0.517 0.509 0.44 0.466 0.479 0.309]';f=@(p,x)

2.350176*p(1)*(1-1/p(2))*(1-(1-x(:,1).^(1/p(2))).^p(2)).^2.*(x(:,1).^ (-1/p(2))-1).^(-p(2)).*x(:,1).^(-1/p(2)-0.5).*x(:,2);p0=[8 0.5]';

opt=optimset('TolFun',1e-3,'TolX',1e-3);%[p R]=nlinfit(x,y,f,p0,opt)2.多项式型的⼀个例⼦

1900-2000年的总⼈⼝情况的曲线拟合

clear all;close all;

%cftool提供了可视化的曲线拟合！

t=[1900 1910 1920 1930 1940 1950 1960 1970 1980 1990 2000]';

y=[75.995 91.972 105.711 123.203 131.669 150.697 179.323 203.212 226.505 249.633 281.4220]';%t太⼤，以t的幂作为基函数会导致设计矩阵尺度太差，列变量⼏乎线性相依。变换为[-1 1]上s=(t-1950)/50;%plot(s,y,'ro');%回归线：y=a+bxmx=mean(s);my=mean(y);sx=std(s);sy=std(y);r=corr(s,y);b=r*sy/sx;a=my-b*mx;rline=a+b.*s;figure;

subplot(3,2,[1 2])plot(s,y,'ro',s,rline,'k');%title('多项式拟合');

set(gca,'XTick',s,'XTickLabel',sprintf('%d|',t));%hold on;n=4;

PreYear=[2010 2015 2020];%预测年份tPreYear=(PreYear-1950)/50;Y=zeros(length(t),n);res=zeros(size(Y));delta=zeros(size(Y));

PrePo=zeros(length(PreYear),n);Predelta=zeros(size(PrePo));for i=1:n

[p S(i)]=polyfit(s,y,i);

[Y(:,i) delta(:,i)]=polyval(p,s,S(i));%拟合的Y

[PrePo(:,i) Predelta(:,i)]=polyval(p,tPreYear,S(i));%预测res(:,i)=y-Y(:,i);%残差end

% plot(s,Y);%2009a⾃动添加不同颜⾊

% legend('data','regression line','1st poly','2nd poly','3rd poly','4th poly',2)% plot(tPreYear,PrePo,'>');% hold off

% plot(Y,res,'o');%残差图r=corr(s,Y).^2 %R^2%拟合误差估计CIYearAdd=[t;PreYear'];tYearAdd=[s;tPreYear'];

CFtit={'⼀阶拟合','⼆阶拟合','三阶拟合','四阶拟合'};for col=1:nsubplot(3,2,col+2);

plot(s,y,'ro',s,Y(:,col),'g-');%原始数据和拟合数据legend('Original','Fitted',2);hold on;

plot(s,Y(:,col)+2*delta(:,col),'r:');%95% CIplot(s,Y(:,col)-2*delta(:,col),'r:');plot(tPreYear,PrePo(:,col),'>');%预测值

plot(tPreYear,PrePo(:,col)+2*Predelta(:,col));%预测95% CIplot(tPreYear,PrePo(:,col)-2*Predelta(:,col));axis([-1.2 1.8 0 400]);

set(gca,'XTick',tYearAdd,'XTickLabel',sprintf('%d|',YearAdd)); title(CFtit{col});hold off;end

figure;%残差图for col=1:nsubplot(2,2,col);plot(Y(:,i),res(:,i),'o');end

⼀个⾮线性的应⽤例⼦(多元情况)要拟合y=a*x1^n1+b*x2^n2+c*x3^n3%注：只是作为应⽤，模型不⼀定正确%x2=x3

y=[1080.94 1083.03 1162.80 1155.61 1092.82 1099.26 1161.06 1258.05 1299.03 1298.30 1440.22 1641.30 1672.211612.73 1658.64 1752.42 1837.99 2099.29 2675.47 2786.33 2881.07]';

x1=[1 1.05 1.1 1.15 1.2 1.25 1.3 1.35 1.4 1.45 1.5 1.55 1.6 1.65 1.7 1.75 1.8 1.85 1.9 1.95 2]';

x2=[1 1.025 1.05 1.075 1.1 1.125 1.15 1.175 1.2 1.225 1.250 1.275 1.3 1.325 1.350 1.375 1.4 1.425 1.45 1.475 1.5]';x3=[1 1.025 1.05 1.075 1.1 1.125 1.15 1.175 1.2 1.225 1.250 1.275 1.3 1.325 1.350 1.375 1.4 1.425 1.45 1.475 1.5]';x=[x1 x2 x3];

f=@(p,x) p(1)*x(:,1).^p(2)+p(3)*x(:,2).^p(4)+p(5)*x(:,3).^p(6);p0=ones(6,1);

p=fminsearch(@(p)sum(y-f(p,x)).^2,p0)res=y-f(p,x);

res2=res.^2 %失败的模型Matlab ⾃定义函数⾃定义函数的途径：M⽂件函数(M file function)在线函数(Inline Function)匿名函数(Anonymous Function)1.M⽂件函数范例

function c=myadd(a,b)

%这⾥可以写函数的使⽤说明，前⾯以%开头%在⼯作区中，help myadd将显⽰此处的说明c=a+b;

%end %⾮必须的

第⼀⾏function告诉Matlab这是⼀个函数，a,b是输⼊，c是输出，myadd是函数名。以m⽂件定义的函数必须保存为函数名的形式，上例中，函数应保存为myadd.m。要使⽤myadd函数，该函数必须在Matlab的搜索路径中。调⽤⽅式：

在Matlab命令符后输⼊a=1;b=2;c=myadd(a,b)

2.在线函数(Inline Function)

通常作为参数传递给另外⼀个函数。⽐如fminsearch，lsqcurvefit等函数需要以函数作为参数。在线函数从字符串表达式创建函数，例如：f=inline('x.^2','x');

创建了函数f(x)=x^2。要计算f(3)，在⼯作区输⼊f(3)即可。f([2 3 4])计算在x=2 3 4时的值f=inline('x+y','x','y')

创建了⼆元函数f(x,y)=x+y，⼯作区输⼊f(2,3)计算2+3，等同于f(f,2,3)。3.匿名函数(Anonymous Function)

匿名函数使⽤函数句柄来表⽰匿名函数，定义形式为函数句柄=@(变量名) 函数表达式

例如：f=@(x) x.^2

定义了函数f(x)=x^2，f(2)计算在x=2处的值。

匿名函数可以调⽤Matlab函数，也可以使⽤⼯作区中存在的变量，例如a=2;

f=@(x) x.^2+a

f(2) %计算时引⽤了变量aa=0;

f(2) %仍然引⽤的是a=2

匿名函数也可以由Matlab的内置函数或M⽂件函数创建，例如f=@sin %f(x)=sin(x)f(pi/2) %sin(pi/2)

functions(f) %查看函数信息

利⽤单元数组可以创建多个函数的句柄，例如f={@sin @cos}f{1}(pi/2) %计算sin(pi/2)f{2}(pi) %计算cos(pi)

函数句柄的另⼀个重要特征是可以⽤来表⽰⼦函数、私有函数和嵌套函数。Matlab 7以后，建议以匿名函数取代在线函数

在创建匿名函数时，Matlab记录了关于函数的信息，当使⽤句柄调⽤该函数的时候，Matlab不再进⾏搜索，⽽是⽴即执⾏该函数，极⼤提⾼了效率。