scipy的插值函数_scipy.interpolate

什么是插值函数？

插值是在直线或曲线上的两点之间找到值的过程。

scipy.interpolate包进行插值

一维插值：内插值interp1d

class scipy.interpolate.interp1d(x, y, kind='linear', axis=- 1, copy=True, bounds_error=None, fill_value=nan, assume_sorted=False)

主要参数：

• x：一维数组，给定数据点集的 x 值。
• y：N 维数组，给定数据点集的 y 值，数组长度必须与 x 相等。
• kind：字符串或整数，可选项，指定使用的样条曲线的种类或插值方法。
  • 可选的字符串：‘linear’, ‘nearest’, ‘nearest-up’, ‘zero’, ‘slinear’, ‘quadratic’, ‘cubic’, ‘previous’, ‘next’；
  • ‘zero’, ‘slinear’, ‘quadratic’, ‘cubic’ 分别表示零次、一次、二次、三次样条插值；
  • ‘previous’, ‘next’ 分别表示只前点插值或后点插值；
  • ‘nearest’ 表示向下舍入， ‘nearest-up’ 表示向上舍入；
  • 默认值为 ‘linear’，即线性插值。

interp1d 允许通过参数 bounds_error、fill_value 设置外推时的边界值，但这并不是进行外推插值。

# 1. 一维插值使用示例
import numpy as np
import matplotlib.pyplot as plt  # 导入 Matplotlib 工具包
from scipy.interpolate import interp1d  # 导入 scipy 中的一维插值工具 interp1d
​
# 已知数据点集 (x,y)
x = [0.0, 2.0, 4.0, 6.0, 8.0, 10.0]  # 已知数据 x
y = [3.1, 2.7, 1.5, 0.1, 1.0, 3.9]   # 已知数据 y
# 由给定数据点集 (x,y) 求插值函数 fx
fx = interp1d(x, y, kind='linear')  # 由已知数据 (x,y) 求出插值函数 fx
# 由插值函数 fx 计算插值点的函数值
xInterp = np.linspace(0,10,100) # 指定需插值的数据点集 xInterp
yInterp = fx(xInterp)  # 调用插值函数 fx，计算 xInterp 的函数值
# 绘图
plt.plot(xInterp, yInterp, label="linear interpolate")
plt.show()

#插值举例比较
# mathmodel24_v1.py
# Demo24 of mathematical modeling algorithm
# Demo of interpolate with Scipy.interpolate
# Copyright 2021 YouCans, XUPT
# Crated：2021-08-01
​
# 2. 一维插值方法(内插)比较
import numpy as np
import matplotlib.pyplot as plt  # 导入 Matplotlib 工具包
from scipy.interpolate import interp1d  # 导入 scipy 中的一维插值工具 interp1d
​
# 生成已知数据点集 (x,y)，需插值的数据点集 xnew
np.random.seed(5)
x = np.linspace(0, 5, 10)  # 生成已知数据点集的 x
y = np.cos(x/10)*2 + 0.5*np.random.rand(10)  # 生成已知数据点集的 y
xnew = np.linspace(0, 5, 100)  # 指定需插值的数据点集 xnew
​
# 使用不同插值方法，由给定数据点集 (x,y) 求插值函数 fx
f1 = interp1d(x, y, kind="linear")  # 线性插值
f2 = interp1d(x, y, kind="zero")  # 零阶样条插值
f3 = interp1d(x, y, kind="slinear")  # 一次样条插值
f4 = interp1d(x, y, kind="quadratic")  # 二次样条插值
f5 = interp1d(x, y, kind="cubic")  # 三次样条插值
f6 = interp1d(x, y, kind="nearest")  # 临近点插值，向下舍入
# f7 = interp1d(x, y, kind="nearest-up")  # 临近点插值，向上舍入
f8 = interp1d(x, y, kind="previous")  # 前点插值
f9 = interp1d(x, y, kind="next")  # 后点插值
​
# 绘图
plt.figure(figsize=(8,6))
plt.suptitle("Data interpolate")  # 全局标题
plt.subplot(221)
plt.plot(x, y, "o",  label="data")  # 已知数据点
plt.plot(xnew, f2(xnew), label="0-order spline")  # 零阶样条插值
plt.plot(xnew, f3(xnew), label="1-order spline")  # 一阶样条插值
plt.legend(loc="lower left")
plt.subplot(222)
plt.plot(x, y, "o",  label="data")  # 已知数据点
plt.plot(xnew, f4(xnew), label="2-order spline")  # 二阶样条插值
plt.plot(xnew, f5(xnew), label="3-order spline")  # 三阶样条插值
plt.legend(loc="lower left")
plt.subplot(223)
plt.plot(x, y, "o",  label="data")  # 已知数据点
plt.plot(xnew, f1(xnew), label="linear")  # 线性插值
plt.plot(xnew, f6(xnew), label="nearest")  # 临近点插值，向下舍入
# plt.plot(xnew, f7(xnew), label="nearest-up")  # 临近点插值，向上舍入
plt.legend(loc="lower left")
plt.subplot(224)
plt.plot(x, y, "o",  label="data")  # 已知数据点
plt.plot(xnew, f8(xnew), label="previous")  # 前点插值
plt.plot(xnew, f9(xnew), label="next")  # 后点插值
plt.legend(loc="lower left")
plt.text(3, 2.4, "youcans-xupt", color='gainsboro')
plt.show()

一维插值方法：外插值UnivariateSpline

class scipy.interpolate.UnivariateSpline(x, y, w=None, bbox=[None, None], k=3, s=None, ext=0, check_finite=False)

主要参数：

• x：一维数组，数值必须递增。
• y：一维数组，数组长度必须与 x 相等。
• w：一维数组，正数，可选项。每个数据点的权重，默认所有点的权重相等。
• k：整数，可选项。样条函数的阶数，​，默认值为 3。
• s：实数，可选项，平滑参数：
  • s=0，数据插值，样条曲线必须通过所有数据点；
  • s>0，数据拟合，满足 ​；
  • 默认不设置 s，则 s=len(w)。

• ext：整数或字符串，可选项。用于控制外推插值的方案：
  • ext= 0 或 "extrapolate"，返回外推值，默认值；
  • ext= 1 或 "zeros"，返回 0；
  • ext= 2 或 "raise"，抛出异常值 ValueError；
  • ext= 3 或 "const"，返回边界值。

# 3. 一维插值方法(外插)
import numpy as np
import matplotlib.pyplot as plt  # 导入 Matplotlib 工具包
from scipy.interpolate import UnivariateSpline  # 导入 scipy 中的一维插值工具 UnivariateSpline
​
# 生成已知数据点集 (x,y)，需插值的数据点集 xnew
x = np.linspace(0, 10, 11)  # 生成已知数据点集的 x
y = np.cos((x)**2/30)*2+2  # 生成已知数据点集的 y
xnew = np.linspace(-0.5, 10.5, 110)  # 指定需插值的数据点集 xnew
​
# 使用 UnivariateSpline 插值工具，由给定数据点集 (x,y) 求插值函数 fSpl
fSpl1 = UnivariateSpline(x, y, s=0)  # 三次样条插值，s=0：插值函数经过所有数据点
y1 = fSpl1(xnew)  # 由插值函数 fSpl1 计算插值点的函数值 y1
​
fSpl2 = UnivariateSpline(x, y)  # 三次样条插值，默认 s= len(w)
y2 = fSpl2(xnew)  # 由插值函数 fSpl2 计算插值点的函数值 y2
​
fSpl2.set_smoothing_factor(0.1)  # 设置光滑因子 sf
y3 = fSpl2(xnew)   # 由插值函数 fSpl2(sf=0.1) 计算插值点的函数值 y3
​
# 绘图
fig, ax = plt.subplots(figsize=(8,6))
plt.plot(x, y, 'ro', ms=5, label="data")
plt.plot(xnew, y1, 'm', label="3rd spline interpolate")
plt.plot(xnew, y2, 'g', label="3rd spline fitting")
plt.plot(xnew, y3, 'b--', label="smoothing factor")
ax.set_title("Data interpolate with extrapolation")
plt.legend(loc="best")
plt.show()

二维插值类 interp2d

 scipy.interpolate.interp2d(x,y,z,kind='linear',copy=True,bounds_error=False,fill_value=None))

• x,y：一维数组，给定数据点集的 x,y 值。
• z：一维数组，给定数据点集对应的函数值 z。
• kind：字符串或整数，可选项，指定使用的样条曲线的种类或插值方法：‘linear’ 表示线性插值，‘cubic’ 表示三次插值，‘quintic’ 表示五次插值。默认值为 ‘linear’，即线性插值。

返回值：

• 类 interp2d() 返回一个函数，其调用方法使用插值来查找新点的值。

# 4. 二维插值方法
import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import interp2d  # 导入 scipy 中的二维插值工具 interp2d
​
# 生成已知数据网格点集 (xx,yy,z)
x = np.linspace(-1, 1.5, 25)  # x 是一维数组
y = np.linspace(-1, 1, 20)  # y 是一维数组
xx, yy = np.meshgrid(x, y)  # 生成网格点的坐标 xx,yy (二维数组)
z = np.sin((xx+yy+xx**2+yy**2))  # 计算数据网格点的值 z=f(xx,yy)
print("shape of original dataset:\n\txx:{},yy:{},z:{}".format(xx.shape,yy.shape,z.shape))
​
# 由给定数据网格点集 (xx,yy,z) 求插值函数 fInterp: xx,yy,z 都是形状相同的二维数组
fInterp = interp2d(xx, yy, z, kind='cubic')  # 三阶样条插值
​
# 由插值函数 fInterp 计算需插值的数据点的函数值
xnew = np.linspace(-1, 1.5, 150)  # xnew 是一维数组
ynew = np.linspace(-1, 1, 100)  # ynew 是一维数组
znew = fInterp(xnew, ynew)  # 计算插值函数 fInterp 在 (xnew, ynew) 所描述网格点集的函数值
print("shape of interpolation dataset:\n\txnew:{},ynew:{},znew:{}".format(xnew.shape,ynew.shape,znew.shape))
​
# 绘图
fig = plt.figure(figsize=(10, 6))
ax1 = plt.subplot(1, 2, 1, projection='3d')
ax1.set_title("2-D original data")
# ax1.plot_wireframe(xx, yy, z, rstride=2, cstride=2, linewidth=1)
surf = ax1.plot_surface(xx, yy, z, rstride=2, cstride=2, cmap=plt.cm.coolwarm)
ax1.set_zlabel('zData')
xxnew, yynew = np.meshgrid(xnew, ynew)  # 将一维数组 xnew, ynew 转换为网格点集（二维数组）
print("\txxnew:{},yynew:{},znew:{}".format(xxnew.shape,yynew.shape,znew.shape))
ax2 = plt.subplot(1, 2, 2, projection='3d')  # 3D 绘图要求 x,y,z 都是 n*m 二维数组
ax2.set_title("2-D interpolation data")
ax2.plot_wireframe(xxnew, yynew, znew, rstride=2, cstride=2,linewidth=1)
surf2 = ax2.plot_surface(xxnew, yynew, znew, rstride=2, cstride=2, cmap=plt.cm.coolwarm)
ax2.set_zlabel('zInterp')
plt.show()

通过设置 interp1d 类的参数kind，可以指定使用的样条曲线的种类或插值方法。

import numpy as np
    import matplotlib.pyplot as plt
    from scipy.interpolate import interp2d  # 导入 scipy 中的二维插值工具 interp2d
    # 生成已知数据网格点集 (xx,yy,z)
    yy, xx = np.mgrid[-2:2:20j, -3:3:30j]  # 生成网格点 30x20 = 600
    z = (1 - 0.5 * xx + xx ** 5 + yy ** 3) * np.exp(-xx ** 2 - 2 * yy ** 2)  # 计算网格点的值 z
    x, y = xx[0, :], yy[:, 0]  # 由数据网格点 xx,yy 转换一维数组 x, y
 
    # 由给定数据点集 (x,y,z) 求插值函数 fInterp： x,y 是一维数组，z 是 len(x)*len(y) 二维数组
    f1 = interp2d(x, y, z, kind='linear')  # 线性插值
    f2 = interp2d(x, y, z, kind='cubic')  # 三阶样条插值
    f3 = interp2d(x, y, z, kind='quintic')  # 五阶样条插值
    # 由插值函数 fInterp 计算需插值的网格点集 ynew,xnew 的函数值
    xnew = np.linspace(-3, 3, 120)  # xnew 是一维数组
    ynew = np.linspace(-2, 2, 80)  # ynew 是一维数组
    z1 = f1(xnew, ynew)  # 根据线性插值函数 f1 计算需插值的网格点集的函数值
    z2 = f2(xnew, ynew)  # 根据三阶样条插值函数 f2 计算需插值的网格点集的函数值
    z3 = f3(xnew, ynew)  # 根据五阶样条插值函数 f3 计算需插值的网格点集的函数值
    print("shape of interpolation dataset:\n\txnew:{},ynew:{},znew:{}".format(xnew.shape, ynew.shape, z1.shape))
    # 绘图
    plt.figure(figsize=(8, 6))
    plt.suptitle("2-D data interpolate")  # 全局标题
    plt.subplot(221)
    plt.pcolor(xx, yy, z, cmap=plt.cm.hsv, shading='auto')
    plt.title("original")
    plt.colorbar()
    plt.subplot(222)
    plt.pcolor(xnew, ynew, z1, cmap=plt.cm.hsv, shading='auto')
    plt.title("linear")
    plt.colorbar()
    plt.subplot(223)
    plt.pcolor(xnew, ynew, z2, cmap=plt.cm.hsv, shading='auto')
    plt.title("cubic")
    plt.colorbar()
    plt.subplot(224)
    plt.pcolor(xnew, ynew, z3, cmap=plt.cm.hsv, shading='auto')
    plt.title("quintic")
    plt.colorbar()
    plt.show()