MATLAB和Python数据插补

拉格朗日插值

拉格朗日插值多项式的构造

在拉格朗日插值中，构造了$n+1$个多项式${\mathbf{L}}_0\left(\mathbf{x}\right),{\mathbf{L}}_1\left(\mathbf{x}\right),\cdots ,{\mathbf{L}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 。 每个多项式${\mathbf{L}}_{\mathbf{j}}\left(\mathbf{x}\right),\left(\mathbf{j}=0,\cdots ,\mathbf{n}\right)$ 的次数为$n$，使得

$${\mathbf{L}}_{\mathbf{j}}\left({\mathbf{x}}_{\mathbf{j}}\right)=\{\begin{array}{l}1,\mathbf{i}=\mathbf{j}\\ 0,\mathbf{i}\ne \mathbf{j}\end{array}$$

然后，插值多项式${\mathbf{P}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 由以下公式定义：

$${\mathbf{P}}_{\mathbf{n}}\left(\mathbf{x}\right)={\mathbf{L}}_0\left(\mathbf{x}\right){\mathbf{y}}_0+{\mathbf{L}}_1\left(\mathbf{x}\right){\mathbf{y}}_1+\cdots +{\mathbf{L}}_{\mathbf{n}}\left(\mathbf{x}\right){\mathbf{y}}_{\mathbf{n}}=\sum _{\mathbf{i}=0}^{\mathbf{n}}{\mathbf{L}}_{\mathbf{i}}\left(\mathbf{x}\right){\mathbf{y}}_{\mathbf{i}}$$

多项式${\mathbf{P}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 的次数为n，因为对于所有$\mathbf{i}=0,\cdots ,\mathbf{n}$，每个${\mathbf{L}}_{\mathbf{i}}\left(\mathbf{x}\right)$ 都是n次多项式。 此外，${\mathbf{P}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 满足：

$${\mathbf{P}}_{\mathbf{n}}\left({\mathbf{x}}_{\mathbf{j}}\right)=\sum _{\mathbf{i}=0}^{\mathbf{n}}{\mathbf{L}}_{\mathbf{i}}\left({\mathbf{x}}_{\mathbf{j}}\right){\mathbf{y}}_{\mathbf{i}}={\mathbf{L}}_{\mathbf{j}}\left({\mathbf{x}}_{\mathbf{j}}\right){\mathbf{y}}_{\mathbf{j}}=1\cdot {\mathbf{y}}_{\mathbf{j}}={\mathbf{y}}_{\mathbf{j}},\mathbf{j}=0,\cdots ,\mathbf{n}$$

这证明${\mathbf{P}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 内插给定数据。

为了构造$n+1$个多项式${\mathbf{L}}_{\mathbf{i}}\left(\mathbf{x}\right)$ ，$\mathbf{i}=0,\cdots ,\mathbf{n}$，重要的是要注意

$$\mathbf{\psi }\left(\mathbf{x}\right)=\left(\mathbf{x}-{\mathbf{x}}_0\right)\left(\mathbf{x}-{\mathbf{x}}_1\right)\cdots \left(\mathbf{x}-{\mathbf{x}}_{\mathbf{n}}\right)$$

的${\mathbf{x}}_0,{\mathbf{x}}_1,\cdots ,{\mathbf{x}}_{\mathbf{n}}$，如果$\tilde{\mathbf{L}}\left(\mathbf{x}\right)$ 定义为

$$\tilde{\mathbf{L}}\left(\mathbf{x}\right)=\frac{\mathbf{\psi }\left(\mathbf{x}\right)}{\mathbf{x}-{\mathbf{x}}_{\mathbf{i}}}=\left(\mathbf{x}-{\mathbf{x}}_0\right)\left(\mathbf{x}-{\mathbf{x}}_1\right)\cdots \left(\mathbf{x}-{\mathbf{x}}_{\mathbf{i}-1}\right)\left(\mathbf{x}-{\mathbf{x}}_{\mathbf{i}+1}\right)\cdots \left(\mathbf{x}-{\mathbf{x}}_{\mathbf{n}}\right)$$

则$x_i$不是$\tilde{\mathbf{L}}\left(\mathbf{x}\right)$ 的根。

现在，可以使用以下公式从$\tilde{\mathbf{L}}\left(\mathbf{x}\right)$ 获得${\mathbf{L}}_{\mathbf{i}}\left(\mathbf{x}\right)$ ：

$${\mathbf{L}}_{\mathbf{i}}\left(\mathbf{x}\right)=\frac{{\tilde{\mathbf{L}}}_{\mathbf{i}}\left(\mathbf{x}\right)}{{\tilde{\mathbf{L}}}_{\mathbf{i}}\left({\mathbf{x}}_{\mathbf{i}}\right)}=\frac{\left(\mathbf{x}-{\mathbf{x}}_0\right)\left(\mathbf{x}-{\mathbf{x}}_1\right)\cdots \left(\mathbf{x}-{\mathbf{x}}_{\mathbf{i}-1}\right)\left(\mathbf{x}-{\mathbf{x}}_{\mathbf{i}+1}\right)\cdots \left(\mathbf{x}-{\mathbf{x}}_{\mathbf{n}}\right)}{\left({\mathbf{x}}_{\mathbf{i}}-{\mathbf{x}}_0\right)\left({\mathbf{x}}_{\mathbf{i}}-{\mathbf{x}}_1\right)\cdots \left({\mathbf{x}}_{\mathbf{i}}-{\mathbf{x}}_{\mathbf{i}-1}\right)\left({\mathbf{x}}_{\mathbf{i}}-{\mathbf{x}}_{\mathbf{i}+1}\right)\cdots \left({\mathbf{x}}_{\mathbf{i}}-{\mathbf{x}}_{\mathbf{n}}\right)}$$

多项式${\mathbf{L}}_{\mathbf{i}}\left(\mathbf{x}\right)$ 也可以写成：

$${\mathbf{L}}_{\mathbf{i}}\left(\mathbf{x}\right)=\prod _{\mathbf{j}=0}^{\mathbf{n}}\frac{\mathbf{x}-{\mathbf{x}}_{\mathbf{j}}}{{\mathbf{x}}_{\mathbf{i}}-{\mathbf{x}}_{\mathbf{j}}}$$

拉格朗日多项式${\mathbf{P}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 可以写成：

$${\mathbf{P}}_{\mathbf{n}}\left(\mathbf{x}\right)=\sum _{\mathbf{i}=0}^{\mathbf{n}}\prod _{\mathbf{j}=0}^{\mathbf{n}}\frac{\mathbf{x}-{\mathbf{x}}_{\mathbf{j}}}{{\mathbf{x}}_{\mathbf{i}}-{\mathbf{x}}_{\mathbf{j}}}{\mathbf{y}}_{\mathbf{i}}$$

下面的MATLAB函数LagrangeInterp.m使用给定的数据点是两个向量x和y，每个向量的长度为n，构造拉格朗日插值多项式${\mathbf{P}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 。 它在$x=t$处评估${\mathbf{P}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 并返回插值结果。

function p = LagrangeInterp(t, x, y)
 n = length(x) ;
 p = 0 ;
 for i = 1 : n
 s = 1 ;
 for j = 1 : n
 if j , i
 s = s*(t-x(j))/(x(i)-x(j)) ;
 else
 continue
 end
 end
 p = p + s*y(i) ;
 end

执行上面的MATLAB代码，我们获得：

>> x = 0:pi/11:pi ;
>> y = sin(x) ;
>> p = LagrangeInterp(pi/6, x, y)
p =
0.499999999914085

Python脚本LagInterp.py使用以下代码实现拉格朗日插值多项式：

执行该代码，显示以下结果：

runfile(’D:/PyFiles/LagInterp.py’, wdir=’D:/PyFiles’)
0.4999999997868132

拉格朗日插值多项式的唯一性

在本节中，将使用代数的基本定理证明拉格朗日插值多项式是唯一的。 它指出，任何具有复数系数的n阶的多项式最多具有n个零。 因此，如果n阶的多项式$\mathbf{p}\left(\mathbf{x}\right)$ 的根数大于n，则它必须是零多项式$\mathbf{p}\left(\mathbf{x}\right)=0,{\forall }_{\mathbf{x}}\in \mathbb{C}$。

为了证明拉格朗日插值多项式的唯一性，假设${\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 和${\mathbf{q}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 是两个n阶多项式，它们对给定数据$\left({\mathbf{x}}_{\mathbf{i}},{\mathbf{y}}_{\mathbf{i}}\right),\mathbf{i}=0,\cdots ,\mathbf{n}$进行插值。 然后，对于所有$\mathbf{i}=0,\cdots ,\mathbf{n}$，${\mathbf{p}}_{\mathbf{n}}\left({\mathbf{x}}_{\mathbf{i}}\right)={\mathbf{q}}_{\mathbf{n}}\left({\mathbf{x}}_{\mathbf{i}}\right)={\mathbf{y}}_{\mathbf{i}}$。对于所有$\mathbf{x}\in \left[{\mathbf{x}}_0,{\mathbf{x}}_{\mathbf{n}}\right]$ ，令$\mathbf{r}\left(\mathbf{x}\right)={\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right)-{\mathbf{q}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 。由于${\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 和${\mathbf{q}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 都均为n阶，因此$\mathbf{r}\left(\mathbf{x}\right)$ 也是。现在，对于$\mathbf{i}=0,\cdots ,\mathbf{n}$，$\mathbf{r}\left({\mathbf{x}}_{\mathbf{i}}\right)={\mathbf{p}}_{\mathbf{n}}\left({\mathbf{x}}_{\mathbf{i}}\right)-{\mathbf{q}}_{\mathbf{n}}\left({\mathbf{x}}_{\mathbf{i}}\right)=0$。 即$\mathbf{r}\left(\mathbf{x}\right)$ 在$\left[{\mathbf{x}}_0,{\mathbf{x}}_{\mathbf{n}}\right]$ 中具有$n+1$个根，而其阶为$n$。 根据代数的基本定理，除非$\mathbf{r}\left(\mathbf{x}\right)$ 是$\left[{\mathbf{x}}_0,{\mathbf{x}}_{\mathbf{n}}\right]$ 中的零多项式，即$\mathbf{r}\left(\mathbf{x}\right)=0,\mathbf{x}\in \left[{\mathbf{x}}_0,{\mathbf{x}}_1\right]$ ，否则这不可能发生。 这证明对于$\mathbf{x}\in \left[{\mathbf{x}}_0,{\mathbf{x}}_1\right],{\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right)={\mathbf{q}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 。

拉格朗日插值误差

给定一个函数$\mathbf{f}:\left[\mathbf{a},\mathbf{b}\right]\to \mathbb{R}$且$\left[\mathbf{a},\mathbf{b}\right]$ 中有$n+1$个点${\mathbf{x}}_0,\cdots ,{\mathbf{x}}_{\mathbf{n}}$。 令${\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 是一个$n$阶多项式，它对数据$\left({\mathbf{x}}_0,\mathbf{f}\left({\mathbf{x}}_0\right)\right),\cdots ,\left({\mathbf{x}}_{\mathbf{n}},\mathbf{f}\left({\mathbf{x}}_{\mathbf{n}}\right)\right)$ 进行插值。 对于所有$\mathbf{j}=0,\cdots ,\mathbf{n}$，${\mathbf{p}}_{\mathbf{n}}\left({\mathbf{x}}_{\mathbf{j}}\right)=\mathbf{f}\left({\mathbf{x}}_{\mathbf{j}}\right)$ 。

根据拉格朗日插值，${\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 由下式给出：

$${\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right)=\sum _{\mathbf{i}=0}^{\mathbf{n}}{\mathbf{L}}_{\mathbf{i}}\left(\mathbf{x}\right)\mathbf{f}\left({\mathbf{x}}_{\mathbf{i}}\right)$$

通过以下方式定义插值误差：

$${\mathbf{E}}_{\mathbf{n}}\left(\mathbf{x}\right)=\mathbf{f}\left(\mathbf{x}\right)-{\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right),\mathbf{x}\in \left[\mathbf{a},\mathbf{b}\right]$$

假设$f$是连续可微的，将推导出函数${\mathbf{E}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 的公式。

以下定理给出了拉格朗日插值的估计。

定理，设${\mathbf{x}}_0,{\mathbf{x}}_1,\cdots ,{\mathbf{x}}_{\mathbf{n}}$为$n+1$个不同的值，$t$为任何实值。 令 a = min ⁡ { x 0 , ⋯   , x n } \boldsymbol{a}=\min \left\{ \boldsymbol{x}_0,\cdots ,\boldsymbol{x}_{\boldsymbol{n}} \right\} a=min{x0​,⋯,xn​} 和 b = max ⁡ { x 0 , ⋯   , x n } \boldsymbol{b}=\max \left\{ \boldsymbol{x}_0,\cdots ,\boldsymbol{x}_{\boldsymbol{n}} \right\} b=max{x0​,⋯,xn​} 。 令$f$为$[a,b]$上的实值函数，可微分$n+1$次。 令${\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 为插值数据$\left({\mathbf{x}}_0,\mathbf{f}\left({\mathbf{x}}_0\right)\right),\cdots ,\left({\mathbf{x}}_{\mathbf{n}},\mathbf{f}\left({\mathbf{x}}_{\mathbf{n}}\right)\right)$ 的拉格朗日多项式。 然后，存在$\mathbf{\xi }\in \left[\mathbf{a},\mathbf{b}\right]$ 使得

$${\mathbf{E}}_{\mathbf{n}}\left(\mathbf{x}\right)=\mathbf{f}\left(\mathbf{x}\right)-{\mathbf{p}}_{\mathbf{n}}\left(\mathbf{x}\right)=\frac{\mathbf{\psi }\left(\mathbf{x}\right)}{\left(\mathbf{n}+1\right)!}{\mathbf{f}}^{\left(\mathbf{n}+1\right)}\left(\mathbf{\xi }\right)$$

证明：

显然，对于$\mathbf{x}={\mathbf{x}}_{\mathbf{i}},{\forall }_{\mathbf{i}}=0,\cdots ,\mathbf{n}$，结果为真。 假设当$\mathbf{x}\ne {\mathbf{x}}_{\mathbf{i}}$时，结果为真。 固定值$\mathbf{t}\in \left[\mathbf{a},\mathbf{b}\right]$ 并定义为：

$$\mathbf{G}\left(\mathbf{x}\right)=\frac{\mathbf{\psi }\left(\mathbf{x}\right)}{\mathbf{\psi }\left(\mathbf{t}\right)}{\mathbf{E}}_{\mathbf{n}}\left(\mathbf{t}\right),\mathbf{x}\in \left[\mathbf{a},\mathbf{b}\right]$$

其中$\mathbf{\psi }\left(\mathbf{x}\right)={\prod }_{\mathbf{j}=0}^{\mathbf{n}}\left(\mathbf{x}-{\mathbf{x}}_{\mathbf{j}}\right)$。

由于$\mathbf{\psi }\left(\mathbf{x}\right)$ 和${\mathbf{E}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 都存在前$n+1$个导数，所以$G(x)$也是如此。因此，$G(x)$在$[a,b]$中具有$n+2$个不同的根。 根据中间值定理，${\mathbf{G}}^{\prime }\left(\mathbf{x}\right)$ 在$\left(\mathbf{a},\mathbf{b}\right)$ 中具有$n+1$个不同的根，而${\mathbf{G}}^{\prime \prime }\left(\mathbf{x}\right)$ 在$\left(\mathbf{a},\mathbf{b}\right)$ 中具有n个不同的根。通过数学归纳法，${\mathbf{G}}^{\left(\mathbf{k}\right)}$在$\left(\mathbf{a},\mathbf{b}\right)$ 中具有$\mathbf{n}-\mathbf{k}+2$个不同的根。 因此，${\mathbf{G}}^{\left(\mathbf{n}+1\right)}$在$\left(\mathbf{a},\mathbf{b}\right)$ 中至少有一个根。 令$\mathbf{\xi }$为$\left(\mathbf{a},\mathbf{b}\right)$ 中${\mathbf{G}}^{\left(\mathbf{n}+1\right)}$的根，即${\mathbf{G}}^{\left(\mathbf{n}+1\right)}\left(\mathbf{\xi }\right)=0$。

根据${\mathbf{E}}_{\mathbf{n}}\left(\mathbf{x}\right)$ 和$\mathbf{\psi }\left(\mathbf{x}\right)$ 的定义，${\mathbf{E}}_{\mathbf{n}}^{\left(\mathbf{n}+1\right)}\left(\mathbf{x}\right)={\mathbf{f}}^{\left(\mathbf{n}+1\right)}\left(\mathbf{x}\right)$ 和${\mathbf{\psi }}^{\left(\mathbf{n}+1\right)}=\left(\mathbf{n}+1\right)!$。 于是，

$${\mathbf{G}}^{\left(\mathbf{n}+1\right)}\left(\mathbf{x}\right)={\mathbf{f}}^{\left(\mathbf{n}+1\right)}\left(\mathbf{x}\right)-\frac{\left(\mathbf{n}+1\right)!}{\mathbf{\psi }\left(\mathbf{t}\right)}{\mathbf{E}}_{\mathbf{n}}\left(\mathbf{t}\right)$$

当$\mathbf{x}=\mathbf{\xi }$，

$${\mathbf{G}}^{\left(\mathbf{n}+1\right)}\left(\mathbf{\xi }\right)={\mathbf{f}}^{\left(\mathbf{n}+1\right)}\left(\mathbf{\xi }\right)-\frac{\left(\mathbf{n}+1\right)!}{\mathbf{\psi }\left(\mathbf{t}\right)}{\mathbf{E}}_{\mathbf{n}}\left(\mathbf{t}\right)=0$$

导出

$${\mathbf{E}}_{\mathbf{n}}\left(\mathbf{t}\right)=\frac{\mathbf{\psi }\left(\mathbf{t}\right)}{\left(\mathbf{n}+1\right)}{\mathbf{f}}^{\left(\mathbf{n}+1\right)}\left(\mathbf{\xi }\right)$$

牛顿插值

MATLAB插值工具

Python数据插值

详情请参考 - 亚图跨际