第6章 数据分析——6.2 数据插值

6.2  数据插值

插值法是一种古老的数学方法。插值问题的数学定义如下：

由实验或测量的方法得到函数y = f(x)在互异点x₀,x₁,…,x_n处的数值y₀,y₁,…,y_n（见图6-2），然后构造一个函数φ(x)作为y = f(x)的近似表达式，即y = f(x)≈φ(x)

使得φ(x₀) = y₀,φ(x₁) = y₁, …,φ(x_n) = y_n。这类问题称为插值问题。y = f(x)称为被插值函数，φ(x)称为插值函数，x₀,x₁,…,x_n称为插值节点。

图6-2  一维插值示意图

插值的任务是由已知的观测点为物理量建立一个简单、连续的解析模型，以便能根据该模型推测该物理量在非观测点处的特性。

插值包括多项式插值、埃尔米特插值、分段插值与样条插值、三角函数插值、辛克插值等几种。插值在数据分析、信号处理、图像处理等诸多领域有着十分重要的应用。

6.2.1  一维插值

当被插值函数y = f(x)为一元函数时，为一维插值。MATLAB使用interp1函数来实现一维插值。interp1函数的调用格式如下。

● Vq=interp1(X,V,Xq,METHOD)：X为自变量的取值范围；V为函数值，或者V为一个向量，其长度必须与X保持一致；Xq为插值点向量或数组；METHOD是字符串变量，用来设定插值方法。

MATLAB提供以下几种插值方法。

● METHOD='nearest'：邻近点插值。插值点函数值估计为与该插值点最近的数据点函数值。

● METHOD='linear'：线性插值。根据相邻数据点的线性函数估计落在该区域内插值数据点的函数值。

● METHOD='spline'：三次样条插值。这种方法在相邻数据点间建立三次多项式函数，根据多项式函数确定插值数据点的函数值。

● METHOD='pchip'或'cubic'：立方插值。通过分段立方埃尔米特插值方法计算插值结果。

● METHOD='v5cubic'：用MATLAB 5版本中的三次样条插值。

上述几种插值方法的比较如下。

● 邻近点插值方法的速度最快，但平滑性最差。

● 线性插值占用的计算机内存比邻近点插值多，运算时间也长；与邻近点插值不同，其结果是连续的，但顶点处的斜率会改变。

● 三次样条插值方法的运算时间最长，但内存的占用较立方插值少，其插值数据和导数都是连续的。在这几种插值方法中，三次样条插值结果的平滑性最好，但如果输入数据不一致或数据点过近，就可能出现很差的插值效果。

选择一种插值方法，需要考虑的因素包括运算时间、占用计算机内存的大小和插值的光滑程度。

下面介绍几种较为常用的一维插值方法。

1．分段线性插值（linear）

分段线性插值的算法是在每个小区间[x_i,x_i+1 ]上采用简单的线性插值。在区间[x_i,x_i+1 ]上的子插值多项式为：

在整个区间[x_i,x_n]上的插值函数为：

其中l_i(x)的定义如下：

例6-16：利用interp1函数对y=sin(x)进行分段线性插值。

在命令行窗口中输入：

clear all;
x = 0 : 2 * pi;
y = sin(x);
xx = 0 : 0.5 : 2 * pi;
yy = interp1(x, y, xx);
plot(x, y, 's', xx, yy)

得到的结果如图6-3所示。

图6-3  分段线性插值结果

2．一维快速傅里叶插值

一维快速傅里叶插值通过函数interpft()来实现。该函数用傅里叶变换把输入数据变换到频域，然后用更多点的傅里叶逆变换变回时域，其结果是对数据进行增采样。函数interpft()的调用格式如下。

● y=interpft(x,n)：对x进行傅里叶变换，然后采用n点傅里叶逆变换变回时域。如果x是一个向量，数据x的长度为m，采样间隔为dx，则数据y的采样间隔是dxm/n。

● y=interpft(x,n,dim)：在dim指定的维度上进行操作。

例6-17：利用一维快速傅里叶插值实现数据增采样。

在命令行窗口中输入：

clear all;
x = 0 : 1.2 : 10;
y = sin(x);
n = 2 * length(x);      % 增采样1倍
yi = interpft(y, n);    % 一维快速傅里叶插值
xi = 0 : 0.6 : 10.4;
hold on;
plot(x, y, 'ro');       % 画图
plot(xi, yi, 'b.-');
title('一维快速傅里叶插值');
legend('原始数据', '插值结果');

得到的结果如图6-4所示。

图6-4  一维快速傅里叶插值结果

3．快速fourier算法

当数据点呈现周期分布时，用上面的几种插值算法效果都不是很好，这时可以使用interpft函数进行插值，此函数使用快速fourier算法做一维插值，其调用格式如下：

y = interpft(x, n)

注意：它返回周期函数在重采样的n个等距点的插值。注意，n必须大于x的长度。

例6-18：采用interpft函数对sin函数插值。

在命令行窗口中输入：

clear all;
x = 0 : 2 * pi;
y = sin(x);
z = interpft(y, 15);
xx = linspace(0, 2 * pi, 15);
plot(x, y, '-o', xx, z, ':o')

得到的结果如图6-5所示。

图6-5  快速fourier算法插值结果

6.2.2  二维插值

当被插值函数y = f(x)为二元函数时，为二维插值。MATLAB使用interp2函数来实现二维插值。interp2函数的调用格式如下。

● Vq=interp2(X,Y,V,Xq,Yq,METHOD)：X、Y、V是具有相同大小的矩阵，V(i,j)是数据点[X(i,j),Y(i,j)]上的函数值；Xq、Yq为待插值数据网格；METHOD是一个字符串变量，表示不同的插值方法。

MATLAB提供以下4种插值方法。

● METHOD='nearest'：领近点插值。将插值点周围4个数据点中离该插值点最近的数据点函数值作为该插值点的函数值的估计值。

● METHOD='linear'：双线性插值，是MATLAB的interp2函数默认使用的插值方法。该方法将插值点周围4个数据点的函数值的线性组合作为插值点的函数值的估计值。

● METHOD='spline'：三次样条插值。该方法的计算效率高，得到的曲面光滑，是用户经常使用的插值方法。

● METHOD='cubic'：双立方插值。该方法利用插值点周围的16个数据点，相对于nearest和linear方法，其需要消耗较多的内存和计算时间，故计算效率不高，但是得到的曲面更加光滑。

例6-19：二维插值示例。

在命令行窗口中输入：

clear all;
[X, Y] = meshgrid(-3 : .25 : 3);        % 产生已知的数据栅格点
Z = peaks(X, Y);            % 计算已知点的函数值
[XI, YI] = meshgrid(-3 : .125 : 3);     % 产生更精密的插值点
ZI = interp2(X, Y, Z, XI, YI);
mesh(X, Y, Z), hold, mesh(XI, YI, ZI + 15)
hold off
axis([-3 3 -3 3 -5 20])

输出结果如图6-6所示。

图6-6  函数二维插值结果