文中内容仅限技术学习与代码实践参考,市场存在不确定性,技术分析需谨慎验证,不构成任何投资建议。
NumPy 快速入门
前提条件
你需要掌握一些 Python 知识。如果需要复习,可以查看 Python 教程。
为了运行示例代码,你需要安装 matplotlib
,除了 NumPy 之外。
学习者画像
这是对 NumPy 中数组的一个快速概览。它展示了如何表示和操作多维( n ≥ 2 n \geq 2 n≥2)数组。特别是,如果你不知道如何对多维数组应用常见函数(而不使用 for 循环),或者你想了解多维数组的轴和形状属性,这篇文章可能会对你有所帮助。
学习目标
阅读完本文后,你应该能够:
- 理解 NumPy 中一维、二维和多维数组的区别;
- 理解如何对多维数组应用一些线性代数运算,而无需使用 for 循环;
- 理解多维数组的轴和形状属性。
基础知识
NumPy 的主要对象是同质多维数组。它是一个由元素(通常是数字)组成的表格,所有元素类型相同,通过非负整数元组进行索引。在 NumPy 中,维度被称为轴。
例如,3D 空间中一个点的坐标数组 [1, 2, 1]
有一个轴。该轴包含 3 个元素,因此我们说它的长度是 3。在下面的示例中,数组有两个轴。第一个轴的长度为 2,第二个轴的长度为 3。
[[1., 0., 0.],
[0., 1., 2.]]
NumPy 的数组类名为 ndarray
,它也被称为 array
。注意,numpy.array
与 Python 标准库中的 array.array
不同,后者仅支持一维数组,并且功能较少。ndarray
对象的一些重要属性如下:
ndarray.ndim
:数组的轴(维度)数量。ndarray.shape
:数组的维度。这是一个整数元组,表示数组在每个维度上的大小。对于一个 n × m n \times m n×m 的矩阵,shape
为(n, m)
。shape
元组的长度就是轴的数量,即ndim
。ndarray.size
:数组中元素的总数。它等于shape
中各元素的乘积。ndarray.dtype
:描述数组中元素类型的对象。可以使用标准 Python 类型创建或指定dtype
。此外,NumPy 还提供了自己的类型,例如numpy.int32
、numpy.int16
和numpy.float64
。ndarray.itemsize
:数组中每个元素的字节大小。例如,类型为float64
的数组的itemsize
为 8(64/8),而类型为complex32
的数组的itemsize
为 4(32/8)。它等同于ndarray.dtype.itemsize
。ndarray.data
:包含数组实际元素的缓冲区。通常,我们不会直接使用这个属性,因为我们可以通过索引功能访问数组中的元素。
示例
>>> import numpy as np
>>> a = np.arange(15).reshape(3, 5)
>>> a
array([[ 0, 1, 2, 3, 4],
[ 5, 6, 7, 8, 9],
[10, 11, 12, 13, 14]])
>>> a.shape
(3, 5)
>>> a.ndim
2
>>> a.dtype.name
'int64'
>>> a.itemsize
8
>>> a.size
15
>>> type(a)
<class 'numpy.ndarray'>
>>> b = np.array([6, 7, 8])
>>> b
array([6, 7, 8])
>>> type(b)
<class 'numpy.ndarray'>
数组创建
创建数组有多种方法。
例如,可以使用 array
函数从普通的 Python 列表或元组创建数组。结果数组的类型将根据序列中元素的类型推断得出。
>>> import numpy as np
>>> a = np.array([2, 3, 4])
>>> a
array([2, 3, 4])
>>> a.dtype
dtype('int64')
>>> b = np.array([1.2, 3.5, 5.1])
>>> b.dtype
dtype('float64')
一个常见的错误是将多个参数传递给 array
,而不是传递一个序列作为参数。
>>> a = np.array(1, 2, 3, 4) # 错误
Traceback (most recent call last):
...
TypeError: array() takes from 1 to 2 positional arguments but 4 were given
>>> a = np.array([1, 2, 3, 4]) # 正确
array
函数会将序列的序列转换为二维数组,将序列的序列的序列转换为三维数组,依此类推。
>>> b = np.array([(1.5, 2, 3), (4, 5, 6)])
>>> b
array([[1.5, 2. , 3. ],
[4. , 5. , 6. ]])
也可以在创建数组时显式指定数组的类型:
>>> c = np.array([[1, 2], [3, 4]], dtype=complex)
>>> c
array([[1.+0.j, 2.+0.j],
[3.+0.j, 4.+0.j]])
通常,数组的元素最初是未知的,但数组的大小是已知的。因此,NumPy 提供了几种函数来创建带有初始占位符内容的数组。这些函数可以减少数组扩展的需要,因为数组扩展是一种代价较高的操作。
zeros
函数创建一个全为零的数组,ones
函数创建一个全为一的数组,而 empty
函数创建一个初始内容随机且取决于内存状态的数组。默认情况下,创建的数组的 dtype
是 float64
,但可以通过关键字参数 dtype
指定。
>>> np.zeros((3, 4))
array([[0., 0., 0., 0.],
[0., 0., 0., 0.],
[0., 0., 0., 0.]])
>>> np.ones((2, 3, 4), dtype=np.int16)
array([[[1, 1, 1, 1],
[1, 1, 1, 1],
[1, 1, 1, 1]],
[[1, 1, 1, 1],
[1, 1, 1, 1],
[1, 1, 1, 1]]], dtype=int16)
>>> np.empty((2, 3))
array([[3.73603959e-262, 6.02658058e-154, 6.55490914e-260], # 可能会变化
[5.30498948e-313, 3.14673309e-307, 1.00000000e+000]])
为了创建数字序列,NumPy 提供了 arange
函数,它类似于 Python 内置的 range
,但返回的是一个数组。
>>> np.arange(10, 30, 5)
array([10, 15, 20, 25])
>>> np.arange(0, 2, 0.3) # 它接受浮点参数
array([0. , 0.3, 0.6, 0.9, 1.2, 1.5, 1.8])
当 arange
使用浮点参数时,由于浮点精度有限,通常无法预测得到的元素数量。因此,通常最好使用 linspace
函数,该函数接受我们想要的元素数量作为参数,而不是步长。
>>> from numpy import pi
>>> np.linspace(0, 2, 9) # 从 0 到 2 的 9 个数字
array([0. , 0.25, 0.5 , 0.75, 1. , 1.25, 1.5 , 1.75, 2. ])
>>> x = np.linspace(0, 2 * pi, 100) # 用于在多个点上评估函数
>>> f = np.sin(x)
打印数组
当你打印一个数组时,NumPy 会以类似于嵌套列表的方式显示它,但布局如下:
- 最后一个轴从左到右打印;
- 倒数第二个轴从上到下打印;
- 其余轴也从上到下打印,每个切片之间用空行分隔。
一维数组被打印成行,二维数组被打印成矩阵,三维数组被打印成矩阵的列表。
>>> a = np.arange(6) # 一维数组
>>> print(a)
[0 1 2 3 4 5]
>>>
>>> b = np.arange(12).reshape(4, 3) # 二维数组
>>> print(b)
[[ 0 1 2]
[ 3 4 5]
[ 6 7 8]
[ 9 10 11]]
>>>
>>> c = np.arange(24).reshape(2, 3, 4) # 三维数组
>>> print(c)
[[[ 0 1 2 3]
[ 4 5 6 7]
[ 8 9 10 11]]
[[12 13 14 15]
[16 17 18 19]
[20 21 22 23]]]
有关 reshape
的更多详细信息,请参见下文。
如果数组太大而无法打印,NumPy 会自动跳过数组的中间部分,只打印角落部分:
>>> print(np.arange(10000))
[ 0 1 2 ... 9997 9998 9999]
>>>
>>> print(np.arange(10000).reshape(100, 100))
[[ 0 1 2 ... 97 98 99]
[ 100 101 102 ... 197 198 199]
[ 200 201 202 ... 297 298 299]
...
[9700 9701 9702 ... 9797 9798 9799]
[9800 9801 9802 ... 9897 9898 9899]
[9900 9901 9902 ... 9997 9998 9999]]
要禁用此行为并强制 NumPy 打印整个数组,可以使用 set_printoptions
更改打印选项。
>>> np.set_printoptions(threshold=sys.maxsize) # 需要导入 sys 模块
基本操作
数组上的算术运算符是逐元素应用的。会创建一个新的数组并填充结果。
>>> a = np.array([20, 30, 40, 50])
>>> b = np.arange(4)
>>> b
array([0, 1, 2, 3])
>>> c = a - b
>>> c
array([20, 29, 38, 47])
>>> b**2
array([0, 1, 4, 9])
>>> 10 * np.sin(a)
array([ 9.12945251, -9.88031624, 7.4511316 , -2.62374854])
>>> a < 35
array([ True, True, False, False])
与许多矩阵语言不同,在 NumPy 数组中,乘法运算符 *
是逐元素运算的。矩阵乘法可以使用 @
运算符(在 Python >=3.5 中)或 dot
函数或方法来完成:
>>> A = np.array([[1, 1],
... [0, 1]])
>>> B = np.array([[2, 0],
... [3, 4]])
>>> A * B # 逐元素乘积
array([[2, 0],
[0, 4]])
>>> A @ B # 矩阵乘积
array([[5, 4],
[3, 4]])
>>> A.dot(B) # 另一种矩阵乘积
array([[5, 4],
[3, 4]])
一些操作,如 +=
和 *=
,会就地修改现有数组,而不是创建一个新的数组。
>>> rg = np.random.default_rng(1) # 创建默认随机数生成器实例
>>> a = np.ones((2, 3), dtype=int)
>>> b = rg.random((2, 3))
>>> a *= 3
>>> a
array([[3, 3, 3],
[3, 3, 3]])
>>> b += a
>>> b
array([[3.51182162, 3.9504637 , 3.14415961],
[3.94864945, 3.31183145, 3.42332645]])
>>> a += b # b 不会自动转换为整数类型
Traceback (most recent call last):
...
numpy._core._exceptions._UFuncOutputCastingError: Cannot cast ufunc 'add' output from dtype('float64') to dtype('int64') with casting rule 'same_kind'
当操作不同类型的数组时,结果数组的类型会对应更通用或更精确的类型(这种行为被称为向上转型,即 upcasting)。
>>> a = np.ones(3, dtype=np.int32)
>>> b = np.linspace(0, pi, 3)
>>> b.dtype.name
'float64'
>>> c = a + b
>>> c
array([1. , 2.57079633, 4.14159265])
>>> c.dtype.name
'float64'
>>> d = np.exp(c * 1j)
>>> d
array([ 0.54030231+0.84147098j, -0.84147098+0.54030231j,
-0.54030231-0.84147098j])
>>> d.dtype.name
'complex128'
许多一元操作,例如计算数组中所有元素的总和,作为 ndarray
类的方法实现。
>>> a = rg.random((2, 3))
>>> a
array([[0.82770259, 0.40919914, 0.54959369],
[0.02755911, 0.75351311, 0.53814331]])
>>> a.sum()
3.1057109529998157
>>> a.min()
0.027559113243068367
>>> a.max()
0.8277025938204418
默认情况下,这些操作会将数组视为一个数字列表,无论其形状如何。然而,通过指定 axis
参数,你可以在数组的指定轴上应用操作:
>>> b = np.arange(12).reshape(3, 4)
>>> b
array([[ 0, 1, 2, 3],
[ 4, 5, 6, 7],
[ 8, 9, 10, 11]])
>>>
>>> b.sum(axis=0) # 每列的和
array([12, 15, 18, 21])
>>>
>>> b.min(axis=1) # 每行的最小值
array([0, 4, 8])
>>>
>>> b.cumsum(axis=1) # 每行的累积和
array([[ 0, 1, 3, 6],
[ 4, 9, 15, 22],
[ 8, 17, 27, 38]])
通用函数(Universal functions)
NumPy 提供了许多熟悉的数学函数,例如 sin
、cos
和 exp
。在 NumPy 中,这些被称为通用函数(ufunc
)。这些函数在 NumPy 中逐元素操作数组,输出也是数组。
>>> B = np.arange(3)
>>> B
array([0, 1, 2])
>>> np.exp(B)
array([1. , 2.71828183, 7.3890561 ])
>>> np.sqrt(B)
array([0. , 1. , 1.41421356])
>>> C = np.array([2., -1., 4.])
>>> np.add(B, C)
array([2., 0., 6.])
索引、切片和迭代
一维数组可以像 Python 中的其他序列(如列表)一样进行索引、切片和迭代。
>>> a = np.arange(10)**3
>>> a
array([ 0, 1, 8, 27, 64, 125, 216, 343, 512, 729])
>>> a[2]
8
>>> a[2:5]
array([ 8, 27, 64])
>>> # 相当于 a[0:6:2] = 1000;
>>> # 从起始位置到第 6 个位置(不包括),每隔一个元素设置为 1000
>>> a[:6:2] = 1000
>>> a
array([1000, 1, 1000, 27, 1000, 125, 216, 343, 512, 729])
>>> a[::-1] # 反转 a
array([ 729, 512, 343, 216, 125, 1000, 27, 1000, 1, 1000])
>>> for i in a:
... print(i**(1 / 3.))
...
9.999999999999998 # 可能会变化
1.0
9.999999999999998
3.0
9.999999999999998
4.999999999999999
5.999999999999999
6.999999999999999
7.999999999999999
8.999999999999998
多维数组可以为每个轴提供一个索引。这些索引用逗号分隔,放在一个元组中:
>>> def f(x, y):
... return 10 * x + y
...
>>> b = np.fromfunction(f, (5, 4), dtype=int)
>>> b
array([[ 0, 1, 2, 3],
[10, 11, 12, 13],
[20, 21, 22, 23],
[30, 31, 32, 33],
[40, 41, 42, 43]])
>>> b[2, 3]
23
>>> b[0:5, 1] # b 的第二列中的每一行
array([ 1, 11, 21, 31, 41])
>>> b[:, 1] # 与上一个例子等效
array([ 1, 11, 21, 31, 41])
>>> b[1:3, :] # b 的第二行和第三行中的每一列
array([[10, 11, 12, 13],
[20, 21, 22, 23]])
如果提供的索引少于轴的数量,则缺少的索引被视为完整的切片 :
。
>>> b[-1] # 最后一行。等同于 b[-1, :]
array([40, 41, 42, 43])
在 b[i]
中的括号内的表达式被视为一个 i
,后面跟着足够多的 :
,以形成一个完整的索引元组。NumPy 还允许你使用点号(...
)来表示足够多的 :
,以形成一个完整的索引元组。例如,如果 x
是一个有 5 个轴的数组,那么:
x[1, 2, ...]
等同于x[1, 2, :, :, :]
;x[..., 3]
等同于x[:, :, :, :, 3]
;x[4, ..., 5, :]
等同于x[4, :, :, 5, :]
。
>>> c = np.array([[[ 0, 1, 2], # 一个 3D 数组(两个堆叠的 2D 数组)
... [ 10, 12, 13]],
... [[100, 101, 102],
... [110, 112, 113]]])
>>> c.shape
(2, 2, 3)
>>> c[1, ...] # 等同于 c[1, :, :] 或 c[1]
array([[100, 101, 102],
[110, 112, 113]])
>>> c[..., 2] # 等同于 c[:, :, 2]
array([[ 2, 13],
[102, 113]])
对多维数组进行迭代是针对第一个轴进行的:
>>> for row in b:
... print(row)
...
[0 1 2 3]
[10 11 12 13]
[20 21 22 23]
[30 31 32 33]
[40 41 42 43]
然而,如果要在数组的每个元素上执行操作,可以使用 flat
属性,它是一个迭代器,可以遍历数组中的所有元素:
>>> for element in b.flat:
... print(element)
...
0
1
2
3
10
11
12
13
20
21
22
23
30
31
32
33
40
41
42
43
形状操作
改变数组的形状
数组的形状由每个轴上的元素数量决定。
>>> a = np.floor(10 * rg.random((3, 4)))
>>> a
array([[3., 7., 3., 4.],
[1., 4., 2., 2.],
[7., 2., 4., 9.]])
>>> a.shape
(3, 4)
可以通过多种命令改变数组的形状。注意,以下三个命令都返回一个修改后的数组,但不会改变原始数组:
>>> a.ravel() # 返回一个展平后的数组
array([3., 7., 3., 4., 1., 4., 2., 2., 7., 2., 4., 9.])
>>> a.reshape(6, 2) # 返回一个修改形状后的数组
array([[3., 7.],
[3., 4.],
[1., 4.],
[2., 2.],
[7., 2.],
[4., 9.]])
>>> a.T # 返回一个转置后的数组
array([[3., 1., 7.],
[7., 4., 2.],
[3., 2., 4.],
[4., 2., 9.]])
>>> a.T.shape
(4, 3)
>>> a.shape
(3, 4)
从 ravel
得到的结果通常是“C 风格”的,即最右边的索引变化最快,因此 a[0, 0]
之后的元素是 a[0, 1]
。如果数组被重新塑形,同样会按照“C 风格”处理。NumPy 通常会按照这种顺序创建数组,因此 ravel
通常不需要复制其参数。但如果数组是通过对另一个数组进行切片或以不寻常的方式创建的,它可能需要被复制。ravel
和 reshape
函数也可以通过可选参数指定使用 FORTRAN 风格的数组,即最左边的索引变化最快。
reshape
函数返回其参数的一个修改形状后的版本,而 ndarray.resize
方法会直接修改数组本身:
>>> a
array([[3., 7., 3., 4.],
[1., 4., 2., 2.],
[7., 2., 4., 9.]])
>>> a.resize((2, 6))
>>> a
array([[3., 7., 3., 4., 1., 4.],
[2., 2., 7., 2., 4., 9.]])
在重塑操作中,如果给定的维度为 -1
,其他维度会自动计算:
>>> a.reshape(3, -1)
array([[3., 7., 3., 4.],
[1., 4., 2., 2.],
[7., 2., 4., 9.]])
将不同数组堆叠在一起
可以沿不同轴将多个数组堆叠在一起:
>>> a = np.floor(10 * rg.random((2, 2)))
>>> a
array([[9., 7.],
[5., 2.]])
>>> b = np.floor(10 * rg.random((2, 2)))
>>> b
array([[1., 9.],
[5., 1.]])
>>> np.vstack((a, b))
array([[9., 7.],
[5., 2.],
[1., 9.],
[5., 1.]])
>>> np.hstack((a, b))
array([[9., 7., 1., 9.],
[5., 2., 5., 1.]])
column_stack
函数可以将一维数组作为列堆叠成二维数组。对于二维数组,它与 hstack
等效:
>>> from numpy import newaxis
>>> np.column_stack((a, b)) # 对于二维数组
array([[9., 7., 1., 9.],
[5., 2., 5., 1.]])
>>> a = np.array([4., 2.])
>>> b = np.array([3., 8.])
>>> np.column_stack((a, b)) # 返回一个二维数组
array([[4., 3.],
[2., 8.]])
>>> np.hstack((a, b)) # 结果不同
array([4., 2., 3., 8.])
>>> a[:, newaxis] # 将 `a` 视为二维列向量
array([[4.],
[2.]])
>>> np.column_stack((a[:, newaxis], b[:, newaxis]))
array([[4., 3.],
[2., 8.]])
>>> np.hstack((a[:, newaxis], b[:, newaxis])) # 结果相同
array([[4., 3.],
[2., 8.]])
对于多于两个维度的数组,hstack
沿第二个轴堆叠,vstack
沿第一个轴堆叠,而 concatenate
允许通过可选参数指定沿哪个轴进行连接。
注意
在复杂情况下,r_
和 c_
对于通过沿一个轴堆叠数字来创建数组很有用。它们允许使用范围字面量 :
。
>>> np.r_[1:4, 0, 4]
array([1, 2, 3, 0, 4])
当与数组作为参数一起使用时,r_
和 c_
的行为类似于 vstack
和 hstack
,但允许通过可选参数指定沿哪个轴进行连接。
将一个数组拆分成多个较小的数组
使用 hsplit
,可以在水平方向(沿着最后一个轴)拆分数组,可以通过指定返回的等形状数组的数量,或者指定拆分的列的位置来实现。
>>> a = np.floor(10 * rg.random((2, 12)))
>>> a
array([[6., 7., 6., 9., 0., 5., 4., 0., 6., 8., 5., 2.],
[8., 5., 5., 7., 1., 8., 6., 7., 1., 8., 1., 0.]])
>>> # 将 `a` 拆分成 3 个子数组
>>> np.hsplit(a, 3)
[array([[6., 7., 6., 9.],
[8., 5., 5., 7.]]), array([[0., 5., 4., 0.],
[1., 8., 6., 7.]]), array([[6., 8., 5., 2.],
[1., 8., 1., 0.]])]
>>> # 在第三列和第四列之后拆分 `a`
>>> np.hsplit(a, (3, 4))
[array([[6., 7., 6.],
[8., 5., 5.]]), array([[9.],
[7.]]), array([[0., 5., 4., 0., 6., 8., 5., 2.],
[1., 8., 6., 7., 1., 8., 1., 0.]])]
vsplit
沿着垂直方向(沿着第一个轴)拆分数组,而 array_split
允许指定沿着哪个轴进行拆分。
拷贝与视图
在操作和处理数组时,它们的数据有时会被拷贝到一个新的数组中,有时则不会。这通常是初学者感到困惑的地方。主要有以下三种情况:
完全不拷贝
简单的赋值操作不会创建对象或其数据的副本。
>>> a = np.array([[ 0, 1, 2, 3],
... [ 4, 5, 6, 7],
... [ 8, 9, 10, 11]])
>>> b = a # 没有创建新的对象
>>> b is a # a 和 b 是同一个 ndarray 对象的两个名称
True
Python 以引用的方式传递可变对象,因此函数调用也不会创建副本。
>>> def f(x):
... print(id(x))
...
>>> id(a) # id 是对象的唯一标识符
148293216 # 可能会变化
>>> f(a)
148293216 # 可能会变化
视图(浅拷贝)
不同的数组对象可以共享相同的数据。view
方法会创建一个新的数组对象,但它们会查看相同的数据。
>>> c = a.view()
>>> c is a
False
>>> c.base is a # c 是 a 数据的视图
True
>>> c.flags.owndata
False
>>>
>>> c = c.reshape((2, 6)) # a 的形状没有改变,重新赋值后的 c 仍然是 a 的视图
>>> a.shape
(3, 4)
>>> c[0, 4] = 1234 # a 的数据会改变
>>> a
array([[ 0, 1, 2, 3],
[1234, 5, 6, 7],
[ 8, 9, 10, 11]])
对数组进行切片也会返回它的视图:
>>> s = a[:, 1:3]
>>> s[:] = 10 # s[:] 是 s 的视图。注意 s = 10 和 s[:] = 10 之间的区别
>>> a
array([[ 0, 10, 10, 3],
[1234, 10, 10, 7],
[ 8, 10, 10, 11]])
深拷贝
copy
方法会创建一个完整的数组及其数据的副本。
>>> d = a.copy() # 创建一个新的数组对象和新的数据
>>> d is a
False
>>> d.base is a # d 不与 a 共享任何内容
False
>>> d[0, 0] = 9999
>>> a
array([[ 0, 10, 10, 3],
[1234, 10, 10, 7],
[ 8, 10, 10, 11]])
在某些情况下,当原始数组不再需要时,建议在切片后调用 copy
。例如,假设 a
是一个巨大的中间结果,而最终结果 b
只包含 a
的一小部分,那么在通过切片构建 b
时应该进行深拷贝:
>>> a = np.arange(int(1e8))
>>> b = a[:100].copy()
>>> del a # 可以释放 `a` 的内存
如果使用 b = a[:100]
而不调用 copy
,那么即使执行了 del a
,a
也会因为被 b
引用而保留在内存中。
函数和方法概览
以下是一些按类别排列的有用的 NumPy 函数和方法名称。完整的列表可以在 Routines and objects by topic 中找到。
数组创建
arange
array
copy
empty
empty_like
eye
fromfile
fromfunction
identity
linspace
logspace
mgrid
ogrid
ones
ones_like
r_
zeros
zeros_like
类型转换
ndarray.astype
atleast_1d
atleast_2d
atleast_3d
mat
数组操作
array_split
column_stack
concatenate
diagonal
dsplit
dstack
hsplit
hstack
ndarray.item
newaxis
ravel
repeat
reshape
resize
squeeze
swapaxes
take
transpose
vsplit
vstack
问题查询
all
any
nonzero
where
排序
argmax
argmin
argsort
max
min
ptp
searchsorted
sort
数组操作
choose
compress
cumprod
cumsum
inner
ndarray.fill
imag
prod
put
putmask
real
sum
基本统计
cov
mean
std
var
基本线性代数
cross
dot
outer
linalg.svd
vdot
广播规则
广播允许通用函数以有意义的方式处理形状不完全相同的输入。
广播的第一条规则是:如果所有输入数组的维度数量不相同,则会在较小数组的形状前面重复添加“1”,直到所有数组的维度数量相同。
第二条广播规则确保在某个维度上大小为 1 的数组会表现得像在该维度上大小为最大形状的数组一样。对于“广播”数组,假设其在该维度上的数组元素值是相同的。
应用广播规则后,所有数组的大小必须匹配。更多详细信息可以在 Broadcasting 中找到。
高级索引与索引技巧
NumPy 提供了比普通 Python 序列更强大的索引功能。除了前面提到的通过整数和切片进行索引之外,数组还可以通过整数数组和布尔数组进行索引。
通过索引数组进行索引
>>> a = np.arange(12)**2 # 前 12 个平方数
>>> i = np.array([1, 1, 3, 8, 5]) # 一个索引数组
>>> a[i] # a 中位置为 i 的元素
array([ 1, 1, 9, 64, 25])
>>>
>>> j = np.array([[3, 4], [9, 7]]) # 一个二维索引数组
>>> a[j] # 结果形状与 j 相同
array([[ 9, 16],
[81, 49]])
当被索引的数组 a
是多维的时,单个索引数组指的是 a
的第一个维度。以下示例通过将标签图像转换为颜色图像,展示了这种行为,使用了调色板。
>>> palette = np.array([[0, 0, 0], # 黑色
... [255, 0, 0], # 红色
... [0, 255, 0], # 绿色
... [0, 0, 255], # 蓝色
... [255, 255, 255]]) # 白色
>>> image = np.array([[0, 1, 2, 0], # 每个值对应调色板中的一种颜色
... [0, 3, 4, 0]])
>>> palette[image] # 生成一个 (2, 4, 3) 的颜色图像
array([[[ 0, 0, 0],
[255, 0, 0],
[ 0, 255, 0],
[ 0, 0, 0]],
[[ 0, 0, 0],
[ 0, 0, 255],
[255, 255, 255],
[ 0, 0, 0]]])
我们还可以为多个维度提供索引数组。每个维度的索引数组必须具有相同的形状。
>>> a = np.arange(12).reshape(3, 4)
>>> a
array([[ 0, 1, 2, 3],
[ 4, 5, 6, 7],
[ 8, 9, 10, 11]])
>>> i = np.array([[0, 1], # 第一维度的索引
... [1, 2]])
>>> j = np.array([[2, 1], # 第二维度的索引
... [3, 3]])
>>>
>>> a[i, j] # i 和 j 必须具有相同的形状
array([[ 2, 5],
[ 7, 11]])
>>>
>>> a[i, 2]
array([[ 2, 6],
[ 6, 10]])
>>>
>>> a[:, j]
array([[[ 2, 1],
[ 3, 3]],
[[ 6, 5],
[ 7, 7]],
[[10, 9],
[11, 11]]])
在 Python 中,arr[i, j]
与 arr[(i, j)]
是完全相同的——因此我们可以将 i
和 j
放在一个元组中,然后用它进行索引。
>>> l = (i, j)
>>> # 等同于 a[i, j]
>>> a[l]
array([[ 2, 5],
[ 7, 11]])
然而,我们不能通过将 i
和 j
放入一个数组中来实现这一点,因为这个数组会被解释为对 a
的第一维度进行索引。
>>> s = np.array([i, j])
>>> # 这不是我们想要的
>>> a[s]
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
IndexError: index 3 is out of bounds for axis 0 with size 3
>>> # 与 `a[i, j]` 相同
>>> a[tuple(s)]
array([[ 2, 5],
[ 7, 11]])
另一个常见的使用索引数组的场景是查找时间序列的最大值:
>>> time = np.linspace(20, 145, 5) # 时间尺度
>>> data = np.sin(np.arange(20)).reshape(5, 4) # 4 个时间依赖序列
>>> time
array([ 20. , 51.25, 82.5 , 113.75, 145. ])
>>> data
array([[ 0. , 0.84147098, 0.90929743, 0.14112001],
[-0.7568025 , -0.95892427, -0.2794155 , 0.6569866 ],
[ 0.98935825, 0.41211849, -0.54402111, -0.99999021],
[-0.53657292, 0.42016704, 0.99060736, 0.65028784],
[-0.28790332, -0.96139749, -0.75098725, 0.14987721]])
>>> # 每个序列的最大值索引
>>> ind = data.argmax(axis=0)
>>> ind
array([2, 0, 3, 1])
>>> # 最大值对应的时间
>>> time_max = time[ind]
>>>
>>> data_max = data[ind, range(data.shape[1])] # => data[ind[0], 0], data[ind[1], 1]...
>>> time_max
array([ 82.5 , 20. , 113.75, 51.25])
>>> data_max
array([0.98935825, 0.84147098, 0.99060736, 0.6569866 ])
>>> np.all(data_max == data.max(axis=0))
True
你还可以使用索引数组作为目标进行赋值:
>>> a = np.arange(5)
>>> a
array([0, 1, 2, 3, 4])
>>> a[[1, 3, 4]] = 0
>>> a
array([0, 0, 2, 0, 0])
然而,当索引列表中包含重复项时,赋值会多次执行,最终保留最后一次赋值的值:
>>> a = np.arange(5)
>>> a[[0, 0, 2]] = [1, 2, 3]
>>> a
array([2, 1, 3, 3, 4])
这种行为是合理的,但要注意,如果你使用 Python 的 +=
构造,它可能不会按你期望的方式工作:
>>> a = np.arange(5)
>>> a[[0, 0, 2]] += 1
>>> a
array([1, 1, 3, 3, 4])
即使 0 在索引列表中出现了两次,第 0 个元素也只被加了 1 次。这是因为 Python 要求 a += 1
等同于 a = a + 1
。
通过布尔数组进行索引
当我们使用整数索引数组时,我们提供了一个要选择的索引列表。布尔索引的方式则不同;我们明确地选择数组中我们想要的元素,以及我们不想要的元素。
布尔索引最自然的方式是使用与原数组形状相同的布尔数组:
>>> a = np.arange(12).reshape(3, 4)
>>> b = a > 4
>>> b # b 是一个与 a 形状相同的布尔数组
array([[False, False, False, False],
[False, True, True, True],
[ True, True, True, True]])
>>> a[b] # 一维数组,包含被选中的元素
array([ 5, 6, 7, 8, 9, 10, 11])
这种特性在赋值操作中非常有用:
>>> a[b] = 0 # a 中所有大于 4 的元素都变成 0
>>> a
array([[0, 1, 2, 3],
[4, 0, 0, 0],
[0, 0, 0, 0]])
你可以查看以下示例,了解如何使用布尔索引生成曼德勃罗集(Mandelbrot set)的图像:
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> def mandelbrot(h, w, maxit=20, r=2):
... """返回一个大小为 (h,w) 的曼德勃罗分形图像。"""
... x = np.linspace(-2.5, 1.5, 4*h+1)
... y = np.linspace(-1.5, 1.5, 3*w+1)
... A, B = np.meshgrid(x, y)
... C = A + B*1j
... z = np.zeros_like(C)
... divtime = maxit + np.zeros(z.shape, dtype=int)
...
... for i in range(maxit):
... z = z**2 + C
... diverge = abs(z) > r # 哪些点发散
... div_now = diverge & (divtime == maxit) # 哪些点现在发散
... divtime[div_now] = i # 记录发散时间
... z[diverge] = r # 避免过度发散
...
... return divtime
>>> plt.clf()
>>> plt.imshow(mandelbrot(400, 400))
布尔索引的第二种方式更类似于整数索引;对于数组的每个维度,我们提供一个一维布尔数组来选择所需的切片:
>>> a = np.arange(12).reshape(3, 4)
>>> b1 = np.array([False, True, True]) # 第一维度的选择
>>> b2 = np.array([True, False, True, False]) # 第二维度的选择
>>>
>>> a[b1, :] # 选择行
array([[ 4, 5, 6, 7],
[ 8, 9, 10, 11]])
>>>
>>> a[b1] # 与上一个例子相同
array([[ 4, 5, 6, 7],
[ 8, 9, 10, 11]])
>>>
>>> a[:, b2] # 选择列
array([[ 0, 2],
[ 4, 6],
[ 8, 10]])
>>>
>>> a[b1, b2] # 一个比较奇怪的操作
array([ 4, 10])
注意,一维布尔数组的长度必须与你要切片的维度(或轴)的长度一致。在上面的例子中,b1
的长度为 3(a
的行数),而 b2
(长度为 4)适用于 a
的第二维度(列)。
ix_()
函数
ix_()
函数可以用来组合不同的向量,从而为每个 n 元组获得结果。例如,如果你想计算所有从每个向量 a、b 和 c 中取出的三元组的 a + b * c
:
>>> a = np.array([2, 3, 4, 5])
>>> b = np.array([8, 5, 4])
>>> c = np.array([5, 4, 6, 8, 3])
>>> ax, bx, cx = np.ix_(a, b, c)
>>> ax
array([[[2]],
[[3]],
[[4]],
[[5]]])
>>> bx
array([[[8],
[5],
[4]]])
>>> cx
array([[[5, 4, 6, 8, 3]]])
>>> ax.shape, bx.shape, cx.shape
((4, 1, 1), (1, 3, 1), (1, 1, 5))
>>> result = ax + bx * cx
>>> result
array([[[42, 34, 50, 66, 26],
[27, 22, 32, 42, 17],
[22, 18, 26, 34, 14]],
[[43, 35, 51, 67, 27],
[28, 23, 33, 43, 18],
[23, 19, 27, 35, 15]],
[[44, 36, 52, 68, 28],
[29, 24, 34, 44, 19],
[24, 20, 28, 36, 16]],
[[45, 37, 53, 69, 29],
[30, 25, 35, 45, 20],
[25, 21, 29, 37, 17]]])
>>> result[3, 2, 4]
17
>>> a[3] + b[2] * c[4]
17
你也可以通过以下方式实现归约操作:
>>> def ufunc_reduce(ufct, *vectors):
... vs = np.ix_(*vectors)
... r = ufct.identity
... for v in vs:
... r = ufct(r, v)
... return r
然后可以这样使用它:
>>> ufunc_reduce(np.add, a, b, c)
array([[[15, 14, 16, 18, 13],
[12, 11, 13, 15, 10],
[11, 10, 12, 14, 9]],
[[16, 15, 17, 19, 14],
[13, 12, 14, 16, 11],
[12, 11, 13, 15, 10]],
[[17, 16, 18, 20, 15],
[14, 13, 15, 17, 12],
[13, 12, 14, 16, 11]],
[[18, 17, 19, 21, 16],
[15, 14, 16, 18, 13],
[14, 13, 15, 17, 12]]])
与普通的 ufunc.reduce
相比,这种版本的归约操作的优势在于它利用了广播规则,从而避免创建一个大小为输出乘以向量数量的参数数组。
通过字符串进行索引
参见 结构化数组。
技巧与提示
这里列出了一些简短且有用的技巧。
“自动”重塑
要改变数组的维度,可以省略其中一个大小,它将被自动推导出来:
>>> a = np.arange(30)
>>> b = a.reshape((2, -1, 3)) # -1 表示“需要多少就是多少”
>>> b.shape
(2, 5, 3)
>>> b
array([[[ 0, 1, 2],
[ 3, 4, 5],
[ 6, 7, 8],
[ 9, 10, 11],
[12, 13, 14]],
[[15, 16, 17],
[18, 19, 20],
[21, 22, 23],
[24, 25, 26],
[27, 28, 29]]])
向量堆叠
如何从一系列等长的行向量构造一个二维数组?在 MATLAB 中这很容易实现:如果 x
和 y
是两个长度相同的向量,你只需要执行 m = [x; y]
。在 NumPy 中,这可以通过 column_stack
、dstack
、hstack
和 vstack
函数实现,具体取决于堆叠的维度。
>>> x = np.arange(0, 10, 2)
>>> y = np.arange(5)
>>> m = np.vstack([x, y])
>>> m
array([[0, 2, 4, 6, 8],
[0, 1, 2, 3, 4]])
>>> xy = np.hstack([x, y])
>>> xy
array([0, 2, 4, 6, 8, 0, 1, 2, 3, 4])
对于超过两个维度的数组,这些函数的逻辑可能会显得有些复杂。
直方图
NumPy 的 histogram
函数应用于一个数组时,会返回一对向量:数组的直方图和一个表示箱体边界的向量。注意:matplotlib
也有一个用于构建直方图的函数(名为 hist
,与 MATLAB 类似),它与 NumPy 中的函数有所不同。主要区别在于,pylab.hist
会自动绘制直方图,而 numpy.histogram
只生成数据。
>>> import numpy as np
>>> rg = np.random.default_rng(1)
>>> import matplotlib.pyplot as plt
>>> # 构建一个均值为 2、标准差为 0.5 的正态分布随机数向量,包含 10000 个样本
>>> mu, sigma = 2, 0.5
>>> v = rg.normal(mu, sigma, 10000)
>>> # 使用 matplotlib 绘制包含 50 个箱体的归一化直方图
>>> plt.hist(v, bins=50, density=True) # 使用 matplotlib 绘图
>>> # 使用 NumPy 计算直方图,然后绘制
>>> (n, bins) = np.histogram(v, bins=50, density=True) # NumPy 版本(不绘图)
>>> plt.plot(.5 * (bins[1:] + bins[:-1]), n)
如果你使用的是 Matplotlib >= 3.4,也可以使用 plt.stairs(n, bins)
。
进一步阅读
风险提示与免责声明
本文内容基于公开信息研究整理,不构成任何形式的投资建议。历史表现不应作为未来收益保证,市场存在不可预见的波动风险。投资者需结合自身财务状况及风险承受能力独立决策,并自行承担交易结果。作者及发布方不对任何依据本文操作导致的损失承担法律责任。市场有风险,投资须谨慎。