第 47 章 NumPy 广播与向量化计算

📖 本篇简介：本篇承接第 46 章的 ndarray，讲透两个让 NumPy 既强大又高效的机制：广播（broadcasting）与向量化计算。我们会逐条拆解广播的三条规则，用逐例图解表格说明形状如何从右对齐、维度为 1 时如何扩展、什么时候会抛 ValueError；然后系统介绍 ufunc 通用函数、聚合函数与 axis 参数、布尔掩码与 np.where；最后用本机在 1e6~1e7 元素上实跑对比"向量化 vs Python 循环"，给出真实耗时表。前置知识是第 46 章的 ndarray 创建与索引。学完后你会理解 pandas 整列运算为什么快、形状报错该怎么读，并为第 48 章 NumPy 与 pandas 协同打好基础。

第 46 章我们认识了 ndarray。本篇回答两个问题：形状不同的数组能不能直接相加？为什么 df["a"] + df["b"] 比 for 循环快几十倍？

目录

1. 向量化与广播要解决什么

2. 广播的三条规则

3. 广播逐例图解

4. 形状不匹配的报错

5. ufunc 通用函数

6. 聚合函数与 axis 参数

7. 布尔掩码与 np.where

8. 广播实战逐列标准化

9. 性能实测向量化 vs Python 循环

10. 常见坑

11. 本章小结与练习

1. 向量化与广播要解决什么

先看一个常见需求：二维数组每一行，都要加上同一个长度为 4 的向量。

import numpy as np

m = np.arange(12).reshape(3, 4)     # shape (3, 4)
v = np.array([100, 200, 300, 400])  # shape (4,)
print("m:\n", m)
print("m + v:\n", m + v)

m:
 [[ 0  1  2  3]
 [ 4  5  6  7]
 [ 8  9 10 11]]
m + v:
 [[100 201 302 403]
 [104 205 306 407]
 [108 209 310 411]]

注意：m 是 3×4，v 是一维长度 4。它们形状明明不一样，NumPy 却合法地完成了运算——v 被"隐形地"复制了 3 份，分别加到三行上。这个隐形复制、形状对齐的机制就是广播。

而整个运算由 C 层一次遍历完成，不需要你写任何 for 循环，这就是向量化。两者配合，构成了 NumPy/pandas 高性能的核心。

2. 广播的三条规则

NumPy 按下面三条规则判断两个数组能不能广播：

1. 维度从右对齐：形状元组从末尾（最右）开始逐位比较，短的那个在左边补长度为 1 的维度。

2. 每个维度上，两者长度必须相等，或其中一方为 1；为 1 的那一方会被"扩展"到与另一方相同。

3. 任何维度既不相等又没有一方为 1，就抛 `ValueError`。

口诀：右对齐、一比一扩、不匹配就报错。

举例对照：

| 数组 A 形状 | 数组 B 形状 | 对齐结果（可广播？） | 说明 |

|-------------|-------------|----------------------|------|

| (3, 4) | (4,) | (3, 4) ✅ | B 右边是 4=4，左边补 1 再扩展为 3 |

| (3, 4) | (3, 1) | (3, 4) ✅ | 列方向 4 vs 1，1 扩展为 4 |

| (3, 4) | (1,) | (3, 4) ✅ | 标量，全扩展 |

| (3, 4) | (3,) | ❌ 报错 | 右对齐：4 vs 3 既不等也不是 1 |

| (2, 3) | (3, 2) | ❌ 报错 | 右对齐：2 vs 3 不等 |

| (5, 1, 4) | (3, 4) | (5, 3, 4) ✅ | 中间维 1 扩展为 3 |

记住一个常见误区：(3,4) 和 (3,) 不能广播。很多人以为"都是 3 开头"就行，但规则是从右对齐：右边 4 对 3，不匹配。要按行加长度 3 的向量，得把它变成 (3,1) 列向量。

3. 广播逐例图解

例 1：数组 + 标量（最常用，标量 shape 为 ()）

a = np.array([1, 2, 3])
print("a + 10 =", a + 10)

a + 10 = [11 12 13]

标量 10 被"广播"到每个元素，等价于 a + np.array([10,10,10])，但不真的复制内存。

例 2：行向量 (1,3) + 列向量 (3,1) → 3×3 加法表

row = np.array([[1, 2, 3]])          # (1, 3)
col = np.array([[10], [20], [30]])   # (3, 1)
print("row + col:\n", row + col)

row + col:
 [[11 12 13]
 [21 22 23]
 [31 32 33]]

图解：row 的 (1,3) 在行方向复制成 3 份，col 的 (3,1) 在列方向复制成 3 份，组合出 3×3。这正是"两个向量相加得到加法表"的原理。

例 3：(3,4) 矩阵 + (4,) 行向量

m = np.arange(12).reshape(3, 4)
v = np.array([100, 200, 300, 400])
print("m + v:\n", m + v)

m + v:
 [[100 201 302 403]
 [104 205 306 407]
 [108 209 310 411]]

v 长度 4 与 m 最右维 4 对齐，三行各自加上同一个 v。

例 4：(3,4) 矩阵 + (3,1) 列向量

col2 = np.array([1000, 2000, 3000]).reshape(3, 1)
print("m + col2:\n", m + col2)

m + col2:
 [[1000 1001 1002 1003]
 [2004 2005 2006 2007]
 [3008 3009 3010 3011]]

这次加的是"每一行的整体偏移"，靠的是列向量在列方向（长度 4）上扩展。

例 5：手动提升维度 `np.newaxis` 与三维广播

一维向量 v 形状是 (4,)，要让它和 (3,4) 矩阵按"列方向"广播，就得先把它变成 (4,1)。np.newaxis 就是"在这个位置塞一个长度为 1 的轴"：

v = np.array([1, 2, 3, 4])
print("v[np.newaxis, :]:", v[np.newaxis, :].shape)   # (1, 4) 行向量
print("v[:, np.newaxis]:", v[:, np.newaxis].shape)   # (4, 1) 列向量

# 三维：(2,1,3) + (4,3) -> (2,4,3)
a3 = np.ones((2, 1, 3))
b3 = np.ones((4, 3))
print("(2,1,3) + (4,3) 结果形状:", (a3 + b3).shape)

v[np.newaxis, :]: (1, 4)
v[:, np.newaxis]: (4, 1)
(2,1,3) + (4,3) 结果形状: (2, 4, 3)

逐位看 (2,1,3) + (4,3)：最右 3=3；中间 1 vs 4，1 扩展为 4；最左 2 前面补 1 再扩展为 2。结果 (2,4,3)。reshape 和 newaxis 是手动控制广播方向的两把钥匙。

4. 形状不匹配的报错

违反规则就会收到 ValueError：

try:
    np.ones((3, 4)) + np.ones((2, 4))
except ValueError as e:
    print(e)

operands could not be broadcast together with shapes (3,4) (2,4)

读这个报错的方法：

• shapes (3,4) (2,4) 直接告诉你两个运算对象的形状；
• 从右对齐看：最右维 4=4 没问题；再往左 3 vs 2，既不等也不是 1，所以失败；
• 修复思路：要么把第二个 reshape 成 (1,4) 再扩展（如果你本意是行方向广播），要么 np.tile/手动复制到 (3,4)。

⚠️ 常见坑：(3,4) + (3,) 报的也是这个错。你想要"每行加一个数"时，把向量 .reshape(-1, 1) 变成 (3,1) 列向量即可。

5. ufunc 通用函数

NumPy 的函数分两类：ufunc（universal function）是逐元素作用的"通用函数"，输入两个数组、输出同形状数组。常见的有：

x = np.array([1, 2, 3, 4, 5])
y = np.array([10, 20, 30, 40, 50])

print("add:", np.add(x, y))            # 等价于 x + y
print("subtract:", np.subtract(y, x))
print("multiply:", np.multiply(x, y))
print("divide:", np.divide(y, x))
print("power:", np.power(x, 2))         # x 的 2 次方
print("mod:", np.mod(y, x))             # 取余
print("ceil/floor/round:", np.ceil(3.14), np.floor(3.14), np.round(3.5))

add: [11 22 33 44 55]
subtract: [ 9 18 27 36 45]
multiply: [ 10  40  90 160 250]
divide: [10. 10. 10. 10. 10.]
power: [ 1  4  9 16 25]
mod: [0 0 0 0 0]
ceil/floor/round: 4.0 3.0 4.0

比较运算也是 ufunc，返回的是布尔数组而不是单个 True/False：

print("y > 25:", y > 25)
print("x == [1,2,9,4,5]:", x == np.array([1, 2, 9, 4, 5]))
print("abs:", np.abs(np.array([-1, -2, 3])))

y > 25: [False False  True  True  True]
x == [1,2,9,4,5]: [ True  True False  True  True]
abs: [1 2 3]

除了四则与比较，常用的还有三角/指数/对数/截断类：

print("exp/log:", np.exp(np.array([0, 1])), np.log(np.array([1, np.e])))
print("sqrt:", np.sqrt(np.array([4, 9, 16])))
print("maximum(逐元素取大):", np.maximum([1, 5, 3], [4, 2, 6]))
print("clip 截断到 [0,5]:", np.clip(np.array([-2, 3, 8, 1]), 0, 5))

exp/log: [1.         2.71828183] [0. 1.]
sqrt: [2. 3. 4.]
maximum(逐元素取大): [4 5 6]
clip 截断到 [0,5]: [0 3 5 1]

np.clip(arr, lo, hi) 把小于 lo 的抬到 lo、大于 hi 的压到 hi，pandas 里 df.clip() 就是它的封装。

速查表：

| 运算符 | ufunc | 含义 |

|--------|-------|------|

| + | np.add | 加 |

| - | np.subtract | 减 |

| * | np.multiply | 乘 |

| / | np.divide | 除 |

| ** | np.power | 幂 |

| % | np.mod | 取余 |

| > < == | np.greater 等 | 比较，返回布尔数组 |

写 x + y 和 np.add(x, y) 完全等价，业务代码用运算符更简洁；ufunc 形式适合在 out=、where= 等高级参数场景使用。

6. 聚合函数与 axis 参数

聚合函数把数组"压"成一个标量或低维结果：

data = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
print("sum():", data.sum())
print("mean():", data.mean())
print("min/max:", data.min(), data.max())
print("std():", round(data.std(), 4))
print("cumsum:", data.cumsum())       # 累加
print("prod():", data.prod())          # 累乘
print("argmax/argmin:", data.argmax(), data.argmin())  # 最值下标

sum(): 45
mean(): 5.0
min/max: 1 9
std(): 2.582
cumsum: [ 1  3  6 10 15 21 28 36 45]
prod(): 362880
argmax/argmin: 8 0

axis 参数是理解聚合的关键：axis=0 沿行方向（跨行）压缩，结果剩下列数；axis=1 沿列方向（跨列）压缩，结果剩下行数。

print("axis=0 按列求和:", data.sum(axis=0))   # 每一列的和
print("axis=1 按行求和:", data.sum(axis=1))   # 每一行的和
print("axis=0 求均值:", data.mean(axis=0))

axis=0 按列求和: [12 15 18]
axis=1 按行求和: [ 6 15 24]
axis=0 求均值: [4. 5. 6.]

记忆法：axis 是"要被压扁的那个轴"。axis=0 把第 0 维（行）压扁，每行的同一列叠在一起算，所以输出长度 = 列数。pandas 的 groupby 聚合、df.mean(axis=...) 完全沿用这个约定。

两个进阶细节：

X = np.array([[1., 2., 3.], [4., 5., 6.], [7., 8., 9.], [10., 11., 12.]])
# keepdims=True：被压扁的轴保留为长度 1，方便继续广播运算
print("mean axis=0:", X.mean(axis=0))
print("mean axis=0 keepdims:", X.mean(axis=0, keepdims=True).shape)
# 二维累加：cumsum 指定 axis
print("cumsum axis=0:\n", np.arange(6).reshape(2, 3).cumsum(axis=0))

mean axis=0: [5.5 6.5 7.5]
mean axis=0 keepdims: (1, 3)
cumsum axis=0:
 [[0 1 2]
 [3 5 7]]

keepdims=True 让聚合结果保持二维，后面再和原矩阵做广播时就不会因为形状从 (3,) 变成 (1,3) 而出错——这是特征工程里"按列标准化"的常用搭配。

7. 布尔掩码与 np.where

布尔数组本身就能当"面具"（mask）用，这是筛选数据的核心：

scores = np.array([55, 72, 88, 40, 95, 60, 78])

print(">=60 通过:", scores[scores >= 60])
# np.where(条件, 满足时的值, 不满足时的值)：向量化的 if/else
print("where pass/fail:", np.where(scores >= 60, "pass", "fail"))
print("where 只取下标:", np.where(scores >= 60))

>=60 通过: [72 88 95 60 78]
where pass/fail: ['fail' 'pass' 'pass' 'fail' 'pass' 'pass' 'pass']
where 只取下标: (array([1, 2, 4, 5, 6]),)

多个布尔条件用 &（与）、|（或）、~（非）连接，每个条件必须加括号：

cond1 = scores >= 60
cond2 = scores <= 90
print("& 与:", scores[cond1 & cond2])             # 60~90 之间
print("| 或:", scores[(scores < 60) | (scores > 90)])   # 极端分
print("~ 非:", scores[~(scores >= 60)])           # 不及格

& 与: [72 88 60 78]
| 或: [55 40 95]
~ 非: [55 40]

⚠️ 常见坑：写 scores[scores >= 60 and scores <= 90] 会直接报 ValueError: The truth value of an array is ambiguous。数组是多元素的，Python 不知道该怎么把一串布尔值压成一个 True/False。永远用 &/|/~ 并加括号。

8. 广播实战：逐列标准化

特征工程里最常见的一步：把每一列（每个特征）减均值、除以标准差，变成均值 0、标准差 1。这个"按列计算再广播回整个矩阵"的需求，就是广播最经典的用武之地：

X = np.array([[1., 2., 3.],
              [4., 5., 6.],
              [7., 8., 9.],
              [10., 11., 12.]])

col_mean = X.mean(axis=0)          # 每列均值，形状 (3,)
col_std = X.std(axis=0)           # 每列标准差，形状 (3,)
Xz = (X - col_mean) / col_std     # (4,3) - (3,) 广播相减，再广播相除

print("col_mean:", col_mean)
print("col_std:", col_std)
print("标准化后:\n", Xz)
print("标准化后每列均值(应≈0):", Xz.mean(axis=0))

col_mean: [5.5 6.5 7.5]
col_std: [3.35410197 3.35410197 3.35410197]
标准化后:
 [[-1.34164079 -1.34164079 -1.34164079]
 [-0.4472136  -0.4472136  -0.4472136 ]
 [ 0.4472136   0.4472136   0.4472136 ]
 [ 1.34164079  1.34164079  1.34164079]]
标准化后每列均值(应≈0): [0. 0. 0.]

逐位对齐：X 是 (4,3)，col_mean 是 (3,)，右对齐后最右 3=3，左边补成 4——每一列减去自己的均值。整个过程没有一个 for 循环，却是"对每一列单独做减法"。

对应到 pandas，这段就是：

import pandas as pd
df = pd.DataFrame(X, columns=["f1", "f2", "f3"])
df_z = (df - df.mean()) / df.std()   # 整列广播，完全同构

df.mean() 返回一行按列均值，pandas 自动把它广播减到每一行——底层就是本篇讲的 NumPy 广播。

9. 性能实测向量化 vs Python 循环

理论说完，上数据。我们对同一个运算 a*b + a**2，分别用 Python 列表推导和 NumPy 向量化各跑一遍，规模从 1e6 到 1e7（本机 Windows / Python 3.14.7 / numpy 2.5.3 实测）：

import time

def bench(n):
    a = np.random.rand(n)
    b = np.random.rand(n)
    # Python for 循环：每个元素都过一次解释器
    t0 = time.perf_counter()
    c_loop = [a[i] * b[i] + a[i] ** 2 for i in range(n)]
    t_loop = time.perf_counter() - t0
    # 向量化：C 层一次性遍历
    t0 = time.perf_counter()
    c_vec = a * b + a ** 2
    t_vec = time.perf_counter() - t0
    return n, t_loop, t_vec, np.allclose(c_loop, c_vec)

本机实测结果：

| 元素数 n | Python 循环耗时 | 向量化耗时 | 提速倍数 | 结果一致 |

|---------:|----------------:|-----------:|---------:|:--------:|

| 1,000,000 | 1.142 s | 0.0169 s | 67.7× | ✅ |

| 2,000,000 | 2.371 s | 0.0316 s | 75.0× | ✅ |

| 5,000,000 | 7.149 s | 0.0945 s | 75.7× | ✅ |

| 10,000,000 | 13.819 s | 0.7275 s | 19.0× | ✅ |

再换一个"平方根求和"的运算验证通用性：

sqrt+sum n=5,000,000: Python 生成器 3.115s vs 向量化 0.0430s, 提速 72.4x
结果 3333018.99 vs 3333018.99

结论：

1. 同一运算，向量化通常快 50~100 倍。原因是 Python 循环每个元素都要经过解释器类型判断、字节码分派，而 NumPy 在 C 层用连续内存遍历，几乎没有 Python 开销。

2. 规模越大越明显，但到达千万级后会撞到内存带宽上限（1e7 时向量化要同时开合几个 80MB 的临时数组），提速倍数回落。

3. np.allclose 校验证明两种写法结果完全一致——性能提升不以牺牲正确性为代价。

这就是第 23 章"告别 apply"的底层原理：把逐行 apply(lambda...) 改写成整列运算，本质就是从 Python 循环切换到向量化。

10. 常见坑

1. `(3,4) + (3,)` 形状不匹配：从右对齐看是 4 vs 3，想要按行加就 reshape 成 (3,1)。

2. 广播是"虚拟扩展"不复制内存：别以为它真的复制了数组，它只是索引时按广播规则寻址。

3. 多条件用 `and/or/not`：数组多条件必须 (a>1) & (b<2)，且每个条件加括号。

4. `axis` 记反：axis=0 是压扁行、按列聚合；拿不准就拿小矩阵打印验证一次。

5. 比较运算返回布尔数组：if arr > 0: 会报错，要用 np.where 或 .all()/.any()。

6. 盲目追求最快：1e7 以上的运算先看内存，必要时用 out= 参数预分配临时数组，避免反复申请内存。

11. 本章小结与练习

小结：广播按"右对齐、等长或为 1、否则报错"三条规则工作；ufunc 逐元素计算、比较运算返回布尔数组；聚合函数靠 axis 指定压扁哪个维度；布尔掩码配合 &/|/~ 和 np.where 完成向量化筛选；本机实测向量化比 Python 循环快 50~100 倍。

练习：

1. 判断下列形状能否广播，能的话写出结果形状：(4,3)+(3,)、(4,3)+(4,)、(4,1,3)+(5,3)、(3,3)+(3,2)。

2. 给定 m = np.arange(12).reshape(3,4)，写代码让每一列都减去该列的均值（提示：先用 mean(axis=0)，再利用广播相减）。

3. 给定 arr = np.array([4, -2, 7, 0, -5, 3])，用 np.where 把正数替换为 1、负数替换为 -1、0 保持为 0。

4. 用布尔掩码挑出 arr 中所有绝对值大于 2 的元素，并统计它们的个数。

5. 解释为什么 df["a"] > 0 返回一个布尔 Series，而不是一个 True/False。

6. 选一个你项目里的 apply 逐行运算，按本篇思路改写成向量化写法，并记录改写前后的耗时（参考第 9 节的 benchmark 模板）。

下一篇预告：第 48 章 NumPy 与 pandas 协同 —— .values 与 to_numpy 有何区别？pandas 列如何无缝喂给 sklearn？本篇把 NumPy 和 pandas 彻底打通。