有效地重塑numpy数组

Question

我正在使用 NumPy 数组。

我有一个2N长度向量D，想将它的一部分重塑为一个N x N数组C。

现在这段代码做了我想要的，但对于更大的来说是一个瓶颈N：

``

import numpy as np
M = 1000
t = np.arange(M)
D = np.sin(t)    # initial vector is a sin() function
N = M / 2
C = np.zeros((N,N))
for a in xrange(N):
    for b in xrange(N):
        C[a,b] = D[N + a - b]

``

一旦C由我继续前进，这样做就可以了，等一些矩阵算法

这个嵌套循环很慢，但由于这个操作本质上是索引的变化，我想我可以使用 NumPy 的内置 reshape ( numpy.reshape) 来加速这部分。

不幸的是，我似乎无法找到转换这些索引的好方法。

有什么帮助加快这部分的速度吗？

Answer 1

您可以使用NumPy broadcasting删除那些嵌套循环 -

C = D[N + np.arange(N)[:,None] - np.arange(N)]

一个也可以np.take用来替换索引，像这样 -

C = np.take(D,N + np.arange(N)[:,None] - np.arange(N))

仔细观察会发现该模式接近于toeplitz和hankel矩阵。因此，使用这些，我们将有另外两种方法来解决它，尽管与广播具有可比的加速。实现看起来像这样 -

from scipy.linalg import toeplitz
from scipy.linalg import hankel

C = toeplitz(D[N:],np.hstack((D[0],D[N-1:0:-1])))
C = hankel(D[1:N+1],D[N:])[:,::-1]

运行时测试

In [230]: M = 1000
     ...: t = np.arange(M)
     ...: D = np.sin(t)    # initial vector is a sin() function
     ...: N = M / 2
     ...: 

In [231]: def org_app(D,N):
     ...:     C = np.zeros((N,N))
     ...:     for a in xrange(N):
     ...:         for b in xrange(N):
     ...:             C[a,b] = D[N + a - b]
     ...:     return C
     ...: 

In [232]: %timeit org_app(D,N)
     ...: %timeit D[N + np.arange(N)[:,None] - np.arange(N)]
     ...: %timeit np.take(D,N + np.arange(N)[:,None] - np.arange(N))
     ...: %timeit toeplitz(D[N:],np.hstack((D[0],D[N-1:0:-1])))
     ...: %timeit hankel(D[1:N+1],D[N:])[:,::-1]
     ...: 
10 loops, best of 3: 83 ms per loop
100 loops, best of 3: 2.82 ms per loop
100 loops, best of 3: 2.84 ms per loop
100 loops, best of 3: 2.95 ms per loop
100 loops, best of 3: 2.93 ms per loop