数学建模美赛（一）：拟合与插值

本博客为备考2022年数学建模美赛，进行相关知识梳理和总结，ppt内容和代码来自课程《数学建模和入门》

一、预测类数学模型——拟合与插值

1.曲线拟合问题

已知一组（二维）数据，即平面上 n 个点（xi, yi),i=1,…n, 寻求一个函数（曲线）y = f(x) , 使 f(x) 在某种准则下与所有数据点最为接近，即曲线拟合得最好。

其中会引入误差平方和等概念来表示曲线拟合的效果，该值越小，则拟合效果越好。

2.拟合与插值的关系

问题：给定一批数据点，需确定满足特定要求的曲线或曲面

解决方案：

• 若要求所求曲线（面）通过所给所有数据点，就是插值问题；

• 若不要求曲线（面）通过所有数据点，而是要求它反映对象整体的变化趋势，这就是数据拟合，又称曲线拟合或曲面拟合。

函数插值与曲线拟合都是要根据一组数据构造一个函数作为近似，由于近似的要求不同，二者的数学方法上完全不同。

二、拟合基本实例

t=[20.5 32.5 51 73 95.7];
r=[765 826 873 942 1032];
%polyfit进行线性拟合，t为自变量，r为因变量，1代表变量维数
aa=polyfit(t,r,1);
a=aa(1)
b=aa(2)
%利用polyval计算t对应的拟合函数值
y=polyval(aa,t);
%plot代表二维绘图，t为自变量，r为因变量，r，k+代表图像标记方式，字母为颜色（r代表red)，
%+代表表现符号
plot(t,r,'k+',t,y,'r')

t=[0.25 0.5 1 1.5 2 3 4 6 8];
c=[19.21 18.15 15.36 14.10 12.89 9.32 7.45 5.24 3.01];
plot(t,c,'r*')
%pause代表停顿，用任意键使程序继续运行
pause
%semilogy代表使y轴为对数坐标轴
semilogy(t,c,'+')
pause 
plot(t,c,'r*')
%hold on指在绘制下一张图时保留前一张图
hold on
c0=23;
k=0.3;
c1=c0*exp(-k*t);
plot(t,c1,'b*')

三、三种插值

三种插值分别为最近临插值、线性插值、样条插值。以下为三种插值的实现和曲线拟合对比图。

x=[1 2 4 7 9 12 13 15 17];
f=[1.5 3.9 6 11.7 12.6 18.8 20.3 20.6 21.1];
axis([0 18 1  22])
xlabel('x')
ylabel('f')
y=1:0.1:17;
%代表绘图1，后续还有绘图2和3分别是产生新的绘图

figure(1)
plot(x,f,'o')
%在图中产生标记点位
gtext('已知数据点')
hold on;pause

%插值函数，总共提供四个参数，前两个参数为提供的数据点横纵坐标，第三个参数为实际曲线的
%所有横指标值，第四个参数为插值类型，nearest指最近临插值。
bb1=interp1(x,f,y,'nearest')
plot(y,bb1)
gtext('nearest')
hold on;pause
a=polyfit(x,f,3)
aa=polyval(a,y)
plot(y,aa)
gtext('曲线拟合')
hold off;pause

figure(2)

plot(x,f,'o')
gtext('已知数据点')
hold on;pause
bb1=interp1(x,f,y,'linear')
plot(y,bb1)
gtext('linear')
hold on;pause
a=polyfit(x,f,3)
aa=polyval(a,y)
plot(y,aa)
gtext('曲线拟合')
hold off;pause

figure(3)

plot(x,f,'o')
gtext('已知数据点')
hold on;pause
bb1=interp1(x,f,y,'spline')
plot(y,bb1)
gtext('spline')
hold on;pause
a=polyfit(x,f,3)
aa=polyval(a,y)
plot(y,aa)
gtext('曲线拟合')
hold off;pause

四、用matlab解拟合问题

1. 线性最小二乘拟合

例子一：

%多项式拟合方法
   x=0:0.1:1;
   y=[-0.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2];
   A=polyfit(x,y,2)
   z=polyval(A,x);
   plot(x,y,'k+',x,z,'r')
%超定方程解法
  x=0:0.1:1;
  y=[-0.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2];
%x'代表x的转置
  R=[(x.^2)' x' ones(11,1)];
  A=R\y'

2. 非线性最小二乘拟合

非线性最小二乘拟合有两种基本函数，以下用两个例子来分别介绍函数的实现。

(1).lsqcurvefit

%编写M文件 curvefun1.m
function f=curvefun1(x,tdata)
    f=x(1)+x(2)*exp(-0.02*x(3)*tdata)   
                   % x(1)=a; x(2)=b,x(3)=k;
%输入命令
clear
  tdata=100:100:1000
  cdata=1e03*[4.54,4.99,5.35,5.65,5.90,6.10,6.26,6.39,6.50,6.59];
  x0=[0.2,0.05,0.05];
  x=lsqcurvefit('curvefun1',x0,tdata,cdata)
  f=curvefun1(x,tdata)  

 figure
 plot(tdata,cdata,'r*')
 hold on
 plot(tdata,f,'bO')
 hold on
 t=100:0.1:1000;
 y=x(1)+x(2)*exp(-0.02*x(3)*t);  
 plot(t,y,'b')

(2).lsqnonlin

%编写M文件 curvefun2.m,cdata,tdata值提前写入function中
function f=curvefun2(x)
tdata=100:100:1000;
cdata=1e03*[4.54,4.99,5.35,5.65,5.90,6.10,6.26,6.39,6.50,6.59];
f=x(1)+x(2)*exp(-0.02*x(3)*tdata)-cdata


%键入命令
clear
x0=[0.2,0.05,0.05];
  x=lsqnonlin('curvefun2',x0)
  f= curvefun2(x)