Python大数据矩阵操作

Question

我想我有一个大数据（N=1e6，dimension=3）的场景。我需要在我的代码中做一些矩阵操作，比如einsum，矩阵求逆等等。为了给大家一个想法，我想做如下的事情。在

import numpy.random as rd

ndata, kdata = 1e6, 1e5

x = rd.normal(0,1,(ndata, kdata,3,3))

y = rd.normal(0,1,(ndata, kdata,3,3))

对于小型ndata，遵循kdata是一种高效、方便的方法

^{pr2}$

由于我有大量的ndata和kdata，上述方法变成了一个内存受限的问题，所以下一个赌注是在ndata和kdata上嵌套for循环的点积，如下所示：

xyloop1 = np.empty((ndata, kdata, 3, 3))

for j in xrange(ndata):

    for k in xrange(kdata):

        xyloop1[j,k] =  np.dot(x[j,k], y[j,k] )

考虑到我所学的，在python中循环是很讨厌的。另外，我想利用numpy的优点，因此认为块矩阵方法更可取，如下所示：

nstep = 200
ndiv  = ndata/nstep   

kstep = 200
kdiv  = kdata/kstep   

xyloop2 = np.empty((ndata, kdata, 3, 3))

for j in xrange(ndiv):

    ji, jf = j*nstep, (j+1)*nstep     

    for k in xrange(kdiv):

        ki, kf = k*kstep, (k+1)*kstep     

        xyloop2[ji:jf,ki:kf] =  einsum('pqrs, pqsu -> pqru', x[ji:jf,ki:kf], y[ji:jf,ki:kf] )

另外，我还需要这些xy或xyloop1或xyloop2来进行进一步的计算。所以每次计算后我都要写和读。考虑到系统I/O的带宽，您认为最佳方法是方法3，因为与方法2相比，方法3意味着更少的I/O和更少的迭代次数？如果你有任何其他想法或需要更多的信息，请告诉我。在

我是新手，所以请对我温柔点：）。任何帮助将不胜感激。顺便说一句，我正试图解决一个大数据的混合建模问题。谢谢！在

Answer 1

很明显，有时einsum是有效的。p、q、r、s为100、50、3、3的示例

示例一：

%timeit tt=np.einsum('pqrs, pqsu->pqru',x,y)
100 loops, best of 3: 3.45 ms per loop

%timeit zz= np.sum(x[:,:,:,None,:]*y[:,:,:,None],axis=-2)
10000 loops, best of 3: 153 µs per loop

例二：

^{pr2}$

Answer 2

虽然我同意这样的评论，即确定知道的唯一方法是为自己分析情况，但是有几个指导原则可以帮助您在第一次尝试时编写高效的numpy代码。以下是一些针对您的问题的建议：

• 创建一个新的numpy数组的开销大约是加法/乘法开销的1000倍，因此方法2应该效率低下，因为每次调用np.dot都会创建一个数组，但只执行27次加法乘法。在
• 如果要在python中有一个慢for循环，尽可能在最左边的轴上执行（对于C顺序数组）。在
• 很难高效地编写非常通用的N维代码，因此我的猜测是，作为一系列更简单的numpy调用将比np.einsum更有效。尝试C = np.sum(A[...,:,None] * B[...,:,:], axis=-2)（尽管这是相当推测性的）。在

所以我会尝试以下方法：

xyloop2 = np.empty((ndata, kdata, 3, 3))

for i in xrange(ndata):
    xyloop2[i] = np.sum(x[i,:,:,:,None] * y[i,:,None,:,:], axis=-2)

与方法2类似，但for循环要简单得多（而且效率更高）。我也把矩阵乘法换成了我认为可能更快的方法。在