在 Scipy 中用 curve_fit 拟合向量函数
Fitting a vector function with curve_fit in Scipy
我想使用 Scipy 的 curve_fit
(或更合适的东西,如果可用的话)来拟合带有向量输出的函数。例如,考虑以下函数:
import numpy as np
def fmodel(x, a, b):
return np.vstack([a*np.sin(b*x), a*x**2 - b*x, a*np.exp(b/x)])
每个组件都有不同的功能,但它们共享我希望适合的参数。理想情况下,我会做这样的事情:
x = np.linspace(1, 20, 50)
a = 0.1
b = 0.5
y = fmodel(x, a, b)
y_noisy = y + 0.2 * np.random.normal(size=y.shape)
from scipy.optimize import curve_fit
popt, pcov = curve_fit(f=fmodel, xdata=x, ydata=y_noisy, p0=[0.3, 0.1])
但是 curve_fit
不适用于带有矢量输出的函数,并抛出错误 Result from function call is not a proper array of floats.
。我所做的是像这样展平输出:
def fmodel_flat(x, a, b):
return fmodel(x[0:len(x)/3], a, b).flatten()
popt, pcov = curve_fit(f=fmodel_flat, xdata=np.tile(x, 3),
ydata=y_noisy.flatten(), p0=[0.3, 0.1])
这行得通。如果不是矢量函数,我实际上也用不同的输入拟合了几个函数,但它们共享模型参数,我可以连接输入和输出。
是否有更合适的方法来将向量函数与 Scipy 或某些附加模块相匹配?对我来说,一个主要的考虑因素是效率——实际的拟合函数要复杂得多,拟合可能需要一些时间,所以如果 curve_fit
的这种使用被破坏并导致运行时间过长,我想知道我应该怎么做改为做。
我认为从效率的角度来看,您所做的非常好。我将尝试查看实现并提出更多量化的东西,但暂时这是我的推理。
你在曲线拟合期间所做的是优化参数 (a,b)
使得
res = sum_i |f(x_i; a,b)-y_i|^2
是最小的。我的意思是你有任意维度的数据点 (x_i,y_i)
,两个参数 (a,b)
和一个近似查询点数据的拟合模型 x_i
.
曲线拟合算法从起始 (a,b)
对开始,将其放入计算上述平方误差的黑匣子中,并尝试提出一个新的 (a',b')
对来产生更小的误差。我的观点是,上面的错误实际上是拟合算法的黑匣子:拟合的配置 space 仅由 (a,b)
参数定义。如果您想象如何实现一个简单的曲线拟合函数,您可以想象您尝试做梯度下降,并将平方误差作为成本函数。
现在黑盒如何计算误差应该与拟合过程无关。很容易看出 x_i
的维数与标量函数无关,因为无论您有 1000 个 1d 查询点还是 3d space 中的 10x10x10 网格都无关紧要。重要的是您有 1000 个点 x_i
需要从模型中计算 f(x_i) ~ y_i
。
唯一需要进一步注意的微妙之处是,在向量值函数的情况下,误差的计算并不微不足道。在我看来,使用向量值函数的 2-范数定义每个 x_i
点的误差是可以的。但是嘿:在这种情况下,点 x_i
处的平方误差是
|f(x_i; a,b)-y_i|^2 == sum_k (f(x_i; a,b)[k]-y_i[k])^2
这意味着累积每个分量的平方误差。这只是意味着您现在正在做的事情是正确的:通过复制您的 x_i
点并分别考虑函数的每个组件,您的平方误差将恰好包含每个误差的 2 范数点.
所以我的观点是您所做的在数学上是正确的,我不希望拟合过程的任何行为取决于 multivariate/vector-valued 函数的处理方式。
如果我能直截了当地推荐我自己的软件包 symfit
,我认为它正是您所需要的。可以在 docs.
中找到有关使用共享参数进行拟合的示例
您的上述具体问题将变为:
from symfit import variables, parameters, Model, Fit, sin, exp
x, y_1, y_2, y_3 = variables('x, y_1, y_2, y_3')
a, b = parameters('a, b')
a.value = 0.3
b.value = 0.1
model = Model({
y_1: a * sin(b * x),
y_2: a * x**2 - b * x,
y_3: a * exp(b / x),
})
xdata = np.linspace(1, 20, 50)
ydata = model(x=xdata, a=0.1, b=0.5)
y_noisy = ydata + 0.2 * np.random.normal(size=(len(model), len(xdata)))
fit = Fit(model, x=xdata, y_1=y_noisy[0], y_2=y_noisy[1], y_3=y_noisy[2])
fit_result = fit.execute()
查看 docs 了解更多信息!
我想使用 Scipy 的 curve_fit
(或更合适的东西,如果可用的话)来拟合带有向量输出的函数。例如,考虑以下函数:
import numpy as np
def fmodel(x, a, b):
return np.vstack([a*np.sin(b*x), a*x**2 - b*x, a*np.exp(b/x)])
每个组件都有不同的功能,但它们共享我希望适合的参数。理想情况下,我会做这样的事情:
x = np.linspace(1, 20, 50)
a = 0.1
b = 0.5
y = fmodel(x, a, b)
y_noisy = y + 0.2 * np.random.normal(size=y.shape)
from scipy.optimize import curve_fit
popt, pcov = curve_fit(f=fmodel, xdata=x, ydata=y_noisy, p0=[0.3, 0.1])
但是 curve_fit
不适用于带有矢量输出的函数,并抛出错误 Result from function call is not a proper array of floats.
。我所做的是像这样展平输出:
def fmodel_flat(x, a, b):
return fmodel(x[0:len(x)/3], a, b).flatten()
popt, pcov = curve_fit(f=fmodel_flat, xdata=np.tile(x, 3),
ydata=y_noisy.flatten(), p0=[0.3, 0.1])
这行得通。如果不是矢量函数,我实际上也用不同的输入拟合了几个函数,但它们共享模型参数,我可以连接输入和输出。
是否有更合适的方法来将向量函数与 Scipy 或某些附加模块相匹配?对我来说,一个主要的考虑因素是效率——实际的拟合函数要复杂得多,拟合可能需要一些时间,所以如果 curve_fit
的这种使用被破坏并导致运行时间过长,我想知道我应该怎么做改为做。
我认为从效率的角度来看,您所做的非常好。我将尝试查看实现并提出更多量化的东西,但暂时这是我的推理。
你在曲线拟合期间所做的是优化参数 (a,b)
使得
res = sum_i |f(x_i; a,b)-y_i|^2
是最小的。我的意思是你有任意维度的数据点 (x_i,y_i)
,两个参数 (a,b)
和一个近似查询点数据的拟合模型 x_i
.
曲线拟合算法从起始 (a,b)
对开始,将其放入计算上述平方误差的黑匣子中,并尝试提出一个新的 (a',b')
对来产生更小的误差。我的观点是,上面的错误实际上是拟合算法的黑匣子:拟合的配置 space 仅由 (a,b)
参数定义。如果您想象如何实现一个简单的曲线拟合函数,您可以想象您尝试做梯度下降,并将平方误差作为成本函数。
现在黑盒如何计算误差应该与拟合过程无关。很容易看出 x_i
的维数与标量函数无关,因为无论您有 1000 个 1d 查询点还是 3d space 中的 10x10x10 网格都无关紧要。重要的是您有 1000 个点 x_i
需要从模型中计算 f(x_i) ~ y_i
。
唯一需要进一步注意的微妙之处是,在向量值函数的情况下,误差的计算并不微不足道。在我看来,使用向量值函数的 2-范数定义每个 x_i
点的误差是可以的。但是嘿:在这种情况下,点 x_i
处的平方误差是
|f(x_i; a,b)-y_i|^2 == sum_k (f(x_i; a,b)[k]-y_i[k])^2
这意味着累积每个分量的平方误差。这只是意味着您现在正在做的事情是正确的:通过复制您的 x_i
点并分别考虑函数的每个组件,您的平方误差将恰好包含每个误差的 2 范数点.
所以我的观点是您所做的在数学上是正确的,我不希望拟合过程的任何行为取决于 multivariate/vector-valued 函数的处理方式。
如果我能直截了当地推荐我自己的软件包 symfit
,我认为它正是您所需要的。可以在 docs.
您的上述具体问题将变为:
from symfit import variables, parameters, Model, Fit, sin, exp
x, y_1, y_2, y_3 = variables('x, y_1, y_2, y_3')
a, b = parameters('a, b')
a.value = 0.3
b.value = 0.1
model = Model({
y_1: a * sin(b * x),
y_2: a * x**2 - b * x,
y_3: a * exp(b / x),
})
xdata = np.linspace(1, 20, 50)
ydata = model(x=xdata, a=0.1, b=0.5)
y_noisy = ydata + 0.2 * np.random.normal(size=(len(model), len(xdata)))
fit = Fit(model, x=xdata, y_1=y_noisy[0], y_2=y_noisy[1], y_3=y_noisy[2])
fit_result = fit.execute()
查看 docs 了解更多信息!