理乐志

chihchao's blog

0%

NumPy学习笔记-2-进阶

发表于 分类于

NumPy的进阶操作，包括常用的数学通用函数、聚合、广播、排序以及花哨的索引等。

1 通用函数

1.1 通用函数的介绍

NumPy数组的计算有时非常快，有时也非常慢，关键是利用向量化操作。NumPy中向量操作是通过通用函数（ufunc）实现的，通用函数的主要目的是对NumPy数组中的值执行更快的重复操作。

NumPy 通用函数的使用方式非常自然，因为它用到了Python原生的算术运算符，标准的加、减、乘、除都可以直接使用：

1
import numpy as np
1
2
3
4
5
6
7
x = np.arange(1, 5)

print('x     =', x)
print('x + 5 =', x + 5)
print('x - 5 =', x - 5)
print('x * 2 =', x * 2)
print('x / 2 =', x / 2)
x     = [1 2 3 4]
x + 5 = [6 7 8 9]
x - 5 = [-4 -3 -2 -1]
x * 2 = [2 4 6 8]
x / 2 = [0.5 1.  1.5 2. ]

其他的一元通用函数还有取负数、取整除、取余数、求幂等等：

1
2
3
4
print('-x     =', - x)
print('x // 2 =', x // 2)  # 返回整除数
print('x % 2  =', x % 2)  # 返回余数
print('x ** 2 =', x ** 2)
-x     = [-1 -2 -3 -4]
x // 2 = [0 1 1 2]
x % 2  = [1 0 1 0]
x ** 2 = [ 1  4  9 16]

还可以将这些运算符组合使用：

1
- (0.5 * x + 1) ** 2
array([-2.25, -4.  , -6.25, -9.  ])

所有这些算术运算符都是NumPy内置函数的简单封装器，例如”+“运算符就是一个add函数的封装器：

1
np.add(x, 5)
array([6, 7, 8, 9])

下表为NumPy 实现的算术运算符：

运算符 对应的通用函数 描述
+ np.add 加法运算
- np.subtract 减法运算
* np.multiply 乘法运算
/ np.divide 除法运算
// np.floor_divide 地板除法运算（取整除）
% np.mod 模运算（取余数）
- np.negative 负数运算

1.2 指数和对数

NumPy中对应的指数函数分别np.power(x,y)、np.exp(x)、np.exp2(x)、np.expm1(x)、np.sqrt(x)等；其中np.expm1(x) 表示exp(x)-1 ，当x值非常小时，计算精度更高。

1
2
3
4
5
6
7
x = [1, 2, 3]

print('x:  ', x)
print('e^x:', np.exp(x))
print('2^x:', np.exp2(x))
print('3^x:', np.power(3, x))
print('x^0.5:', np.sqrt(x))
x:   [1, 2, 3]
e^x: [ 2.71828183  7.3890561  20.08553692]
2^x: [2. 4. 8.]
3^x: [ 3  9 27]
x^0.5: [1.         1.41421356 1.73205081]
1
2
3
4
x = 0.0000001

print('expm1(x):', np.expm1(x))
print('exp(x)-1:', np.exp(x)-1)
expm1(x): 1.0000000500000016e-07
exp(x)-1: 1.0000000494336803e-07

对数函数分别有：np.ln(x)、np.log2(x)、np.log10(x)、np.log1p(x) 等；其中，np.log1p(x) 表示np.log(1+x) ，当x值非常小时，计算精度更高。

1
2
3
4
5
6
x = [1, 2, 4, 10]

print('x:       ', x)
print('ln(x):   ', np.log(x))
print('log2(x): ', np.log2(x))
print('log10(x):', np.log10(x))
x:        [1, 2, 4, 10]
ln(x):    [0.         0.69314718 1.38629436 2.30258509]
log2(x):  [0.         1.         2.         3.32192809]
log10(x): [0.         0.30103    0.60205999 1.        ]

需要注意的是，NumPy中并没有np.log(x,base)函数，因为np.log的第二个参数不是base而是out array，如果只是执行普通的log操作，可以选择使用np.math.log(x, base) 或使用python自带的math模块的log函数。

1
np.math.log(81, 3)
4.0

1.3 角度转换与三角函数

1
2
3
4
a = np.degrees(np.pi)
b = np.radians(180)

print('a:', a, '\nb:', b)
a: 180.0 
b: 3.141592653589793
1
2
3
4
5
6
7
8
theta = np.linspace(0, np.pi, 3)

sin_theta = np.sin(theta)
cos_theta = np.cos(theta)
tan_theta = np.tan(theta)

print('theta:', theta, '\nsin_theta:', sin_theta,
      '\ncos_theta:', cos_theta, '\ntan_theta:', tan_theta)
theta: [0.         1.57079633 3.14159265] 
sin_theta: [0.0000000e+00 1.0000000e+00 1.2246468e-16] 
cos_theta: [ 1.000000e+00  6.123234e-17 -1.000000e+00] 
tan_theta: [ 0.00000000e+00  1.63312394e+16 -1.22464680e-16]

以上结果因为是在机器精度内计算，所以有些本应该是0的值并没有精确到0。

1
2
3
4
5
6
7
8
values = [-1, 0, 1]

arcsin_values = np.arcsin(values)
arccos_values = np.arccos(values)
arctan_values = np.arctan(values)

print('values:', values, '\narcsin_values:', arcsin_values,
      '\narccos_values:', arccos_values, '\narctan_values', arctan_values)
values: [-1, 0, 1] 
arcsin_values: [-1.57079633  0.          1.57079633] 
arccos_values: [3.14159265 1.57079633 0.        ] 
arctan_values [-0.78539816  0.          0.78539816]

1.4 其他

NumPy中还有其他函数运算和一些常量值。

  • 一些常量值
1
2
3
4
print('Pi:', np.pi)
print('e:', np.e)
print('NaN:', np.nan)
print('inf:', np.inf)
Pi: 3.141592653589793
e: 2.718281828459045
NaN: nan
inf: inf
  • 求绝对值

求绝对值的运算操作，只是没有对应的运算符而已。求绝对值的函数为：np.absolute() ，可简写为np.abs() 。此函数如果用来处理复数，则返回复数的模。

1
2
print(np.absolute(x))
print(np.abs(x))
[ 1  2  4 10]
[ 1  2  4 10]
1
2
y = np.array([3-4j, 4-3j, 2, 1j])
np.abs(y)
array([5., 5., 2., 1.])
  • 取整

np.ceil(x) 向上取整，即返回大于或者等于 x 的最小整数；np.floor(x) 向下取整，即返回小于或等于 x 的最大整数。

1
2
print(np.ceil(3.14))
print(np.floor(3.14))
4.0
3.0
  • 分离浮点数

np.modf(x) 返回浮点数x的整数和小数部分。

1
np.modf(3.14)
(0.14000000000000012, 3.0)

2 高级通用函数特性

2.1 指定输出

在进行大量运算时，有时候指定一个用于存放运算结果的数组是非常有用的。不同于创建临时数组，你可以用这个特性将计算结果直接写入你期望的存储位置。所有的通用函数都可以通过 out 参数来指定计算结果的存放位置：

1
2
3
4
x = np.arange(5)
y = np.zeros(5)
np.multiply(x, 10, out=y)
print('x:', x, '\ny:', y)
x: [0 1 2 3 4] 
y: [ 0. 10. 20. 30. 40.]

这个特性也可以被用作数组视图，例如将计算结果写入特定数组的每隔一个元素的位置：

1
2
3
y = np.zeros(10)
np.power(2, x, out=y[::2])
print('y:', y)
y: [ 1.  0.  2.  0.  4.  0.  8.  0. 16.  0.]

如果这里写的是 y[::2]=2**x ，那么结果将是创建一个临时数组，该数组存放的是 2**x 的结果，并且接下来会将这些值复制到y数组中。对于上述例子中比较小的计算量来说，着两种方式的差别并不大。但是对于较大的数组，通过慎重使用 out 参数将能够有效节约内存。

2.2 聚合

二元通用函数有些非常有趣的聚合功能，这些聚合可以直接在对象上计算。例如，如果我们希望用一个特定的运算 reduce 一个数组，那么可以用任何通用函数的 reduce 方法。一个 reduce 方法会对给定的元素和操作重复执行，直至得到单个的结果。

例如，对 add 通用函数调用 reduce 方法会返回数组中所有元素的和：

1
2
x = np.arange(1, 6)
np.add.reduce(x)
15

同样，对 multiply 通用函数调用 reduce 方法会返回数组中所有元素的乘积：

1
np.multiply.reduce(x)
120

如果需要存储每次计算的中间结果，可以使用 accumulate ：

1
np.add.accumulate(x)
array([ 1,  3,  6, 10, 15], dtype=int32)
1
np.multiply.accumulate(x)
array([  1,   2,   6,  24, 120], dtype=int32)

注意 ：在一些特殊情况中，NumPy 提供了专用的函数（np.sum、 np.prod、np.cumsum、np.cumprod），他们也可以实现以上 reduce 的功能：

1
2
3
4
print(np.sum(x))
print(np.cumsum(x))
print(np.prod(x))
print(np.cumprod(x))
15
[ 1  3  6 10 15]
120
[  1   2   6  24 120]

2.3 外积

最后，任何通用函数都可以用 outer 方法获取连个不同输入数组所有元素对的函数运算结果。这意味着你可以用一行代码实现一个乘法表：

1
2
x = np.arange(1, 10)
np.multiply.outer(x, x)
array([[ 1,  2,  3,  4,  5,  6,  7,  8,  9],
       [ 2,  4,  6,  8, 10, 12, 14, 16, 18],
       [ 3,  6,  9, 12, 15, 18, 21, 24, 27],
       [ 4,  8, 12, 16, 20, 24, 28, 32, 36],
       [ 5, 10, 15, 20, 25, 30, 35, 40, 45],
       [ 6, 12, 18, 24, 30, 36, 42, 48, 54],
       [ 7, 14, 21, 28, 35, 42, 49, 56, 63],
       [ 8, 16, 24, 32, 40, 48, 56, 64, 72],
       [ 9, 18, 27, 36, 45, 54, 63, 72, 81]])

通用函数另外一个非常有用的特性是它能操作不同大小和形状的数组，一组这样的操作被称为广播（broadcasting）。有关通用函数更多的信息（包括可用的通用函数的完整列表）可以在 NumPy 和 SciPy 文档的网站中找到。

3 聚合：最小值、最大值和其他值

在介绍 NumPy 时讲到，NumPy计算要比Python内置函数计算快的多，下面以实例验证这一点：

1
2
3
4
big_array = np.random.rand(1000000)

%timeit sum(big_array)
%timeit np.sum(big_array)
131 ms ± 2.39 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
712 µs ± 94.2 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

由结果可以看出，在计算大数据时，NumPy 比 Python 内置函数快了整整一个量级！NumPy 中内置聚合函数还有很多：

1
2
print(np.min(big_array))
print(np.max(big_array))
1.959273376650472e-06
0.9999965430780281

对于 min、max、sum 和其他 NumPy 聚合，一种更简洁的语法形式是对数组对象直接调用这些方法：

1
2
3
print(big_array.min())
print(big_array.max())
print(big_array.sum())
1.959273376650472e-06
0.9999965430780281
500131.080477

在实际应用中经常会遇到沿着一行或一列聚合的情况，此时可通过指定轴（axis）来进行：

1
2
a = np.random.random((3, 4))
print('a:\n', a)
a:
 [[0.45288685 0.52669344 0.68856631 0.80570042]
 [0.88351024 0.90855213 0.02373143 0.41029715]
 [0.99508181 0.56756179 0.93662538 0.12835504]]
1
2
3
4
5
6
7
8
# 每一列的最小值
print(a.min(axis=0))

# 每一行的最小值
print(a.min(axis=1))

# 整个数组的最小值
print(a.min())
[0.45288685 0.52669344 0.02373143 0.12835504]
[0.45288685 0.02373143 0.12835504]
0.023731431665718228

如果数组中有缺失值（NaN），那么处理时就需要采用安全处理策略（NaN-safe），即计算时忽略所有的缺失值，否则会影响输出结果：

1
2
a[0, 0] = np.nan
print('a:\n', a)
a:
 [[       nan 0.52669344 0.68856631 0.80570042]
 [0.88351024 0.90855213 0.02373143 0.41029715]
 [0.99508181 0.56756179 0.93662538 0.12835504]]
1
2
3
# 缺失值的存在会影响函数的输出结果

print(a.min())
nan

在操作时改为NaN安全版本，只需要在函数前加上nan即可，如：np.min() 变为 np.nanmin() 等。需要注意的是，采用安全版本时，只能使用 np.nan... 形式的函数，不能再使用数组对象直接调用方法的形式，否则会引起报错。

1
2
3
print(np.nanmin(a, axis=0))
print(np.nanmin(a, axis=1))
print(np.nanmin(a))
[0.88351024 0.52669344 0.02373143 0.12835504]
[0.52669344 0.02373143 0.12835504]
0.023731431665718228

NumPy 中提供了很多聚合函数，详见下表：

函数名称 NaN安全版本 描述
np.sum np.nansum 计算元素的和
np.prod np.nanprod 计算元素的积
np.mean np.nanmean 计算元素的平均值
np.std np.nanstd 计算元素的标准差
np.var np.nanvar 计算元素的方差
np.min np.nanmin 找出最小值
np.max np.nanmax 找出最大值
np.argmin np.nanargmin 找出最小值的索引
np.argmax np.nanargmax 找出最大值的索引
np.median np.nanmedian 计算元素的中位数
np.percentile np.nanpercentile 计算基于元素排序的统计值
np.any N/A 验证任何一个元素是否为真
np.all N/A 验证所有元素是否为真

下面以一组身高数据练习常用聚合函数：

1
2
heights = np.random.randint(150, 190, size=100)
print(heights)
[161 160 183 162 167 150 171 157 177 182 168 166 170 172 151 154 169 171
 187 174 171 161 186 189 169 183 175 152 167 152 185 172 187 186 165 172
 155 151 169 184 151 164 162 179 188 159 156 176 161 155 174 166 168 180
 175 166 185 187 180 162 152 185 151 155 152 157 178 169 169 151 161 173
 166 183 168 161 172 165 155 181 176 172 152 182 175 167 181 177 175 170
 171 164 161 160 178 152 153 185 176 164]
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
# 身高平均值
print("Mean height:", heights.mean())

# 身高标准差
print("Standard deviation:", heights.std())

# 身高最低值
print("Minimun height:", heights.min())

# 身高最高值
print("Maximun height:", heights.max())

# 身高值的25%分位数
print("25th percentile:", np.percentile(heights, 25))

# 身高值的中位数
print("Median:", np.median(heights))

# 身高值的75%分位数
print("75th percentile:", np.percentile(heights, 75))
Mean height: 168.72
Standard deviation: 11.048149166263098
Minimun height: 150
Maximun height: 189
25th percentile: 161.0
Median: 169.0
75th percentile: 177.0
1
2
3
4
5
6
7
8
9
10
import seaborn
import matplotlib.pyplot as plt
%matplotlib inline

seaborn.set()

plt.hist(heights)
plt.title("Height Distribution")
plt.xlabel("height(cm)")
plt.ylabel("number")

4 数组的计算：广播

4.1 广播的介绍

前面介绍了 NumPy 如何通过通用函数的向量化操作来减少缓慢的Python循环、另外一种向量化操作的方法是利用NumPy的广播功能。广播可以简单理解为用于不同大小数组的二进制通用函数（加、减、乘）的一组规则。广播允许这些二进制操作可以用于不同大小的数组。下面用三个例子进行展示NumPy的功能。

  • 例一

将一个标量（可以认为是一个零维数组）和一个一维数组相加：

1
2
3
4
5
a = np.arange(3)
print('a:\n', a, '\na_ndim:', a.ndim, '\na_shape:', a.shape)

b = a + 5
print('\nb:\n', b, '\nb_ndim:', b.ndim, '\nb_shape:', b.shape)
a:
 [0 1 2] 
a_ndim: 1 
a_shape: (3,)

b:
 [5 6 7] 
b_ndim: 1 
b_shape: (3,)

我们可以认为这个操作是将数值5扩展或重复至数组 [5,5,5] ，然后执行加法，NumPy广播功能的好处是，这种对值的重复实际上并没有发生，但是这是一种很好用的理解广播的模型。

我们同样可以将这个原理扩展到更高维度的数组。以下将一个一维数组和二维数组相加的结果。

  • 例二
1
2
3
4
5
6
7
8
a = np.arange(3)
print('a:\n', a, '\na_ndim:', a.ndim, '\na_shape:', a.shape)

b = np.ones((3, 3))
print('\nb:\n', b, '\nb_ndim:', b.ndim, '\nb_shape:', b.shape)

c = a + b
print('\nc:\n', c, '\nc_ndim:', c.ndim, '\nc_shape:', c.shape)
a:
 [0 1 2] 
a_ndim: 1 
a_shape: (3,)

b:
 [[1. 1. 1.]
 [1. 1. 1.]
 [1. 1. 1.]] 
b_ndim: 2 
b_shape: (3, 3)

c:
 [[1. 2. 3.]
 [1. 2. 3.]
 [1. 2. 3.]] 
c_ndim: 2 
c_shape: (3, 3)

这里这个一维数组就被扩展或者广播了。它被扩展到匹配b数组的形状。以上这些例子理解起来相对容易，更复杂的情况会涉及对两个数组的同时广播。

  • 例三
1
2
3
4
5
6
7
8
a = np.arange(3)
print('a:\n', a, '\na_ndim:', a.ndim, '\na_shape:', a.shape)

b = np.arange(3).reshape((3, 1))
print('\nb:\n', b, '\nb_ndim:', b.ndim, '\nb_shape:', b.shape)

c = a + b
print('\nc:\n', c, '\nc_ndim:', c.ndim, '\nc_shape:', c.shape)
a:
 [0 1 2] 
a_ndim: 1 
a_shape: (3,)

b:
 [[0]
 [1]
 [2]] 
b_ndim: 2 
b_shape: (3, 1)

c:
 [[0 1 2]
 [1 2 3]
 [2 3 4]] 
c_ndim: 2 
c_shape: (3, 3)

正如此前将一个值扩展或广播以匹配另外一个数组的形状，这里将 a 和 b 都进行了扩展来匹配一个公共的形状，最终的结果是一个二维数组。下图可以更直观的展示以上三个例子的广播过程：

上图中浅色盒子表示广播的值。同样需要注意的是，这个额外的内存并没有在实际操作中进行分配，但是这样的想象方式更方便我们从概念上理解。

4.2 广播的规则

NumPy 的广播遵循一组严格的规则，设定这组规则是为了决定两个数组间的操作：

  • 规则1：如果两个数组的维度数不同，那么小维度数组的形状将会在最左边补1。

  • 规则2：如果两个数组的形状在任何一个维度上都不匹配，那么数组的形状会沿着维度为 1 的维度扩展以匹配另外一个数组的形状。

  • 规则3：如果两个数组的形状在任何一个维度上都不匹配并且没有任何一个维度等于1，那么会引发异常。

为了更清楚的理解这些规则，下面按照这些规则来回顾一下上面的示例：

在例一中，a 是一维的，标量5可以理解为零维的，按照规则1，标量5维度数更小，所以在其左边补1，此时标量5也变为一维的，此时两个维度相同，最终的形状都是一维。

在例二中，根据规则1，a的维度更小，a的形状左边补1，变为 (1,3)；b的形状仍为 (3,3) ；根据规则2，第一个维度不匹配，因此扩展这个维度以匹配数组，此时a 的形状变为 (3,3) 。a 和 b 最终的形状都为 (3,3) ，相加后的形状也为 (3,3) 。

在例三中，a为1维，b为2维，a 和 b 的形状分别为 (3,) 和 (3, 1)。 根据规则1，a的形状左边补1，变为 (1,3)；根据规则2，需要更新这两个数组的维度来相互匹配，a的形状更新为 (3,3)，b的形状也更新为 (3,3)。a 和 b 的形状相同，最终相加后的形状也是 (3,3) 。

上面三个例子都是符合规则的运算，下面看一个因不符合规则而引发错误的例子：

1
2
3
4
5
a = np.arange(3)
print('a:\n', a, '\na_ndim:', a.ndim, '\na_shape:', a.shape)

b = np.ones((3, 2))
print('\nb:\n', b, '\nb_ndim:', b.ndim, '\nb_shape:', b.shape)
a:
 [0 1 2] 
a_ndim: 1 
a_shape: (3,)

b:
 [[1. 1.]
 [1. 1.]
 [1. 1.]] 
b_ndim: 2 
b_shape: (3, 2)

在此例中，根据规则1，a 的形状左边补1，变为 (1,3)；b的形状仍为 (3,2) ；根据规则2，a 的形状扩展为 (3,3)；而此时，两个数组的第一个维度是匹配的，但第二个维度是不匹配的，而且没有任何一个维度等于1。因此，根据规则3，a + b 相加会引发异常：

1
a + b
---------------------------------------------------------------------------

ValueError                                Traceback (most recent call last)

<ipython-input-41-bd58363a63fc> in <module>
----> 1 a + b


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

在上面的例子中主要以相加作为示例，实际上，这些广播规则对于任意二进制通用函数都是适用的。比如 logaddexp(a, b) 函数，该函数比 log(exp(a)+exp(b)) 计算的更为准确：

1
2
3
4
5
6
7
8
a = np.arange(3)[:, np.newaxis]
print('a:\n', a, '\na_ndim:', a.ndim, '\na_shape:', a.shape)

b = np.ones((3, 2))
print('\nb:\n', b, '\nb_ndim:', b.ndim, '\nb_shape:', b.shape)

c = np.logaddexp(a, b)
print('\nc:\n', c, '\nc_ndim:', c.ndim, '\nc_shape:', c.shape)
a:
 [[0]
 [1]
 [2]] 
a_ndim: 2 
a_shape: (3, 1)

b:
 [[1. 1.]
 [1. 1.]
 [1. 1.]] 
b_ndim: 2 
b_shape: (3, 2)

c:
 [[1.31326169 1.31326169]
 [1.69314718 1.69314718]
 [2.31326169 2.31326169]] 
c_ndim: 2 
c_shape: (3, 2)

4.3 广播的实际应用

广播是NumPy的核心，下面通过几个简单的示例来展示广播功能的作用。

  • 数组的归一化
1
2
3
4
5
a = np.random.random((10, 3))
print('a:\n', a)

a_mean = a.mean(axis=0)
print('\na_mean:\n', a_mean)
a:
 [[0.4943422  0.3015842  0.66961799]
 [0.36327177 0.81985191 0.10642248]
 [0.70012129 0.63598483 0.97612182]
 [0.1715966  0.86723366 0.36162404]
 [0.52428464 0.57423658 0.50148675]
 [0.81345869 0.51409075 0.94525601]
 [0.03809255 0.58624396 0.58671427]
 [0.01312058 0.93794146 0.64981756]
 [0.68367647 0.70328816 0.87909459]
 [0.63869101 0.91279874 0.13474207]]

a_mean: 
 [0.44406558 0.68532543 0.58108976]

现在通过a数组的元素中减去这个均值实现归一化（该操作是一个广播过程）。为了进一步核查处理过程是否正确，可以查看归一化的数组的均值是否接近0：

1
2
a_centered = a - a_mean
a_centered.mean(axis=0)
array([-6.66133815e-17,  2.22044605e-17, -3.33066907e-17])

在机器精度范围内，该值为0。

  • 画一个二维函数

广播另外一个非常有用的地方在于，它能基于二维函数显示图像。下面先定义一个函数 z = f(x, y) ，可以用广播沿数值区间计算该函数，并使用 Matplotlib 画出这个二维数组。

1
2
3
4
x = np.linspace(0, 5, 50)
y = np.linspace(0, 5, 50)[:, np.newaxis]

z = np.sin(x)**10 + np.cos(10+y*x) * np.cos(x)
1
2
3
4
import matplotlib.pyplot as plt

plt.imshow(z, origin='lower', extent=[0, 5, 0, 5], cmap='viridis')
plt.colorbar()

5 比较、掩码和布尔逻辑

下面学习用布尔掩码来查看和操作NumPy数组中的值。当你想基于某些准则来抽取、修改、计数或对一个数组中的值进行其他操作时，掩码就可以排上用场了。例如你可能希望统计数组中又多少值大于某一个给定值，或者删除所有超出某些限值的异常点。在NumPy中，布尔掩码通常是完成这类任务的最高效方式。

5.1 比较操作

前面介绍了通用函数，并且特别关注了算术运算符。我们看到加、减、乘、除和其他一些运算符实现了数组的逐元素操作。NumPy还实现了如小于、大于等的逐元素比较的通用函数。这些运算的结果是一个布尔数据类型的数组。一共又六种标准的比较操作。

运算符 对应的通用函数 描述
== np.equal 等于
!= np.not_equal 不等于
< np.less 小于
<= np.less_equal 小于等于
> np.greater 大于
>= np.greater_equal 大于等于 
1
2
a = np.arange(1, 6)
print('a:', a)
a: [1 2 3 4 5]
1
2
3
print('a>3:', a > 3)
print('a<3:', a < 3)
print("(2*a)==(a**2):", (2*a) == (a**2))
a>3: [False False False  True  True]
a<3: [ True  True False False False]
(2*a)==(a**2): [False  True False False False]

和算术运算通用函数一样，这些比较运算通用函数也可以用于任意形状、大小的数组。下面是一个二维数组的示例：

1
2
3
4
np.random.seed(0)
a = np.random.randint(10, size=(3, 4))
print('a:\n', a)
print('a<6:\n', a < 6)
a:
 [[5 0 3 3]
 [7 9 3 5]
 [2 4 7 6]]
a<6:
 [[ True  True  True  True]
 [False False  True  True]
 [ True  True False False]]

5.2 操作布尔数组

5.2.1 统计记录的个数

如果你需要统计布尔数组中 True 记录的个数，有两种方式可以实现，一种是直接使用 np.count_nonzero 函数；另外一种是利用 np.sum ，在这种方式中 False 会被解释成0，True会被解释成1。

1
2
3
4
# 有多少个值小于6

print(np.count_nonzero(a < 6))
print(np.sum(a < 6))
8
8

sum() 的好处是，和其他 NumPy 聚合函数一样，这个求和也可以沿着行或列进行：

1
2
3
4
5
# 每一行有多少值小于6
print(np.sum(a < 6, axis=1))

# 每一列有多少值小于6
print(np.sum(a < 6, axis=0))
[4 2 2]
[2 2 2 2]

如果要快速检查任意或者所有值是否为 True ，可以用 np.any() 或 np.all() ，这两个函数也可以沿着特定的坐标轴运用：

1
2
3
4
5
6
7
8
9
10
11
# 有没有值大于8
print('Is there a value greater than 8?', np.any(a > 8))

# 有没有值小于8
print('Is there a value less than 8?', np.any(a < 8))

# 是否所有值都等于8
print('Whether all values are equal to 6?', np.all(a == 8))

# 是否每一行的值都小于8
print('Whether the value of each row is less than 10?', np.all(a < 8, axis=1))
Is there a value greater than 8? True
Is there a value less than 8? True
Whether all values are equal to 6? False
Whether the value of each row is less than 10? [ True False  True]

5.2.2 布尔运算符

如果需要满足多个条件，可以使用Python的 逐位逻辑运算符 （或| 、且&、非~ 、相异^ ）来实现。同标准的算术运算符一样，NumPy 用通用函数重载了这些逻辑运算符，这样可以实现数组的逐位运算（通常是布尔运算）。

运算符 对应的通用函数 描述
& np.bitwise_and 按位与：如果两个相应位都为1，则该位的结果为1，否则为0
| np.bitwise_or 按位或：只要对应的两个二进位有一个为1时，结果位就为1
^ np.bitwise_xor 按位异或：对应的两个二进位相异时，结果为1
~ np.bitwise_not 按位取反：对数据的每个二进制位取反，即把1变为0，把0变为1 
1
np.sum((a > 3) & (a < 8))
6
1
np.sum(~ ((a <= 3) | (a >= 6)))
3

5.3 将布尔数组作为掩码

在前面的小节中，我们看到了如何直接对布尔数组进行聚合计算。一种更强大的模式是使用布尔数组作为掩码，通过该掩码选择数据的子数据集。如前面的介绍过的方法，利用比较运算符可以得到一个布尔数组：

1
2
3
4
np.random.seed(0)
a = np.random.randint(10, size=(3, 4))
print('a:\n', a)
print('a<5:\n', a < 5)
a:
 [[5 0 3 3]
 [7 9 3 5]
 [2 4 7 6]]
a<5:
 [[False  True  True  True]
 [False False  True False]
 [ True  True False False]]

现在为了将这些值从数组中选出，可以进行简单的索引，即掩码操作：

1
a[a < 5]
array([0, 3, 3, 3, 2, 4])

现在返回的是一个一维数组，它包含了所有满足条件的值。换句话说，所有这些值是掩码数组对应位置为 True 的值。接下来，就可以对这些值做任何操作了。

1
2
x = (a > 3) & (a < 8)
a[x]
array([5, 7, 5, 4, 7, 6])
1
2
y = ~ ((a <= 3) | (a >= 6))
a[y]
array([5, 5, 4])

5.4 关键字与逻辑运算符

关键字（and、or、not）与 逻辑操作运算符（&、|、/）是有区别的，什么时候该使用哪一种经常容易让人产生混淆。关键字判断的整个对象是真是假，而逻辑运算符是指每个对象中的比特位。

当使用and 或 or 时，就等于让Python将这个对象当作整个布尔实体。在Python中，所有非零的整数都会被当作是True：

1
bool(42), bool(0)
(True, False)
1
bool(42 and 0)
False
1
bool(42 or 0)
True

当你对整数使用 & 和 | 时，表达式操作的是元素的比特，将 and 或 or 应用于组成该数字的每个比特。& 和 | 运算时，对应的二进制比特位进行比较以得到最终结果：

1
bin(42), bin(59)
('0b101010', '0b111011')
1
bin(42 & 59)
'0b101010'
1
bin(42 | 59)
'0b111011'

当你在 NumPy中有一个布尔数组时，该数组可以被当作是由比特字符组成的、其中1=True、0=False。这样的数组可以用上面介绍的方式进行 & 和 | 的操作：

1
2
3
A = np.array([1, 0, 1, 0, 1, 0], dtype=bool)
B = np.array([1, 1, 1, 0, 1, 1], dtype=bool)
A | B
array([ True,  True,  True, False,  True,  True])

如果使用 or 来计算这两个数组，Python会计算整个数组对象的真或假，这会导致程序出错：

1
A or B
---------------------------------------------------------------------------

ValueError                                Traceback (most recent call last)

<ipython-input-66-ea2c97d9d9ee> in <module>
----> 1 A or B


ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()

同样，对给定数组进行逻辑运算时，也应该使用 & 或 | ，而不是 and 或 or：

1
2
a = np.arange(10)
(a > 4) & (a < 8)
array([False, False, False, False, False,  True,  True,  True, False,
       False])

如果试图计算整个数组的真或假，程序也同样会给出ValueError的错误：

1
(a > 4) and (a < 8)
---------------------------------------------------------------------------

ValueError                                Traceback (most recent call last)

<ipython-input-68-9b76ab08d34b> in <module>
----> 1 (a > 4) and (a < 8)


ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()

因此，and 和 or 对整个对象执行单个布尔运算，而 & 和 | 对一个对象的内容（单个比特或字节）执行多个布尔运算。对于NumPy布尔数组，后者是常用的操作。

6 花哨的索引

下面学习另外一种数组索引，也称作花哨的索引（fancy indexing）。花哨的索引与前面那些简单的索引非常类似，但是传递的是索引数组，而不是单个标量。花哨的索引让我们能够快速获得并修改复杂的数组值的子数据集。

6.1 探索花哨的索引

花哨的索引在概念上非常简单，它意味着传递一个索引数组来一次性获得多个数组元素。

1
2
3
4
rand = np.random.RandomState(42)

a = rand.randint(100, size=10)
print('a:', a)
a: [51 92 14 71 60 20 82 86 74 74]
1
[a[3], a[7], a[2]]
[71, 86, 14]
1
2
ind = [3, 7, 2]
a[ind]
array([71, 86, 14])

利用花哨的索引，结果的形状与索引数组的形状一致，而不是与被索引数组的形状一致：

1
2
3
ind = np.array([[3, 7],
                [4, 5]])
a[ind]
array([[71, 86],
       [60, 20]])

花哨的索引也对多个维度适用。和标准的索引一样，第一个索引指的是行，第二个索引指的是列：

1
2
b = np.arange(12).reshape((3, 4))
b
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])
1
2
3
4
row = np.array([0, 1, 2])
col = np.array([2, 1, 3])

b[row, col]
array([ 2,  5, 11])

在上面这个例子中，第一个值是a[0,2]，第二个值是a[1,1]，第三个值是a[2,3]。在花哨的索引中，索引值的配对遵循广播的规则。因此，当我们将一个列向量和一个行向量组合在一个索引中时，会得到一个二维的结果：

1
b[row[:, np.newaxis], col]
array([[ 2,  1,  3],
       [ 6,  5,  7],
       [10,  9, 11]])

这里特别需要记住的是，花哨的索引返回的值反映的是广播后的索引数组的形状，而不是被索引的数组的形状。

6.2 组合索引

花哨的索引可以和其他索引方案结合起来形成更强大的索引操作。

  • 花哨的索引和简单的索引组合
1
b[2, [2, 0, 1]]
array([10,  8,  9])
  • 花哨的索引和掩码组合
1
2
mask = np.array([1, 0, 1, 0], dtype=bool)
b[row[:, np.newaxis], mask]
array([[ 0,  2],
       [ 4,  6],
       [ 8, 10]])

6.3 用花哨的索引修改值

1
2
3
4
5
6
a = np.arange(10)
print('a:', a)

i = np.array([2, 1, 8, 4])
a[i] = 99
print('a:', a)
a: [0 1 2 3 4 5 6 7 8 9]
a: [ 0 99 99  3 99  5  6  7 99  9]

可以用任何的赋值操作来实现，如：

1
2
a[i] += 1
print('a:', a)
a: [  0 100 100   3 100   5   6   7 100   9]

不过需要注意，操作中重复的索引会导致一些出乎意料的结果产生：

1
2
3
4
5
a = np.zeros(10)
print('a:', a)

a[[0, 0]] = [4, 6]
print('a:', a)
a: [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]
a: [6. 0. 0. 0. 0. 0. 0. 0. 0. 0.]

上面这个操作首先赋值 a[0] = 4 ，然后赋值 a[0] = 6，因此结果 a[0] 的值为6 。以上还算合理，但是设想以下操作：

1
2
3
4
5
6
print('a:', a)

i = [2, 3, 3, 4, 4, 4]
a[i] += 1

print('a:', a)
a: [6. 0. 0. 0. 0. 0. 0. 0. 0. 0.]
a: [6. 0. 1. 1. 1. 0. 0. 0. 0. 0.]

在这个例子中，我们可能会期望 a[3] 的值为2， a[4] 的值为3，因为这是索引值重复的次数。但结果却并不是我们所预想。从概念的角度理解，这是因为 a[i] += 1 是 a[i] = a[i] +1 的简写。 a[i] +1 计算后，这个结果被赋值给了a相应的索引值。记住这个原理后，我们发现数组并没有发生多次累加，而是发生了赋值，显然这不是我们希望的结果。

如果你希望累加，该怎么做呢？可以借助通用函数中的 at() 方法来实现：

1
2
3
4
5
a = np.zeros(10)
print('a:', a)

np.add.at(a, i, 1)
print('a:', a)
a: [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]
a: [0. 0. 1. 2. 3. 0. 0. 0. 0. 0.]

at() 函数在这里对给定的操作、给定的索引（这里是i）以及给定的值（这里是1）执行的是就地操作。另一个可以实现该功能的类似方法是通用函数中的 reduceat() 函数。

6.4 示例：选择随机点

花哨的索引的一个常见用途是从一个矩阵中选择行的子集。例如，我们有一个 N × D 的矩阵，表示在D个维度的N个点。以下是一个二维正态分布的点组成的数组：

1
2
3
4
5
mean = [0, 0]
cov = [[1, 2], [2, 5]]
rand = np.random.RandomState(42)
a = rand.multivariate_normal(mean, cov, 100)
a.shape
(100, 2)
1
2
3
4
5
import matplotlib.pyplot as plt
import seaborn
seaborn.set()

plt.scatter(a[:, 0], a[:, 1])

我们利用花哨的索引随机选取20个点——选择20个随机的、不重复的索引值，并利用这些索引值选取到原始数组对应的值：

1
2
indices = np.random.choice(a.shape[0], 20, replace=False)
indices
array([26, 48, 10, 66, 96,  7, 91,  2, 90, 33, 22, 78, 75, 70, 82, 92, 16, 13, 21, 45])
1
2
selection = a[indices]
selection.shape
(20, 2)

现在看看哪些点被选中了，将选中的点在图上用大圆圈标示出来：

1
2
3
plt.scatter(a[:, 0], a[:, 1], alpha=0.3)
plt.scatter(selection[:, 0], selection[:, 1],
            facecolor='none', edgecolors='b', s=200)

这种方法通常用于快速分割数据，即需要分割训练、测试数据集以验证统计模型时，以及在解答统计问题时的抽样方法中使用。

7 数组的排序

7.1 NumPy中的快速排序

如果想在不修改原始输入数组的基础上返回一个排好序的数组，可以使用np.sort；如果希望用排好序的数组替代原始数组，可以使用数组的 sort方法：

1
2
3
a = np.array([2, 1, 4, 3, 5])
print(np.sort(a))
print(a)
[1 2 3 4 5]
[2 1 4 3 5]
1
2
a.sort()
print(a)
[1 2 3 4 5]

另外一个相关的函数是 argsort ，该函数返回的是原始数组排好序的索引值：

1
2
3
a = np.array([2, 1, 4, 3, 5])
i = np.argsort(a)
print(i)
[1 0 3 2 4]

以上结果的第一个元素是数组中最小元素的索引值，第二个值给出的次小元素的索引值，以此类推。这些索引值可以被用于（通过花哨的索引）创建有序的数组：

1
a[i]
array([1, 2, 3, 4, 5])

NumPy 排序算法的一个有用的功能是通过axis参数，沿着多维数组的行或列进行排序：

1
2
3
rand = np.random.RandomState(42)
a = rand.randint(0, 10, (4, 6))
print('a:\n', a)
a:
 [[6 3 7 4 6 9]
 [2 6 7 4 3 7]
 [7 2 5 4 1 7]
 [5 1 4 0 9 5]]
1
2
# 对a的每一列排序
np.sort(a, axis=0)
array([[2, 1, 4, 0, 1, 5],
       [5, 2, 5, 4, 3, 7],
       [6, 3, 7, 4, 6, 7],
       [7, 6, 7, 4, 9, 9]])
1
2
# 对a的每一行排序
np.sort(a, axis=1)
array([[3, 4, 6, 6, 7, 9],
       [2, 3, 4, 6, 7, 7],
       [1, 2, 4, 5, 7, 7],
       [0, 1, 4, 5, 5, 9]])

需要注意的是，这种处理方式是将行或列当作独立的数组，任何行或列的值之间的关系将会丢失。

7.2 部分排序：分隔

有时候我们不希望对整个数组进行排序，仅仅希望找到数组中第K小的值，NumPy的np.partition函数提供了该功能：

np.partition(a, kth, axis=-1, kind='introselect', order=None) a为数组，kth为整数或由整数组成的元组。该函数返回一个新数组，所有小于第k值的值会被排在K值之前，所有大于或等于第K值的值会排在后面。

1
2
b = np.array([7, 2, 3, 1, 6, 5, 4])
np.partition(b, 3)
array([2, 1, 3, 4, 6, 5, 7])

输出结果中前三个值是数组中最小的三个值，剩下的位置是原始数组剩下的值。在这两个分隔区间中，元素都是任意排列的。

与排序类似，也可以沿着多维数组任意的轴进行分离：

1
np.partition(a, 2, axis=1)
array([[3, 4, 6, 7, 6, 9],
       [2, 3, 4, 7, 6, 7],
       [1, 2, 4, 5, 7, 7],
       [0, 1, 4, 5, 9, 5]])

正如 np.argsort 函数计算的是排序的索引值，也有一个 np.argpartition函数计算的分隔的索引值：

1
2
3
print(b)
print(np.partition(b, 3))
print(np.argpartition(b, 3))
[7 2 3 1 6 5 4]
[2 1 3 4 6 5 7]
[1 3 2 6 4 5 0]

8 NumPy的结构化数组

假定现在有关于一些人的分类数据，我们需要存储这些数据，那么一种可行的方法是将他们存在三个单独的数组中：

1
2
3
name = ['Alice', 'Bob', 'Cathy', 'Doug']
age = [25, 45, 37, 19]
weight = [55.0, 85.5, 68.0, 61.5]

这种方法的一个缺点是没有任何信息告诉我们这三个数组是相关联的。如果用一种单一结构来存储所有的数据，那么看起来会更自然。NumPy可以用结构化数组实现这种存储，这些结构化数组是复合数据类型。

1
2
3
4
data = np.zeros(4, dtype={'names': ('name', 'age', 'weight'),
                          'formats': ('U10', 'i4', 'f8')})
print(data)
print(data.dtype)
[('', 0, 0.) ('', 0, 0.) ('', 0, 0.) ('', 0, 0.)]
[('name', '<U10'), ('age', '<i4'), ('weight', '<f8')]
1
2
3
4
data['name'] = name
data['age'] = age
data['weight'] = weight
print(data)
[('Alice', 25, 55. ) ('Bob', 45, 85.5) ('Cathy', 37, 68. )
 ('Doug', 19, 61.5)]

结构化数组的方便之处在于，你可以通过索引或名称查看相应的值：

1
2
3
4
5
6
7
8
# 获取所有的名字
print(data['name'])

# 获取数据第一行
print(data[0])

# 获取最后一行的名字
print(data[-1]['name'])
['Alice' 'Bob' 'Cathy' 'Doug']
('Alice', 25, 55.)
Doug

利用布尔掩码，还可以做一些更复杂的操作，如按照年龄进行筛选：

1
2
# 获取年龄小于30岁的人的名字
data[data['age'] < 30]['name']
array(['Alice', 'Doug'], dtype='<U10')

请注意，如果你希望实现比上面更复杂的操作，那么你应该考虑使用 Pandas 包，Pandas 提供了一个 DataFrame 对象，该结构是构建于 NumPy 数组之上的，提供了很多有用的数据操作功能。如果你每天都需要使用结构化数据，那么 Pandas 包是更好的选择。