插值算法简介
数据分析是在大数据时代下不可获取的一环,合理、全面地分析数据,能够使得决策者在决策时作出最为明智的决定。在数据分析过程中,常常可以使用插值算法来根据已知的数据估算出未知的数据,从而模拟产生一些新的值来满足要求。
一维插值
在许多插值问题中,我们常常研究的是一维插值:
设函数 y=f(x)y=f(x)y=f(x) 在区间 [a,b][a, b][a,b] 上有定义,且已知在点 a≤x0<x1<⋯<xn≤ba \leq x_0 < x_1 < \cdots < x_n \leq ba≤x0<x1<⋯<xn≤b 上的值分别为:y0,y1,⋯,yn,y_0, y_1, \cdots, y_n,y0,y1,⋯,yn,
若存在一简单函数P(x)P(x)P(x),使 P(xi)=yi(i=0,1,2,⋯,xn)P(x_i) = y_i (i = 0, 1, 2, \cdots, x_n)P(xi)=yi(i=0,1,2,⋯,xn)则称 P(x)P(x)P(x) 为f(x)f(x)f(x) 的插值函数,点 x0,x1,⋯,xnx_0, x_1, \cdots, x_nx0,x1,⋯,xn 称为插值节点,包含插值节点的区间 [a,b][a, b][a,b] 称为插值区间,求插值函数 P(x)P(x)P(x) 的方法称为插值法。
主要插值法
- 若 P(x)P(x)P(x) 是次数不超过 nnn 的代数多项式,即P(x)=a0+a1x+⋯+anxnP(x)=a_0+a_1x+\cdots+a_nx^nP(x)=a0+a1x+⋯+anxn
- 分段插值:即一段一段的进行插值,得到的插值函数比较复杂,但是准确度较高。
- 三角插值:主要会使用傅里叶变换。
插值算法
-
拉格朗日插值法
∙\bullet∙ 两个点:(x0,y0)(x1,y1)(x_0, y_0)(x_1, y_1)(x0,y0)(x1,y1)
可以设出插值函数为:f(x)=x−x0x0−x1y0+x−x0x1−x0y1f(x)=\frac{x-x_0}{x_0-x_1}y_0 + \frac{x-x_0}{x_1-x_0}y_1f(x)=x0−x1x−x0y0+x1−x0x−x0y1∙\bullet∙ 三个点:(x0,y0)(x1,y1)(x2,y2)(x_0, y_0)(x_1, y_1)(x_2, y_2)(x0,y0)(x1,y1)(x2,y2)
可以设出插值函数为:f(x)=(x−x1)(x−x2)(x0−x1)(x0−x2)y0+(x−x0)(x−x2)(x1−x0)(x1−x2)y1+(x−x0)(x−x1)(x2−x0)(x2−x1)y2f(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_2f(x)=(x0−x1)(x0−x2)(x−x1)(x−x2)y0+(x1−x0)(x1−x2)(x−x0)(x−x2)y1+(x2−x0)(x2−x1)(x−x0)(x−x1)y2
以此类推,可以得出4个、5个……点的拉格朗日插值函数。然而,拉格朗日插值法在平常的插值问题分析中并不常用,原因是会产生龙格现象(即多项式的次数越高,函数两端的波动便会越大,越不准确)。
-
分段线性以及分段二次插值
顾名思义,就是将函数分为一段一段的,每一段都是用线性函数或者二次函数来进行估计。 -
牛顿插值
f(x)=f(x0)+f[x0,x1](x−x0)+f[x0,x1,x2](x−x0)(x−x1)+⋯+f[x0,x1,⋯,xn−2,xn−1](x−x0)(x−x1)⋯(x−xn−3)(x−xn−2)+f[x0,x1,⋯,xn−1,xn](x−x0)(x−x1)⋯(x−xn−2)(x−xn−1)\begin{aligned}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) + \cdots \\ &+f[x_0, x_1, \cdots, x_{n-2}, x_{n-1}](x-x_0)(x-x_1)\cdots(x-x_{n-3})(x-x_{n-2}) \\ &+f[x_0, x_1, \cdots, x_{n-1}, x_n](x-x_0)(x-x_1)\cdots(x-x_{n-2})(x-x_{n-1}) \end{aligned}f(x)=f(x0)+f[x0,x1](x−x0)+f[x0,x1,x2](x−x0)(x−x1)+⋯+f[x0,x1,⋯,xn−2,xn−1](x−x0)(x−x1)⋯(x−xn−3)(x−xn−2)+f[x0,x1,⋯,xn−1,xn](x−x0)(x−x1)⋯(x−xn−2)(x−xn−1)
其中 f[x0,x1]f[x_0, x_1]f[x0,x1] 是指 f(x)f(x)f(x) 关于点 x0,x1x_0, x_1x0,x1 的差商。一阶差商:f[x0,xk]=f(xk)−f(x0)xk−x0f[x_0, x_k] = \frac{f(x_k)-f(x_0)}{x_k-x_0}f[x0,xk]=xk−x0f(xk)−f(x0)
二阶差商:f[x0,x1,x2]=f[x1,x2]−f[x0,x1]x2−x0f[x_0, x_1, x_2] = \frac{f[x_1, x_2]-f[x_0, x_1]}{x_2 - x_0}f[x0,x1,x2]=x2−x0f[x1,x2]−f[x0,x1]
⋯\cdots⋯
k阶差商:f[x0,x1,⋯,xk]=f[x1,⋯,xk−1,xk]−f[x0,x1,⋯,xk−1]xk−x0f[x_0, x_1, \cdots, x_k] = \frac{f[x_1, \cdots, x_{k-1}, x_{k}]-f[x_0, x_1, \cdots, x_{k-1}]}{x_k-x_0}f[x0,x1,⋯,xk]=xk−x0f[x1,⋯,xk−1,xk]−f[x0,x1,⋯,xk−1]但是,牛顿插值法也会存在龙格现象的问题。
-
最为常用的两种插值方法:三次埃尔米特插值以及三次样条插值
4.1 三次埃尔米特插值
∙\bullet∙ 埃尔米特插值原理
设函数f(x)f(x)f(x)在区间[a,b][a, b][a,b]上有 n + 1 个互异节点 a=x0<x1<x2<⋯<xn=ba=x_0 < x_1 < x_2 < \cdots < x_n=ba=x0<x1<x2<⋯<xn=b,定义在[a,b][a, b][a,b]上函数f(x)f(x)f(x)在节点上满足:
f(xi)=yi,f′(xi)=yi′(i=0,1,2,⋯,n)(2n+2个条件)f(x_i)=y_i, f'(x_i)=y_i'(i=0, 1, 2, \cdots, n)(2n+2\text{个条件})f(xi)=yi,f′(xi)=yi′(i=0,1,2,⋯,n)(2n+2个条件)
可唯一确定一个次数不超过2n+12n+12n+1的多项式H2n+1(x)=H(x)H_{2n+1}(x)=H(x)H2n+1(x)=H(x)满足:
H(xj)=yj,H′(xj)=mj(j=0,1,⋯,n)H(x_j)=y_j, H'(x_j)=m_j(j=0, 1, \cdots, n)H(xj)=yj,H′(xj)=mj(j=0,1,⋯,n)其余项为:R(x)=f(x)−H(x)=f2n+2(ξ)(2n+2)!ω2n+2(x)R(x)=f(x)-H(x)=\frac{f^{2n+2}(\xi)}{(2n+2)!}\omega_{2n+2}(x)R(x)=f(x)−H(x)=(2n+2)!f2n+2(ξ)ω2n+2(x)
三次埃尔米特插值可在 MatlabMatlabMatlab 中调用 phipphipphip 函数进行直接求取,详见后面的例题。4.2 三次样条插值
∙\bullet∙ 三次样条插值条件
设 y=f(x)y=f(x)y=f(x) 在点 x0,x1,⋯,xnx_0, x_1, \cdots, x_nx0,x1,⋯,xn 的值为 y0,y1,⋯,yny_0, y_1, \cdots, y_ny0,y1,⋯,yn,若函数 S(x)S(x)S(x) 满足下列条件:
△\triangle△ S(xi)=f(xi)=yi,i=0,1,2,⋯,nS(x_i) = f(x_i) = y_i, i = 0, 1, 2,\cdots,nS(xi)=f(xi)=yi,i=0,1,2,⋯,n;
△\triangle△ 在每个子区间 [xix_ixi,xi+1x_{i+1}xi+1](i=0,1,2,⋯,n−1i=0,1,2,\cdots,n-1i=0,1,2,⋯,n−1)上S(x)S(x)S(x)是三次多项式;
△\triangle△ S(x)S(x)S(x)在[a,b]上二阶连续可微,则称S(x)S(x)S(x)为函数f(x)f(x)f(x)的三次样条插值函数;同样,三次样条插值也可在 MatlabMatlabMatlab 中调用 splinesplinespline 函数进行直接求取,详见后面的例题。
例题:淡水养殖中池塘水体质量的评估
-
原始数据
-
三次埃尔米特插值
MatlabMatlabMatlab代码为:
[n,m]=size(Pool_1);x=1:2:15;new_x=1:15;res_2=zeros(n,size(new_x,2));%pchip埃尔米特
数据处理结果:
可以看到,原始数据只有15周内奇数周的数据,经过插值过后,获得了完整的15周池塘水体质量的数据。
-
三次样条插值
MatlabMatlabMatlab代码为:
[n,m]=size(Pool_1);
x=1:2:15;
new_x=1:15;
res_1=zeros(n,size(new_x,2));%spline 得到的y
数据处理结果为:
同样也获得了15周内完整的池塘水体质量数据。
注意,这样的插值方法还可以对未来进行短期的预测。方法不同,预测结果同样会有些差异。
如果想要深入了解插值算法,推荐
刘春凤教授:中国大学MOOC数值计算方法
这本书。
有什么好的建议,请一定告诉我哦~~~