一、插值
一维插值问题可描述为：已知函数在x0，x1，x2，…，xn处的值为y0，y1，…，yn，求简单函数p(x)，使p(x)=yi。
💥实际应用中不应使用七次以上的插值多项式（Runge现象）。
避免Runge现象常用方法是：将插值区间分成若干小区间，在小区间内用低次（二次、三次）插值，即分段低次插值，如样条函数插值。
1）Lagrange插值法
👇在MATLAB里面定义lagrange函数👇
function y=lagrange(x0,y0,x)
n=length(x0);
m=length(x);
for i=1:m
z=x(i);
s=0.0;
for k=1:n
p=1.0;
for j=1:n
if j~=k
p=p*(z-x0(j))/(x0(k)-x0(j));
end
end
s=p*y0(k)+s;
end
y(i)=s;
end
2）Newton插值法
MATLAB插值
一维插值interp1
yi = interp1 ( x , y , xi , ’ method ’ )
x,y为插值点，xi,yi为被插值和插值结果，x,y和xi,yi通常为向量；
’ method '表示插值方法
①‘nearest’——最邻近插值 ②’linear’‘——线性插值（默认） ③ ‘spline’——三次样条插值 ④’cubic’——立方插值
✨代码演示①
x=0:2:24;
y=[12 9 9 10 18 24 28 27 25 20 18 15 13];
x1=13;
y1=interp1(x,y,x1,'spline')
xi=0:1/3600:24;
yi=interp1(x,y,xi,'spline');
plot(x,y, '*',xi, yi)
✨代码演示②
新建脚本后输入如下代码
function plane
x0=[0 3 5 7 9 11 12 13 14 15];
y0=[0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6];
x=0:0.1:15;
y1=lagrange(x0,y0,x);
y2=interp1(x0,y0,x);
y3=interp1(x0,y0,x,'spline');
subplot(3,1,1);
plot(x0,y0,'k+',x,y1,'r');
grid;
title('lagrange');
subplot(3,1,2);
plot(x0,y0,'k+',x,y2,'r');
grid;
title('piecewise linear');
subplot(3,1,3);
plot(x0,y0,'k+',x,y3, 'r');
结果如下
二维插值interp2
zi = interp2 ( x , y , z, xi , yi , ’ method ’ )
x,y,z为插值点，z为被插值函数在(x,y)处的值；
xi,yi为被插值点，zi为输出的插值结果（插值函数在(xi,yi)处的值）
x,y为向量，xi,yi为向量或矩阵，z和zi为矩阵
①‘nearest’——最邻近插值 ②’linear’‘——双线性插值（默认） ③ ‘spline’——双三次样条插值 ④’cubic’——双立方插值
✨代码演示①
x=1:5;
y=1:3;
temps=[82 81 80 82 84;79 63 61 65 81;84 84 82 85 86];
figure(1);
mesh(x,y,temps);
xi=1:0.2:5;
yi=1:0.2:3;
zi=interp2(x,y,temps,xi,yi','cubic');
figure(2);
mesh(xi,yi,zi);
figure(3);
contour(xi,yi,zi,20,'r');
[i,j]=find(zi==min(min(zi)));
x=xi(j),y=yi(i),zmin=zi(i,j)
[i,j]=find(zi==max(max(zi)));
x=xi(j),y=yi(i),zmax=zi(i,j)
①plot3空间曲线②mesh空间曲面网格图③surf空间曲面表面图④contour等高线
📢contour(x,y,z,n)：作出由点(x,y,z)插值而成曲面的n条等高线。
📢meshc和surfc可在曲面下方画等高线；meshz和surfz画垂帘图。
✨应用：题例
在某山区测得一些地点的高程如下表。平面区域为0<x<5600, 0<y<4800
试用Matlab中的最邻近插值、双线性插值和双三次插值3种方法作出该山区的地貌图和等高线图，并求出最高和最低点。
✍MATLAB代码：
>> x=0:400:5600;
y=0:400:4800;
z=[370 470 550 600 670 690 670 620 580 450 400 300 100 150 250;...
510 620 730 800 850 870 850 780 720 650 500 200 300 350 320;...
650 760 880 970 1020 1050 1020 830 900 700 300 500 550 480 350;...
740 880 1080 1130 1250 1280 1230 1040 900 500 700 780 750 650 550;...
830 980 1180 1320 1450 1420 1400 1300 700 900 850 840 380 780 750;...
880 1060 1230 1390 1500 1500 1400 900 1100 1060 950 870 900 930 950;...
910 1090 1270 1500 1200 1100 1350 1450 1200 1150 1010 880 1000 1050 1100;...
950 1190 1370 1500 1200 1100 1550 1600 1550 1380 1070 900 1050 1150 1200;...
1430 1430 1460 1500 1550 1600 1550 1600 1600 1600 1550 1500 1500 1550 1550;...
1420 1430 1450 1480 1500 1550 1510 1430 1300 1200 980 850 750 550 500;...
1380 1410 1430 1450 1470 1320 1280 1200 1080 940 780 620 460 370 350;...
1370 1390 1410 1430 1440 1140 1110 1050 950 820 690 540 380 300 210;...
1350 1370 1390 1400 1410 960 940 880 800 690 570 430 290 210 150];
图①
figure(1);
meshz (x,y,z);
xlabel('X'), ylabel ('Y'), zlabel ('Z');
图②
[xi,yi]=meshgrid(0:50:5600,0:50:4800);
figure(2);
z1i=interp2(x,y,z,xi,yi,'nearest');
surfc(xi,yi,z1i);
xlabel('X'), ylabel ('Y'), zlabel ('Z');
图③
figure (3);
z2i=interp2(x,y,z,xi,yi);
surfc(xi,yi,z2i) ;
xlabel('X'), ylabel ('Y'), zlabel ('Z');
图④
figure(4);
z3i=interp2(x,y,z,xi,yi,'cubic');
surfc(xi,yi,z3i);
xlabel('X'), ylabel ('Y'), zlabel ('Z');
图⑤
figure(5);
subplot(1,3,1),contour(xi,yi,z1i,10,'r');
subplot(1,3,2),contour(xi,yi,z2i,10,'r');
subplot(1,3,3),contour(xi,yi,z3i,10,'r');
散乱点插值griddata
griddata ( x , y , z , xi , yi , ’ method ’ )——(x,y)为散乱点
✨应用：题例
在某海域测得一些点(x,y)处的水深z如下表，船的吃水深度为5英尺，在矩形区域(75,200)*(-50,150)内的哪些地方船要避免进入。
✍MATLAB代码：
x=[129 140 103.5 88 185.5 195 105.5 157.5 107.5 77 81 162 162 117.5];
y=[7.5 141.5 23 147 22.5 137.5 85.5 -6.5 -81 3 56.5 -66.5 84 -33.5];
z=[-4 -8 -6 -8 -6 -8 -8 -9 -9 -8 -8 -9 -4 -9];
[xi,yi]=meshgrid(75:0.5:200,-70:0.5:150);
zi=griddata(x,y,z,xi,yi,'cubic');
figure(1);
meshz(xi,yi,zi);
xlabel('X'),ylabel('Y'),zlabel('Z');
figure(2), contour(xi,yi,zi,[-5 -5], 'b');
grid;
hold on;
plot(x,y,'+');
xlabel('X'),ylabel('Y');
在等高线图中，-5等高线圈内易搁浅，圈外不易搁浅。
二、拟合
曲线拟合需要解决两个问题：
（1）线型选择
（2）线型参数计算
线型拟合参数计算采用最小二乘法；非线性则使用Gauss-Newton迭代法。
MATLAB拟合
1）多项式拟合
格式：[ a , S ] = polyfit ( x , y , n )
x,y是被拟合数据的自变量和因变量；n为拟合多项式的次数
a为拟合多项式系数构成的向量
S为分析拟合效果所需的指标（可省略）
✨代码演示
x=1:12;
y=[5, 8, 9, 15, 25, 29, 31, 30, 22, 25, 27, 24];
a=polyfit(x,y,9)
xp=1:0.1:12;
yp=polyval(a,xp);
plot(x,y,'.k' ,xp,yp,'r');
返回结果
a =
0.0001 -0.0044 0.1107 -1.5423 13.1415 -70.6122 237.0006 -470.8685 492.7583 -195.0000
2）非线性拟合
格式：[ b , r ] = polyfit ( x , y , fun , b0 , option )
x,y是被拟合数据的自变量和因变量；fun为拟合函；b0为拟合参数的初始迭代值；option为拟合选项
b拟合参数；r为拟合残差
✨代码演示
x=1:16;
y=[4.00 6.40 8.00 8.80 9.22 9.50 9.70 9.86 10.00 10.20 10.32 10.42 10.50 10.55 10.58 10.60];
y1=@(b,t)b(1)*exp(-t/b(2))+b(3)*exp(-t/b(4))+b(5);
b0=[-1 1 -1 1 1];
a=nlinfit(x,y,y1,b0)
xp=1:0.1:16;
yp=y1(a,xp);
plot(x,y,'.k',xp,yp,'r');
返回结果
a =
-9.3518 1.5112 -2.7661 11.2383 11.3240
MATLAB拟合工具箱
在命令窗口键入cftool启动拟合工具箱
❗❗拟合问题与插值问题的区别
1️⃣插值函数过已知点，拟合函数不一定过已知点。（本质区别）
2️⃣插值主要求函数值，拟合主要求的是函数关系（从而进行预测等进一步分析）
1212
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
