您将如何使用 numpy 有效地矢量化这种操作？

Question

输入数据

生成 n 个给定大小（此处为 3x2）的矩阵。我也选择了 n=25，但我让 n 强调我们所拥有的是一堆矩阵这一事实。

import numpy as np
n = 25
data = np.random.rand(n, 3, 2)

这只是一个格式示例：我无法更改它。或者，如果我这样做了，则必须考虑这一更改的计算成本。

当前实现

我想以原子方式实现的是：

output = []
for datum in data: # This outputs on (3x2) matrix after the other
    d0 = datum[0]
    dr = datum[1:]
    output.append(dr-d0)

或者，以更快的方式：

output = [dr-d0 for (dr, d0) in zip(datum[:,0], datum[:,1:])]

问题

这太慢了，而且：

output = datum[:,1:] - datum[:,0]

不起作用，因为在这种情况下，减法运算的行为没有得到很好的定义。另外，这种切片效率不高。

Cython/Nuitka/PyPy 等是可能的解决方案，但如果可能的话，我现在想坚持使用原始 Numpy。也许某种函数可以非常快速地应用于 numpy 数组的外部循环的元素，而无需 python 的开销......

np.vectorize 函数不起作用：

def get_diff(mat):
    return mat[1:] - mat[0]

所以我呼吁你们，Numpy 的大祭司，Python 的仆人来启迪我可怜的灵魂！

编辑：

XY 问题

（我不知道它有名字）

我真正想做的是确定很多单纯形（读作“四面体”）的内容（读作“体积”）。 AFAIK 最简单、最有效的方法是计算：

np.linalg.det(mat[:1]-mat[0])

然后让我重新表述我的问题：如何使用普通 python 和 numpy 有效地计算任何维度为 k 的单纯形集合的内容？

Answer 1

我建议data[:,1:] - data[:,0,None]。 None 创建一个新轴（官方应该使用 np.newaxis，这让你在做什么很清楚），然后减法将按照你想要的方式运行。

纠正我认为您的列表理解中的错误：

def loop(data):
    output = []
    for datum in data: # This outputs on (3x2) matrix after the other
        d0 = datum[0]
        dr = datum[1:]
        output.append(dr-d0)
    return output

def listcomp(data):
    output = [dr-d0 for (d0, dr) in zip(data[:,0], data[:,1:])]
    return output

def sub(data):
    output = data[:,1:] - data[:,0,None]
    return output

我们有

>>> import numpy as np
>>> n = 25
>>> data = np.random.rand(n, 3, 2)
>>> res_loop = loop(data)
>>> res_listcomp = listcomp(data)
>>> res_sub = sub(data)
>>> np.allclose(res_loop, res_listcomp)
True
>>> np.allclose(res_loop, res_sub)
True
>>> 
>>> %timeit loop(data)
10000 loops, best of 3: 184 µs per loop
>>> %timeit listcomp(data)
10000 loops, best of 3: 158 µs per loop
>>> %timeit sub(data)
100000 loops, best of 3: 12.8 µs per loop