相关疑难解决方法(0)

我有一个k*n矩阵X和一个k*k矩阵A.对于每一列X,我想计算标量

X[:, i].T.dot(A).dot(X[:, i])

(或者,数学上Xi' * A * Xi).

目前,我有一个for循环:

out = np.empty((n,))
for i in xrange(n):
    out[i] = X[:, i].T.dot(A).dot(X[:, i])

但是因为n很大,我想尽可能快地做到这一点(即使用一些NumPy函数而不是循环).

由于我np.dot被OpenBlas和Openmpi加速,我想知道是否有可能写出双倍数额

for i in range(N):
     for j in range(N):
         B[k,l] += A[i,j,k,l] * X[i,j]

作为内在产品.就在我正在使用的那一刻

B = np.einsum("ijkl,ij->kl",A,X)

但不幸的是它很慢,只使用一个处理器.有任何想法吗？

编辑:我用一个简单的例子对给出的答案进行了基准测试,似乎它们都处于同一数量级:

A = np.random.random([200,200,100,100])
X = np.random.random([200,200])
def B1():
    return es("ijkl,ij->kl",A,X) 
def B2():
    return np.tensordot(A, X, [[0,1], [0, 1]])
def B3():
    shp = A.shape
    return np.dot(X.ravel(),A.reshape(shp[0]*shp[1],1)).reshape(shp[2],shp[3])

%timeit B1()
%timeit B2()
%timeit B3()

1 loops, best of 3: 300 ms per loop
10 loops, best of 3: 149 ms per loop
10 loops, best of 3: 150 ms per loop …

给定两个维度为 (2,2,2) 的 3-D 数组：

A = [[[ 0,  0],
    [92, 92]],

   [[ 0, 92],
    [ 0, 92]]]

B = [[[ 0,  0],
    [92,  0]],

   [[ 0, 92],
    [92, 92]]]

你如何有效地找到 A 和 B 中每个向量的欧几里德距离？

我尝试过 for 循环，但这些循环很慢，我正在按 (>>2, >>2, 2) 的顺序处理 3-D 数组。

最终我想要一个形式的矩阵：

C = [[d1, d2],
     [d3, d4]]

编辑：

我尝试了以下循环，但最大的问题是丢失了我想要保留的尺寸。但距离是正确的。

[numpy.sqrt((A[row, col][0] - B[row, col][0])**2 + (B[row, col][1] -A[row, col][1])**2) for row in range(2) for col in range(2)]

我有一个关于如何尽可能快地计算numpy距离的问题,

def getR1(VVm,VVs,HHm,HHs):
    t0=time.time()
    R=VVs.flatten()[numpy.newaxis,:]-VVm.flatten()[:,numpy.newaxis]
    R*=R
    R1=HHs.flatten()[numpy.newaxis,:]-HHm.flatten()[:,numpy.newaxis]
    R1*=R1
    R+=R1
    del R1
    print "R1\t",time.time()-t0, R.shape, #11.7576191425 (108225, 10500) 
    print numpy.max(R) #4176.26290975
    # uses 17.5Gb ram
    return R


def getR2(VVm,VVs,HHm,HHs):
    t0=time.time()
    precomputed_flat = numpy.column_stack((VVs.flatten(), HHs.flatten()))
    measured_flat = numpy.column_stack((VVm.flatten(), HHm.flatten()))
    deltas = precomputed_flat[None,:,:] - measured_flat[:, None, :]
    #print time.time()-t0, deltas.shape # 5.861109972 (108225, 10500, 2)
    R = numpy.einsum('ijk,ijk->ij', deltas, deltas)
    print "R2\t",time.time()-t0,R.shape, #14.5291359425 (108225, 10500)
    print numpy.max(R) #4176.26290975
    # uses 26Gb ram
    return R


def getR3(VVm,VVs,HHm,HHs):
    from numpy.core.umath_tests import inner1d
    t0=time.time()
    precomputed_flat = …

原始问题

我将大小为n的行P与大小为n×m的矩阵O的每列相关联.我精心设计了以下代码:

import numpy as np

def ColumnWiseCorrcoef(O, P):
    n = P.size
    DO = O - (np.sum(O, 0) / np.double(n))
    DP = P - (np.sum(P) / np.double(n))
    return np.dot(DP, DO) / np.sqrt(np.sum(DO ** 2, 0) * np.sum(DP ** 2))

它比天真的方法更有效:

def ColumnWiseCorrcoefNaive(O, P):
    return np.corrcoef(P,O.T)[0,1:O[0].size+1]

以下是我在Intel核心上使用numpy-1.7.1-MKL的时间:

O = np.reshape(np.random.rand(100000), (1000,100))
P = np.random.rand(1000)

%timeit -n 1000 A = ColumnWiseCorrcoef(O, P)
1000 loops, best of 3: 787 us per loop
%timeit -n 1000 B = ColumnWiseCorrcoefNaive(O, P)
1000 loops, best of …

给定一个4D数组M: (m, n, r, r),我如何求和所有m * n内部矩阵(形状(r, r))以获得新的形状矩阵(r * r)？

例如,

    M [[[[ 4,  1],
         [ 2,  1]],

        [[ 8,  2],
         [ 4,  2]]],

       [[[ 8,  2],
         [ 4,  2]],

        [[ 12, 3],
         [ 6,  3]]]]

我希望结果应该是

      [[32, 8],
       [16, 8]]

我有一个大的矢量场,其中场很大(例如512 ^ 3;但不一定是正方形),矢量是2D或3D(例如形状是[512,512,512,2]或[512,512, 512,3]).

计算向量的平方幅度的标量场的最快方法是什么？

我可以绕过每个方向,即

import numpy as np
shp = [256,256,256,3]                       # Shape of vector field
vf = np.arange(3*(256**3)).reshape(shp)     # Create vector field
sf = np.zeros(shp[:3])                      # Create scalar field for result

for ii in range(shp[0]):
    for jj in range(shp[1]):
        for kk in range(shp[2]):
            sf[ii,jj,kk] = np.dot( vf[ii,jj,kk,:] , vf[ii,jj,kk,:] )

但那是相当缓慢的,有什么更快的吗？