numpy: 列向点积

Question

给定一个二维的 numpy 数组，我需要计算每一列与自身的点积，并把结果存储在一个一维数组里。下面的代码可以实现这个功能：

In [45]: A = np.array([[1,2,3,4],[5,6,7,8]])

In [46]: np.array([np.dot(A[:,i], A[:,i]) for i in xrange(A.shape[1])])
Out[46]: array([26, 40, 58, 80])

有没有简单的方法可以避免使用 Python 的循环？上面的做法虽然可以，但如果有 numpy 的现成方法可以用，那我想用这个。

补充说明 实际上，这个矩阵有很多行，但列相对较少。因此，我不太想创建比 O(A.shape[1]) 更大的临时数组。而且我也不能直接修改 A。

Answer 1

在线性代数中，行 i 和行 j 的点积就是矩阵 AA^T 中第 i,j 个位置的值。同样，列 i 和列 j 的点积就是矩阵 (A^T)A 中第 i,j 个位置的值。

所以，如果你想计算矩阵 A 中每个列向量和它自己的点积，可以用 ColDot = np.dot(np.transpose(A), A).diagonal() 这个代码。另一方面，如果你想计算每个行向量和它自己的点积，可以用 RowDot = np.dot(A, np.transpose(A)).diagonal() 这个代码。

这两行代码都会返回一个数组。

Answer 2

你可以计算所有元素的平方，然后按列求和，方法是使用下面的代码：

np.sum(np.square(A),0);

（我对sum函数的第二个参数不是很确定，这个参数用来指定求和的方向，而我现在没有安装numpy。也许你需要自己试试看 :) ...）

编辑

看了DSM的帖子，似乎你应该使用axis=0。使用square函数可能比用A*A更高效一些。

Answer 3

这样怎么样：

>>> A = np.array([[1,2,3,4],[5,6,7,8]])
>>> (A*A).sum(axis=0)
array([26, 40, 58, 80])

---

编辑：嗯，好吧，你不想要中间的大对象。也许可以试试：

>>> from numpy.core.umath_tests import inner1d
>>> A = np.array([[1,2,3,4],[5,6,7,8]])
>>> inner1d(A.T, A.T)
array([26, 40, 58, 80])

这似乎会快一点。这个方法应该能在后台完成你想要的，因为 A.T 是一个视图（如果我没理解错的话，它不会自己复制），而 inner1d 似乎会按照需要的方式循环。

非常迟来的更新：另一个选择是使用 np.einsum：

>>> A = np.array([[1,2,3,4],[5,6,7,8]])
>>> np.einsum('ij,ij->j', A, A)
array([26, 40, 58, 80])
>>> timeit np.einsum('ij,ij->j', A, A)
100000 loops, best of 3: 3.65 us per loop
>>> timeit inner1d(A.T, A.T)
100000 loops, best of 3: 5.02 us per loop
>>> A = np.random.randint(0, 100, (2, 100000))
>>> timeit np.einsum('ij,ij->j', A, A)
1000 loops, best of 3: 363 us per loop
>>> timeit inner1d(A.T, A.T)
1000 loops, best of 3: 848 us per loop
>>> (np.einsum('ij,ij->j', A, A) == inner1d(A.T, A.T)).all()
True