目录
- NumPy在内部将数据存储在连续的内存块上,这与其他的Python内建数据结构是不同的。NumPy的算法库是用C语言写的,所以在操作数据内存时,不需要任何类型检查或者其他管理操作。NumPy数组使用的内存量也小于其他Python内建序列。
- NumPy可以针对全量数组进行复杂计算而不需要写Python循环。
# 假设一个NumPy数组包含100万个整数,还有一个同样数据内容的Python列表
>>> import numpy as np
>>> my_arr = np.arange(1000000)
>>> my_list = list(range(1000000))
1. NumPy ndarray:多维数组对象
>>> import numpy as np
# 生成随机数组
>>> data = np.random.randn(2, 3)
>>> data
array([[-0.2047, 0.4789, -0.5194],
[-0.5557, 1.9658, 1.3934]])
>>> data * 10
array([[ -2.0471, 4.7894, -5.1944],
[ -5.5573, 19.6578, 13.9341]])
>>> data + data
array([[-0.4094, 0.9579, -1.0389],
[-1.1115, 3.9316, 2.7868]])
一个ndarray是一个通用的多维同类数据容器,也就是说,它包含的每一个元素均为相同类型。每一个数组都有一个shape属性,用来表征数组每一维度的数量;每一个数组都有一个dtype属性,用来描述数组的数据类型:
>>> data.shape
(2, 3)
>>> data.dtype
dtype('float64')
1.1 生成ndarray
生成数组最简单的方式就是使用array函数。array函数接收任意的序列型对象(当然也包括其他的数组),生成一个新的包含传递数据的NumPy数组。例如,列表的转换就是一个好例子:
>>> data1 = [6, 7.5, 8, 0, 1]
>>> arr1 = np.array(data1)
>>> arr1
array([ 6. , 7.5, 8. , 0. , 1. ])
嵌套序列,例如同等长度的列表,将会自动转换成多维数组:
>>> data2 = [[1, 2, 3, 4], [5, 6, 7, 8]]
>>> arr2 = np.array(data2)
>>> arr2
array([[1, 2, 3, 4],
[5, 6, 7, 8]])
因为data2是一个包含列表的列表,所以Numpy数组arr2形成了二维数组。我们可以通过检查ndim和shape属性来确认这一点:
>>> arr2.ndim
2
>>> arr2.shape
(2, 4)
除非显式地指定,否则np.array会自动推断生成数组的数据类型。数据类型被存储在一个特殊的元数据dtype中。例如,之前的两个例子:
>>> arr1.dtype
dtype('float64')
>>> arr2.dtype
dtype('int32')
除了np.array,还有很多其他函数可以创建新数组。例如,给定长度及形状后,zeros可以一次性创造全0数组,ones可以一次性创造全1数组。empty则可以创建一个没有初始化数值的数组。想要创建高维数组,则需要为shape传递一个元组:
>>> np.zeros(10)
array([ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.])
>>> np.zeros((3, 6))
array([[ 0., 0., 0., 0., 0., 0.],
[ 0., 0., 0., 0., 0., 0.],
[ 0., 0., 0., 0., 0., 0.]])
>>> np.empty((2, 3, 2))
array([[[ 0., 0.],
[ 0., 0.],
[ 0., 0.]],
[[ 0., 0.],
[ 0., 0.],
[ 0., 0.]]])
想要使用np.empty来生成一个全0数组,并不安全,有些时候它可能会返回未初始化的垃圾数值。
arange是Python内建函数range的数组版:
>>> np.arange(15)
array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14])
由于NumPy专注于数值计算,如果没有特别指明的话,默认的数据类型是float64(浮点型)。
数组生成函数
1.2 ndarray的数据类型
数据类型,即dytpe,是一个特殊的对象,它包含了ndarray需要为某一种类型数据所申明的内存块信息(也称为元数据,即表示数据的数据)
>>> arr1 = np.array([1, 2, 3], dtype=np.float64)
>>> arr2 = np.array([1, 2, 3], dtype=np.int32)
>>> arr1.dtype
dtype('float64')
>>> arr2.dtype
dtype('int32')
dtype是NumPy能够与其他系统数据灵活交互的原因。通常,其他系统提供一个硬盘或内存与数据的对应关系,使得利用C或Fortran等底层语言读写数据变得十分方便。数据的dtype通常都是按照一个方式命名:类型名,比如float和int,后面再接上表明每个元素位数的数字。一个标准的双精度浮点值(Python中数据类型为float),将使用8字节或64位。因此,这个类型在NumPy中称为float64。
NumPy数据类型
可以使用astype方法显式地转换数组的数据类型:
>>> arr = np.array([1, 2, 3, 4, 5])
>>> arr.dtype
dtype('int32')
# 整数被转换成了浮点数
>>> float_arr = arr.astype(np.float64)
>>> float_arr.dtype
dtype('float64')
# 把浮点数转换成整数,则小数点后的部分将被消除
>>> arr = np.array([3.7, -1.2, -2.6, 0.5, 12.9, 10.1])
>>> arr
array([ 3.7, -1.2, -2.6, 0.5, 12.9, 10.1])
>>> arr.astype(np.int32)
array([ 3, -1, -2, 0, 12, 10], dtype=int32)
如果你有一个数组,里面的元素都是表达数字含义的字符串,也可以通过astype将字符串转换为数字:
>>> numeric_strings = np.array(['1.25', '-9.6', '42'], dtype=np.string_)
>>> numeric_strings.astype(float)
array([ 1.25, -9.6 , 42. ])
在NumPy中,当使用numpy.string_类型作字符串数据要小心,因为NumPy会修正它的大小或删除输入且不发出警告。pandas在处理非数值数据时有更直观的开箱型操作。
如果因为某些原因导致转换类型失败(比如字符串无法转换为float64位时),将会抛出一个ValueError。
>>> int_array = np.arange(10)
>>> calibers = np.array([.22, .270, .357, .380, .44, .50], dtype=np.float64)
>>> int_array.astype(calibers.dtype)
array([ 0., 1., 2., 3., 4., 5., 6., 7., 8., 9.])
# 也可以使用类型代码来传入数据类型
>>> empty_uint32 = np.empty(8, dtype='u4')
>>> empty_uint32
array([ 0, 1075314688, 0, 1075707904, 0,
1075838976, 0, 1072693248], dtype=uint32)
使用astype时总是生成一个新的数组,即使你传入的dtype与之前一样。
1.3 NumPy数组算术
数组之所以重要是因为它允许你进行批量操作而无须任何for循环。NumPy用户称这种特性为向量化。任何在两个等尺寸数组之间的算术操作都应用了逐元素操作的方式:
>>> arr = np.array([[1., 2., 3.], [4., 5., 6.]])
>>> arr
array([[ 1., 2., 3.],
[ 4., 5., 6.]])
>>> arr * arr
array([[ 1., 4., 9.],
[ 16., 25., 36.]])
>>> arr - arr
array([[ 0., 0., 0.],
[ 0., 0., 0.]])
带有标量计算的算术操作,会把计算参数传递给数组的每一个元素:
>>> 1 / arr
array([[ 1. , 0.5 , 0.3333],
[ 0.25 , 0.2 , 0.1667]])
>>> arr ** 0.5
array([[ 1. , 1.4142, 1.7321],
[ 2. , 2.2361, 2.4495]])
同尺寸数组之间的比较,会产生一个布尔值数组:
>>> arr2 = np.array([[0., 4., 1.], [7., 2., 12.]])
>>> arr2
array([[ 0., 4., 1.],
[ 7., 2., 12.]])
>>> arr2 > arr
array([[False, True, False],
[ True, False, True]], dtype=bool)
1.4 基础索引与切片
>>> arr = np.arange(10)
>>> arr
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
>>> arr[5]
5
>>> arr[5:8]
array([5, 6, 7])
>>> arr[5:8] = 12
>>> arr
array([ 0, 1, 2, 3, 4, 12, 12, 12, 8, 9])
>>> arr_slice = arr[5:8]
>>> arr_slice
array([12, 12, 12])
# 改变arr_slice
>>> arr_slice[1] = 12345
>>> arr
array([ 0, 1, 2, 3, 4, 12, 12345, 12, 8, 9])
不写切片值的[:]将会引用数组的所有值
>>> arr_slice[:] = 64
>>> arr
array([ 0, 1, 2, 3, 4, 64, 64, 64, 8, 9])
如果你还是想要一份数组切片的拷贝而不是一份视图的话,你就必须显式地复制这个数组,例如arr[5:8].copy()
对更高维度的数组,你会有更多选择。在一个二维数组中,每个索引值对应的元素不再是一个值,而是一个一维数组:
>>> arr2d = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
>>> arr2d[2]
array([7, 8, 9])
单个元素可以通过递归的方式获得。但是要多写点代码,你可以通过传递一个索引的逗号分隔列表去选择单个元素,以下两种方式效果一样:
>>> arr2d[0][2]
3
>>> arr2d[0, 2]
3
可以将0轴看作“行”,将1轴看作“列”
在多维数组中,你可以省略后续索引值,返回的对象将是降低一个维度的数组。因此在一个2×2×3的数组arr3d中:
>>> arr3d = np.array([[[1, 2, 3], [4, 5, 6]], [[7, 8, 9], [10, 11, 12]]])
>>> arr3d
array([[[ 1, 2, 3],
[ 4, 5, 6]],
[[ 7, 8, 9],
[10, 11, 12]]])
arr3d[0]是一个2×3的数组:
>>> arr3d[0]
array([[1, 2, 3],
[4, 5, 6]])
标量和数组都可以传递给arr3d[0]:
>>> old_values = arr3d[0].copy()
>>> arr3d[0] = 42
>>> arr3d
array([[[42, 42, 42],
[42, 42, 42]],
[[ 7, 8, 9],
[10, 11, 12]]])
>>> arr3d[0] = old_values
>>> arr3d
array([[[ 1, 2, 3],
[ 4, 5, 6]],
[[ 7, 8, 9],
[10, 11, 12]]])
arr3d[1,0]返回的是一个一维数组:
>>> arr3d[1, 0]
array([7, 8, 9])
上面的表达式可以分解为下面两步
>>> x = arr3d[1]
>>> x
array([[ 7, 8, 9],
[10, 11, 12]])
>>> x[0]
array([7, 8, 9])
以上的数组子集选择中,返回的数组都是视图
1.4.1 数组的切片索引
与Python列表的一维对象类似,数组可以通过类似的语法进行切片:
>>> arr
array([ 0, 1, 2, 3, 4, 64, 64, 64, 8, 9])
>>> arr[1:6]
array([ 1, 2, 3, 4, 64])
再回想下前面的二维数组,arr2d,对数组进行切片略有不同:
>>> arr2d
array([[1, 2, 3],
[4, 5, 6],
[7, 8, 9]])
# 数组沿着轴0进行了切片。表达式arrzd[:2]的含义为选择arr2d的前两“行”。
>>> arr2d[:2]
array([[1, 2, 3],
[4, 5, 6]])
>>> arr2d[:2, 1:]
array([[2, 3],
[5, 6]])
# 选择第二行但是只选择前两列
>>> arr2d[1, :2]
array([4, 5])
# 选择第三列,但是只选择前两行
>>> arr2d[:2, 2]
array([3, 6])
# 单独一个冒号表示选择整个轴上的数组
>>> arr2d[:, :1]
array([[1],
[4],
[7]])
>>> arr2d[:2, 1:] = 0
>>> arr2d
array([[1, 0, 0],
[4, 0, 0],
[7, 8, 9]])
二维数组的切片
1.5 布尔索引
使用numpy.random中的randn函数来生成一些随机正态分布的数据:
>>> names = np.array(['Bob', 'Joe', 'Will', 'Bob', 'Will', 'Joe', 'Joe'])
>>> data = np.random.randn(7, 4)
>>> names
array(['Bob', 'Joe', 'Will', 'Bob', 'Will', 'Joe', 'Joe'],
dtype='<U4')
>>> data
array([[ 0.0929, 0.2817, 0.769 , 1.2464],
[ 1.0072, -1.2962, 0.275 , 0.2289],
[ 1.3529, 0.8864, -2.0016, -0.3718],
[ 1.669 , -0.4386, -0.5397, 0.477 ],
[ 3.2489, -1.0212, -0.5771, 0.1241],
[ 0.3026, 0.5238, 0.0009, 1.3438],
[-0.7135, -0.8312, -2.3702, -1.8608]])
假设每个人名都和data数组中的一行相对应,并且我们想要选中所有'Bob'对应的行。与数学操作类似,数组的比较操作(比如==)也是可以向量化的。因此,比较names数组和字符串'Bob'会产生一个布尔值数组:
>>> names == 'Bob'
array([ True, False, False, True, False, False, False], dtype=bool)
在索引数组时可以传入布尔值数组:
>>> data[names == 'Bob']
array([[ 0.0929, 0.2817, 0.769 , 1.2464],
[ 1.669 , -0.4386, -0.5397, 0.477 ]])
布尔值数组的长度必须和数组轴索引长度一致。
>>> data[names == 'Bob', 2:]
array([[ 0.769 , 1.2464],
[-0.5397, 0.477 ]])
>>> data[names == 'Bob', 3]
array([ 1.2464, 0.477 ])
为了选择除了'Bob'以外的其他数据,你可以使用!=或在条件表达式前使用~对条件取反:
>>> names != 'Bob'
array([False, True, True, False, True, True, True], dtype=bool)
>>> data[~(names == 'Bob')]
array([[ 1.0072, -1.2962, 0.275 , 0.2289],
[ 1.3529, 0.8864, -2.0016, -0.3718],
[ 3.2489, -1.0212, -0.5771, 0.1241],
[ 0.3026, 0.5238, 0.0009, 1.3438],
[-0.7135, -0.8312, -2.3702, -1.8608]])
~符号可以在你想要对一个通用条件进行取反时使用:
>>> cond = names == 'Bob'
>>> data[~cond]
array([[ 1.0072, -1.2962, 0.275 , 0.2289],
[ 1.3529, 0.8864, -2.0016, -0.3718],
[ 3.2489, -1.0212, -0.5771, 0.1241],
[ 0.3026, 0.5238, 0.0009, 1.3438],
[-0.7135, -0.8312, -2.3702, -1.8608]])
当要选择三个名字中的两个时,可以对多个布尔值条件进行联合,需要使用数学操作符如&(and)和|(or)
>>> mask = (names == 'Bob') | (names == 'Will')
>>> mask
array([ True, False, True, True, True, False, False], dtype=bool)
>>> data[mask]
array([[ 0.0929, 0.2817, 0.769 , 1.2464],
[ 1.3529, 0.8864, -2.0016, -0.3718],
[ 1.669 , -0.4386, -0.5397, 0.477 ],
[ 3.2489, -1.0212, -0.5771, 0.1241]])
使用布尔值索引选择数据时,总是生成数据的拷贝,即使返回的数组并没有任何变化。
Python的关键字and和or对布尔值数组并没有用,请使用&(and)和|(or)来代替。
将data中所有的负值设置为0
>>> data[data < 0] = 0
>>> data
array([[ 0.0929, 0.2817, 0.769 , 1.2464],
[ 1.0072, 0. , 0.275 , 0.2289],
[ 1.3529, 0.8864, 0. , 0. ],
[ 1.669 , 0. , 0. , 0.477 ],
[ 3.2489, 0. , 0. , 0.1241],
[ 0.3026, 0.5238, 0.0009, 1.3438],
[ 0. , 0. , 0. , 0. ]])
利用一维布尔值数组对每一行设置数值也是非常简单的
>>> data[names != 'Joe'] = 7
>>> data
array([[ 7. , 7. , 7. , 7. ],
[ 1.0072, 0. , 0.275 , 0.2289],
[ 7. , 7. , 7. , 7. ],
[ 7. , 7. , 7. , 7. ],
[ 7. , 7. , 7. , 7. ],
[ 0.3026, 0.5238, 0.0009, 1.3438],
[ 0. , 0. , 0. , 0. ]])
1.6 神奇索引
神奇索引是NumPy中的术语,用于描述使用整数数组进行数据索引。
假设我们有一个8×4的数组:
>>> arr = np.empty((8, 4))
>>> for i in range(8):
>>> arr[i] = i
>>> arr
array([[ 0., 0., 0., 0.],
[ 1., 1., 1., 1.],
[ 2., 2., 2., 2.],
[ 3., 3., 3., 3.],
[ 4., 4., 4., 4.],
[ 5., 5., 5., 5.],
[ 6., 6., 6., 6.],
[ 7., 7., 7., 7.]])
为了选出一个符合特定顺序的子集,你可以简单地通过传递一个包含指明所需顺序的列表或数组来完成:
>>> arr[[4, 3, 0, 6]]
array([[ 4., 4., 4., 4.],
[ 3., 3., 3., 3.],
[ 0., 0., 0., 0.],
[ 6., 6., 6., 6.]])
如果使用负的索引,将从尾部进行选择
>>> arr[[-3, -5, -7]]
array([[ 5., 5., 5., 5.],
[ 3., 3., 3., 3.],
[ 1., 1., 1., 1.]])
传递多个索引数组时情况有些许不同,这样会根据每个索引元组对应的元素选出一个一维数组:
>>> arr = np.arange(32).reshape((8, 4))
>>> arr
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, 30, 31]])
>>> arr[[1, 5, 7, 2], [0, 3, 1, 2]]
array([ 4, 23, 29, 10])
在上述例子中,元素(1,0)、(5,3)、(7,1)和(2,2)被选中。如果不考虑数组的维数(本例中是二维),神奇索引的结果总是一维的。
在本例中,神奇索引的行为和一些用户所设想的并不相同。通常情况下,我们所设想的结果是通过选择矩阵中行列的子集所形成的矩形区域。下面是实现我们想法的一种方式:
>>> arr[[1, 5, 7, 2]][:, [0, 3, 1, 2]]
array([[ 4, 7, 5, 6],
[20, 23, 21, 22],
[28, 31, 29, 30],
[ 8, 11, 9, 10]])
神奇索引与切片不同,它总是将数据复制到一个新的数组中。
1.7 数组转置和换轴
转置是一种特殊的数据重组形式,可以返回底层数据的视图而不需要复制任何内容。数组拥有transpose方法,也有特殊的T属性:
>>> arr = np.arange(15).reshape((3, 5))
>>> arr
array([[ 0, 1, 2, 3, 4],
[ 5, 6, 7, 8, 9],
[10, 11, 12, 13, 14]])
>>> arr.T
array([[ 0, 5, 10],
[ 1, 6, 11],
[ 2, 7, 12],
[ 3, 8, 13],
[ 4, 9, 14]])
当进行矩阵计算时,你可能会经常进行一些特定操作,比如,当计算矩阵内积会使用np.dot:
>>> arr = np.random.randn(6, 3)
>>> arr
array([[-0.8608, 0.5601, -1.2659],
[ 0.1198, -1.0635, 0.3329],
[-2.3594, -0.1995, -1.542 ],
[-0.9707, -1.307 , 0.2863],
[ 0.378 , -0.7539, 0.3313],
[ 1.3497, 0.0699, 0.2467]])
>>> np.dot(arr.T, arr)
array([[ 9.2291, 0.9394, 4.948 ],
[ 0.9394, 3.7662, -1.3622],
[ 4.948 , -1.3622, 4.3437]])
对于更高维度的数组,transpose方法可以接收包含轴编号的元组,用于置换轴
# 轴已经被重新排序,使得原先的第二个轴变为第一个,原先的第一个轴变成了第二个,最后一个轴并没有改变。
>>> arr = np.arange(16).reshape((2, 2, 4))
>>> arr
array([[[ 0, 1, 2, 3],
[ 4, 5, 6, 7]],
[[ 8, 9, 10, 11],
[12, 13, 14, 15]]])
>>> arr.transpose((1, 0, 2))
array([[[ 0, 1, 2, 3],
[ 8, 9, 10, 11]],
[[ 4, 5, 6, 7],
[12, 13, 14, 15]]])
ndarray有一个swapaxes方法,该方法接收一对轴编号作为参数,并对轴进行调整用于重组数据
>>> arr
array([[[ 0, 1, 2, 3],
[ 4, 5, 6, 7]],
[[ 8, 9, 10, 11],
[12, 13, 14, 15]]])
>>> arr.swapaxes(1, 2)
array([[[ 0, 4],
[ 1, 5],
[ 2, 6],
[ 3, 7]],
[[ 8, 12],
[ 9, 13],
[10, 14],
[11, 15]]])
swapaxes返回的是数据的视图,而没有对数据进行复制
2. 通用函数:快速的逐元素数组函数
通用函数,也可以称为ufunc,是一种在ndarray数据中进行逐元素操作的函数。某些简单函数接收一个或多个标量数值,并产生一个或多个标量结果,而通用函数就是对这些简单函数的向量化封装。
有很多ufunc是简单的逐元素转换,比如sqrt或exp函数:
>>> arr = np.arange(10)
>>> arr
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
>>> np.sqrt(arr)
array([ 0. , 1. , 1.4142, 1.7321, 2. , 2.2361, 2.4495,
2.6458, 2.8284, 3. ])
>>> np.exp(arr)
array([ 1. , 2.7183, 7.3891, 20.0855, 54.5982,
148.4132, 403.4288, 1096.6332, 2980.958 , 8103.0839])
这些是所谓的一元通用函数。还有一些通用函数,比如add或maximum则会接收两个数组并返回一个数组作为结果,因此称为二元通用函数:
>>> x = np.random.randn(8)
>>> y = np.random.randn(8)
>>> x
array([-0.0119, 1.0048, 1.3272, -0.9193, -1.5491, 0.0222, 0.7584,
-0.6605])
>>> y
array([ 0.8626, -0.01 , 0.05 , 0.6702, 0.853 , -0.9559, -0.0235,
-2.3042])
# 逐个元素将x和y中的最大值计算出来
>>> np.maximum(x, y)
array([ 0.8626, 1.0048, 1.3272, 0.6702, 0.853 , 0.0222, 0.7584,
-0.6605])
也有一些通用函数返回多个数组。比如modf,是Python内建函数divmod的向量化版本。它返回了一个浮点值数组的小数部分和整数部分:
>>> arr = np.random.randn(7) * 5
>>> arr
array([-3.2623, -6.0915, -6.663 , 5.3731, 3.6182, 3.45 , 5.0077])
>>> remainder, whole_part = np.modf(arr)
>>> remainder
array([-0.2623, -0.0915, -0.663 , 0.3731, 0.6182, 0.45 , 0.0077])
>>> whole_part
array([-3., -6., -6., 5., 3., 3., 5.])
通用函数接收一个可选参数out,允许对数组按位置操作:
>>> arr
array([-3.2623, -6.0915, -6.663 , 5.3731, 3.6182, 3.45 , 5.0077])
>>> np.sqrt(arr)
array([ nan, nan, nan, 2.318, 1.9022, 1.8574, 2.2378])
>>> np.sqrt(arr, arr)
array([ nan, nan, nan, 2.318, 1.9022, 1.8574, 2.2378])
>>> arr
array([ nan, nan, nan, 2.318, 1.9022, 1.8574, 2.2378])
一元通用函数
二元通用函数
3. 使用数组进行面向数组编程
使用NumPy数组可以使你利用简单的数组表达式完成多种数据操作任务,而无须写些大量循环。这种利用数组表达式来替代显式循环的方法,称为向量化。通常,向量化的数组操作会比纯Python的等价实现在速度上快一到两个数量级(甚至更多),这对所有种类的数值计算产生了最大的影响。
假设我们想要对一些网格数据来计算函数sqrt(x^2+y^2)的值。np.meshgrid函数接收两个一维数组,并根据两个数组的所有(x,y)对生成一个二维矩阵:
>>> points = np.arange(-5, 5, 0.01) # 1000 equally spaced points
>>> xs, ys = np.meshgrid(points, points)
>>> ys
array([[-5. , -5. , -5. , ..., -5. , -5. , -5. ],
[-4.99, -4.99, -4.99, ..., -4.99, -4.99, -4.99],
[-4.98, -4.98, -4.98, ..., -4.98, -4.98, -4.98],
...,
[ 4.97, 4.97, 4.97, ..., 4.97, 4.97, 4.97],
[ 4.98, 4.98, 4.98, ..., 4.98, 4.98, 4.98],
[ 4.99, 4.99, 4.99, ..., 4.99, 4.99, 4.99]])
用matplotlib来生成这个二维数组的可视化
import matplotlib.pyplot as plt
plt.imshow(z, cmap=plt.cm.gray);
plt.colorbar()
plt.title("Image plot of $\sqrt{x^2 + y^2}$ for a grid of values")
3.1 将条件逻辑作为数组操作
numpy.where函数是三元表达式x if condition else y的向量化版本。假设我们有一个布尔值数组和两个数值数组:
>>> xarr = np.array([1.1, 1.2, 1.3, 1.4, 1.5])
>>> yarr = np.array([2.1, 2.2, 2.3, 2.4, 2.5])
>>> cond = np.array([True, False, True, True, False])
假设cond中的元素为True时,我们取xarr中的对应元素值,否则取yarr中的元素。我们可以通过列表推导式来完成,像下列代码这样:
>>> result = [(x if c else y)
>>> for x, y, c in zip(xarr, yarr, cond)]
>>> result
[1.1000000000000001, 2.2000000000000002, 1.3, 1.3999999999999999, 2.5]
这样会产生多个问题。首先,如果数组很大的话,速度会很慢(因为所有的工作都是通过解释器解释Python代码完成)。其次,当数组是多维时,就无法奏效了。而使用np.where时,就可以非常简单地完成:
>>> result = np.where(cond, xarr, yarr)
>>> result
array([ 1.1, 2.2, 1.3, 1.4, 2.5])
np.where的第二个和第三个参数并不需要是数组,它们可以是标量。where在数据分析中的一个典型用法是根据一个数组来生成一个新的数组。假设你有一个随机生成的矩阵数据,并且你想将其中的正值都替换为2,将所有的负值替换为-2,使用np.where会很容易实现:
>>> arr = np.random.randn(4, 4)
>>> arr
array([[-0.5031, -0.6223, -0.9212, -0.7262],
[ 0.2229, 0.0513, -1.1577, 0.8167],
[ 0.4336, 1.0107, 1.8249, -0.9975],
[ 0.8506, -0.1316, 0.9124, 0.1882]])
>>> arr > 0
array([[False, False, False, False],
[ True, True, False, True],
[ True, True, True, False],
[ True, False, True, True]], dtype=bool)
>>> np.where(arr > 0, 2, -2)
array([[-2, -2, -2, -2],
[ 2, 2, -2, 2],
[ 2, 2, 2, -2],
[ 2, -2, 2, 2]])
可以使用np.where将标量和数组联合,例如,我可以像下面的代码那样将arr中的所有正值替换为常数2:
>>> np.where(arr > 0, 2, arr) # 仅将正值设为2
array([[-0.5031, -0.6223, -0.9212, -0.7262],
[ 2. , 2. , -1.1577, 2. ],
[ 2. , 2. , 2. , -0.9975],
[ 2. , -0.1316, 2. , 2. ]])
传递给np.where的数组既可以是同等大小的数组,也可以是标量。
3.2 数学和统计方法
许多关于计算整个数组统计值或关于轴向数据的数学函数,可以作为数组类型的方法被调用。你可以使用聚合函数(通常也叫缩减函数),比如sum、mean和std(标准差),既可以直接调用数组实例的方法,也可以使用顶层的NumPy函数。
生成了一些正态分布的随机数,并计算了部分聚合统计数据:
>>> arr = np.random.randn(5, 4)
>>> arr
array([[ 2.1695, -0.1149, 2.0037, 0.0296],
[ 0.7953, 0.1181, -0.7485, 0.585 ],
[ 0.1527, -1.5657, -0.5625, -0.0327],
[-0.929 , -0.4826, -0.0363, 1.0954],
[ 0.9809, -0.5895, 1.5817, -0.5287]])
>>> arr.mean()
0.19607051119998253
>>> np.mean(arr)
0.19607051119998253
>>> arr.sum()
3.9214102239996507
像mean、sum等函数可以接收一个可选参数axis,这个参数可以用于计算给定轴向上的统计值,形成一个下降一维度的数组:
>>> arr.mean(axis=1)
array([ 1.022 , 0.1875, -0.502 , -0.0881, 0.3611])
>>> arr.sum(axis=0)
array([ 3.1693, -2.6345, 2.2381, 1.1486])
arr.mean(1)表示“计算每一列的平均值”,而arr.sum(0)表示“计算行轴向的累和”。
其他的方法,例如cumsum和cumprod并不会聚合,它们会产生一个中间结果:
>>> arr = np.array([0, 1, 2, 3, 4, 5, 6, 7])
>>> arr.cumsum()
array([ 0, 1, 3, 6, 10, 15, 21, 28])
在多维数组中,像cumsum这样的累积函数返回相同长度的数组,但是可以在指定轴向上根据较低维度的切片进行部分聚合:
>>> arr = np.array([[0, 1, 2], [3, 4, 5], [6, 7, 8]])
>>> arr
array([[0, 1, 2],
[3, 4, 5],
[6, 7, 8]])
>>> arr.cumsum(axis=0)
array([[ 0, 1, 2],
[ 3, 5, 7],
[ 9, 12, 15]])
>>> arr.cumprod(axis=1)
array([[ 0, 0, 0],
[ 3, 12, 60],
[ 6, 42, 336]])
基础数组统计方法
3.3 布尔值数组的方法
在前面的方法,布尔值会被强制为1(True)和0(False)。因此,sum通常可以用于计算布尔值数组中的True的个数:
>>> arr = np.random.randn(100)
>>> (arr > 0).sum() # 正值的个数
42
对于布尔值数组,有两个非常有用的方法any和all。any检查数组中是否至少有一个True,而all检查是否每个值都是True:
>>> bools = np.array([False, False, True, False])
>>> bools.any()
True
>>> bools.all()
False
这些方法也可适用于非布尔值数组,所有的非0元素都会按True处理。
3.4 排序
和Python的内建列表类型相似,NumPy数组可以使用sort方法按位置排序:
>>> arr = np.random.randn(6)
>>> arr
array([ 0.6095, -0.4938, 1.24 , -0.1357, 1.43 , -0.8469])
>>> arr.sort()
>>> arr
array([-0.8469, -0.4938, -0.1357, 0.6095, 1.24 , 1.43 ])
你可以在多维数组中根据传递的axis值,沿着轴向对每个一维数据段进行排序:
>>> arr = np.random.randn(5, 3)
>>> arr
array([[ 0.6033, 1.2636, -0.2555],
[-0.4457, 0.4684, -0.9616],
[-1.8245, 0.6254, 1.0229],
[ 1.1074, 0.0909, -0.3501],
[ 0.218 , -0.8948, -1.7415]])
>>> arr.sort(1)
>>> arr
array([[-0.2555, 0.6033, 1.2636],
[-0.9616, -0.4457, 0.4684],
[-1.8245, 0.6254, 1.0229],
[-0.3501, 0.0909, 1.1074],
[-1.7415, -0.8948, 0.218 ]])
顶层的np.sort方法返回的是已经排序好的数组拷贝,而不是对原数组按位置排序。下面的例子计算的是一个数组的分位数,并选出分位数所对应的值,这是一种应急的方式:
>>> large_arr = np.random.randn(1000)
>>> large_arr.sort()
>>> large_arr[int(0.05 * len(large_arr))] # 5% quantile
-1.5311513550102103
3.5 唯一值与其他集合逻辑
NumPy包含一些针对一维ndarray的基础集合操作。常用的一个方法是np.unique,返回的是数组中唯一值排序后形成的数组:
>>> names =np.array(['Bob', 'Joe', 'Will', 'Bob', 'Will', 'Joe', 'Joe'])
>>> np.unique(names)
array(['Bob', 'Joe', 'Will'],
dtype='<U4')
>>> ints = np.array([3, 3, 3, 2, 2, 1, 1, 4, 4])
>>> np.unique(ints)
array([1, 2, 3, 4])
将np.unique和纯Python实现相比较:
>>> sorted(set(names))
['Bob', 'Joe', 'Will']
另一个函数,np.in1d,可以检查一个数组中的值是否在另外一个数组中,并返回一个布尔值数组:
>>> values = np.array([6, 0, 0, 3, 2, 5, 6])
>>> np.in1d(values, [2, 3, 6])
array([ True, False, False, True, True, False, True], dtype=bool)
数组的集合操作
4. 使用进行文件输入和输出
NumPy可以在硬盘中将数据以文本或二进制文件的形式进行存入硬盘或由硬盘载入。
np.save和np.load是高效存取硬盘数据的两大工具函数。数组在默认情况下是以未压缩的格式进行存储的,后缀名是.npy:
>>> arr = np.arange(10)
>>> np.save('some_array', arr)
如果文件存放路径中没写.npy时,后缀名会被自动加上。硬盘上的数组可以使用np.load进行载入:
>>> np.load('some_array.npy')
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
你可以使用np.savez并将数组作为参数传递给该函数,用于在未压缩文件中保存多个数组:
>>> np.savez('array_archive.npz', a=arr, b=arr)
当载入一个.npy文件的时候,你会获得一个字典型的对象,并通过该对象很方便地载入单个数组:
>>> arch = np.load('array_archive.npz')
>>> arch['b']
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
如果你的数据已经压缩好了,你可能会想要使用numpy.savez_compressed将数据存入已经压缩的文件:
# numpy.savez_compressed将数据存入已经压缩的文件:
>>> np.savez_compressed('arrays_compressed.npz', a=arr, b=arr)
5. 线性代数
线性代数,比如矩阵乘法、分解、行列式等方阵数学,是所有数组类库的重要组成部分。和Matlab等其他语言相比,NumPy的线性代数中所不同的是*是矩阵的逐元素乘积,而不是矩阵的点乘积。因此NumPy的数组方法和numpy命名空间中都有一个函数dot,用于矩阵的操作:
>>> x = np.array([[1., 2., 3.], [4., 5., 6.]])
>>> y = np.array([[6., 23.], [-1, 7], [8, 9]])
>>> x
array([[ 1., 2., 3.],
[ 4., 5., 6.]])
>>> y
array([[ 6., 23.],
[ -1., 7.],
[ 8., 9.]])
>>> x.dot(y)
array([[ 28., 64.],
[ 67., 181.]])
# x.dot(y)等价于np.dot(x,y)
>>> np.dot(x, y)
array([[ 28., 64.],
[ 67., 181.]])
一个二维数组和一个长度合适的一维数组之间的矩阵乘积,其结果是一个一维数组
>>> np.dot(x, np.ones(3))
array([ 6., 15.])
特殊符号@也作为中缀操作符,用于点乘矩阵操作:
>>> x @ np.ones(3)
array([ 6., 15.])
numpy.linalg拥有一个矩阵分解的标准函数集,以及其他常用函数,例如求逆和行列式求解。这些函数都是通过在MATLAB和R等其他语言使用的相同的行业标准线性代数库来实现的,例如BLAS、LAPACK或英特尔专有的MKL(数学核心库)(是否使用MKL取决于使用NumPy的版本):
from numpy.linalg import inv, qr
>>> X = np.random.randn(5, 5)
>>> mat = X.T.dot(X)
>>> inv(mat)
array([[ 933.1189, 871.8258, -1417.6902, -1460.4005, 1782.1391],
[ 871.8258, 815.3929, -1325.9965, -1365.9242, 1666.9347],
[-1417.6902, -1325.9965, 2158.4424, 2222.0191, -2711.6822],
[-1460.4005, -1365.9242, 2222.0191, 2289.0575, -2793.422 ],
[ 1782.1391, 1666.9347, -2711.6822, -2793.422 , 3409.5128]])
>>> mat.dot(inv(mat))
array([[ 1., 0., -0., -0., -0.],
[-0., 1., 0., 0., 0.],
[ 0., 0., 1., 0., 0.],
[-0., 0., 0., 1., -0.],
[-0., 0., 0., 0., 1.]])
>>> q, r = qr(mat)
>>> r
array([[-1.6914, 4.38 , 0.1757, 0.4075, -0.7838],
[ 0. , -2.6436, 0.1939, -3.072 , -1.0702],
[ 0. , 0. , -0.8138, 1.5414, 0.6155],
[ 0. , 0. , 0. , -2.6445, -2.1669],
[ 0. , 0. , 0. , 0. , 0.0002]])
表达式X.T.dot(X)计算的是X和它的转置矩阵X.T的点乘积。
常用numpy.linalg函数
6. 伪随机数生成
numpy.random模块填补了Python内建的random模块的不足,可以高效地生成多种概率分布下的完整样本值数组。例如,你可以使用normal来获得一个4×4的正态分布样本数组:
>>> samples = np.random.normal(size=(4, 4))
>>> samples
array([[ 0.5732, 0.1933, 0.4429, 1.2796],
[ 0.575 , 0.4339, -0.7658, -1.237 ],
[-0.5367, 1.8545, -0.92 , -0.1082],
[ 0.1525, 0.9435, -1.0953, -0.144 ]])
然而Python内建的random模块一次只能生成一个值。你可以从下面的示例中看到,numpy.random在生成大型样本时比纯Python的方式快了一个数量级:
>>> from random import normalvariate
>>> N = 1000000
>>> %timeit samples = [normalvariate(0, 1) for _ in range(N)]
1.77 s +- 126 ms per loop (mean +- std. dev. of 7 runs, 1 loop each)
>>> %timeit np.random.normal(size=N)
61.7 ms +- 1.32 ms per loop (mean +- std. dev. of 7 runs, 10 loops each)
我们可以称这些为伪随机数,因为它们是由具有确定性行为的算法根据随机数生成器中的随机数种子生成的。你可以通过np.random.seed更改NumPy的随机数种子:
>>> np.random.seed(1234)
numpy.random中的数据生成函数公用了一个全局的随机数种子。为了避免全局状态,你可以使用numpy.random.RandomState生成一个随机数生成器,使数据独立于其他的随机数状态:
>>> rng = np.random.RandomState(1234)
>>> rng.randn(10)
array([ 0.4714, -1.191 , 1.4327, -0.3127, -0.7206, 0.8872, 0.8596,
-0.6365, 0.0157, -2.2427])
numpy.random中的部分函数列表
7. 随机漫步
首先,让我们考虑一个简单的随机漫步,从0开始,步进为1和-1,且两种步进发生的概率相等。
以下是使用内建random模块利用纯Python实现的一个1,000步的随机漫步:
# 1,000步的随机漫步:
>>> import random
>>> position = 0
>>> walk = [position]
>>> steps = 1000
>>> for i in range(steps):
>>> step = 1 if random.randint(0, 1) else -1
>>> position += step
>>> walk.append(position)
# 随机漫步的前100步可视化
>>> plt.plot(walk[:100])
使用np.random模块一次性抽取1,000次投掷硬币的结果,每次投掷的结果为1或-1,然后计算累积值:
>>> nsteps = 1000
>>> draws = np.random.randint(0, 2, size=nsteps)
>>> steps = np.where(draws > 0, 1, -1)
>>> walk = steps.cumsum()
由此我们开始从漫步轨道上提取一些统计数据,比如最大值、最小值等:
>>> walk.min()
-3
>>> walk.max()
31
更复杂的统计是第一次穿越时间,即随机漫步的某一步达到了某个特定值。这里假设我们想要知道漫步中是何时连续朝某个方向连续走了10步。np.abs(walk)>=10给我们一个布尔值数组,用于表明漫步是否连续在同一方向走了十步,但是我们想要的是第一次走了10步或-10步的位置。我们可以使用argmax来计算,该函数可以返回布尔值数组中最大值的第一个位置(True就是最大值):
>>> (np.abs(walk) >= 10).argmax()
37
这里使用argmax效率并不高,因为它总是完整地扫描整个数组。在这个特殊的示例中,一旦True被发现,我们就知道最大值了。
7.1 一次性模拟多次随机漫步
模拟多次随机漫步,比如说5,000步。可以稍微地修改下之前的代码来生成所有的随机步。如果传入一个2个元素的元组,numpy.random中的函数可以生成一个二维的抽取数组,并且我们可以一次性地跨行算出全部5,000个随机步的累积和:
>>> nwalks = 5000
>>> nsteps = 1000
>>> draws = np.random.randint(0, 2, size=(nwalks, nsteps)) # 0或1
>>> steps = np.where(draws > 0, 1, -1)
>>> walks = steps.cumsum(1)
>>> walks
array([[ 1, 0, 1, ..., 8, 7, 8],
[ 1, 0, -1, ..., 34, 33, 32],
[ 1, 0, -1, ..., 4, 5, 4],
...,
[ 1, 2, 1, ..., 24, 25, 26],
[ 1, 2, 3, ..., 14, 13, 14],
[ -1, -2, -3, ..., -24, -23, -22]])
现在我们可以计算出这些随机步的最大值和最小值了:
>>> walks.max()
138
>>> walks.min()
-133
让我们在这些随机步中计算出30或-30的最小穿越时间。这有点棘手,因为不是所有的5,000个都达到了30。我们可以用any方法来检查:
>>> hits30 = (np.abs(walks) >= 30).any(1)
>>> hits30
array([False, True, False, ..., False, True, False], dtype=bool)
>>> hits30.sum() # 达到30或-30的数字
3410
我们可以使用布尔值数组来选出绝对步数超过30的步所在的行,并使用argmax从轴向1上获取穿越时间:
>>> crossing_times = (np.abs(walks[hits30]) >= 30).argmax(1)
>>> crossing_times.mean()
498.88973607038122
利用其他分布而不是等概率的掷硬币实验来随机漫步也是很容易的。你只需要使用一个不同的随机数生成函数,比如normal,再根据特定的均值和标准差即可生成正态分布下的随机步:
>>> steps = np.random.normal(loc=0, scale = 0.25, size=(nwalks, nsteps))