加速Numpy阵列切片

我需要在两个维度上切分中等大小的2d Numpy阵列。作为实例，加速Numpy阵列切片

import numpy as np 
X = np.random.normal(loc=0, scale=1, size=(3000, 100))

从该阵列中，我需要选择大量的行和相当小数量的列，例如

row_idx = np.random.random_integers(0, 2999, 2500) 
col_idx = np.random.random_integers(0, 99, 10)

现在，我这样做是通过以下命令：

X.take(col_idx, axis=1).take(row_idx, axis=0)

这需要我的电脑上大约115μs。问题是我需要每次运行数百万次执行此步骤。

你是否看到有机会加速这个速度？

编辑（附加信息）：我有一个矩阵X是nxk。这n行包含1xk向量。有三组：活动组（V），左组（L）和右组（R）。此外，还有系数v0和v。我需要计算此数量：http://goo.gl/KNoSl3（对不起，我无法发布图像）。问题中的公式选择左（右）集中的所有X行以及活动集中的所有列。

编辑2

我发现了另一个小的改进。

X.take(col_idx, axis=1, mode='clip').take(row_idx, axis=0, mode='clip')

有点快（我的机器上大概是25％）。

来源

2014-03-13 BayerSe

'take（）'方法需要将所选行复制一列。你应该调整你的算法来制作这个不需要的东西。我们无法告诉你如何在没有进一步的背景下做到这一点。 –

你的指数多久改变一次？ – Daniel

对于几十个观察值，行索引保持不变（更确切地说：我有k个变量，这些变量分为活动集和非活动集，我需要检查哪个变量是非活动集最合适 - 即行索引保持不变不变，只要我检查非活动集中的变量） – BayerSe

你可以使用二维看中索引：

X[row_idx,col_idx[:,None]]

但是，使用你的方法需要1毫秒〜我的机器上，VS〜300US。

除非您有关于row_idx和col_idx中的值的其他信息，否则似乎您的方法是您可以做的最好的方法。

来源

2014-03-13 14:37:06 perimosocordiae

让我们做一些事情，我们做一个指数的一维数组，满足我们的n维网格的条件。

def make_multi_index(arr, *inds): 
    tmp = np.meshgrid(*inds, indexing='ij') 
    idx = np.vstack([x.ravel() for x in tmp]) 
    return np.ravel_multi_index(idx, X.shape)

使用您的测试阵列和原来的情况下，以供参考：

%timeit X.take(col_idx, axis=1).take(row_idx, axis=0) 
10000 loops, best of 3: 95.4 µs per loop

让我们使用这个功能来构建指数，追究他们，然后用取来回报您所需的输出：

inds = make_multi_index(X, row_idx, col_idx) 
tmp = np.take(X,inds).reshape(row_idx.shape[0], col_idx.shape[0]) 

np.allclose(tmp, X.take(col_idx, axis=1).take(row_idx, axis=0)) 
Out[128]: True

因此，建立我们的指数，并保持它们似乎工作，现在的时机：

%timeit make_multi_index(X, row_idx, col_idx) 
1000 loops, best of 3: 356 µs per loop 

%timeit np.take(X,inds).reshape(row_idx.shape[0], col_idx.shape[0]) 
10000 loops, best of 3: 59.9 µs per loop

因此，它发生的并不是非常好 - 这可能会得到更好的维度，你想从中取得更好。无论如何，如果您保留这些索引超过10-15次迭代，它可以帮助一些或如果您添加一个额外的维度，并同时采取所有非活动的数据集。

来源

2014-03-13 15:25:29 Daniel

这听起来很有希望。我看我是否可以利用这种方法，也许有一些方法可以最小化索引更改的数量。谢谢！ – BayerSe

加速Numpy阵列切片

回答

相关问题