数学建模入门笔记(3) 插值与拟合

插值与拟合

插值和拟合的区别：拟合不要求过每一个已知点，而插值要求过每一个已知点，因而插值可以看作过每一个点的拟合。

插值适用于补全缺失值，因为使用一般拟合就有可能使已知值偏移，不符合需求。据说PS用某种样条插值，放大的时候最大程度的保留连续性，因此显得不是那么模糊.

数学建模笔记1.3——插值 - Infty的文章 - 知乎
https://zhuanlan.zhihu.com/p/390028714

1 插值

1.1 分段线性插值

从几何上解释就是用直线把已知点相连形成折线

数学表示

分段线性插值函数记为 $I_n(x)$
$$I_n(x)=\sum _{i=0}^n{y}_i{l}_i(x)$$
其中基函数记为 $l_i(x)$
$$l_i(x)=\{\begin{array}{l}\frac{x-x_{i-1}}{x_i-x_{i-1}},x\in \left[x_{i-1},x_i\right],i\ne 0,\\ \frac{x-x_{i+1}}{x_i-x_{i+1}},x\in \left[x_i,x_{i+1}\right],i\ne n,\\ 0,\text{ 其他 }\end{array}$$
$I_n(x)$ 有良好的收敛性, 即对于 $x\in [a,b]$, 有
$$\underset{n\to \infty }{\lim }{I}_n(x)=f(x)$$

特点

计算量小，用 $I_n(x)$ 计算 $x$ 点的插值时，只用到 $x$ 左右的两个节点， 计算量与节点个数 $n$ 无关。但 $n$ 越大插值误差越小。

不够光滑。

1.2 拉格朗日插值

数学表示

拉格朗日 (Lagrange) 插值的基函数为
$$\begin{array}{rl}{l}_i(x)& =\frac{\left(x-x_0\right)\cdots \left(x-x_{i-1}\right)\left(x-x_{i+1}\right)\cdots \left(x-x_n\right)}{\left(x_i-x_0\right)\cdots \left(x_i-x_{i-1}\right)\left(x_i-x_{i+1}\right)\cdots \left(x_i-x_n\right)}\\ & =\prod _{\begin{array}{c}j=0\\ j\ne i\end{array}}^n\frac{x-x_j}{x_i-x_j},i=0,1,\cdots ,n\end{array}$$
$l_i(x)$ 是 $n$ 次多项式, 满足
$$l_i\left(x_j\right)=\{\begin{array}{l}0,j\ne i,\\ 1,j=i。\end{array}$$

拉格朗日插值函数
$$L_n(x)=\sum _{i=0}^n{y}_i{l}_i(x)=\sum _{i=0}^n{y}_i\left(\prod _{\begin{array}{c}j=0\\ j\ne i\end{array}}^n\frac{x-x_j}{x_i-x_j}\right)$$

特点

$(1)$ 增加或改变已知节点则需要重新计算生成拉格朗日插值函数，计算量大。

$(2)$ 会出现龙格现象（Runge）：
指的是对于某些函数，使用均匀节点构造高次多项式差值时，在插值区间的边缘的误差可能很大的现象。它是由Runge在研究多项式差值的误差时发现的，这一发现很重要，因为它表明，并不是插值多项式的阶数越高，效果就会越好。

龙格现象(Runge Phenomenon) - sslchi的文章 - 知乎
https://zhuanlan.zhihu.com/p/138747068

1.3 牛顿插值

数学表示

差商

自变量之差与因变量之差之比叫差商
定义: 函数 $y=f(x)$ 在区间 $\left[x_i,x_{i+1}\right]$ 上的平均变化率
$$f\left[x_i,x_{i+1}\right]=\frac{f\left(x_{i+1}\right)-f\left(x_i\right)}{x_{i+1}-x_i}$$

称为 $f(x)$ 关于 $x_i,x_{i+1}$ 的一阶差商，并记为 $f\left[x_i,x_{i+1}\right]$
二阶差商:
$$f\left[x_i,x_{i+1},x_{i+2}\right]=\frac{f\left[x_{i+1},x_{i+2}\right]-f\left[x_i,x_{i+1}\right]}{x_{i+2}-x_i}$$
$\mathrm{m}$ 阶差商:
$$f\left[x_0,x_1,\cdots ,x_m\right]=\frac{f\left[x_1,x_2,\cdots ,x_m\right]-f\left[x_0,x_1,\cdots ,x_{m-1}\right]}{x_m-x_0}$$

牛顿插值公式

$$\begin{array}{c}{a}_0=f\left(x_0\right)\\ a_1=f\left[x_0,x_1\right]\\ a_2=f\left[x_0,x_1,x_2\right]\end{array}$$

其中一般式:
$$a_k=f\left[x_0,x_1,\cdots ,x_k\right](k=0,1,\cdots ,n)$$

将求得系数代入多项式中即可得到n次牛顿插值公式
$$N_n(x)=f\left(x_0\right)+f\left[x_0,x_1\right]\left(x-x_0\right)+\cdots +f\left[x_0,x_1,\cdots ,x_n\right]\left(x-x_0\right)\left(x-x_1\right)\cdots \left(x-x_n\right)$$

其余项为
$$\begin{array}{c}{R}_n(x)=f\left[x_0,x_1,\cdots ,x_n,x\right]\left(x-x_0\right)\left(x-x_1\right)\cdots \left(x-x_n\right)\\ =f\left[x_0,x_1,\cdots ,x_n,x\right]\prod _{i=0}^n\left(x-x_i\right)=\frac{f^{(n+1)}(\xi )}{(n+1)!}\prod _{i=0}^n\left(x-x_i\right)\\ f\left[x_0,x_1,\cdots ,x_n\right]=\frac{f^{(n+1)}(\xi )}{(n+1)!}\end{array}$$

特点

$(1)$ 添加新节点计算量小

$(2)$ 仍会出现龙格现象

http://t.csdnimg.cn/ksi89

牛顿插值的几何解释是怎么样的？ - 马同学的回答 - 知乎
https://www.zhihu.com/question/22320408/answer/141973314

1.4 三次样条插值

样条插值

即取插值函数为样条函数, 称为样条插值。例如分段线性插值是一次样条插值。

三次样条插值的数学表示

定义

即已知函数 $y=f(x)$ 在区间 $[a,b]$ 上的 $n+1$ 个节点
$$a=x_0<x_1<\cdots <x_{n-1}<x_n=b$$
上的值 $y_i=f\left(x_i\right)(i=0,1,\cdots ,n)$, 求插值函数 $S(x)$, 使得
$(1)$ $S\left(x_i\right)=y_i(i=0,1,\cdots ,n)$ 。
$(2)$ 在每个小区间 $\left[x_i,x_{i+1}\right](i=0,1,\cdots ,n-1)$ 上 $S(x)$ 是三次多项式, 记为 $S_i(x)$ 。
$(3)$ $S(x)$ 在 $[a,b]$ 上二阶连续可微。

函数 $S(x)$ 称为 $f(x)$ 的三次样条插值函数。

固定条件

由条件(2), 不妨记
S ( x ) = { S i ( x ) , x ∈ [ x i , x i + 1 ] , i = 0 , 1 , ⋯   , n − 1 } , S i ( x ) = a i x 3 + b i x 2 + c i x + d i , \begin{aligned} & S(x)=\left\{S_i(x), x \in\left[x_i, x_{i+1}\right], i=0,1, \cdots, n-1\right\}, \\ & S_i(x)=a_i x^3+b_i x^2+c_i x+d_i, \end{aligned} ​S(x)={Si​(x),x∈[xi​,xi+1​],i=0,1,⋯,n−1},Si​(x)=ai​x3+bi​x2+ci​x+di​,​

式中: $a_i,b_i,c_i,d_i$ 为待定系数, 共 $4n$ 个。
由条件 (3), 有
$$\{\begin{array}{l}{S}_i\left(x_{i+1}\right)=S_{i+1}\left(x_{i+1}\right),\\ S_i^{\prime }\left(x_{i+1}\right)=S_{i+1}^{\prime }\left(x_{i+1}\right),i=0,1,\cdots ,n-2\\ S^{\prime \prime }{}_i\left(x_{i+1}\right)=S^{\prime \prime }{}_{i+1}\left(x_{i+1}\right),\end{array}$$

容易看出, 式 (5.1) 和式 (5.2) 共含有 $4n-2$ 个方程, 为确定 $S(x)$ 的 $4n$ 个待定参数, 尚需再给出两个边界条件。

边界条件

常用的三次样条函数的边界条件有 3 种类型:
$(1)$ $S^{\prime }(a)=y_0^{\prime },S^{\prime }(b)=y_n^{\prime }$ 。 由这种边界条件建立的样条插值函数称为 $f(x)$ 的完备三次样条插值函数。
特别地, $y^{\prime }{}_0=y_n^{\prime }=0$ 时, 样条曲线在端点处呈水平状态。
如果 $f^{\prime }(x)$ 不知道, 我们可以要求 $S^{\prime }(x)$ 与 $f^{\prime }(x)$ 在端点处近似相等。这时以 $x_0,x_1$, $x_2,x_3$ 为节点作一个三次 Newton 插值多项式 $N_a(x)$, 以 $x_n,x_{n-1},x_{n-2},x_{n-3}$ 作一个三次 Newton 插值多项式 $N_b(x)$,要求
$$S^{\prime }(a)=N_a^{\prime }(a),S^{\prime }(b)=N_b^{\prime }(b)。$$

由这种边界条件建立的三次样条称为 $f(x)$ 的 Lagrange 三次样条插值函数。
$(2)$ $S^{\prime \prime }(a)=y_0^{\prime \prime },S^{\prime \prime }(b)=y_n^{\prime \prime }$ 。特别地, $y^{\prime \prime }{}_0=y_n^{\prime \prime }=0$ 时, 称为自然边界条件。
$(3)$ $S^{\prime }(a+0)=S^{\prime }(b-0),S^{\prime \prime }(a+0)=S^{\prime \prime }(b-0)$, 此条件称为周期条件。

2 拟合

曲线拟合问题的提法是, 已知一组 (二维) 数据, 即平面上的 $n$ 个点 $\left(x_i,y_i\right),i=1$, $2,\cdots ,n,x_i$ 互不相同, 寻求一个函数 (曲线) $y=f(x)$, 使 $f(x)$ 在某种准则下与所有数据点最为接近,即曲线拟合得最好。

2.1 线性最小二乘法

线性最小二乘法是解决曲线拟合最常用的方法。

基本思路

令
$$f(x)=a_1{r}_1(x)+a_2{r}_2(x)+\cdots +a_m{r}_m(x)$$

其中:

$r_k(x)$ 为事先选定的一组线性无关的函数; $a_k$ 为待定系数 $(k=1,2,\cdots ,m;m<n)$ 。
拟合准则是使 $y_i(i=1,2,\cdots ,n)$ 与 $f\left(x_i\right)$ 的距离 ${\delta }_i$ 的平方和最小, 称为最小二乘准则。

步骤

1. 函数 $r_k(x)$ 的选取

面对一组数据 $\left(x_i,y_i\right),i=1,2,\cdots ,n$, 用线性最小二乘法作曲线拟合时, 首要的也是关键的一步是恰当地选取 $r_1(x),\cdots ,r_m(x)$ 。

如果通过机理分析, 能够知道 $y$ 与 $x$ 之间的函数关系, 则 $r_1(x),\cdots ,r_m(x)$ 容易确定。

若无法知道 $y$ 与 $x$ 之间的关系, 通常可以将数据 $\left(x_i,y_i\right),i=1,2,\cdots ,n$ 作图, 直观地判断应该用什么样的曲线去作拟合。

常用的曲线有：
(1) 直线 $y=a_{1}x+a_2$ 。
(2) 多项式 $y=a_1{x}^m+\cdots +a_{m}x+a_{m+1}$ (一般 $m=2,3$, 不宜太高 $)。$
(3) 双曲线 (一支) $y=\frac{a_1}{x}+a_2$ 。
(4) 指数曲线 $y=a_1{\mathrm{e}}^{a_{2}x}$ 。
对于指数曲线, 拟合前需作变量代换, 化为对 $a_1,a_2$ 的线性函数。

已知一组数据, 用什么样的曲线拟合最好, 可以在直观判断的基础上, 选几种曲线分别拟合,然后比较,看哪条曲线的最小二乘指标 $J$ 最小。

2. 系数 $a_k$ 的确定

记
$$J\left(a_1,\cdots ,a_m\right)=\Vert \mathbf{RA}-\mathbf{Y}{\Vert }_2^2=\sum _{i=1}^n{\delta }_i^2=\sum _{i=1}^n{\left[f\left(x_i\right)-y_i\right]}^2,$$

为求 $a_1,\cdots ,a_m$ 使 $J$ 达到最小, 只需利用极值的必要条件 $\frac{\partial J}{\partial a_j}=0(j=1,\cdots ,m)$, 得到关于 $a_1,\cdots ,a_m$ 的线性方程组
$$\sum _{i=1}^n{r}_j\left(x_i\right)\left[\sum _{k=1}^m{a}_k{r}_k\left(x_i\right)-y_i\right]=0,j=1,\cdots ,m,$$

即
$$\sum _{k=1}^m{a}_k\left[\sum _{i=1}^n{r}_j\left(x_i\right)r_k\left(x_i\right)\right]=\sum _{i=1}^n{r}_j\left(x_i\right)y_i,j=1,\cdots ,m,$$

记
$$\begin{array}{rl}\mathbf{R}& ={\left[\begin{array}{ccc}{r}_1\left(x_1\right)& \cdots & r_m\left(x_1\right)\\ ⋮& ⋮& ⋮\\ r_1\left(x_n\right)& \cdots & r_m\left(x_n\right)\end{array}\right]}_{n\times m},\\ \mathbf{A}& ={\left[a_1,\cdots ,a_m\right]}^{\mathrm{T}},\mathbf{Y}={\left[y_1,\cdots ,y_n\right]}^{\mathrm{T}},\end{array}$$

方程组式 (5.5) 可表为
$${\mathbf{R}}^{\mathrm{T}}\mathbf{R}\mathbf{A}={\mathbf{R}}^{\mathrm{T}}\mathbf{Y}$$

当 { r 1 ( x ) , ⋯   , r m ( x ) } \left\{r_1(x), \cdots, r_m(x)\right\} {r1​(x),⋯,rm​(x)} 线性无关时, $\mathbf{R}$ 列满秩, ${\mathbf{R}}^{\mathrm{T}}\mathbf{R}$ 可逆, 于是方程组式 (5.6) 有唯一解
$$\mathbf{A}={\left({\mathbf{R}}^{\mathrm{T}}\mathbf{R}\right)}^{-1}{\mathbf{R}}^{\mathrm{T}}\mathbf{Y}$$

3. Matlab解法

Matlab 中的线性最小二乘的标准型为
$$\underset{A}{\min }\Vert \mathbf{RA}-\mathbf{Y}{\Vert }_2^2$$

命令为 $\mathbf{A}=\mathbf{R}\backslash \mathbf{Y}$

2.2 最小二乘优化

优化问题中目标函数由若干函数的平方和组成且求其最小，则属于最小二乘优化问题

可以用Matlab优化工具箱解决

2.3 曲线拟合与函数逼近

拟合： 离散$\Rightarrow$ 连续函数，符合最小二乘准则

逼近： 复杂连续函数$\Rightarrow$ 简单连续函数，符合最小平方逼近准则

最小平方逼近准则原理类似最小二乘准则，可看作连续的版本