对于那些可以阅读 Latex 的人,这就是我正在尝试计算的:
$$k_{xyi} =\sum_{j}\left (\left ( x_{i}-x_{j}\right )^{2}+\left ( y_{i}-y_{j}\右)^{2}\右)$$
其中 x 和 y 是矩阵 A 的行。
对于仅计算机语言的人来说,这将翻译为: k(x,y,i) = sum_j( (xi - xj)^2 + (yi - yj)^2 ) 其中 x 和 y 是矩阵 A 的行。
所以 k 是一个 3d 矩阵。
这可以仅通过 API 调用来完成吗? (没有 for 循环)
这里是测试启动:
import numpy as np
A = np.random.rand(4,4)
k = np.empty((4,4,4))
for ix in range(4):
for iy in range(4):
x = A[ix,]
y = A[iy,]
sx = np.power(x - x[:,np.newaxis],2)
sy = np.power(y - y[:,np.newaxis],2)
k[ix,iy] = (sx + sy).sum(axis=1).T
现在,对于编码大师来说,请将两个 for 循环替换为 numpy API 调用。
更新: 忘了说了,我需要一个节省RAM空间的方法,我的A矩阵通常是20-30千平方。因此,如果您的答案不会创建巨大的临时多维数组,那就太好了。
最佳答案
我会改变你的 latex ,使其看起来更像下面的东西 - 在我看来,它不会那么令人困惑:
据此,我认为表达式中的最后一行实际上应该是:
k[ix,iy] = (sx + sy).sum(axis=-1)
如果是这样,您可以按如下方式计算上述表达式:
Axij = (A[:, None, :] - A[..., None])**2
k = np.sum(Axij[:, None, :, :] + Axij, axis=-1)
上面首先扩展了一个内存密集型 4D 数组。如果您担心内存问题,可以通过引入新的 for 循环来跳过此步骤:
k = np.empty((4,4,4))
Axij = (A[:, None, :] - A[..., None])**2
for xi in range(A.shape[0]):
k[xi] = np.sum(Axij[xi, None, :, :] + Axij, axis=-1)
这会慢一些,但不会像你想象的那么慢,因为你仍然在 numpy 中进行大量操作。您可能可以跳过 3D Axij
中间部分,但这样做同样会导致性能损失。
如果您的矩阵在边缘上确实有 20k,您的 3D 输出将是 64TB。你不会在 numpy 甚至内存中执行此操作(除非你有一个大规模的分布式内存系统)。
关于python - numpy 行对行差异平方和,无需 for 循环(仅 api 调用),我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/27189980/