Python中的scipy BSpline拟合

Question

这是我第一次使用BSpline，我想将曲线拟合到我的数据点上。我尝试过使用Univariate Spline，并尝试使用splev和splrep，但我真的想学习如何使用BSpline来完成这个任务。

看起来我的拟合结果非常粗糙，曲线甚至没有经过数据点。

arraymagU = linspace(U_timeband.min(),U_timeband.max(),300) #array for my x data points
UfunctionBS = BSpline(U_timeband,U_magband,k=4,extrapolate=False)
arraymagU2 = UfunctionBS(arraymagU)

plt.plot(arraymagU,arraymagU2)

U_timeband是我的x坐标，U_magband仅仅是我的y坐标。我认为k=4表示一个三次函数拟合？我已经尝试过改变这个值，但好像并没有使结果更好。

它产生了这样的结果：

如何让这个结果更好、一致呢？ 我认为我可能需要定义断点，但我不确定怎么做。

- Sophia Medallon

Answer 1

splrep 返回一个元组 (t,c,k) 包含节向量，B样条系数和样条阶数。 这些可以输入到 interpolate.BSpline 中创建一个 BSpline 对象：

import numpy as np
import scipy.interpolate as interpolate
import matplotlib.pyplot as plt

x = np.array([ 0. ,  1.2,  1.9,  3.2,  4. ,  6.5])
y = np.array([ 0. ,  2.3,  3. ,  4.3,  2.9,  3.1])

t, c, k = interpolate.splrep(x, y, s=0, k=4)
print('''\
t: {}
c: {}
k: {}
'''.format(t, c, k))
N = 100
xmin, xmax = x.min(), x.max()
xx = np.linspace(xmin, xmax, N)
spline = interpolate.BSpline(t, c, k, extrapolate=False)

plt.plot(x, y, 'bo', label='Original points')
plt.plot(xx, spline(xx), 'r', label='BSpline')
plt.grid()
plt.legend(loc='best')
plt.show()

- unutbu

---

谢谢。我已经使用了UnivariateSpline，它的拟合效果更好，但我相信我的导师坚持要使用BSpline进行光曲线拟合（这是来自超新星的光度学数据）。那么如何知道最佳的样条拟合呢？使用UnivariateSpline的样条还是interp1d的样条？ - Sophia Medallon

Answer 2

BSpline允许您在已知系数的情况下构造B样条。如果您想要适应这些系数，您将需要使用类似splrep的方法。另一种选择是在线性回归上使用BSpline.basis_elemement，但对于您的用例来说，使用splrep几乎肯定更好。

通常需要给出结点向量，但幸运的是这并不太复杂。接受的答案(使用s=0进行精确匹配)基本上将节点设置为输入坐标的内部点，但对于带有噪声数据的情况，它可能会过度拟合并且仍然很“棱角分明”：

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import BSpline, splrep, splev

# Generate
np.random.seed(0)
n = 300
ts = np.sort(np.random.uniform(0, 5, size=n))
ys = np.sin(ts) + 0.1*np.random.randn(n)

# Fit
tck = splrep(ts, ys, t=ts[2:-2], k=3)
# Alternative:
# tck = splrep(ts, ys, s=0, k=3)
ys_interp = splev(ts, tck)

# Display
plt.figure(figsize=(12, 6))
plt.plot(ts, ys, '.c')
plt.plot(ts, ys_interp, '-m')
plt.show()

通常更好的方法是将节点定义为输入坐标的分位数，并选择合理数量（我发现对于简单形状，5-10个节点效果很好）：

# Fit
n_interior_knots = 5
qs = np.linspace(0, 1, n_interior_knots+2)[1:-1]
knots = np.quantile(ts, qs)
tck = splrep(ts, ys, t=knots, k=3)
ys_smooth = splev(ts, tck)

# Alternative if one really wants to use BSpline: 
# ys_smooth = BSpline(*tck)(ts)

# Display
plt.figure(figsize=(12, 6))
plt.plot(ts, ys, '.c')
plt.plot(ts, ys_smooth, '-m')
plt.show()

- Thrastylon

Answer 3

请注意，要使用interpolate.BSpline，您需要scipy 0.19或更高版本（请参阅：https://docs.scipy.org/doc/scipy/reference/release.0.19.0.html#scipy-interpolate-improvements）。

- user3097732