当前位置:   article > 正文

Matlab 线性拟合、一维、多维度非线性拟合、多项式拟合_matlab线性拟合

matlab线性拟合

一、线性拟合

  线性拟合

我随便设定一个函数然后通过解方程计算出对应的系数

假设我的函数原型是

y=a*sin(0.1*x.^2+x)+b* squre(x+1)+c*x+d

  1. clc;
  2. clear;
  3. x=0:0.2:10;
  4. % 我们这里假设 a=3.2 b=0.7 c=5.0 d是一个随机
  5. y=3.2*sin(0.1*x.^2+x)+0.7*sqrt(x+1)+5*x +rand(size(x));
  6. plot(x,y,'*');
  7. hold on ;
  8. y1=sin(0.1*x.^2+x);
  9. y2=sqrt(x+1);
  10. y3=x;
  11. y4=rand(size(x));
  12. X=[y1;y2;y3;y4];%将各自的俩带入
  13. P=X'\y' % 通过解方程计算出4个系数
  14. yn=P(1)*y1+P(2)*y2+P(3)*y3+P(4)*y4; % 得到一个新的函数 计算得出的拟合Y的值
  15. plot(x,yn,'r');
  16. legend('原始数据点','红色曲线拟合值','Location','southoutside','Orientation','horizontal')

 拟合系数:

 

  1. clear;
  2. clc;
  3. close all
  4. t=0:0.001:2*pi;%原函数
  5. YS=sin(t);%基函数
  6. N=21;
  7. Yo=[];
  8. for k=1:N
  9. Yn=sawtooth(k*(t+pi/2),0.5);
  10. Yo=[Yo,Yn'];
  11. end
  12. YS=YS';%拟合
  13. a = Yo\YS;
  14. %绘图
  15. figure()
  16. for k=1:N
  17. clf
  18. plot(t,Yo(:,1:k)*a(1:k,:),t,YS,'LineWidth',1)
  19. ylim([-1.3,1.3])
  20. xlim([0,6.3])
  21. pause(0.1)
  22. end

二、非线性拟合

利用matlab实现非线性拟合(三维、高维、参数方程)_matlab多元非线性拟合_hyhhyh21的博客-CSDN博客

上面的这位是真正的大佬,我们都是照猫画虎的学习。

1、一维

简单的一维的拟合:

思路: 将非线性-》线性:

eg:

 将其两边都取对数

用线性的方式计算出a b 

逆变换 ,画出曲线:

  1. clear
  2. clc
  3. close all
  4. % 假设函数 为 y=a* exp(-bx)
  5. x=0:0.1:5;
  6. % 我们这里假设 a=2.4 b=1.2
  7. a=2.4;
  8. b=1.2;
  9. y=a*exp(-b*x);%
  10. y=y+1.3*y.*rand(size(y)); % 增加噪声
  11. plot(x,y,'.');
  12. hold on;
  13. %Lg_y=Lg_a+b*(-x) 变成了ax+b 的形式 ,然而我们的最终的目的是通过x 来计算出a 和 b
  14. % 对等式的两边取对数
  15. lg_y=log(y);
  16. y1=ones(size(x));
  17. y2=-x;
  18. % 同理和上面计算线性的一杨
  19. X=[y1;y2];
  20. P =X'\lg_y'
  21. % 画出拟合后的曲线
  22. a_fit=exp(P(1));
  23. b_fit=P(2);
  24. x2=0:0.01:10;
  25. plot(x2,a_fit*exp(-b_fit*x2),'-','linewidth',1.5,'color','r')

Matlab 中的非线性拟合方法

1、fit 方法

fit是最常用,最经典的方法

 

  1. ft = fittype( 'a*x+b*sin(c*x).*exp(d*x)+e', 'independent', 'x', 'dependent', 'y' );; %函数的表达式,
  2. OP1 = fitoptions( 'Method', 'NonlinearLeastSquares' );% 非线性拟合方法
  3. OP1.StartPoint = 5*rand(1,5);%初始值,越靠近真实值越好
  4. OP1.Lower = [-2, 0, 2, 0, 0];%参数的最小边界
  5. OP1.Upper = [1, 3, 5, 2, 3];%参数的最大边界
  6. % 开始拟合
  7. fitobject = fit(x',y',ft,OP1);
  8. Fit_P=ones(size(P));
  9. Fit_P(1)=fitobject.a;
  10. Fit_P(2)=fitobject.b;
  11. Fit_P(3)=fitobject.c;
  12. Fit_P(4)=fitobject.d;
  13. Fit_P(5)=fitobject.e;

2、nlinfit()函数 Levenberg-Marquard

L-M 非线性迭代
 

  1. % 2 用nlinfit()函数 Levenberg-Marquardt 
  2. % 定义一个函数
  3. Func=@(P,x)( P(1)*x+P(2)*sin(P(3)*x).*exp(P(4)*x)+P(5));% 也就是说定义一个函数模型
  4. OP2 = statset('nlinfit');
  5. %  x,y  modelfun是函数模型  beta表示的是初始值 ,我这里写成最进行的那个参数   OP2  拟合的方法
  6. beta=[-0.17 2.1  3.0  0.25  2.0];% 初始值
  7. Fit_P2 = nlinfit(x,y,Func,beta,OP2);
  8. %拟合
  9. fit_y2 = Fit_P2(1)*x1+Fit_P2(2)*sin(Fit_P2(3)*x1).*exp(Fit_P2(4)*x1)+Fit_P2(5);;
  10. subplot(3,2,2)
  11. hold on;
  12. plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
  13. plot(x1,fit_y2,'-','linewidth',1.5,'color','r');
  14. box on
  15. %ylim(y_lim)
  16. title('nlinfit函数')

3、信赖域法(trust region reflective)

信赖域法(trust region reflective)是通过Hessian 矩阵,逐步试探邻域内的最小化,来求解问题的。相比较之前的那些雅克比相关的方法,信赖域法会占用更多内存和速度,所以适用于中小规模的矩阵。

  1. % 3 信赖区间 IsqNonLin()
  2. func2=@(P)(P(1)*x+P(2)*sin(P(3)*x).*exp(P(4)*x)+P(5) -y);
  3. % lsqnonlin方法
  4. % 'Algorithm','trust-region-reflective' 算法是trust-region-reflective
  5. % MaxFunctionEvaluations MaxFunctionEvaluations可以理解为试探的次数,
  6. % 比如算法在一个点的四周试探了三个邻近点的值,然后确定下一步要往其中的某个点走,
  7. % 这个时候FunctionEvaluations对应3次,即试探了3次,而Iteration是一次,即走了一步,完成了一步迭代
  8. % MaxIterations 最大迭代次数
  9. OP3=optimoptions(@lsqnonlin,'Algorithm','trust-region-reflective','MaxFunctionEvaluations',1e4,'MaxIterations',1e3);
  10. %[-2, 0, 2, 0, 0];%参数的最小边界
  11. %[1, 3, 5, 2, 3];%参数的最大边界
  12. lower=[-2, 0, 2, 0, 0];
  13. up=[1, 3, 5, 2, 3];
  14. % 计算出系数
  15. [Fit_P3,~]=lsqnonlin(func2,beta,lower,up,OP3);
  16. fit_y3=Fit_P3(1)*x1+Fit_P3(2)*sin(Fit_P3(3)*x1).*exp(Fit_P3(4)*x1)+Fit_P3(5);
  17. subplot(3,2,3);
  18. hold on
  19. plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
  20. plot(x1,fit_y3,'-','linewidth',1.5,'color','r')
  21. hold off
  22. box on
  23. %ylim(y_lim)
  24. title('lsqnonlin函数');
  25. %4 lsqcurvefit()函数 trust-region-reflective
  26. modelfun2 = @(p,x)(p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5)) ;
  27. OP4=optimoptions('lsqcurvefit','Algorithm','trust-region-reflective','MaxFunctionEvaluations',1e4,'MaxIterations',1e3);
  28. %[-2, 0, 2, 0, 0];%参数的最小边界
  29. %[1, 3, 5, 2, 3];%参数的最大边界
  30. lower=[-2, 0, 2, 0, 0];
  31. up=[1, 3, 5, 2, 3];
  32. % 计算出系数
  33. %[p,~] = lsqcurvefit(modelfun,p0,x,y,[-2,0,2,0,0],[1,3,5,3,3],OP4);
  34. [Fit_P4,~]=lsqcurvefit(modelfun2,beta,x,y,lower,up,OP4);
  35. fit_y4=Fit_P4(1)*x1+Fit_P4(2)*sin(Fit_P4(3)*x1).*exp(Fit_P4(4)*x1)+Fit_P4(5);
  36. subplot(3,2,4);
  37. hold on
  38. plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
  39. plot(x1,fit_y4,'-','linewidth',1.5,'color','r')
  40. hold off
  41. box on
  42. %ylim(y_lim)
  43. title('lsqcurvefit函数');

4、fsolve()函数

默认的算法为trust-region-dogleg,属于信赖域法。

5、粒子群法

 所有代码:

  1. clear;
  2. close all;
  3. clc;
  4. % 自定义一个非线性的函数 y=a*x+b*sin(c*x).*exp(d*x)+e 那将函数
  5. x = 0:0.05:10;
  6. P=[-0.2 2.4 3.4 0.3 1.7];
  7. y = P(1)*x+P(2)*sin(P(3)*x).*exp(P(4)*x)+P(5);
  8. y=y+0.5*randn(size(x)); % 添加噪声
  9. figure();
  10. % 1 .fit 函数开始拟合
  11. ft = fittype( 'a*x+b*sin(c*x).*exp(d*x)+e', 'independent', 'x', 'dependent', 'y' );; %函数的表达式,
  12. OP1 = fitoptions( 'Method', 'NonlinearLeastSquares' );% 非线性拟合方法
  13. %OP1.StartPoint = 5*rand(1,5);%初始值,越靠近真实值越好
  14. OP1.StartPoint = [-0.17 2.1 3.0 0.25 2.0];
  15. OP1.Lower = [-2, 0, 2, 0, 0];%参数的最小边界
  16. OP1.Upper = [1, 3, 5, 2, 3];%参数的最大边界
  17. % 开始拟合
  18. fitobject = fit(x',y',ft,OP1);
  19. Fit_P=ones(size(P));
  20. Fit_P(1)=fitobject.a;
  21. Fit_P(2)=fitobject.b;
  22. Fit_P(3)=fitobject.c;
  23. Fit_P(4)=fitobject.d;
  24. Fit_P(5)=fitobject.e;
  25. %plot(x,y,'.');
  26. % 开始计算拟合后的y
  27. x1 = 0:0.01:10;
  28. fit_y1 = Fit_P(1)*x1+Fit_P(2)*sin(Fit_P(3)*x1).*exp(Fit_P(4)*x1)+Fit_P(5);
  29. subplot(3,2,1)
  30. hold on
  31. plot(x,y,'LineStyle','none','MarkerSize',10,'Marker','.','color','k');
  32. plot(x1,fit_y1,'-','linewidth',1.5,'color','r');
  33. % 开始计算拟合后的y
  34. fit_y1 = Fit_P(1)*x+Fit_P(2)*sin(Fit_P(3)*x).*exp(Fit_P(4)*x)+Fit_P(5);
  35. subplot(3,2,1)
  36. plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
  37. plot(x,fit_y1,'-','linewidth',1.5,'color','r');
  38. hold on;
  39. title('经典fit函数');
  40. box on;
  41. % 2 用nlinfit()函数 Levenberg-Marquardt
  42. % 定义一个函数
  43. Func=@(P,x)( P(1)*x+P(2)*sin(P(3)*x).*exp(P(4)*x)+P(5));% 也就是说定义一个函数模型
  44. OP2 = statset('nlinfit');%
  45. % x,y modelfun是函数模型 beta表示的是初始值 ,我这里写成最进行的那个参数 OP2 拟合的方法
  46. beta=[-0.17 2.1 3.0 0.25 2.0];% 初始值
  47. Fit_P2 = nlinfit(x,y,Func,beta,OP2);
  48. %拟合
  49. fit_y2 = Fit_P2(1)*x1+Fit_P2(2)*sin(Fit_P2(3)*x1).*exp(Fit_P2(4)*x1)+Fit_P2(5);
  50. subplot(3,2,2)
  51. hold on;
  52. plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
  53. plot(x1,fit_y2,'-','linewidth',1.5,'color','r');
  54. box on
  55. %ylim(y_lim)
  56. title('nlinfit函数')
  57. % 3 信赖区间 IsqNonLin()
  58. func2=@(P)(P(1)*x+P(2)*sin(P(3)*x).*exp(P(4)*x)+P(5) -y);
  59. % lsqnonlin方法
  60. % 'Algorithm','trust-region-reflective' 算法是trust-region-reflective
  61. % MaxFunctionEvaluations MaxFunctionEvaluations可以理解为试探的次数,
  62. % 比如算法在一个点的四周试探了三个邻近点的值,然后确定下一步要往其中的某个点走,
  63. % 这个时候FunctionEvaluations对应3次,即试探了3次,而Iteration是一次,即走了一步,完成了一步迭代
  64. % MaxIterations 最大迭代次数
  65. OP3=optimoptions(@lsqnonlin,'Algorithm','trust-region-reflective','MaxFunctionEvaluations',1e4,'MaxIterations',1e3);
  66. %[-2, 0, 2, 0, 0];%参数的最小边界
  67. %[1, 3, 5, 2, 3];%参数的最大边界
  68. lower=[-2, 0, 2, 0, 0];
  69. up=[1, 3, 5, 2, 3];
  70. % 计算出系数
  71. [Fit_P3,~]=lsqnonlin(func2,beta,lower,up,OP3);
  72. fit_y3=Fit_P3(1)*x1+Fit_P3(2)*sin(Fit_P3(3)*x1).*exp(Fit_P3(4)*x1)+Fit_P3(5);
  73. subplot(3,2,3);
  74. hold on
  75. plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
  76. plot(x1,fit_y3,'-','linewidth',1.5,'color','r')
  77. hold off
  78. box on
  79. %ylim(y_lim)
  80. title('lsqnonlin函数');
  81. %4 lsqcurvefit()函数 trust-region-reflective
  82. modelfun2 = @(p,x)(p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5)) ;
  83. OP4=optimoptions('lsqcurvefit','Algorithm','trust-region-reflective','MaxFunctionEvaluations',1e4,'MaxIterations',1e3);
  84. %[-2, 0, 2, 0, 0];%参数的最小边界
  85. %[1, 3, 5, 2, 3];%参数的最大边界
  86. lower=[-2, 0, 2, 0, 0];
  87. up=[1, 3, 5, 2, 3];
  88. % 计算出系数
  89. %[p,~] = lsqcurvefit(modelfun,p0,x,y,[-2,0,2,0,0],[1,3,5,3,3],OP4);
  90. [Fit_P4,~]=lsqcurvefit(modelfun2,beta,x,y,lower,up,OP4);
  91. fit_y4=Fit_P4(1)*x1+Fit_P4(2)*sin(Fit_P4(3)*x1).*exp(Fit_P4(4)*x1)+Fit_P4(5);
  92. subplot(3,2,4);
  93. hold on
  94. plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
  95. plot(x1,fit_y4,'-','linewidth',1.5,'color','r')
  96. hold off
  97. box on
  98. %ylim(y_lim)
  99. title('lsqcurvefit函数');
  100. %% 5 fsolve()函数 %默认算法trust-region-dogleg
  101. modelfun3 = @(p)(p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5) -y);
  102. p0 = 5*rand(1,5);
  103. OP5 = optimoptions('fsolve','Display','off');
  104. Fit_P = fsolve(modelfun3,beta,OP5);
  105. fit_y5 = Fit_P(1)*x1+Fit_P(2)*sin(Fit_P(3)*x1).*exp(Fit_P(4)*x1)+Fit_P(5);
  106. subplot(3,2,5)
  107. hold on
  108. plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
  109. plot(x1,fit_y5,'-','linewidth',1.5,'color','r')
  110. hold off
  111. box on
  112. title('fsolve函数')
  113. %% 6 粒子群PSO算法
  114. fun6 = @(p) ( norm(p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5) -y) );%这里需要计算误差的平方和
  115. OP6 = optimoptions('particleswarm','InertiaRange',[0.4,1.2],'SwarmSize',100);
  116. [p,~,~,~] = particleswarm(fun6,5,[-5,-5,-5,-5],[5,5,5,5],OP6);%区间可以稍微放大一些,不怕
  117. y6 = p(1)*x+p(2)*sin(p(3)*x).*exp(p(4)*x)+p(5);
  118. subplot(3,2,6)
  119. hold on
  120. plot(x,y,'LineStyle','none','MarkerSize',15,'Marker','.','color','k')
  121. plot(x,y6,'-','linewidth',1.5,'color','r')
  122. hold off
  123. box on
  124. ylim(y_lim)
  125. title('PSO算法')

三、多项式曲线

 Matlab:

  1. >> x=linspace(0,4*pi,150);
  2. y=cos(x)+10*rand(1);
  3. plot(x,y,'.');
  4. hold on;
  5. [p,s]=polyfit(x,y,9);% 拟合为7阶的函数
  6. x1=linspace(0,4*pi,150);
  7. y1=polyval(p,x1);
  8. plot(x1,y1,'color','r');
  9. p
  10. p =
  11. 0.0000 0.0000 -0.0004 0.0081 -0.0783 0.3753 -0.7660 0.3815 -0.4104 10.6154
  12. % 方程变换
  13. >> x=linspace(0,4*pi,150);
  14. y=2*exp(-(x-1).^2/1.^2)+0.1*rand(1);
  15. plot(x,y,'.')
  16. >> x=linspace(0,4*pi,50);
  17. y=2*exp(-(x-1).^2/1.^2)+0.1*rand(1);
  18. plot(x,y,'.')
  19. >> x=linspace(0,4*pi,50);
  20. y=2*exp(-(x-1).^2/1.^2)+0.1*rand(1);
  21. plot(x,y,'.');
  22. hold on;
  23. [p,s]=polyfit(x,y,9);% 拟合为7阶的函数
  24. x1=linspace(0,4*pi,50);
  25. y1=polyval(p,x1);
  26. plot(x1,y1,'color','r')

 从上图我们可以看出9阶的拟合效果要比7阶的好很多,那么我们用c++实现的时候也就按照9阶的来。

声明:本文内容由网友自发贡献,不代表【wpsshop博客】立场,版权归原作者所有,本站不承担相应法律责任。如您发现有侵权的内容,请联系我们。转载请注明出处:https://www.wpsshop.cn/w/凡人多烦事01/article/detail/104497
推荐阅读
相关标签
  

闽ICP备14008679号