devbox
weixin_40725706
这个屌丝很懒，什么也没留下！
关注作者
热门标签
热门文章
当前位置:   article > 正文

2024美赛数学建模常用数学建模模型之——最小二乘法_数学建模最小二乘多项式拟合方法

作者：weixin_40725706 | 2024-04-28 17:32:08

赞

踩

数学建模最小二乘多项式拟合方法

2.1 线性最小二乘法

       曲线拟合问题的提法 是，已知一组（二维）数据，即平面上的 n 个点 ( x i , y i ) ， i = 1,2, L , n ， x i 互不相同，寻求一个函数（曲线） y = f ( x ) ，使 f ( x ) 在某种准则下 与所有数据点最为接近，即曲线拟合得最好。
   线性最小二乘法 是解决曲线拟合最常用的方法，基本思路是，令

2.2 最小二乘法的 Matlab 实现

2.2.1 解方程组方法
在上面的记号下， 
  1. x=[19 25 31 38 44]'; 
  2. y=[19.0 32.3 49.0 73.3 97.8]'; 
  3. r=[ones(5,1),x.^2]; 
  4. ab=r\y 
  5. x0=19:0.1:44; 
  6. y0=ab(1)+ab(2)*x0.^2; 
  7. plot(x,y,'o',x0,y0,'r')
2.2.2 多项式拟合方法
如果取 { r 1 ( x ), L , r m + 1 ( x )} = {1, x , L , x m } ，即用 m 次多项式拟合给定数据， Matlab 中有现成的函数
                                                a=polyfit(x0,y0,m)
其中输入参数 x0,y0 为要拟合的数据， m 为拟合多项式的次数，输出参数 a 为拟合多项
式 y=a m x m + …+a 1 x+a 0 系数 a=[ a m , …, a 1 , a 0 ]。
多项式在 x 处的值 y 可用下面的函数计算
                                                y=polyval(a,x)。
例 5 某乡镇企业 1990-1996 年的生产利润如表 5
试预测 1997 年和 1998 年的利润。
作已知数据的的散点图， 
  1. x0=[1990 1991 1992 1993 1994 1995 1996]; 
  2. y0=[70 122 144 152 174 196 202]; 
  3. plot(x0,y0,'*')
发现该乡镇企业的年生产利润几乎直线上升。因此，我们可以用 y = a 1 x + a 0 作为
拟合函数来预测该乡镇企业未来的年利润。编写程序如下： 
  1. x0=[1990 1991 1992 1993 1994 1995 1996]; 
  2. y0=[70 122 144 152 174 196 202]; 
  3. a=polyfit(x0,y0,1) 
  4. y97=polyval(a,1997) 
  5. y98=polyval(a,1998)

2.3 最小二乘优化

        在无约束最优化问题中，有些重要的特殊情形，比如目标函数由若干个函数的平方和构成。这类函数一般可以写成：

lsqlin 函数

Matlab 中的函数为：
x=lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0) 
  1. x=[19 25 31 38 44]'; 
  2. y=[19.0 32.3 49.0 73.3 97.8]'; 
  3. r=[ones(5,1),x.^2]; 
  4. ab=lsqlin(r,y) 
  5. x0=19:0.1:44; 
  6. y0=ab(1)+ab(2)*x0.^2; 
  7. plot(x,y,'o',x0,y0,'r')

lsqcurvefit 函数

解 该问题即解最优化问题： 
  1. function f=fun1(x,tdata);
  2. f=x(1)+x(2)*exp(-0.02*x(3)*tdata); %其中 x(1)=a，x(2)=b，x(3)=k
调用函数 lsqcurvefit，编写程序如下： 
  1. td=100:100:1000; 
  2. cd=[4.54 4.99 5.35 5.65 5.90 6.10 6.26 6.39 6.50 6.59]; 
  3. x0=[0.2 0.05 0.05]; 
  4. x=lsqcurvefit(@fun1,x0,td,cd)

lsqnonlin 函数 

  1. function f=fun2(x); 
  2. td=100:100:1000; 
  3. cd=[4.54 4.99 5.35 5.65 5.90 6.10 6.26 6.39 6.50 6.59]; 
  4. f=x(1)+x(2)*exp(-0.02*x(3)*td)-cd;
调用函数 lsqnonlin，编写程序如下： 
  1. x0=[0.2 0.05 0.05]; %初始值是任意取的 
  2. x=lsqnonlin(@fun2,x0)

lsqnonneg 函数

编写程序如下： 
  1. c=[0.0372 0.2869;0.6861 0.7071;0.6233 0.6245;0.6344 0.6170]; 
  2. d=[0.8587;0.1781;0.0747;0.8405]; 
  3. x=lsqnonneg(c,d)

曲线拟合的用户图形界面求法

Matlab 工具箱提供了命令 cftool ，该命令给出了一维数据拟合的交互式环境。具体
执行步骤如下：
（ 1 ）把数据导入到工作空间；
（ 2 ）运行 cftool ，打开用户图形界面窗口；
（ 3 ）对数据进行预处理；
（ 4 ）选择适当的模型进行拟合；
（ 5 ）生成一些相关的统计量，并进行预测。
可以通过帮助（运行 doc cftool ）熟悉该命令的使用细节。

2.4 曲线拟合与函数逼近

前面讲的曲线拟合是已知一组离散数据{( x i , y i), i =1,L, n}，选择一个较简单的函数 f ( x)，如多项式，在一定准则如最小二乘准则下，最接近这些数据。
解 编写程序如下： 
  1. syms x 
  2. base=[1,x^2,x^4]; 
  3. y1=base.'*base 
  4. y2=cos(x)*base.' 
  5. r1=int(y1,-pi/2,pi/2) 
  6. r2=int(y2,-pi/2,pi/2) 
  7. a=r1\r2 
  8. xishu1=double(a) 
  9. digits(8),xishu2=vpa(a)
求得 xishu1=0.9996 -0.4964 0.0372，即所求的最佳平方逼近多项式为

2.5 黄河小浪底调水调沙问题

  1. clc,clear 
  2. load data.txt %data.txt 按照原始数据格式把水流量和排沙量排成 4 行，12 列
  3. liu=data([1,3],:); 
  4. liu=liu';liu=liu(:); 
  5. sha=data([2,4],:); 
  6. sha=sha';sha=sha(:); 
  7. y=sha.*liu;y=y'; 
  8. i=1:24; 
  9. t=(12*i-4)*3600; 
  10. t1=t(1);t2=t(end); 
  11. pp=csape(t,y); 
  12. xsh=pp.coefs %求得插值多项式的系数矩阵，每一行是一个区间上多项式的系数。
  13. TL=quadl(@(tt)ppval(pp,tt),t1,t2)
也可以利用 3 次 B 样条函数进行插值，求得总的排沙量也为 1.844 × 10 9 t ，，计算
的 Matlab 程序如下： 
  1. clc,clear 
  2. load data.txt %data.txt 按照原始数据格式把水流量和排沙量排成 4 行，12 列
  3. liu=data([1,3],:); 
  4. liu=liu';liu=liu(:); 
  5. sha=data([2,4],:); 
  6. sha=sha';sha=sha(:); 
  7. y=sha.*liu;y=y'; 
  8. i=1:24; 
  9. t=(12*i-4)*3600; 
  10. t1=t(1);t2=t(end); 
  11. pp=spapi(4,t,y) %三次 B 样条
  12. pp2=fn2fm(pp,'pp') %把 B 样条函数转化为 pp 格式
  13. TL=quadl(@(tt)fnval(pp,tt),t1,t2)
对于问题（ 2），研究排沙量与水量的关系，从试验数据可以看出，开始排沙量是随
着水流量的增加而增长，而后是随着水流量的减少而减少。显然，变化规律并非是线性
的关系，为此，把问题分为两部分，从开始水流量增加到最大值 2720m 3 /s （即增长的过
程）为第一阶段，从水流量的最大值到结束为第二阶段，分别来研究水流量与排沙量的
关系。 
  1. load data.txt 
  2. liu=data([1,3],:); 
  3. liu=liu';liu=liu(:); 
  4. sha=data([2,4],:); 
  5. sha=sha';sha=sha(:); 
  6. y=sha.*liu; 
  7. subplot(1,2,1), plot(liu(1:11),y(1:11),'*') 
  8. subplot(1,2,2), plot(liu(12:24),y(12:24),'*')
从散点图可以看出，第一阶段基本上是线性关系，第二阶段准备依次用二次、三次、
四次曲线来拟合，看哪一个模型的剩余标准差小就选取哪一个模型。最后求得第一阶段
排沙量 y 与水流量 v 之间的预测模型为 
  1. clc, clear 
  2. load data.txt %data.txt 按照原始数据格式把水流量和排沙量排成 4 行，12 列
  3. liu=data([1,3],:); liu=liu'; liu=liu(:); 
  4. sha=data([2,4],:); sha=sha'; sha=sha(:); 
  5. y=sha.*liu; 
  6. %以下是第一阶段的拟合
  7. format long e 
  8. nihe1_1=polyfit(liu(1:11),y(1:11),1) %拟合一次多项式，系数排列从高次幂到低次幂
  9. nihe1_2=polyfit(liu(1:11),y(1:11),2)
  10. yhat1_1=polyval(nihe1_1,liu(1:11)); %求预测值
  11. yhat1_2=polyval(nihe1_2,liu(1:11)); 
  12. %以下求误差平凡和与剩余标准差
  13. cha1_1=sum((y(1:11)-yhat1_1).^2); rmse1_1=sqrt(cha1_1/9) 
  14. cha1_2=sum((y(1:11)-yhat1_2).^2); rmse1_2=sqrt(cha1_2/8) 
  15. %以下是第二阶段的拟合
  16. for j=1:3 
  17.  str1=char(['nihe2_' int2str(j) '=polyfit(liu(12:24),y(12:24),' int2str(j+1) ')']); 
  18.  eval(str1) 
  19.  str2=char(['yhat2_' int2str(j) '=polyval(nihe2_' int2str(j) ',liu(12:24));']); 
  20.  eval(str2) 
  21.  str3=char(['cha2_' int2str(j) '=sum((y(12:24)-yhat2_' int2str(j) ').^2);'... 
  22.  'rmse2_' int2str(j) '=sqrt(cha2_' int2str(j) '/(11-j))']); 
  23.  eval(str3) 
  24. end 
  25. format

更多美赛资料，请点击下方名片获取：

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

Copyright © 2003-2013 www.wpsshop.cn 版权所有，并保留所有权利。

  

闽ICP备14008679号