【问题标题】:Fitting a 3D array of data to a 1D function with numpy or scipy使用 numpy 或 scipy 将 3D 数据数组拟合到 1D 函数
【发布时间】:2013-02-12 05:29:25
【问题描述】:

4我目前正在尝试将大量数据拟合到正弦函数中。在我只有一组数据(一维数组)的情况下,scipy.optimize.curve_fit() 工作正常。但是,据我所知,如果函数本身只是一维的,则它不允许更高维的数据输入。我不想使用 for 循环遍历数组,因为这在 python 中运行速度非常慢。

到目前为止,我的代码应该类似于:

from scipy import optimize
import numpy as np    
def f(x,p1,p2,p3,p4): return p1 + p2*np.sin(2*np.pi*p3*x + p4)      #fit function

def fit(data,guess):
   n = data.shape[0] 
   leng = np.arange(n)
   param, pcov = optimize.curve_fit(f,leng,data,guess)
   return param, pcov

其中 data 是一个三维数组 (shape=(x,y,z)),我想将每一行 data[:,a,b] 拟合到函数中,param 是一个 (4,y,z) 形状的数组作为输出。 当然,对于多维数据,这会导致

ValueError: operands could not be broadcast together with shapes (2100,2100) (5)

也许有一个简单的解决方案,但我不知道该怎么做。有什么建议吗?

搜索我的问题的答案非常困难,因为使用这些关键字的大多数主题都与高维函数的拟合有关。

【问题讨论】:

  • 不要担心 for 循环很小。无论如何,我很确定curve_fitting 将是您代码中较慢的部分。如果您怀疑循环是瓶颈,请分析代码!
  • 嗯,这个想法是如果曲线拟合可以使用整个数组而不是运行函数 y*z 次,它会更快。这就是我说 for 循环很慢时的意思。
  • 或许你可以通过FFT得到sin的参数。
  • 为什么曲线拟合应该快得多?大部分时间可能会花在拟合例程上,而不是循环数据上。
  • 我预计它会更快,因为这是我在比较 for 循环方法以遍历数据数组与直接将整个数组输入函数时使用其他例程所经历的(最多一个因素100)。我承认这些是相当简单的计算,也许使用像 curve_fit 这样更耗时的方法节省的时间可能不那么重要,但我仍然希望这个过程能稍微加速。我也尝试过使用 FFT,但由于我的采样点很少,结果不太令人满意。

标签: python arrays numpy scipy


【解决方案1】:

使用np.apply_along_axis() 可以解决您的问题。只需这样做:

func1d = lambda y, *args: optimize.curve_fit(f, xdata=x, ydata=y, *args)[0] #<-- [0] to get only popt
param = np.apply_along_axis( func1d, axis=2, arr=data )

请看下面的例子:

from scipy import optimize
import numpy as np
def f(x,p1,p2,p3,p4):
    return p1 + p2*np.sin(2*np.pi*p3*x + p4)
sx = 50  # size x
sy = 200 # size y
sz = 100 # size z
# creating the reference parameters
tmp = np.empty((4,sy,sz))
tmp[0,:,:] = (1.2-0.8) * np.random.random_sample((sy,sz)) + 0.8
tmp[1,:,:] = (1.2-0.8) * np.random.random_sample((sy,sz)) + 0.8
tmp[2,:,:] = np.ones((sy,sz))
tmp[3,:,:] = np.ones((sy,sz))*np.pi/4
param_ref = np.empty((4,sy,sz,sx))     # param_ref in this shape will allow an
for i in range(sx):                    # one-shot evaluation of f() to create 
    param_ref[:,:,:,i] = tmp           # the data sample
# creating the data sample
x = np.linspace(0,2*np.pi)
factor = (1.1-0.9)*np.random.random_sample((sy,sz,sx))+0.9
data = f(x, *param_ref) * factor       # the one-shot evalution is here
# finding the adjusted parameters
func1d = lambda y, *args: optimize.curve_fit(f, xdata=x, ydata=y, *args)[0] #<-- [0] to get only popt
param = np.apply_along_axis( func1d, axis=2, arr=data )

【讨论】:

  • 感谢您的回答,这确实有效。不幸的是,与使用 for 循环的方法相比,它并没有加快这个过程。我猜 David Zwicker 关于需要大量计算时间的 fit 函数是正确的。
  • 是的,你是对的,显然 apply_along_axisvectorizeapply_over_axes 与 Python for 循环相比没有提供任何性能提升,但它们确实简化了代码......
猜你喜欢
  • 2013-01-08
  • 2021-04-24
  • 2020-04-22
  • 2021-10-18
  • 2022-08-03
  • 2021-08-15
  • 1970-01-01
  • 2018-09-11
  • 2020-10-03
相关资源
最近更新 更多