插值

一维数据插值

当我们有一些离散的数据点，但需要在这些数据点之间进行估计或预测时，插值法就可以帮助我们。它能够通过已知的数据点来填补空缺的数据点，使我们能够得到更完整的数据。想象一下你手上有几个散落的点，但你想要画一条平滑的曲线来穿过这些点。插值法就像是使用橡皮筋，你把它放在这些点上，然后让它自然地弯曲，最终形成一条平滑的曲线。这条曲线能够逼近原始数据点，同时也能够在数据点之间进行预测。插值法就是通过已知的数据点来预测或填补未知的数据点，让我们能够获得更完整、更平滑的数据集。

插值方法有很多种，下面列举几种常见的插值方法：

线性插值：线性插值是最简单的插值方法之一。它假设两个数据点之间的变化是线性的，并通过连接这些数据点来估计中间的未知数据点。
分段线性插值：将数据点划分为若干段，每段使用线性插值来估计未知数据点。
拉格朗日插值：使用拉格朗日多项式来逼近数据点。它通过在每个数据点上构建一个多项式，并将这些多项式相加来获得插值函数。拉格朗日插值可以准确地经过所有的已知数据点，但在计算复杂度上可能比较高。

实际使用中，存在一个问题。当使用一些简单的插值方法（比如上述的线性插值或拉格朗日插值）时，如果函数在某些地方很快地变化，或者有剧烈的震荡，那么插值结果可能不准确，还会出现明显的振荡效应。

这种现象被称为龙格现象，是插值法存在的主要问题。举例说明，假设我们想要使用插值方法来逼近一个函数，在两个已知的数据点之间进行插值，如果这个函数在两个数据点之间有一个非常陡峭的峰或谷，那么直接使用简单的插值方法很难直接预测，可能会导致插值结果出现误差和振荡效应。这种情况在越高次的插值中越容易出现。存在两个问题：插值多项式次数越高，精度未必显著提升，而误差可能显著增大。

可以采用分段低次插值。将整个函数分成若干小区间，每个区间内采用低次插值方法（比如线性插值或二次插值）来逼近函数，可以在一定程度上避免高次插值带来的龙格现象。但是，这一方法也不能全面反映被插入函数的特征。解决方案是埃尔米特插值，它保证该点的导数也相等，模拟了形态特征。

综合多种方案的优缺点，我们一般采用以下方法：

分段三次埃尔米特插值 ：不仅考虑数据点之间的线性变化，还考虑了斜率的变化。通过使用三次多项式来逼近每一段数据，可以得到更平滑的插值结果。
三次样条插值：更好的光滑插值方法。它将数据点之间的曲线定义为一系列的三次多项式段，并要求这些段在连接点上具有连续的一阶和二阶导数。这样可以得到更平滑、更自然的插值曲线。

至于具体的数学证明，很显然，我是完全不会的。matlab内置了相关函数，我们可以直接调用。

这里我们举一个例子，对正弦函数进行分段三次埃尔米特插值和三次样条插值，并绘制出原始数据和插值结果的曲线。

% 分段三次埃尔米特插值
x = -pi:pi; y = sin(x);
new_x = -pi:0.1:pi;
p = pchip(x,y,new_x);
figure(1); % 在同一个脚本文件里面，要想画多个图，需要给每个图编号
plot(x, y, 'o', new_x, p, 'r-')
  
% plot函数用法:
% plot(x1,y1,x2,y2)
% 线方式： - 实线  :点线    -. 虚点    - - 波折线
% 点方式： . 圆点  +加号  * 星号  x x形  o 小圆
% 颜色： y黄； r红； g绿； b蓝； w白； k黑； m紫； c青
  
% 三次样条插值和分段三次埃尔米特插值的对比
x = -pi:pi;
y = sin(x);
new_x = -pi:0.1:pi;
p1 = pchip(x,y,new_x);   %分段三次埃尔米特插值
p2 = spline(x,y,new_x);  %三次样条插值
figure(2);
plot(x,y,'o',new_x,p1,'r-',new_x,p2,'b-')
legend('样本点','三次埃尔米特插值','三次样条插值','Location','SouthEast')   %标注显示在东南方向
% LEGEND(string1,string2,string3, …)
% 分别将字符串1、字符串2、字符串3……标注到图中，每个字符串对应的图标为画图时的图标。
% ‘Location’用来指定标注显示的位置

所以我们现在的任务是整理数据，把有用的部分从完整的数据集里取出，单独组合在一起，这样有利于更方便地训练我们的模型。

n维数据插值

上面是一维数据的插值，同样的道理，我们也可以对n维数据进行处理。

在一维插值中，我们只有一个自变量和一个因变量，例如在一条直线上对数据进行插值。而在二维插值中，我们有两个自变量和一个因变量，例如在一个平面上对数据进行插值。而在更高维度的情况下，我们可能有更多的自变量。在多维数据中，我们需要使用多个坐标轴来描述每个数据点的位置。例如，在三维数据中，我们可以使用(x, y, z)来表示每个数据点的位置，其中x、y、z分别表示沿着三个不同方向的坐标值。通过在不同方向上进行插值计算，我们可以推断出缺失数据点的值，进而对整个多维空间进行预测或填充。

对于n维数据插值，优先使用三次样条插值。我也不知道怎么证明，反正用这个效果就是好。

% n维数据的插值
x = -pi:pi; y = sin(x);
new_x = -pi:0.1:pi;
p = interpn (x, y, new_x, 'spline');
% 等价于 p = spline(x, y, new_x)，优先使用三次样条插值。
figure(3);
plot(x, y, 'o', new_x, p, 'r-')

插值法一般不用于预测。如果真的要用它预测，你多半也就完蛋了。

拟合

拟合的实现

拟合的思想就像是给你一些散落的数据点，你需要找到一条曲线或者函数，尽可能地贴近这些数据点，同时还要能够反映出数据的整体趋势。

可以想象一下，面前有一堆螺丝钉，你想要知道它们的整体形状。拟合就像是你拿一个橡皮筋，和插值不同，我们不直接放在螺丝钉上，而是使得它能够尽量贴合这些螺丝钉的位置。你会不断调整线的形状，直到你觉得它与螺丝钉的位置最接近，并且能够反映出整体的形状。如果使用不同的函数形式来进行拟合。通过寻找最适合的函数形式和参数，我们就可以找到一个拟合曲线或函数，使其与已知数据点尽可能接近，并且能够描述数据的总体趋势。尽管不能保证它通过每一个点，但是误差足够小就可以满足要求，得出一个确定的曲线。即使有一些数据点没有给出，也可以根据拟合曲线来预测这些缺失的数据点的位置。

因此我们可以认为，当数据点比较多、分布比较连续，且需要反映数据整体趋势时，拟合是比较合适的选择。

下面给出代码，使用最小二乘法找出x与y之间的拟合曲线：

clear;clc
load  data1
plot(x,y,'o')
% 给x和y轴加上标签
xlabel('x的值')
ylabel('y的值')
n = size(x,1);
k = (n*sum(x.*y)-sum(x)*sum(y))/(n*sum(x.*x)-sum(x)*sum(x))
b = (sum(x.*x)*sum(y)-sum(x)*sum(x.*y))/(n*sum(x.*x)-sum(x)*sum(x)) %最小二乘法
hold on % 继续在之前的图形上来画图形
grid on % 显示网格线

% 画出y=kx+b的函数图像 plot(x,y)
% 传统的画法：模拟生成x和y的序列，比如要画出[0,5]上的图形

xx = 2.5: 0.1 :7  % 间隔设置的越小画出来的图形越准确
yy = k * xx + b  % k和b都是上面已经算出来的已知值
plot(xx,yy,'-')

% 匿名函数的基本用法。
% handle = @(arglist) anonymous_function
% 其中handle为调用匿名函数时使用的名字。
% arglist为匿名函数的输入参数，可以是一个，也可以是多个，用逗号分隔。
% anonymous_function为匿名函数的表达式。
%  z=@(x,y) x^2+y^2; 
%  z(1,2) 
% % ans =  5
% fplot(f,xinterval) 将匿名函数f在指定区间xinterval绘图。xinterval =  [xmin xmax] 表示定义域的范围

f=@(x) k*x + b;
fplot(f,[2.5,7]);
legend('样本数据','拟合函数','location','SouthEast')

评价拟合质量

拟合优度 R^2，也称为决定系数，是用来衡量回归模型对观测数据的解释能力的统计指标，只适用于线性函数。它的取值范围在0到1之间，越接近1，说明误差平方和越接近0，也就意味着拟合效果越好。

R^2 的计算方式是通过回归平方和 (SSR) 除以总体平方和 (SST) 得到的。具体公式为：R^2 = SSR / SST。

回归平方和（SSR）表示因变量的变异程度中可以被回归模型解释的部分的总和。计算方式是将每个预测值与因变量均值的差的平方求和。

误差平方和（SSE）表示因变量的变异程度中不能被回归模型解释的部分的总和。计算方式是将每个观测值与对应的预测值的差的平方求和，是另外一个重要指标，越小越好。

总体平方和（SST）表示因变量（被解释变量）的总方差。计算方式是将每个观测值与因变量均值的差的平方求和。

具体原理不赘述，只需要套公式实现代码就可以了。

代码如下：

y_hat = k*x+b; % y的拟合值
SSR = sum((y_hat-mean(y)).^2)  % 回归平方和
SSE = sum((y_hat-y).^2) % 误差平方和
SST = sum((y-mean(y)).^2) % 总体平方和
SST-SSE-SSR   % 5.6843e-14  =   5.6843*10^-14   matlab浮点数计算的一个误差
R_2 = SSR / SST

曲线拟合工具箱

进入曲线拟合工具箱界面“Curve Fitting tool”

点击“Data”按钮，弹出“Data”窗口；
利用X data和Y data的下拉菜单读入数据x,y，可修改数据集名“Data set name”，然后点击“Create data set”按钮，退出“Data”窗口，返回工具箱界面，这时会自动画出数据集的曲线图；
点击“Fitting”按钮，弹出“Fitting”窗口；
点击“New fit”按钮，可修改拟合项目名称“Fit name”，通过“Data set”下拉菜单选择数据集，然后通过下拉菜单“Type of fit”选择拟合曲线的类型;
类型设置完成后，点击“Apply”按钮，就可以在Results框中得到拟合结果。

工具箱提供的拟合类型有：

Custom Equations：用户自定义的函数类型
Exponential：指数逼近，有2种类型， aexp(bx) 、 aexp(bx) + cexp(dx)
Fourier：傅立叶逼近，有7种类型，基础型是 a0 + a1cos(xw) + b1sin(xw)
Gaussian：高斯逼近，有8种类型，基础型是 a1*exp(-((x-b1)/c1)^2)
Interpolant：插值逼近，有4种类型，不赘述。
Polynomial：多形式逼近，有9种类型。
Power：幂逼近，有2种类型，ax^b 、ax^b + c
Rational：有理数逼近
Smoothing Spline：平滑逼近
Sum of Sin Functions：正弦曲线逼近，有8种类型，基础型是 a1sin(b1x + c1)
Weibull：只有一种，abx^(b-1)exp(-ax^b)

操作结束后，会在工具箱窗口中显示拟合曲线。如果拟合效果不好，还可以在“Fitting”窗口点击“New fit”按钮，按照上述步骤4、5进行新的拟合。

左上角“文件”，选择“generate code”导出代码。这里给出样例，数据是乱编的。我这辈子都不想再看见Fourier了。

clear;clc
year = 1790:10:2000;
population = [3.9,5.3,7.2,9.6,12.9,17.1,23.2,31.4,38.6,50.2,62.9,76.0,92.0,106.5,123.2,131.7,150.7,179.3,204.0,226.5,251.4,281.4];
plot(year,population,'o')
  
[fitresult, gof] = createFit(year, population) % 这个就是从工具箱导出的函数
  
t = 2001:2030;
xm = 342.4;  
r =  0.02735; % 这些数据是工具箱拟合后得到的
predictions = xm./(1+(xm./3.9-1).*exp(-r.*(t-1790)));  % 计算预测值（注意这里要写成点乘和点除,这样可以保证按照对应元素进行计算）
figure(3)
plot(year,population,'o',t,predictions,'.')  % 绘制预测结果图