python - Curve_fit 到 apply_along_axis。如何加快速度？

Question

我有一些大型数据集，我想将其拟合到单指数时间衰减。

数据由在不同时间获取的多个 4D 数据集组成，因此拟合应沿着第五维运行(通过数据集)。

我当前使用的代码如下:

import numpy as np
import scipy.optimize as opt

[... load 4D datasets ....]
data = (dataset1, dataset2, dataset3)
times = (10, 20, 30)

def monoexponential(t, M0, t_const):
    return M0*np.exp(-t/t_const)

# Starting guesses to initiate  descent.
M0_init = 80.0
t_const_init = 50.0
init_guess = (M0_init, t_const_init)

def fit(vector):
    try:
        nlfit, nlpcov = opt.curve_fit(monoexponential, times, vector,
                                      p0=init_guess,
                                      sigma=None,
                                      check_finite=False,
                                      maxfev=100, ftol=0.5, xtol=1,
                                      bounds=([0, 2000], [0, 800]))
        M0, t_const = nlfit
    except:
        t_const = 0

    return t_const

# Concatenate datasets in data into a single 5D array.
concat5D = np.concatenate([block[..., np.newaxis] for block in data],
                     axis=len(data[0].shape))

# And apply the curve fitting along the last dimension.
decay_map = np.apply_along_axis(fit, len(concat5D.shape) - 1, concat5D)

代码工作正常，但需要很长时间(例如，对于dataset1.shape == (100,100,50,500))。我读过一些其他主题，提到 apply_along_axis 非常慢，所以我猜这就是罪魁祸首。不幸的是，我真的不知道这里可以使用什么作为替代方案(除了显式的 for 循环？)。

有人知道我可以做什么来避免 apply_along_axis 并加速 curve_fit 被多次调用吗？

Answer 1

因此，您要对一维数组应用 fit 操作 100*100*50*500 次(示例中包含 3 个值，现实生活中更多？)？

apply_along_axis 会迭代输入数组的所有维度(一维除外)。无需同时在多个轴上进行编译或执行此fit操作。

如果没有 apply_along_axis，最简单的方法是将数组 reshape 为二维数组，将 (100,100,50,500) 压缩到一维 (250...,)，然后对其进行迭代。然后 reshape 结果。

我认为在最后一个轴上连接数据集可能比在第一个轴上连接数据集慢，但时间显示并非如此。

np.stack 是 concatenate 的新版本，可以轻松在任何位置添加新轴。

In [319]: x=np.ones((2,3,4,5),int)
In [320]: d=[x,x,x,x,x,x]

In [321]: np.stack(d,axis=0).shape   # same as np.array(d)
Out[321]: (6, 2, 3, 4, 5)

In [322]: np.stack(d,axis=-1).shape
Out[322]: (2, 3, 4, 5, 6)

对于更大的列表(使用简单的 sum 函数):

In [295]: d1=[x]*1000       # make a big list

In [296]: timeit np.apply_along_axis(sum,-1,np.stack(d1,-1)).shape
10 loops, best of 3: 39.7 ms per loop

In [297]: timeit np.apply_along_axis(sum,0,np.stack(d1,0)).shape
10 loops, best of 3: 39.2 ms per loop

使用数组 reshape 时间的显式循环大约相同

In [312]: %%timeit 
   .....: d2=np.stack(d1,-1)
   .....: d2=d2.reshape(-1,1000)
   .....: res=np.stack([sum(i) for i in d2],0).reshape(d1[0].shape)
   .....: 
10 loops, best of 3: 39.1 ms per loop

但是像 sum 这样的函数可以在整个数组上工作，而且速度更快

In [315]: timeit np.stack(d1,-1).sum(-1).shape
100 loops, best of 3: 3.52 ms per loop

因此改变堆叠和迭代方法不会对速度产生太大影响。但改变“配合”使其可以在多个维度上工作可能会有很大帮助。我对 optimize.fit 的了解不够，不知道这是否可行。

====================

我刚刚深入研究了 apply_along_axis 的代码。它基本上构建了一个类似于 ind=(0,1,slice(None),2,1) 的索引，并执行 func(arr[ind]) ，并且然后递增它，像带有进位的长算术一样排序。因此，它只是系统地逐步遍历所有元素，同时保持一个轴为 : 切片。