插值与拟合

插值算法

数模比赛中，常常需要根据已知函数点进行数据、模型的处理，but现有数据少得可怜，这时就要模拟产生一些靠谱的值来满足需求，这就是插值的作用。

一维插值

即平面中的点，只有$x$和$y$

$y=f(x)$在区间$[a,b]$上有定义，且在已知点$a\le x_0\le x_1\le ...\le b$上的值分别为：$y_0,y_1,...,y_n$

若存在一个简单函数$P(x)$，使$P(x_i)=y_i$ $(i-0,1,2,...,n)$ $(1,1)$

则称$P(x)$为$f(x)$的插值函数，点$x_0,x_1,...,x_n$称为插值节点，包含插值节点的区间$[a,b]$称为插值区间，求插值函数$P(x)$的方法称之为插值法。

$P(x)$不唯一：

多项式插值：若$P(x)$是次数不超过$n$的代数多项式，即$P(x)=a_0+a_{1}x+...+a_n{x}^n$，称为多项式插值

分段插值：若$P(x)$为分段多项式，就称为分段插值

三角插值：若$P(x)$为三角多项式，就称为三角插值

多项式插值

定理:

设有$n+1$个互不相同的节点$(x_i,y_i)$ $(i=0,\mathrm{1...},n)$,则存在唯一的多项式$L(x)=a_0+a_{1}x+a_2{x}^2+...+a_n{x}^n$使得$L_n(x_j)=y_j$ $(j=0,1,\mathrm{2...},n)$

证明:

构造方程组：$x_i=\{\begin{array}{l}{a}_0+a_1{x}_0+a_2{x}_0^2+...+a_n{x}_0^n=y_0\\ a_0+a_1{x}_1+a_2{x}_1^2+...+a_n{x}_1^n=y_1\\ ....\\ a_0+a_1{x}_n+a_2{x}_n^2+...+a_n{x}_n^n=y_n\end{array}$

线性代数解方程： A = \left[ \matrix{ 1 & x_0 & \cdots & x_0^n \\ 1 & x_1 & \cdots & x_1^n \\ \dots \\ 1 & x_n & \cdots & x_n^n } \right] $$\begin{gathered} X=\left[\text{matrix}{a}_0 \\ {a}_1 \\ ... \\ {a}_n\right] \end{gathered}$$ $$\begin{gathered} Y=\left[\text{matrix}{y}_0 \\ {y}_1 \\ ... \\ {y}_n\right] \end{gathered}$$

方程组的矩阵形式如下：$AX=Y$

由于$|A|={\prod }_{i=1}^n{\prod }_{j=0}^{n-1}(x_i=x_j)\ne 0$所以方程组有唯一解

从而$L_n(x)=a_0+a_{x}x+a_2{x}^2+...+a_n{x}^n$唯一存在

ps：

• 只要$n+1$个节点互异，满足上述插值条件的多项式唯一存在
• 如果不限制多项式次数，插值多项式并不唯一

拉格朗日插值法

• 两个点$(x_0,y_0)$&$(x_1,y_1)$：$f(x)=\frac{x-x_1}{x_0-x_1}{y}_0+\frac{x-x_0}{x_1-x_0}{y}_1$

• 三个点：$(x_0,y_0)$&$(x_1,y_1)$&$(x_2,y_2)$:

  $$f(x)=\frac{(x-x_1)(x-x_2)}{(x_0-x_1)(x_0-x_2)}{y}_0+\frac{(x-x_0)(x-x_2)}{(x_1-x_0)(x_1-x_2)}{y}_1+\frac{(x-x_0)(x-x_1)}{(x_2-x_0)(x_2-x_1)}{y}_2$$

• n个点：$L_n(x)={\sum }_{k=0}^n{y}_k\frac{{\omega }_{n+1}(x)}{(x-x_k){\omega }_{n+1}^{{}^{\prime }}(x_k)}$

  其中 ${\omega }_{n+1}(x)=(x-x_0)(x-x_1)...(x-x_n)$ & ${\omega }_{n+1}^{{}^{\prime }}(x_k)=(x_k-x_0)...(x_k-x_{k-1})(x_k-xk+1)...(x_k-x_n)$

存在问题：龙格现象

高次插值会产生龙格现象，即在两端波动较大，产生明显的震荡。

在不熟悉曲线运动趋势的情况下，不要轻易的使用高次插值

如何提高插值精度呢 ⬇

分段插值

• 线性插值：找所求点最近的两点构成直线在此点对应函数值
• 二次插值：取距$x$最近的三点$x_{i-1},x_i,x_{i+1}$进行二次插值，即在$[x_{i-1},x_{i+1}]$上取$f(x)\approx L_2(x)=su{m}_{k=i-1}^{i+1}[y_k{\prod }_{j=i-1,j\ne k}^{i+1}\frac{x_i-x_j}{x_k-x_i}]$

牛顿插值法

$$\begin{array}{rl}f(x)=& f(x_0)+f[x_0,x_1](x-x_0)\\ & +f[x_0,x_1,x_2](x-x_0)(x-x_1)+...\\ & +f[x_0,x_1,...,x_n](x-x_0)(x-x_1)...(x-x_n)\end{array}$$

差商：$f[x_0,x_k]=\frac{f(x_k)-f(x_0)}{x_k-x_0}$为函数$f(x)$的一阶差商

$f[x_0,x_1,x_2]=\frac{f[x_1,x_2]-f[x_0,x_1]}{x_2-x_0}$为函数的二阶差商

同理，$f[x_0,...,x_k]=\frac{f[x_1,...,x_k]-f[x_0,...x_{k-1}]}{x_k-x_0}$为函数的k阶差商

与拉格朗日插值法相比，此方法存在继承性。但是也存在龙格现象。

继承性:

每次插值只和前n项有关，只要添加新项，就可以产生新函数。

埃尔米特插值

一般插值问题的要求：$\phi (x_i)=y_i$

较高要求：${\phi }^{{}^{\prime }}(x_i)=y_i^{{}^{\prime }}$

满足节点函数值相等，且对应导数甚至高阶导数也相等的插值多项式就是埃尔米特插值多项式，它更能反映被插值函数的性态

简略版原理：

设函数$f(x)$在区间[a,b]有$n+1$个互异的节点 $a=x_0<x_1<...<x_n<b$,定义在[a,b]上函数$f(x)$在节点上满足：

​ $f(x_i)=y_i,f^{{}^{\prime }}(x_i)=y_i^{{}^{\prime }}(i=0,1,2,...,n)$可唯一确定一个次数不超过$2n+1$的多项式$H_{2n+1}(x)=H(x)$满足：

$$H(x_j)=y_j,h^{{}^{\prime }}(x_j)=m_j(j=0,1,...,n)$$

其余项为：

$$R(x)=f(x)-H(x)=\frac{f^{2n+2}(\xi )}{（2n+2)!}{\omega }_{2n+2}(x)$$

在实际应用中，往往使用分段三次Hermite插值多项式

% 分段三次埃尔米特插值
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青

样条插值

Matlab有内置函数：

p = spline(x,y,new_x);  %三次样条插值

设$y=f(x)$在点$x_0,x_1,...x_n$的值为$y_0,y_1,...,y_n$若函数$S(x)$满足下列条件：

• $s(x_i)=f(x_i)=y_i$, $i=0,1,2,...,n$

• 在每个子区间$[x_i,x_{i+1}]$上$S(x)$是三次多项式

• $S(x)$在$[a,b]$上二阶连续可微。

三次样条插值和分段三次埃尔米特插值的对比

 
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维数据插值(了解)

p = interpn(x1,x2,...,xn,y,new_x1,new_x2,...,new_xn,method)

method分为：

• ‘linear’线性插值(默认)
• ‘cubic’三次插值
• ‘spline’三次样条插值(最精准)
• ‘nearest’最邻近插值算法

不恰当的例子：

% 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-')

拟合算法

与插值问题不同，在拟合问题中不需要曲线一定经过给定的点。拟合问题的目标是寻求一个函数(曲线)，使得曲线在某种准则下与所有数据点最为接近，即曲线拟合最好(最小化损失函数)

插值与拟合的区别

样本点过多由于次数过高，造成龙格现象。

尽管我们可以通过分段来避免此现象，但我们更多的时候更倾向于得到一个确定的曲线，即使其不能经过每个样本点，只要误差足够小即可。$拟合的结果是得到一个确定的曲线$

最小二乘法

线性代数讲过了QAQ

$$\begin{gathered} \text{DeclareMathOperator}\ast \mathrm{argmin}argmin\hat{k},\hat{b}=\underset{k,b}{\mathrm{argmin}}(\sum _{i=1}^n(y_i-\widehat{y_i}{)}^2) \\ \hat{k}=\frac{n\sum _{i=1}^n{x}_i{y}_i-\sum _{i=1}^n{y}_i\sum _{i=1}^n{x}_i}{n\sum _{i=1}^n{x}_i^2-\sum _{i=1}^n{x}_i\sum _{i=1}^n{x}_i} \\ \hat{b}=\frac{\sum _{i=1}^n{x}_i^2\sum _{i=1}^n{y}_i-\sum _{i=1}^n{x}_i\sum _{i=1}^n{x}_i{y}_i}{n\sum _{i=1}^n{x}_i^2-\sum _{i=1}^n{x}_i\sum _{i=1}^n{x}_i} \end{gathered}$$

评价拟合好坏

拟合优度(可决系数)$R^2$

总体平方和SST(Total sum of squares)：$SST={\sum }_{i=1}^n(y_i-\overline{y}{)}^2$

误差平方和SSE(The sum of squares due to error)：$SSE={\sum }_{i=1}^n(y_i-\widehat{y_i}{)}^2$

回归平方和SSR(Sum of squares of the regression)：$SSR={\sum }_{i=1}^n(\widehat{y_i}-\overline{y}{)}^2$

证明SST=SSE+SSR

$$\begin{gathered} {\sum }_{i=1}^n(y_i-\overline{y}{)}^2=SST \\ ={\sum }_{i=1}^n(y_i-\widehat{y_i}+\widehat{y_i}-\overline{y}{)}^2 \\ ={\sum }_{i=1}^n(y_i-\widehat{y_i}{)}^2+{\sum }_{i=1}^n(\widehat{y_i}-\overline{y}{)}^2+2{\sum }_{i=1}^n(y_i-\widehat{y_i})(\widehat{y_i}-\overline{y}) \end{gathered}$$

$$\begin{gathered} 2\sum _{i=1}^n(y_i-\widehat{y_i})(\widehat{y_i}-\overline{y}) \\ =\sum _{i=1}^n\widehat{y_i}(y_i-\widehat{y_i})-\overline{y}\sum _{i=1}^n(y_i-\widehat{y_i}) \\ =\sum _{i=1}^n(k{x}_i+b)(y_i-\widehat{y_i})-0=0 \end{gathered}$$

拟合优度：$0\le R^2=\frac{SSR}{SST}=\frac{SST-SSE}{SST}=1-\frac{SSE}{SST}\le 1$

$R^2$越接近1，说明误差平方和越接近0，拟合得越好(仅能评价线性函数)，与其他函数比较时直接看SSE即可。