Python SciPy interpolate.make_lsq_spline用法及代码示例

本文简要介绍 python 语言中 scipy.interpolate.make_lsq_spline 的用法。

用法:

scipy.interpolate.make_lsq_spline(x, y, t, k=3, w=None, axis=0, check_finite=True)#

计算基于 LSQ(最小二乘)的拟合 B-spline(的系数)。

结果是线性组合

\[S(x) = \sum_j c_j B_j(x; t)\]

B-spline 基本元素中的 \(B_j(x; t)\) ，它最小化

\[\sum_{j} \left( w_j \times (S(x_j) - y_j) \right)^2\]

参数 ：：

x： 数组, 形状 (m,)

横坐标。

y： 数组, 形状 (m, ...)

纵坐标。

t： 数组，形状(n + k + 1，)。

结。结点和数据点必须满足Schoenberg-Whitney 条件。

k： 整数，可选

B-spline学位。默认为立方体，k = 3 。

w： 数组，形状 (m,)，可选

样条拟合的权重。必须是正的。如果 None ，则权重全部相等。默认为 None 。

axis： 整数，可选

插值轴。默认为零。

check_finite： 布尔型，可选

是否检查输入数组是否仅包含有限数。禁用可能会提高性能，但如果输入确实包含无穷大或 NaN，则可能会导致问题(崩溃、非终止)。默认为真。

返回 ：：

b： 度为 k 且带有结 t 的 BSpline 对象。

注意：

数据点的数量必须大于样条线阶数 k 。

结 t 必须满足 Schoenberg-Whitney 条件，即，必须存在数据点 x[j] 的子集，使得 t[j] < x[j] < t[j+k+1] 对应于 j=0, 1,...,n-k-2 。

例子：

生成一些嘈杂的数据：

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> rng = np.random.default_rng()
>>> x = np.linspace(-3, 3, 50)
>>> y = np.exp(-x**2) + 0.1 * rng.standard_normal(50)

现在拟合带有预定义内部结的平滑三次样条。在这里，我们通过添加边界结来使结向量 (k+1)-regular：

>>> from scipy.interpolate import make_lsq_spline, BSpline
>>> t = [-1, 0, 1]
>>> k = 3
>>> t = np.r_[(x[0],)*(k+1),
...           t,
...           (x[-1],)*(k+1)]
>>> spl = make_lsq_spline(x, y, t, k)

为了比较，我们还为同一组数据构建了一个插值样条曲线：

>>> from scipy.interpolate import make_interp_spline
>>> spl_i = make_interp_spline(x, y)

绘制两者：

>>> xs = np.linspace(-3, 3, 100)
>>> plt.plot(x, y, 'ro', ms=5)
>>> plt.plot(xs, spl(xs), 'g-', lw=3, label='LSQ spline')
>>> plt.plot(xs, spl_i(xs), 'b-', lw=3, alpha=0.7, label='interp spline')
>>> plt.legend(loc='best')
>>> plt.show()

NaN 处理：如果输入数组包含nan值，结果没有用，因为底层样条拟合例程无法处理nan.一种解决方法是对not-a-number 数据点使用零权重：

>>> y[8] = np.nan
>>> w = np.isnan(y)
>>> y[w] = 0.
>>> tck = make_lsq_spline(x, y, t, w=~w)

请注意，需要用数值替换 nan(只要相应的权重为零，精确值无关紧要。)