Numpy矢量化错误

Question

我正在尝试对以下等式进行矢量化。我让它以非矢量化形式工作

print(K)
print(N)
print(gamma.shape)
print(means.shape)
print(img.shape)
print(N_k.shape)

产量

2
1000
(1000, 2)
(2,)
(1000)
(2,)



self.variances = np.sum(gamma * (img  - self.means) ** 2) / N_k

这给了我错误的答案。答案的形状正确，但值相差甚远

然而

self.variances = [ np.sum([gamma[n][k] * (img[n] - self.means[k]) ** 2 for n in range(0, N)]) / N_k[k] for k in range(0, K) ]

有效。我是 numpy 的新手，所以我确定我犯了一个愚蠢的错误，但在这种情况下我找不到它

谢谢

Answer 1

使用

self.variances = np.sum(gamma * (img[:, None] - self.means) ** 2, axis=0) / N_k

而不是

self.variances = np.sum(gamma * (img - self.means) ** 2) / N_k

---

import numpy as np

N, K = 10, 20
gamma = np.random.random((N, K))
means = np.random.random(K)
N_k = np.random.random(K)
img = np.random.random(N)

expected = np.array([ np.sum([gamma[n][k] * (img[n] - means[k]) ** 2 for n in range(0, N)]) / N_k[k] for k in range(0, K) ])

result = np.sum(gamma * (img[:, None]  - means) ** 2, axis=0) / N_k

assert np.allclose(result, expected)

---

注意img - self.means 减去相应的值elementwise。 img 的形状为(N,) 和self.means 的形状为(K,)，因此如果N == K 可以从另一个中减去一个，但如果N != K 则引发ValueError，因为那时形状不兼容。由于您没有收到错误，N 必须等于 K。

由于您要计算一个二维数组值，首先通过添加一个新轴：img[:, None]，使img 成为一个形状为(N, 1) 的二维数组。

然后利用 broadcasting所以 当您减去一维数组self.means 时，它将广播到二维 兼容形状(N, K) 的数组。默认情况下广播会添加新轴on 左边。因此形状(K,) 的数组可以自动广播到一个 形状数组(N, K)。 （这也解释了为什么我们需要使用img[:, None] 在右侧显式添加新轴。）

现在(img[:, None] - self.means) 将是形状为(N, K) 的二维数组。

---

还请注意，在调用np.sum 时指定axis=0 很重要，因为求和是在长度N 的第一个轴上完成的。由于 Python 使用基于 0 的索引，因此第一个轴对应于 axis=0。

如果您不指定轴，则默认情况下np.sum 对所有轴求和。