Numpy广播切片数组和向量

Question

给定三个numpy数组：一个多维数组x，一个带有尾随单例维度的向量y，一个没有尾随单例维度的向量z。

x = np.zeros((M,N))
y = np.zeros((M,1))
z = np.zeros((M,))

广播操作的行为取决于向量表示和上下文：

x[:,0] = y      # error cannot broadcast from shape (M,1) into shape (M)
x[:,0] = z      # OK

x[:,0] += y     # error non-broadcastable output with shape (M) doesn't match 
                # broadcast shape(M,M)
x[:,0] += z     # OK

x - y           # OK
x - z           # error cannot broadcast from shape (M,N) into shape (M)

我意识到我可以做以下事情：

x - z[:,None]   # OK

但我不明白这种显式表示法有什么用处。它并没有提高可读性。我不明白为什么表达式 x - y 是可以的，但 x - z 却是模棱两可的。

为什么Numpy对带或不带尾随单例维度的向量进行不同处理？

编辑：文档 指出：当两个维度相等或其中一个维度为1时，它们是兼容的，但是 y 和 z 都是功能上的 M x 1 向量，因为一个 M x 0 向量不包含任何元素。

- hbristow

---

一个关键点是 x[:,0] 是1维的。 x[:,[0]] 是 (M,1) 的。 z[:]=y 可能会产生相同的错误。 - hpaulj

z 可以在不复制的情况下被重塑为 (M,1) 或 (1,M)。但从功能上讲，它更接近于 (1,M)，因为会自动在开头进行扩展。 - hpaulj

Answer 1

约定俗成的规则是广播会在数组形状的开头插入单例维度。这使得在数组的最后几个维度上执行操作更加方便，所以 (x.T - z).T 应该可以工作。

如果它自动决定了哪个轴与 z 匹配，那么诸如 x - z 的操作仅在 N == M 时产生错误，从而使代码难以测试。因此，这种约定允许一些便利性，同时又能对某些错误具有健壮性。

如果你不喜欢 z[:, None] 的缩写，也许你会发现 z[:, np.newaxis] 更清晰。

要使类似于 x[:,0] = y 这样的分配工作正常，可以使用 x[:,0:1] = y 代替。

- joeln

---

你能举一个在哪里添加单例维度有用的例子吗？我想执行的常见操作是从矩阵中减去列均值，但最直观的方法 x -= np.mean(x, axis=1) 实际上并不能像我们希望的那样工作。 - hbristow

Answer 2

使用Numpy的矩阵接口而不是数组接口可以产生所需的广播行为：

x = np.asmatrix(np.zeros((M,N)))
y = np.asmatrix(np.zeros((M,1)))

x[:,0] = y              # OK
x[:,0] = y[:,0]         # OK
x[:,0] = y[:,0:1]       # OK
x[:,0] += y             # OK
x - y                   # OK
x - np.mean(x, axis=0)  # OK
x - np.mean(x, axis=1)  # OK

- hbristow

---

另一方面，引入了各种其他问题。矩阵几乎总是难以找到错误的根源。例如，* 运算符执行矩阵乘法而不是逐元素操作。如果将其提供给期望数组的函数，则会出现问题。 - Davidmh

Answer 3

将 (M,1) 和 (M,) 区别对待的一个好处是可以让您指定要对齐哪些维度以及要广播哪些维度。

假设您有：

a = np.arange(4)
b = np.arange(16).reshape(4,4)
# i.e a = array([0, 1, 2, 3])
# i.e b = array([[ 0,  1,  2,  3],
#   [ 4,  5,  6,  7],
#   [ 8,  9, 10, 11],
#   [12, 13, 14, 15]])

当你执行 c = a + b 时，a 和 b 将在 axis=1 对齐，并且 a 将沿着 axis=0 广播：

array([[0, 1, 2, 3],
   [0, 1, 2, 3],
   [0, 1, 2, 3],
   [0, 1, 2, 3]])

但是如果您想在axis=0中对a和b进行对齐，并在axis=1中进行广播呢？

array([[0, 0, 0, 0],
   [1, 1, 1, 1],
   [2, 2, 2, 2],
   [3, 3, 3, 3]])

(M,1) 和 (M,) 的区别使您能够指定要对齐和广播的维度。

（例如，如果将(M,1)和(M,)视为相同，则如何告诉numpy您想要在axis=1上进行广播？）

- aha

---

我的错误在于从代数的角度而不是通用容器的角度来处理numpy数组。在您上面强调的情况下，当a.shape == (M,1)时，c = a.T + b和c = a + b会以一种直观的方式进行广播。这是一种哲学上的差异，我正在尝试通过询问为什么在numpy的设计中做出了特定的决策来理解它。感谢您的帮助！ - hbristow

在二维操作中，a.T + b 可以指定广播的首选项，但是在高维操作中应用它会非常困难。我想NumPy的方法更加通用。 - aha