所以我阅读了关于curve_fit
here的文档。它包含以下示例:
import numpy as np
import scipy.optimize as so
def func(x, a,b,c ):
return a * np.exp(-b * x) + c
a,b,c = 2.5, 1.3, 0.5
nx = 500
noiseAlpha = 0.5
#
xdata = np.linspace(0, 4, nx)
y = func(xdata, a,b,c)
ydata = y + noiseAlpha * np.random.normal(size=len(xdata))
如果我现在调用curve_fit
,它将近似于导数,因为我没有提供任何东西。我们来计时:
实际上,curve_fit
调用leastsq
(doc here),它接受一个Dfun
参数来计算雅可比。所以我这样做了:
def myDfun( abc, xdata, ydata, f ) :
a,b,c = abc
ebx = np.exp(-b * xdata)
res = np.vstack( ( ebx, a * -xdata * ebx, np.ones(len(xdata)) ) ).T
return res
我再次计时:
%%timeit
popt, pcov = so.curve_fit(func, xdata, ydata, Dfun=myDfun )
1000 loops, best of 3: 858 µs per loop
呃?使用雅可比比比近似它慢吗?我做错什么了吗?在
这不是一个真正的答案,但我觉得这取决于问题的大小。对于小尺寸(n=500),在计算jacobian(使用提供的jac)上花费的额外时间最终可能不会得到回报。在
n=500,有刺拳:
没有:
^{pr2}$n=5000,使用刺针:
没有:
n=50000,注射:
没有:
另外,您可能需要考虑重写jacobian函数,例如去掉昂贵的
.T()
可以带来大约15%的速度提升:相关问题 更多 >
编程相关推荐