3.7.1 一维插值
一维数据的插值运算可以通过interpld( )完成。其调用形式如下，它实际上不是函数而是一个类：
interp1d(x, y, kind='linear', ...)
其中，x和y参数是一系列已知的数据点，kind参数是插值类型，可以是字符串或整数，它给出插值的B样条曲线的阶数，可以有如下候选值：
●'zero'、'nearest'：阶梯插值，相当于0阶B样条曲线。
●'slinear'、'linear'：线性插值，用一条直线连接所有的取样点，相当于一阶B样条曲线，'slinear'使用扩展库中的相关函数计算，而'linear'则直接使用Python编写的函数进行运算，其结果一样。
●'quadratic'、'cubic'：二阶和三阶B样条曲线，更高阶的曲线可以直接使用整数值指定。
interpld对象可以计算x的取值范围之内任意点的函数值。它可以像函数一样直接调用，和NumPy的ufunc函数一样能对数组中的每个元素进行计算，并返回一个新的数组。
下面的程序演示了kind参数以及与其对应的插值曲线，结果如图3-38所示。程序中我们使用循环对相同的数据进行4种不同阶数的插值运算。❶首先使用数据点创建一个interpld对象f，通过kind参数指定其阶数。❷调用f( )计算出一系列的插值结果。本例中，决定插值曲线的数据点一共有11个，插值之后的曲线数据点有101个。
图3-38 interpld的各阶插值
高次interp1d( )插值的运算量很大，因此对于点数较多的数据，建议使用后面介绍的UnivariateSpline( )。
from scipy import interpolate
x = np.linspace(0, 10, 11)
y = np.sin(x)
xnew = np.linspace(0, 10, 101)
pl.plot(x,y,'ro')
for kind in ['nearest', 'zero', 'slinear', 'quadratic']:
f = interpolate.interp1d(x,y,kind=kind) ❶
ynew = f(xnew) ❷
pl.plot(xnew, ynew, label=str(kind))
pl.legend(loc='lower right')
1．外推和Spline拟合
上节所介绍的interpld类要求其参数x是一个递增的序列，并且只能在x的取值范围之内进行内插计算，不能用它进行外推运算，即计算x的取值范围之外的数据点。UnivariateSpline类的插值运算比interpld更高级，它支持外推和拟合运算，其调用形式如下：
UnivariateSpline(x, y, w=None, bbox=[None, None], k=3, s=None)
●x、y是保存数据点的X-Y坐标的数组，其中x必须是递增序列。
●w是为每个数据点指定的权重值。
●k为样条曲线的阶数。
●s是平滑系数，它使得最终生成的样条曲线满足条件：
。即当s>0时，样条曲线并不一定通过各个数据点。为了让曲线通过所有数据点，必须将s参数设置为0。此外，还可以使用InterpolatedUnivariateSpline类，它与UnivariateSpline的唯一区别就是它通过所有的数据点，相当于将s设置为0。
下面的程序演示了使用UnivariateSpline对数据进行插值、外推以及样条曲线拟合：
❶如图3-39（上）所示，UnivariateSpline能够进行外推运算，虽然输入数据中没有X轴大于10的点，但是它能计算出X轴在0到12的插值结果。在X轴大于10的部分，样条曲线仍然呈现出和正弦波类似的形状，越远离输入数据范围，误差会越大，因此外推的范围是有限的。由于s参数为0，因此插值曲线经过所有的数据点。
图3-39 使用UnivariateSpline进行插值：外推（上）和数据拟合（下）
❷图3-39（下）则显示了s参数不为零时的结果，对于带噪声的输入数据，选择合适的s参数能够使得样条曲线接近无噪声时的波形，可以把它看作使用样条曲线对数据进行拟合运算。
x1 = np.linspace(0, 10, 20)
y1 = np.sin(x1)
sx1 = np.linspace(0, 12, 100)
sy1 = interpolate.UnivariateSpline(x1, y1, s=0)(sx1) ❶
x2 = np.linspace(0, 20, 200)
y2 = np.sin(x2) + np.random.standard_normal(len(x2))*
0.2
sx2 = np.linspace(0, 20, 2000)
spline2 = interpolate.UnivariateSpline(x2, y2, s=8) ❷
sy2 = spline2(sx2)
pl.figure(figsize=(8, 5))
pl.subplot(211)
pl.plot(x1, y1, ".", label=u"数据点")
pl.plot(sx1, sy1, label=u"spline曲线")
pl.legend( )
pl.subplot(212)
pl.plot(x2, y2, ".", label=u"数据点")
pl.plot(sx2, sy2, linewidth=2, label=u"spline曲线")
pl.plot(x2, np.sin(x2), label=u"无噪声曲线")
pl.legend( )
当曲线为三阶曲线时，UnivariateSpline.roots( )可以用于计算曲线与y=0的交点横坐标。下面显示了图3-39（下）中的曲线与X轴的6个交点的横坐标：
print np.array_str( spline2.roots( ), precision=3 )
[ 3.288 6.329 9.296 12.578 15.75 18.805]
如果需要计算曲线与任意横线的交点，可以事先将曲线的拟合数据在Y轴方向上进行平移。但若要计算与多条y=c横线的交点，则需要对同样的数据进行平移和拟合多次。
❶下面的root_at( )通过直接修改拟合曲线的参数，实现曲线的平移，从而可以计算与任意横线的交点。❷将roots_at( )动态添加为UnivariateSpline类的方法。❸对多条横线求交点，并进行绘图，其结果如图3-40所示。
图3-40 计算Spline与水平线的交点
def roots_at(self, v): ❶
coeff = self.get_coeffs( )
coeff -= v
try:
root = self.roots( )
return root
finally:
coeff += v
interpolate.UnivariateSpline.roots_at = roots_at ❷
pl.plot(sx2, sy2, linewidth=2, label=u"spline曲线")
ax = pl.gca( )
for level in [0.5, 0.75, -0.5, -0.75]:
ax.axhline(level, ls=":", color="k")
xr = spline2.roots_at(level) ❸
pl.plot(xr, spline2(xr), "ro")
2．参数插值
前面介绍的插值函数都需要X轴的数据是按照递增顺序排列的，就像一般的y=f(x)函数曲线一样。数学上还有一种参数曲线，它使用参数t和两个函数x=f(t)、y=g(t)来定义二维平面上的一条曲线，例如圆形、心形等曲线都是参数曲线。参数曲线的插值可以通过splprep( )和splev( )实现，这组函数支持高维空间的曲线的插值，这里以二维曲线为例介绍其用法。
❶首先调用splprep( )，其第一个参数为一组一维数组，每个数组是各点在对应轴上的坐标。s参数为平滑系数，与UnivariateSpline的含义相同。splprep( )返回两个对象，其中tck是一个元组，它包含了插值曲线的所有信息。t是自动计算出的参数曲线的参数数组。
❷调用splev( )进行插值运算，其第一个参数为一个新的参数数组，这里将t的取值范围等分200份，第二个参数为splprep( )返回的第一个对象。实际上，参数数组t是正规化之后的各个线段长度的累计，因此t的范围为0到1。
其结果如图3-41所示，图中比较了平滑系数为0和1e-4时的插值曲线。当平滑系数为0时，插值曲线通过所有的数据点。
图3-41 使用参数插值连接二维平面上的点
x = [ 4.913, 4.913, 4.918, 4.938, 4.955, 4.949, 4.911,
4.848, 4.864, 4.893, 4.935, 4.981, 5.01 , 5.021]
y = [ 5.2785, 5.2875, 5.291 , 5.289 , 5.28 , 5.26 , 5.245 ,
5.245 , 5.2615, 5.278 , 5.2775, 5.261 , 5.245 , 5.241]
pl.plot(x, y, "o")
for s in (0, 1e-4):
tck, t = interpolate.splprep([x, y], s=s) ❶
xi, yi = interpolate.splev(np.linspace(t[0], t[-1], 200), tck) ❷
pl.plot(xi, yi, lw=2, label=u"s=%g" % s)
pl.legend( )
3．单调插值
前面介绍的几种插值方法不能保证数据点的单调性，即曲线的最值可能出现在数据点之外的地方。PchipInterpolator类（别名pchip）使用单调三次插值，能够保证曲线的所有最值都出现在数据点之上。下面的程序用pchip( )对数据点进行插值，并绘制其一阶导数曲线，由图3-42的导数曲线可知，所有最值点处的导数都为0。
图3-42 单调插值能保证两个点之间的曲线为单调递增或递减
x = [0, 1, 2, 3, 4, 5]
y = [1, 2, 1.5, 2.5, 3, 2.5]
xs = np.linspace(x[0], x[-1], 100)
curve = interpolate.pchip(x, y)
ys = curve(xs)
dys = curve.derivative(xs)
pl.plot(xs, ys, label=u"pchip")
pl.plot(xs, dys, label=u"一阶导数")
pl.plot(x, y, "o")
pl.legend(loc="best")
pl.grid( )
pl.margins(0.1, 0.1)