精华内容
下载资源
问答
  • 此程序可以用于曲线拟合。包括线性拟合以及非线性拟合
  • Labview非线性拟合

    2017-09-16 16:33:27
    Labview用LM算法非线性拟合,可以实现一组数据的非线性拟合,而不用知道具体的方程系数,我都觉得很烦了,能两句话说完,要那么多废话干嘛
  • 非线性拟合工具

    热门讨论 2011-12-05 16:56:12
    当前非线性拟合和多元拟合的工具较少,这是针对常用的拟合算法,开发的一款数据拟合为主的软件。包括线性拟合的各种算法,非线性拟合的各种算法,以及多元拟合的各种算法。其中提供了很多非线性方程的模型,以满足...
  • 常用的高次非线性拟合matlab源代码,内含一次二次和非线性拟合的教程,还有部分神经网络的源代码
  • 传统的BP神经网络分类和拟合精度不高,很...本程序运用遗传算法初始化BP神经网络的参数,使网络的非线性拟合和分类精度更高。对于想要进行非线性拟合和分类却无法建立数学模型的小伙伴,通过该方法也是一个较好的办法。
  • 使用C#的Mathnet类库实现最小二乘法非线性拟合 作者:linbor tinka
  • 使用C#的Mathnet类库实现非线性拟合 作者:linbor tinka
  • 利用matlab实现非线性拟合0 前言1 线性拟合1.1 多项式拟合1.2 线性拟合2 一维非线性拟合2.1 简单的非线性拟合2.2 matlab中Curve Fitting App2.3 matlab中非线性拟合的实现2.3.1 fit()函数2.3.2 nlinfit()函数2.3.3 ...

    2021.06 更新。补充了非线性拟合中,关于最小二乘定义的问题,放在了最后一章。
    2021.12更新。应评论区建议,补充了3.3章绘图用的代码。

    本文首发于 matlab爱好者 微信公众号,欢迎关注。

    惯例声明:本人没有相关的工程应用经验,只是纯粹对相关算法感兴趣才写此博客。所以如果有错误,欢迎在评论区指正,不胜感激。本文主要关注于算法的实现,对于实际应用等问题本人没有任何经验,所以也不再涉及。

    本人在学习matlab匿名函数时,作为练习有感而发,写下下面内容,非专业,难免有误。
    在这里插入图片描述

    0 前言

    一般而言,通过已有的数据点去推导其它数据点,常见的方法有插值和拟合。插值适用性较广,尤其是线性插值或样条插值已被广泛的应用。但是通过已知的函数去拟合数据,是连接理论与实验重要的桥梁,这一点是插值无法替代的。

    日常学习工作中,经常会遇到下面这种问题:想要用某个具体函数去拟合自己的数据,明明知道这个函数的具体形式,却不知道其中的参数怎么选取。本文就简单介绍一下matlab环境下,如何进行非线性拟合。

    1 线性拟合

    1.1 多项式拟合

    多项式拟合就是利用下面形式的方程去拟合数据:
    y = a + b x + c x 2 + d x 3 + . . . y=a +bx+cx^2+dx^3+... y=a+bx+cx2+dx3+...
    matlab中可以用polyfit()函数进行多项式拟合。下面举一个小例子:

    对于已有的数据点,我们采用4阶多项式拟合。其中已知函数的的表达式为

    y=0.03 x^4 - 0.5 x^3 + 2 x^2 - 4
    

    在此基础上添加了一些噪声点。拟合曲线依然采用4阶进行拟合,结果如下。

    在这里插入图片描述
    可以看到拟合曲线与理论曲线基本一致,说明这种方法能够较好的拟合出原始数据的趋势。源代码如下:

    x=0:0.5:10;
    y=0.03*x.^4-0.5*x.^3+2.0*x.^2-0*x-4+6*(rand(size(x))-0.5);
    p=polyfit(x,y,4);
    x2=0:0.05:10;
    y2=polyval(p,x2);
    figure();
    subplot(1,2,1)
    hold on
    plot(x,y,'linewidth',1.5,'MarkerSize',15,'Marker','.','color','k')
    plot(x,0.03*x.^4-0.5*x.^3+2.0*x.^2-0*x-4,'linewidth',1,'color','b')
    hold off
    legend('原始数据点','理论曲线','Location','southoutside','Orientation','horizontal')
    legend('boxoff')
    box on
    subplot(1,2,2)
    hold on
    plot(x2,y2,'-','linewidth',1.5,'color','r')
    plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
    hold off
    box on
    legend('拟合曲线','数据点','Location','southoutside','Orientation','horizontal')
    legend('boxoff')
    

    对于多项式拟合,不是阶数越多越好。比如还是这个上面这个例子,阶数越多,曲线越能够穿过每一个点,但是曲线的形状与理论曲线偏离越大。所以选择最适合的才是最好的。
    在这里插入图片描述

    1.2 线性拟合

    线性拟合就是能够把拟合函数写成下面这种形式的:
    y = p 1 f 1 ( x ) + p 2 f 2 ( x ) + p 3 f 3 ( x ) + . . . y=p_1 f_1(x) +p_2 f_2(x)+p_3 f_{3}(x)+... y=p1f1(x)+p2f2(x)+p3f3(x)+...

    其中f(x)是关于x的函数,其表达式是已知的。p是常数系数,这个是未知的。

    对于这种形式的拟合,matlab内部有一个及其强悍的函数,可以自动输出p的解,并且满足最小二乘。这个函数就是\。没错,就是斜杠。这个符号通常用于求解方程AX=B的情况,我们用X=A\B可以求出未知数X。我们利用当A行和列不等时,输出X的最小二乘这个特性,就可以求出相应的最佳拟合。

    还是举个例子

    在这里插入图片描述
    虽然看上去很非线性,但是x和y都是已知的,a、b、c是未知的,而且是线性的,所以可以被非常简单的拟合出来。假设a=2.5,b=0.5,c=-1,加入随机扰动。拟合的最终效果为:

    在这里插入图片描述
    最终得到的拟合参数为:a=2.47,b=0.47,c=-0.66。考虑到随机噪声的影响,与原始数据相差不大,源代码如下:

    %线性拟合
    x=0:0.5:10;
    a=2.5;
    b=0.5;
    c=-1;
    %原函数
    y=a*sin(0.2*x.^2+x)+b*sqrt(x+1)+c+0.5*rand(size(x));
    
    %拟合部分
    yn1=sin(0.2*x.^2+x);
    yn2=sqrt(x+1);
    yn3=ones(size(x));
    X=[yn1;yn2;yn3];
    %拟合,得到系数
    Pn=X'\y';
    yn=Pn(1)*yn1+Pn(2)*yn2+Pn(3)*yn3;
    
    %绘图
    figure()
    subplot(1,2,1)
    plot(x,y,'linewidth',1.5,'MarkerSize',15,'Marker','.','color','k')
    legend('原始数据点','Location','southoutside','Orientation','horizontal')
    legend('boxoff')
    subplot(1,2,2)
    hold on
    x2=0:0.01:10;
    plot(x2,Pn(1)*sin(0.2*x2.^2+x2)+Pn(2)*sqrt(x2+1)+Pn(3),'-','linewidth',1.5,'color','r')
    plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
    hold off
    box on
    legend('拟合曲线','数据点','Location','southoutside','Orientation','horizontal')
    legend('boxoff')
    

    事实上,其实常用的拟合方式中,有很多都是线性拟合,比如多项式拟合,傅里叶拟合等。对于多项式拟合,只需要把拟合的部分替换成x,x.^2,x.^3这种形式。对于傅里叶级数,只需要把拟合的部分替换为sin(x),sin(2*x),sin(3*x)这种形式就行。

    下面展示一个利用线性拟合,进行不同频率的三角波级数拟合正弦函数的例子:
    在这里插入图片描述

    clear;
    clc;
    close all
    t=0:0.001:2*pi;%原函数
    YS=sin(t);%基函数
    N=21;
    Yo=[];
    for k=1:N
        Yn=sawtooth(k*(t+pi/2),0.5);
        Yo=[Yo,Yn'];
    end
    YS=YS';%拟合
    a = Yo\YS;
    %绘图
    figure()
    for k=1:N
        clf
        plot(t,Yo(:,1:k)*a(1:k,:),t,YS,'LineWidth',1)
        ylim([-1.3,1.3])
        xlim([0,6.3])
        pause(0.1)
    end
    

    2 一维非线性拟合

    2.1 简单的非线性拟合

    前面介绍了线性的拟合方法。如果一个非线性方程,可以化为上面线性方程中公式中给出的样子,那么我们也可以套用线性的方法去求解。

    比如下面的方程:
    y = a ∗ exp ⁡ ( − b x ) y=a*\exp(-bx) y=aexp(bx)
    经过取对数变换,可以变为下面等效的线性形式:
    lg ⁡ ( y ) = lg ⁡ ( a ) + b ∗ ( − x ) \lg(y)=\lg(a)+b*(-x) lg(y)=lg(a)+b(x)
    这样求出来之后,再反变换回去,就可以得到原方程的系数了。
    在这里插入图片描述
    代码如下:

    clear
    clc
    close all
    
    %方法1
    x=0:0.5:10;
    a=2.5;
    b=0.5;
    %原函数
    y=a*exp(-b*x);
    y=y+0.5*y.*rand(size(y));%加噪声
    %变形函数
    %Lg_y=Lg_a+b*(-x)
    Lg_y=log(y);
    %拟合部分
    yn1=ones(size(x));
    yn2=-x;
    X=[yn1;yn2];
    %拟合,得到系数
    Pn=X'\Lg_y';
    %反变换
    a_fit=exp(Pn(1));
    b_fit=Pn(2);
    %绘图
    figure()
    hold on
    x2=0:0.01:10;
    plot(x2,a_fit*exp(-b_fit*x2),'-','linewidth',1.5,'color','r')
    plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
    hold off
    

    2.2 matlab中Curve Fitting App

    matlab自带了一个Curve Fitting App,它在matlab集成的App里面。
    在这里插入图片描述
    界面左边输入数据,右侧选择方法。
    在这里插入图片描述
    界面里常用的拟合方式都有,而且直接展示拟合效果,非常方便。非常适合鼠标直接拖拖拽拽点点点的操作方式。

    2.3 matlab中非线性拟合的实现

    2.3.1 fit()函数

    matlab中,fit()函数是一个比较通用的用于函数拟合的函数。它的优点就是非常的全面,可以用于各种种类的拟合。上面的App里,很多拟合种类都是间接调用了fit函数来实现的拟合。

    对于非线性拟合,可以使用fit()函数中的Nonlinear Least Squares方法。其大概原理为,首先确定一个初始的点,计算该点的雅可比矩阵J,来近似线性化的判断该点周围的趋势,并将这个点向更小的方向移动。

    因此,这个方法的一个缺点在于,对于初始点的选取非常敏感,最终结果只能在初始点附近的局部最小值点上,而不能保证全局最小值。

    2.3.2 nlinfit()函数

    相比于前面的fit()函数,nlinfit()函数是matlab专门的非线性拟合函数。对于非稳健估计,采用的是Levenberg-Marquardt(LM)方法,也叫阻尼最小二乘法。对于稳健估计,采用的是Iteratively Reweighted Least Squares方法,也就是在Least Squares基础上,对每一个拟合点的权重进行调整的一种方法。这两者方法也都是基于雅克比矩阵的方法。

    2.3.3 lsqnonlin()函数和lsqcurvefit()函数

    lsqnonlin()也是matlab中自带的一个非线性拟合函数。它给出了两种计算非线性拟合的方法,一种是比较经典之前介绍过的LM方法,一种是信赖域方法。信赖域法(trust region reflective)是通过Hessian 矩阵,逐步试探邻域内的最小化,来求解问题的。相比较之前的那些雅克比相关的方法,信赖域法会占用更多内存和速度,所以适用于中小规模的矩阵。

    lsqcurvefit()函数和lsqnonlin()内容上相似,只是引用格式上有所不同。

    2.3.4 fsolve()函数

    这也是一个求解非线性方程的函数,可以求解方程组或者矩阵形式,功能非常强大。默认的算法为trust-region-dogleg,俗称狗腿法,属于信赖域法。这里用到的功能比较基础,所以也不过多介绍了。

    2.3.5 粒子群算法

    说了那么多,发现逐渐从如何非线性拟合,陷入到了最优化的深坑里。而且前面的那么多方法,很多都解决不了陷入到局部最优解的问题。实际上,这种问题如果进入了最优化领域,很多智能算法也可以被考虑进来。所以我也把粒子群PSO算法加入到了里面,尝试将结果收敛到全局最优解。

    2.3.6 不同算法的对比效果

    第一个例子是 y=a.*exp(-b\*(x-c).^2)+d,一个简单的高斯函数形式的非线性方程,其参数给定为:

    abcd
    3.82.14.4-1.3

    在已知函数形式,求解这四个参数条件下,6种不同的函数非拟合效果如下:
    在这里插入图片描述
    可以看到,这几种方法都能够较好的拟合出想要的结果。

    第二个例子是一个指数增长的正弦函数,在很多线性系统中都可以测量到这种信号。函数的形式为:

    y=a*x+b*sin(c*x).*exp(d*x)+e 。其给定的参数为:

    abcde
    -0.32.14.40.31.7

    这个函数的拟合具有一定难度,拟合过程中会遇到非常多的局部解。

    在这里插入图片描述
    由于初始点比较随机,所以最终结果每次都会不一样。
    其中前面的几种方法对于初始值的敏感度比较高,如果初始值选的比较接近原始解,也是可以得到较好的结果。其中nlinfit函数经常会报错,容错率较低。而PSO算法经常能够收敛到最优解(虽然不是每次都可以,偶尔也会陷入局部解)。

    代码如下:

    clear
    clc
    %函数大比拼
    close all
    
    %初始设置
    x = 0:0.1:10;
    a = -0.3;
    b = 2.1;
    c = 4.4;
    d = 0.3;
    e = 1.7;
    
    y = a*x+b*sin(c*x).*exp(d*x)+e;
    noise = 0.05*abs(y-1).*randn(size(x));
    y = y+noise;%加噪声函数
    
    figure();%plot(x,y)
    y_lim = [-40,40];
    
    %% 1 fit()函数 Least Squares
    ft = fittype( 'a*x+b*sin(c*x).*exp(d*x)+e', 'independent', 'x', 'dependent', 'y' );
    OP1 = fitoptions( 'Method', 'NonlinearLeastSquares' );
    OP1.StartPoint = 5*rand(1,5);%初始值,越靠近真实值越好
    OP1.Lower = [-2,0,2,0,0];%参数的最小边界
    OP1.Upper = [1,3,5,3,3];%参数的最大边界
    %拟合
    fitobject = fit(x',y',ft,OP1); 
    a2=fitobject.a;
    b2=fitobject.b;
    c2=fitobject.c;
    d2=fitobject.d;
    e2=fitobject.e;
    %展示结果
    y1 = a2*x+b2*sin(c2*x).*exp(d2*x)+e2;
    subplot(3,2,1)
    hold on
    plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
    plot(x,y1,'-','linewidth',1.5,'color','r')
    hold off
    box on
    ylim(y_lim)
    title('fit函数')
    
    %% 2 nlinfit()函数 Levenberg-Marquardt %容易报错
    modelfun = @(p,x)( p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5) );
    
    OP2 = statset('nlinfit');
    %opts.RobustWgtFun = 'bisquare';
    p0 = 5*rand(1,5);
    p = nlinfit(x,y,modelfun,p0,OP2);
    %拟合
    y2 = p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5);
    subplot(3,2,2)
    hold on
    plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
    plot(x,y2,'-','linewidth',1.5,'color','r')
    hold off
    box on
    ylim(y_lim)
    title('nlinfit函数')
    
    %% 3 lsqnonlin()函数 trust-region-reflective
    modelfun = @(p)(p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5) -y) ;%这里和nlinfit()函数定义不一样
    p0 = 5*rand(1,5);
    OP3 = optimoptions(@lsqnonlin,'Algorithm','trust-region-reflective','MaxFunctionEvaluations',1e4,'MaxIterations',1e3);
    [p,~] = lsqnonlin(modelfun,p0,[-2,0,2,0,0],[1,3,5,3,3],OP3);
    y3 = p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5);
    subplot(3,2,3)
    hold on
    plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
    plot(x,y3,'-','linewidth',1.5,'color','r')
    hold off
    box on
    ylim(y_lim)
    title('lsqnonlin函数')
    
    %% 4 lsqcurvefit()函数 trust-region-reflective
    modelfun = @(p,x)(p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5)) ;%这里和其它函数的定义也不一样
    p0 = 5*rand(1,5);
    OP4 = optimoptions('lsqcurvefit','Algorithm','trust-region-reflective','MaxFunctionEvaluations',1e4,'MaxIterations',1e3);
    [p,~] = lsqcurvefit(modelfun,p0,x,y,[-2,0,2,0,0],[1,3,5,3,3],OP4);
    y4 = p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5);
    subplot(3,2,4)
    hold on
    plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
    plot(x,y4,'-','linewidth',1.5,'color','r')
    hold off
    box on
    ylim(y_lim)
    title('lsqcurvefit函数')
    
    %% 5 fsolve()函数 %默认算法trust-region-dogleg
    modelfun = @(p)(p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5) -y);
    p0 = 5*rand(1,5);
    OP5 = optimoptions('fsolve','Display','off');
    p = fsolve(modelfun,p0,OP5);
    y5 = p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5);
    subplot(3,2,5)
    hold on
    plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
    plot(x,y5,'-','linewidth',1.5,'color','r')
    hold off
    box on
    ylim(y_lim)
    title('fsolve函数')
    
    %% 6 粒子群PSO算法
    fun = @(p) ( norm(p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5) -y) );%这里需要计算误差的平方和
    OP6 = optimoptions('particleswarm','InertiaRange',[0.4,1.2],'SwarmSize',100);
    [p,~,~,~]  = particleswarm(fun,5,[-5,-5,-5,-5],[5,5,5,5],OP6);%区间可以稍微放大一些,不怕
    y6 = p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5);
    subplot(3,2,6)
    hold on
    plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
    plot(x,y6,'-','linewidth',1.5,'color','r')
    hold off
    box on
    ylim(y_lim)
    title('PSO算法')
    

    下图展示了PSO求解过程中逐渐收敛到全局最优解的过程。
    在这里插入图片描述

    3 高维非线性方程组拟合

    3.1 一般形式高维方程或方程组的拟合

    之前的文章中的数据具有一 一对应的特点,所以严格来讲并不是普遍的二维拟合。对于一些复杂的二维函数,比如椭圆,可能原本的拟合就会失效。

    对于这种一般性质的方程或方程组,将原本输入方程整理为f(x,y,z,…)=0的形式。衡量拟合程度的优化函数,就直接取函数f(xi,yi,zi,…)的值即可。

    下面演示最终的两个例子:

    第一个是三维直线,采用两平面式描述。

    Ax+By+Cz-1=0
    Dx+Ey+Fz-1=0
    

    总共2个方程,维度为3维,第一个方程有3个参数,第二个方程也有3个参数。离散点已知的条件下,三维直线的两平面表达式不唯一。

    最终拟合效果如下:

    在这里插入图片描述

    第二个是二维椭圆,方程为:

    x^2+Axy+By^2+Cx+Dy+E=0
    

    总共1个方程,维度为2维。方程共有5个参数。

    最终拟合效果如下:
    在这里插入图片描述

    clear
    clc
    close all
    %% 演示1
    %1 导入数据(这里用的是人工生成的数据)
    %三维直线拟合,函数表示
    %1.0*x+1.9*y+3.0*z=1;
    %1.2*x-0.4*y-1.7*z=1;
    x=0:0.04:1;
    %求解出y和z
    % [ 1.0, 3.0    [ y        [1.0
    %  -0.4,-1.7] *   z] = 1 -  1.2] x
    y=zeros(size(x));z=y;
    for k=1:length(x)
        R=[1.9,3.0;-0.4,-1.7]\[1-1.0*x(k);1-1.2*x(k)];
        [y(k),z(k)]=Value(R);
    end
    rand_n=randperm(length(x));
    %生成随机的初始输入数据
    x1=x(rand_n)+0.05*randn(size(x));
    y1=y(rand_n)+0.05*randn(size(x));
    z1=z(rand_n)+0.05*randn(size(x));
    %2 整理格式,将数据和要拟合的函数进行格式整理
    %准备数据
    XX={x1,y1,z1};
    %准备函数
    F1=@(p1,XX1) p1(1)*XX1{1}+p1(2)*XX1{2}+p1(3)*XX1{3}-1;
    F2=@(p2,XX2) p2(1)*XX2{1}+p2(2)*XX2{2}+p2(3)*XX2{3}-1;
    FF={F1,F2};
    %3 生成最终优化函数,带入到优化方程中求解
    fun=@(p) Dis(p,{3,3},XX,FF);%p参数,{3,3}第一个方程3个参数,第二个方程3个参数。XX离散点。FF函数表达式。
    OP=optimoptions('particleswarm','FunctionTolerance',0);
    [p,fval,~,~]=particleswarm(fun,6,[-5,-5,-5,-5,-5,-5],[5,5,5,5,5,5],OP);
    %4 对比拟合效果
    figure()
    x2=0:0.01:1;
    y2=zeros(size(x2));
    z2=y2;
    for k=1:length(x2)
        R=[p(2),p(3);p(5),p(6)]\[1-p(1)*x2(k);1-p(4)*x2(k)];
        [y2(k),z2(k)]=Value(R);
    end
    %系数可能不同。因为直线的两平面表示不唯一
    hold on
    plot3(x2,y2,z2)
    plot3(x1,y1,z1,'*');
    hold off
    view(3)
    
    %% 演示2
    %1 导入数据(这里用的是人工生成的数据)
    %二维椭圆拟合
    th=0:0.15:2*pi;
    a=3.2;%椭圆轴1
    b=4.8;%椭圆轴2
    x0=-1.9;
    y0=-4.1;
    beta=1.1;%椭圆旋转角度
    %绘制椭圆
    x=a*cos(th);
    y=b*sin(th);
    %旋转beta角度
    x_r=x*cos(beta)-y*sin(beta);
    y_r=x*sin(beta)+y*cos(beta);
    rand_n=randperm(length(x));
    %生成随机的初始输入数据
    x1=x_r(rand_n)+0.1*randn(size(x));
    y1=y_r(rand_n)+0.1*randn(size(x));
    %2 整理格式,将数据和要拟合的函数进行格式整理
    %准备数据
    XX={x1,y1};
    %准备函数
    F1=@(p1,XX1) XX1{1}.*XX1{1}+p1(1)*XX1{1}.*XX1{2}+p1(2)*XX1{2}.*XX1{2}+p1(3)*XX1{1}+p1(4)*XX1{2}+p1(5);
    FF={F1};
    %3 生成最终优化函数,带入到优化方程中求解
    fun=@(p) Dis(p,{5},XX,FF);
    OP=optimoptions('particleswarm','FunctionTolerance',0);
    [p,fval,~,~]=particleswarm(fun,5,[-20,-20,-20,-20,-20],[20,20,20,20,20],OP);
    %4 对比拟合效果
    figure()
    hold on
    plot(x1,y1,'*')
    Fun_Plot=@(xp,yp) xp.*xp+p(1)*xp.*yp+p(2)*yp.*yp+p(3)*xp+p(4)*yp+p(5);
    fimplicit(Fun_Plot,[-6 6 -6 6],'LineWidth',1)
    hold off
    
    %% 用到的函数
    function varargout=Value(f)
    %多元素赋值,例子:
    %[a,b,c]=Value([1,2,3]);%多变量赋值
    %[xy(1),xy(2),t]=Value([2,5,3]);%复合赋值
    %[b,a]=Value([a,b]);%元素交换
    %来源:hyhhhyh21,matlab爱好者
    N=numel(f);
    varargout=cell(1,N);
    if nargout~=N
        error('输入输出数量不一致')
    end
    for k=1:N
        varargout{k}=f(k);
    end
    end
    
    function Sum_N=Dis(p,p_num,XX,FF)
    %用函数值评价拟合程度
    %输入:要拟合的参数p
    %输入:p_num cell格式,每个方程的参数数量
    %输入:XX 数据,以cell形式储存
    %输入:FF 拟合函数,以cell形式储存
    N_F=numel(FF);%要联立的方程数量
    L=length(XX{1});%离散点的数量
    N_L=numel(XX);%离散点维度
    Sum_N=0;
    XXm=cell(1,N_L);
    %直接计算函数的值
    p_n=1;%参数的索引
    for k=1:N_F
        %对每一个方程进行计算
        FF_k=FF{k};%方程
        F_p=p_num{k};%该方程用到的参数数量
        for m=1:L
            for n=1:N_L
                XXm{n}=XX{n}(m);
            end
            Distance=FF_k(p(p_n:(p_n+F_p-1)),XXm);
            Sum_N=Sum_N+Distance^2;
        end
        p_n=p_n+F_p;%更新函数索引
    end
    end
    

    3.2 一般形式参数方程的拟合

    高维参数方程的拟合比较困难。因为原本方程中x、y、z的坐标点都是已知的。但是参数方程中,x、y、z的坐标点已知,但是与参数u、v往往未知。所以相当于原本的方程中引入了额外的未知信息。

    但是基本思路和普通方程是一样的。离散点距离假想方程的距离,需要用该点距方程上每个点的距离中的最小值,来近似判断。

    依然是展示两个例子。

    第一个是计算三维李萨如图形。参数方程为:

    x=sin(A*u)
    y=cos(B*u)
    z=sin(C*u)
    

    方程为三维参数方程,只有一个参数u。第一个方程有1个未知量A,第二个方程有1个未知量B,第三方程有1个未知量C。

    最终拟合效果如下。由于李萨如图形中,只要频率的比例固定,图案就会固定。所以最终ABC的值不唯一,但是它们的比例肯定唯一。

    在这里插入图片描述
    第二个例子是一个三维旋转曲面。参数方程为:

    x= A*u.*sin(v+B)
    y=-C*u.*cos(v+D)
    z=v
    

    方程为三维参数方程,有2个参数u、v。第一个方程有2个未知量AB,第二个方程有2个未知量CD,第三方程有0个未知量。

    最终拟合效果如下:

    在这里插入图片描述
    这两个例子的演示代码如下。这里参数方程的Dis_P()由于频繁的计算点与点之间的距离,所以运算速度比第一章单纯函数的Dis()较慢。

    clear
    clc
    close all
    %% 演示1
    %三维李萨如图形拟合
    %1 导入数据(这里用的是人工生成的数据)
    t=0:0.01:4*pi;
    x=sin(4*t);
    y=cos(5*t);
    z=sin(6*t);
    
    rand_n=randperm(length(t));
    x1=x(rand_n)+0.02*randn(size(t));
    y1=y(rand_n)+0.02*randn(size(t));
    z1=z(rand_n)+0.02*randn(size(t));
    %2 整理格式,将数据和要拟合的函数进行格式整理
    %输入参数方程
    XX={x1,y1,z1};%离散点输入
    F1=@(p1,uu1) sin(p1(1)*uu1{1});
    F2=@(p2,uu1) cos(p2(1)*uu1{1});
    F3=@(p3,uu1) sin(p3(1)*uu1{1});
    FF={F1,F2,F3};%方程输入
    u1=0:0.05:13;%设置参数的最大最小范围以及精度,能达到绘图精度即可
    uu={u1};%参数输入
    %3 生成最终优化函数,带入到优化方程中求解
    fun=@(p) Dis_P(p,{1,1,1},uu,XX,FF);%使得DisP函数输出的Sum_N最小
    %p参数,{1,1,1}代表有3个方程,每个方程的代求参数p个数为1。uu为参数输入。XX为离散点输入。FF为方程输入。
    OP=optimoptions('particleswarm','FunctionTolerance',0,'InertiaRange',[0.3,1.2],'MaxStallIterations',200);
    [pp,fval,~,~]=particleswarm(fun,3,[1,1,1],[10,10,10],OP);
    %4 对比拟合效果
    figure()
    hold on
    tt=0:0.01:4*pi;
    %pp=pp/pp(1)*4;%这里不一定必须是4,5,6,只需要比例为4:5:6就行
    [a2,b2,c2]=Value(pp);
    x2=sin(a2*tt);
    y2=cos(b2*tt);
    z2=sin(c2*tt);
    
    plot3(x2,y2,z2);
    plot3(x1,y1,z1,'*');
    hold off
    view(3)
    
    %% 演示2
    %三维螺旋面拟合
    %1 导入数据(这里用的是人工生成的数据)
    F1=@(p1,uu1) p1(1).*uu1{1}.*sin(uu1{2}+p1(2));
    F2=@(p2,uu1) -p2(1).*uu1{1}.*cos(uu1{2}+p2(2));
    F3=@(p3,uu1) uu1{2};
    
    u_rand=rand_AB([-4,4],100,1);
    v_rand=rand_AB([-5,5],100,1);
    x=F1([2,0.1],{u_rand,v_rand});
    y=F2([1,0.1],{u_rand,v_rand});
    z=F3([],{u_rand,v_rand});
    
    rand_n=randperm(length(x));
    x1=x(rand_n)+0.01*randn(size(x));
    y1=y(rand_n)+0.01*randn(size(x));
    z1=z(rand_n)+0.01*randn(size(x));
    
    %2 整理格式,将数据和要拟合的函数进行格式整理
    %输入参数方程
    XX={x1,y1,z1};%离散点输入
    FF={F1,F2,F3};%方程输入
    u1=-4:0.1:4;%设置参数的最大最小范围以及精度,能达到绘图精度即可
    v1=-5:0.1:5;%设置参数的最大最小范围以及精度,能达到绘图精度即可
    uu={u1,v1};%参数输入
    %3 生成最终优化函数,带入到优化方程中求解
    fun=@(p) Dis_P(p,{2,2,0},uu,XX,FF);%使得DisP函数输出的Sum_N最小
    OP=optimoptions('particleswarm','FunctionTolerance',0,'InertiaRange',[0.3,1.2],'MaxStallIterations',200);
    [pp,fval,~,~]=particleswarm(fun,4,[0.1,0,0.1,0],[10,2*pi,10,2*pi],OP);
    
    %4 对比拟合效果
    figure()
    hold on
    plot3(x1,y1,z1,'*');
    funx = @(u,v) pp(1)*u.*sin(v+pp(2));
    funy = @(u,v) -pp(3)*u.*cos(v+pp(4));
    funz = @(u,v) v;
    fsurf(funx,funy,funz,[-4 4 -5 5],'LineStyle','none','FaceAlpha',0.5)
    xlim([-8,8]);
    ylim([-8,8]);
    zlim([-5,5]);
    colormap(hsv)
    camlight
    plot3([0,8],[0,0],[0,0],'linewidth',2,'color','k')
    plot3([0,0],[0,8],[0,0],'linewidth',2,'color','k')
    plot3([0,0],[0,0],[0,5],'linewidth',2,'color','k')
    hold off
    view(3)
    
    
    function Sum_N=Dis_P(p,p_num,uu,XX,FF)
    %用距离曲线的距离评价拟合程度(参数方程)
    %输入:p 要拟合的参数
    %输入:p_num 数组,每个方程的参数数量
    %输入:uu 参数方程中的参数,以cell形式储存
    %输入:XX 数据,以cell形式储存
    %输入:FF 拟合函数,以cell形式储存
    
    N_F=numel(FF);%要联立的方程数量
    L=length(XX{1});%离散点的数量
    N_L=numel(XX);%拟合参数p的数量
    N_u=numel(uu);%参数方程中参数的数量
    if N_u>1
        uu_new=ndgrid_h(uu{:});
        uu=uu_new;
    end
    %循环每个点,求到直线的距离
    %在假定参数p的情况下,计算假定函数
    Point_NF=cell(N_F,1);
    p_n=1;%参数的索引
    for k=1:N_F
        F_p=p_num{k};%该方程用到的参数数量
        Point_NF{k}=FF{k}(p(p_n:(p_n+F_p-1)),uu);%计算假定函数各个点坐标
        p_n=p_n+F_p;%更新函数索引
    end
    Sum_N=0;
    for k=1:L
        %分别求每个假定函数的点,与真实离散点之间距离的平方和
        Point_distance2=0;
        for m=1:N_F
            %读取真实点坐标
            Point_distance2=Point_distance2+(Point_NF{m}-XX{m}(k)).^2;
        end
        Min_distance2=min(Point_distance2);%求出最小距离,即为点与假定函数之间的距离
        Sum_N=Sum_N+Min_distance2;
    end
    end
    
    
    function varargout=Value(f)
    %多元素赋值,例子:
    %[a,b,c]=Value([1,2,3]);%多变量赋值
    %[xy(1),xy(2),t]=Value([2,5,3]);%复合赋值
    %[b,a]=Value([a,b]);%元素交换
    %来源:hyhhhyh21,matlab爱好者
    N=numel(f);
    varargout=cell(1,N);
    if nargout~=N
        error('输入输出数量不一致')
    end
    for k=1:N
        varargout{k}=f(k);
    end
    end
    
    function X=rand_AB(AB,varargin)
    %生成区间[A,B]之间的随机分布
    %除了AB之外,其余输入与rand相同
    [A,B]=Value(AB);
    X=rand(varargin{1:end});
    X=A+X*(B-A);
    end
    
    function S=ndgrid_h(varargin)
    %来源于matlab自带的源代码。
    %用法和ndgrid用法一样,但是将输出更改为cell
    N=nargin;
    S=cell(1,N);
    if N==1
        S{1}=varargin;
    else 
        j = 1:N;
        siz = cellfun(@numel,varargin);
        if N == 2 % Optimized Case for 2 dimensions
            x = full(varargin{1}(:));
            y = full(varargin{2}(:)).';
            S{1} = repmat(x,size(y));
            S{2} = repmat(y,size(x));
        else
            for i=1:N
                x = full(varargin{j(i)});
                s = ones(1,N);
                s(i) = numel(x);
                x = reshape(x,s);
                s = siz;
                s(i) = 1;
                S{i} = repmat(x,s);
            end
        end
    end
    S2=cell(1,N);
    for k=1:N
        S2{k}=S{k}(:);
    end
    S=S2;
    end
    

    主要函数就是Dis和Dis_P,准备工作就是把方程和离散点整理成可以输入的形式。优化用到的函数就是PSO(particleswarm),需要更改未知参数数量和范围就可以。

    3.3 补充

    这一部分的拟合优度或者置信区间就没办法拿出来说了,毕竟有试凑的嫌疑。适合单纯的验证。

    另外这里和最小二乘求出来的结果会有所不同。这主要是定义的问题。其中最小二乘指的是距离y方向上最小的二乘距离。而这一章节中定义的最小二乘,是指的点到拟合曲线距离的最小二乘(有点类似于主成分分析)。虽然一般情况下两者并不会有太大差别,但是如果有误差的话,还请注意这一点,

    如下图,这里是一个椭圆分布的离散点。
    在这里插入图片描述
    如果将方程的形式定义为y=ax+b,要求距离y最小,则最终结果为红色的线。如果将方程定义为ax+by+c=0,要求距离直线的距离最小,则最终拟合结果为黄色的线。这里无所谓对错,只是要求的最小定义不同。

    单独在后面补充了这么一个章节,希望在应用中注意。

    上面的图,绘图代码如下:

    %构建原始数据
    t1=rand(3000,1);t2=rand(3000,1);
    x=2*sqrt(t2).*cos(2*pi*t1);
    y=sqrt(t2).*sin(2*pi*t1);
    XY=[x,y];
    XY2=XY*[1,1;-1,1];%变形
    x2=XY2(:,1);y2=XY2(:,2);
    x=x2;y=y2;
    %% 演示1
    %线性最小二乘拟合,以|y-yi|^2作为衡量指标
    yn1=x;
    yn2=ones(size(x));
    X=[yn1,yn2];
    %拟合,得到系数
    Pn=X\y;
    yn=Pn(1)*yn1+Pn(2)*yn2;
    %% 演示2
    %二维直线非线性拟合,以|ax+by+c|^2作为衡量指标
    %准备数据
    XX={x,y};
    %准备函数
    F1=@(p1,XX1) p1(1)*XX1{1}-p1(2)*XX1{2}+p1(3);%ax-by+c=0
    FF={F1};
    %3 生成最终优化函数,带入到优化方程中求解
    fun=@(p) Dis(p,{3},XX,FF);%p参数,{2}第一个方程2个参数。XX离散点。FF函数表达式。
    % Sum_N=Dis([1,0],{2},XX,FF)
    OP=optimoptions('particleswarm','FunctionTolerance',0);
    [p,fval,~,~]=particleswarm(fun,3,[1,1,-1],[5,5,1],OP);
    %% 对比
    figure()
    hold on
    plot(x,y,'.')
    plot(x,yn,'linewidth',4)
    plot(x,(p(1)*x+p(3))/p(2),'linewidth',4)
    hold off
    legend({'原始数据','最小二乘','非线性拟合'})
    

    参考文献

    信赖域法 https://zhuanlan.zhihu.com/p/205555114

    展开全文
  • python完成非线性拟合

    2021-12-26 22:48:05
    在之前的博客"使用python来完成数据的线性拟合"当中,介绍了基于python,使用三种方法完成线性拟合的理论和代码实现。同样经常会碰到样本分布呈现非线性关系的情况,那么如何拟合出来呢?本文侧重对数据已经有建模,...

            在之前的博客"使用python来完成数据的线性拟合"当中,介绍了基于python,使用三种方法完成线性拟合的理论和代码实现。同样经常会碰到样本分布呈现非线性关系的情况,那么如何拟合出来呢?本文侧重对数据已经有建模,但是准确的关系需要得以确定的情况。

            如果想直接求出拟合系数,而不清楚原本模型的话,直接利用theta = np.polyfit(X, Y_noise, deg=4)得到y=a*x^4+b*x^3+c*x^2+d方程的各个权重系数,theta=[a,b,c,d]。这里deg=4表明多项式的最高次方为4.

           非线性拟合,也可以利用线性拟合的方法,只不过需要将非线性模型进行整理,做成线性模型来处理。举例:比如我有一个模型,它的曲线方程形式为 y=a*x^4+b*x^2+c. 很明显,这是一个偶函数,关于y轴对称。输入为x,输出为y。权重参数a、b、c需要通过已知的大量样本点(x,y)来求得。

           将问题转换为线性拟合的问题,即已知x,那么x0=x^4,  x1=x^2,  x2=x^0=1的x0,x1,x2就能够直接得到了。更近一步,模型 y=a*x^4+b*x^2+c其实就可以表示为 y=a*x0+b*x1+c*x2,而后者恰恰就是线性模型了,因为这里a,b,c为未知权重系数,而且都是一次幂,x0,x1,x2都是已知。

           所以,就可以直接套用线性拟合的方法来求取此本质为非线性模型,被我们转换为线性模型来计算的最佳参数a,b,c了。

           代码实现功能如下,

    1. 完成了利用预设模型y= a * x ** 4 + b * x**2 + c (已假设a=10,b=6,c=20)来制造样本数据。

     

    2. 然后利用不同方法来拟合得到a,b,c。

    3.可视化拟合结果。

    import numpy as np
    from sklearn.linear_model import LinearRegression
    from matplotlib import pyplot as plt
    
    # 生成数据
    SAMPLE_NUM = 1000
    X = np.linspace(-2, 2, SAMPLE_NUM)
    a = 10
    b = 6
    c = 20
    Y = list(map(lambda x: a * x ** 4 + b * x ** 2 + c, X))
    Y_noise = [np.random.randn() * 5 + y for y in Y]
    plt.figure()
    plt.legend()
    plt.xlabel("x")
    plt.ylabel("y")
    plt.scatter(X, Y_noise, c='blue', s=10, linewidths=0.3)
    plt.show()
    
    # 直接求解,每个x阶次都要有对应的权重参数
    # 使用最小二乘法做多项式拟合
    theta = np.polyfit(X, Y_noise, deg=4)
    print(theta)
    
    # 如果曲线方程形式是已知的,要确定准确的方程参数。那么就可以用如下方式:
    x0 = np.array(list(map(lambda x: x ** 4, X)))
    x1 = np.array(list(map(lambda x: x ** 2, X)))
    x2 = np.ones(len(X))
    
    # shape=(SAMPLE_NUM,3)
    A = np.stack((x0, x1, x2), axis=1)
    b = np.array(Y_noise).reshape((SAMPLE_NUM, 1))
    
    print("方法列表如下:"
                  "1.最小二乘法 least square method "
                  "2.常规方程法 Normal Equation "
                  "3.线性回归法 Linear regression")
    method = int(input("method="))
    
    if method == 1:
        theta, _, _, _ = np.linalg.lstsq(A, b, rcond=None)
        theta = theta.flatten()
        a_ = theta[0]
        b_ = theta[1]
        c_ = theta[2]
        print("拟合结果为: y={:.4f}*x^4+{:.4f}*x^2+{:.4f}".format(a_, b_, c_))
        Y_predict = list(map(lambda x: a_ * x ** 4 + b_ * x ** 2 + c_, X))
        plt.scatter(X, Y_noise, s=10, linewidths=0.3)
        plt.plot(X, Y_predict, c='red')
        plt.title("method {}: y={:.4f}*x^4+{:.4f}*x^2+{:.4f}".format(method, a_, b_, c_))
        plt.show()
    elif method == 2:
        AT = A.T
        A1 = np.matmul(AT, A)
        A2 = np.linalg.inv(A1)
        A3 = np.matmul(A2, AT)
        A4 = np.matmul(A3, b)
        A4 = A4.flatten()
        print("A4=", A4)
        a_ = A4[0]
        b_ = A4[1]
        c_ = A4[2]
        print("拟合结果为: y={:.4f}*x^4+{:.4f}*x^2+{:.4f}".format(a_, b_, c_))
        Y_predict = list(map(lambda x: a_ * x ** 4 + b_ * x ** 2 + c_, X))
        plt.scatter(X, Y_noise, s=10, linewidths=0.3)
        plt.plot(X, Y_predict, c='red')
        plt.title("method {}: y={:.4f}*x^4+{:.4f}*x^2+{:.4f}".format(method, a_, b_, c_))
        plt.show()
    elif method == 3:
        # 利用线性回归构建模型,拟合数据
        model = LinearRegression()
        X_normalized = A
        Y_noise_normalized = np.array(Y_noise).reshape((SAMPLE_NUM, 1))  #
        model.fit(X_normalized, Y_noise_normalized)
        # 利用已经拟合到的模型进行预测
        Y_predict = model.predict(X_normalized)
        a_ = model.coef_.flatten()[0]
        b_ = model.coef_.flatten()[1]
        c_ = model.intercept_[0]
        print(model.coef_)
        print("拟合结果为: y={:.4f}*x^4+{:.4f}*x^2+{:.4f}".format(a_, b_, c_))
        plt.scatter(X, Y_noise, s=10, linewidths=0.3)
        plt.plot(X, Y_predict, c='red')
        plt.title("method {}: y={:.4f}*x^4+{:.4f}*x^2+{:.4f}".format(method, a_, b_, c_))
        plt.show()
    else:
        print("请重新输入")

    展开全文
  • Matlab线性拟合和非线性拟合

    万次阅读 多人点赞 2019-01-23 09:25:36
    线性拟合 已知如下图像的x,y坐标,x = [1.0, 1.5, 2.0, 2.5, 3.0],y = [0.9, 1.7, 2.2, 2.6, 3.0],如何用一条直线去拟合下列散点? 代码: x = [1.0, 1.5, 2.0, 2.5, 3.0]'; y = [0.9, 1.7, 2.2, 2.6, 3.0]'...

    线性拟合

    已知如下图像的x,y坐标,x = [1.0, 1.5, 2.0, 2.5, 3.0],y = [0.9, 1.7, 2.2, 2.6, 3.0],如何用一条直线去拟合下列散点?

    代码:

    x = [1.0, 1.5, 2.0, 2.5, 3.0]';
    y = [0.9, 1.7, 2.2, 2.6, 3.0]';
    a = polyfit(x,y,1)  % a会返回两个值,[斜率,x=0时y的值]
    xi = 1:0.1:3;
    yi = polyval(a,xi);
    plot(x,y,'o',xi,yi);

    拟合后的图像如下:

    非线性拟合

    普通非线性拟合:已知如下图像的x,y坐标,x = [1.0, 1.5, 2.0, 2.5, 3.0],y = [0.9, 1.7, 2.2, 2.6, 3.0],如何用一条曲线a*x+b*sin(x)+c去拟合下列散点?

    代码:

    x = [1.0, 1.5, 2.0, 2.5, 3.0]';
    y = [0.9, 1.7, 2.2, 2.6, 3.0]';
    p = fittype('a*x+b*sin(x)+c');
    f = fit(x,y,p)   % f会返回a,b,c的值
    plot(f,x,y);

    指数非线性拟合:如何拟合1790-1900年美国人口指数增长模型

    思路:

    实现代码:

    t = 1790:10:1900;
    p = [3.9 5.3 7.2 9.6 ...
        12.9 17.1 23.2 31.4 ...
        38.6 50.2 62.9 76.0];
    Y = log(p); % Y = log(p) 返回数组p中每个元素的自然对数ln(x)
    X = t;
    a = polyfit(X,Y,1)
    x0 = exp(a(2)); r = a(1);
    ti = 1790:1900;
    pti= x0*exp(r*ti);
    plot(t,p,'o',ti,pti,'m')
    xlabel('Year')
    ylabel('Population')

    对数形式的非线性拟合:

    %% 对数形式非线性回归
    x = [1.5, 4.5, 7.5,10.5,13.5,16.5,19.5,22.5,25.5];
    y = [7.0,4.8,3.6,3.1,2.7,2.5,2.4,2.3,2.2];
    plot(x, y, '*', 'linewidth', 1) % 这里的linewidth指的是散点大小
    m1 = @(b,x) b(1) + b(2)*log(x);
    nonlinfit1 = fitnlm(x,y,m1,[0.01;0.01])
    b = nonlinfit1.Coefficients.Estimate;
    Y1 = b(1,1) + b(2,1)*log(x);
    hold on 
    plot(x, Y1, '--k', 'linewidth',2)

    展开全文
  • MATLAB多元非线性拟合

    千次阅读 2021-04-20 10:48:21
    使用“regress” 线性的不行,用二次函数。 format long A=[... 0.2 13.6 8503 251 27.4 7.7 9.9 3658 314 13.9 5.8 10.8 7307 433 26.8 7.70 9.70 6717 257 23.8 7.5 9.8 7609 280 21.7 ...

    辛苦几天收集的资料:

    1.使用“nlinfit”

    x1=1150,1000,900,850,700,625,550,475,3350,3500,5900,5800,5700,4600,4625,4725,11650,11200,11200]';

    x2=[175,100,25,0,75,100,150,200,50,600,500,225,100,1225,1600,2000,1200,1000,1550

    ]';

    x=[x1,x2];

    y=[1.44E-02,1.80E-02,6.08E-02,5.59E-02,3.42E-02,7.74E-03,1.17E-03,6.16E-03,1.91E-04,1.91E-04,1.02E-03,2.83E-03,9.52E-05,3.77E-04,2.70E-04,1.87E-04,3.98E-04,4.04E-04,4.02E-04]';

    beta0=[0.1 0.1 1 1];

    myfun=@(a,x)4030.0./pi./4.2./(a(1).*x(:,1).^a(2).*a(3).*x(:,1).^a(4)).*exp(-(x(:,2).^2./2./(a(1).*x(:,1).^a(2)).^2+30.0.^2./2./(a(3).*x(:,1).^a(4)).^2));

    [a,b,c,d,res]=nlinfit(x,y,myfun,beta0);a,res

    plot3(x1,x2,y,'o',x1,x2,myfun(a,x))

    值的选取没有定法,与实际问题的模型有关。

    2.使用“regress”

    线性的不行,用二次函数。

    format long

    A=[...

    0.2 13.6 8503 251 27.4 7.7 9.9 3658 314 13.9 5.8 10.8 7307 433 26.8 7.70 9.70 6717 257 23.8 7.5 9.8 7609 280 21.7 5.6 11.3 4271 533 14.6 6.2 7.6 52169 48 225 3.23 9.16 16516 80 44.1 0.33 11.3 17366 85 54.1 0.14 9.5 14245 91 56.6 5.5 9.7 18184 3 31.6 2.3 8.9 33612 250 114.9 3.3 4.6 73927 5 166 1.9 9.7 32175 150 107.5 0.6 9.9 33088 242 142.3 0.22 11.7 18620 567 60.4 1.88 11.76 27885 267 71.6 2.78 10.9 21780 76 58.7]

    x=A(:,1:4),Y=A(:,5)

    x11=x(:,1).*x(:,1);

    x12=x(:,1).*x(:,2);

    x13=x(:,1).*x(:,3);

    x14=x(:,1).*x(:,4);

    x22=x(:,2).*x(:,2);

    x23=x(:,2).*x(:,3);

    x24=x(:,2).*x(:,4);

    x33=x(:,3).*x(:,3);

    x34=x(:,3).*x(:,4);

    x44=x(:,4).*x(:,4);

    X=[x(:,:),x11,x12,x13,x14,x22,x23,x24,x33,x34,x44]

    [B,BINT,R] = REGRESS(Y,[ones(length(Y),1),X]) 结果:

    B =

    1.0e+003 *

    -1.426098928217992

    -0.004076772421011

    0.255534919787513

    0.000012942581436

    0.000845938681439

    0.000607150442496

    -0.000574488595437

    0.000000405451807

    -0.000042626483419

    -0.011775830339062

    -0.000000876232149

    0.000008150156703

    -0.000000000013441

    -0.000000013991054

    -0.000000969496753

    R =

    3.122573422039807

    0.447341267999400

    -7.343326306615449

    2.107836742251767

    -6.239492394117182

    9.044235126157025

    2.238791755625499

    4.285551199892858

    -2.231536057549363

    -1.979307925154075

    3.503835830046878

    1.414933242530537

    -1.426757776398972

    -12.052007973319576

    14.597045597468522

    -5.024612350970848

    -1.747668123505179

    -2.717435276394376

    B就是系数,R就是预测值与实际值的差值。

    3. 使用“lsqcurvefit”

    clear

    clc

    x=[40 50 60 70 80 90 100 110 120 135 150];

    y=[0.0096 0.0145 0.0194 0.0348 0.0501 0.0751 0.1000 0.1497 0.1993 0.2496 0.2999];

    z=[0.2400 0.2865 0.3330 0.3600 0.3870 0.4010 0.4150 0.4390 0.4630 0.4875 0.5120];

    X0=[1 1 1 1 1 1];

    %只要这样写就可以了

    f=@(p,x)( p(1) + p(2)*x(1,:) +

    p(3)*x(2,:) + p(4)*x(1,:).^2 + p(5)*x(1,:).*x(2,:) +

    p(6)*x(2,:).^2);

    p=lsqcurvefit(f,X0,[x;y],z)

    展开全文
  • 对于某些非线性函数,如指数函数y=e^(ax+b),也可以对函数转化后,求得精确的拟合结果,如上述指数函数可转化为x=(ln y)/a -b/a,同样可以求得具有全局最优拟合误差的拟合函数。上述函数都可以用MATLAB的regress函数...
  • 基于MATLAB洛伦兹线型非线性拟合算法实现.pdf
  • 如何用matlab进行多元非线性拟合

    千次阅读 2021-04-20 10:48:20
    % 使用最小二乘拟合: % opt指定拟合选项(注意查看命令窗口提示的优化终止条件,如对结果不满意考虑适当修改) % b0为初值(要慎重选择,不同初值得到的结果可能不同) opt = optimset('MaxFunEvals',50000,'MaxIter...
  • python里面多元非线性回归有哪些方法SciPy 里面的子函数库optimize, 一般情况下可用curve_fit函数直接拟合或者leastsq做最小二乘第九句:简单的事重复做,你就是专家;重复的事用心做,你就是赢家。Python怎么实现...
  • MATLAB进行非线性拟合

    千次阅读 2020-04-20 16:59:39
    matlab进行非线性拟合常用最小二乘法,适用于:已经求解出函数,但含有未知数,不过已经收集到了一系列数据 1.lsqcurvefit 格式:[x, resnorm,r,flag]=lsqcurvefit(fun, c0,xdata,ydata) c0为初始解向量;xdata,...
  • 基于matlab的神经网络非线性函数拟合问题
  • 本期视频时长约106分钟,通过两个具体的应用案例,深入地解读了MATLAB非线性拟合函数lsqcurvefit的具体用法;生动地讲解了局部最优解的基本概念;对拟合效果的数学评价进行了MATLAB代码实现;视频的最后,还对拟合...
  • #资源达人分享计划#
  • 文章目录[MATLAB 在科学计算中的应用] 使用MATLAB 进行非线性拟合前言引述MATLAB 曲线拟合函数简述一二维数据非线性拟合一维数据拟合例子二维数据拟合例子高维数据非线性拟合lsqcurvefitnlinfit 函数数据拟合工具箱...
  • 线性拟合 % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % % 注意:在进行线性拟合之前,必须要具备两个条件 % (1)数据点 % (2)要进行拟合的...
  • origin 2016非线性拟合

    2016-07-01 08:54:09
    origin 2016官方文件的非线性拟合方法介绍
  • 使用Origin非线性拟合曲线案例,记载了使用origin进行非线性拟合曲线的方法。
  • 多元非线性拟合是非常困难的事情,我选用matlab进行拟合,以下为我找到的三种方法以及具体数据。1.使用“nlinfit”x1=1150,1000,900,850,700,625,550,475,3350,3500,5900,5800,5700,4600,4625,4725,11...
  • 此篇,我们来介绍一些其他的拟合
  • Matlab 线性拟合 非线性拟合

    千次阅读 2019-01-24 12:04:24
    Matlab 线性拟合 非线性拟合
  • Matlab 线性拟合 & 非线性拟合

    万次阅读 2018-07-07 11:45:53
    使用Matlab进行拟合是图像处理中线条变换的一个重点内容,本文将详解Matlab中的直线拟合和曲线拟合用法。 关键函数: fittype Fit type for cu
  • 非线性函数拟合

    2018-04-20 16:05:58
    神经网络非线性函数拟合神经网络非线性函数拟合神经网络非线性函数拟合神经网络非线性函数拟合
  • 非线性拟合/GAM

    千次阅读 2019-02-27 19:51:46
    前两种方法都属于参数方法,下面介绍一种结合了参和参数模型的方法:lowess或者又叫做loess,通过对局部进行线性拟合(非常数值拟合),来解决NW的缺点。 lowess 算法: 1)给定数据集 D = ( x i , y i ) D=...
  • 比如z=f(x,y),给出(x,y,z)的多个数据点,然后怎么拟合出f函数,给出源码就更好了</p>

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 80,320
精华内容 32,128
关键字:

非线性拟合