Numpy基础:数组和矢量计算

Numpy是高性能科学计算和数学分析的基础包,是大多数高级工具的构建基础。具备以下功能:

  • ndarray,一个具有矢量算术运算和复杂广播能力的快速且节省空间的多维数组

  • 对数组数据进行快速运算的标准数学函数(无需编写循环)

  • 读写磁盘数据的工具

  • 线性代数、随机数生成以及傅里叶变换

  • 集成由C、C++、Fortran等语言编写代码的工具

标准的NumPy约定导入方法:

1
import numpy as np

NumPy的ndarray:一种多维数组对象

NumPy最重要的特点是其N维数组对象(即ndarray),能够对整块数据执行一些数学运算,其语法与基本数据类型之间的运算一样。

1
2
data = np.random.random((2,3))
data
array([[ 0.08455582,  0.69068026,  0.0977948 ],
       [ 0.42529949,  0.9647802 ,  0.09677386]])
1
data * 10
array([[ 0.84555825,  6.90680264,  0.97794796],
       [ 4.25299487,  9.647802  ,  0.96773858]])
1
data + data
array([[ 0.16911165,  1.38136053,  0.19558959],
       [ 0.85059897,  1.9295604 ,  0.19354772]])

特点:

  1. ndarray是一个通用的同构数据容器,即其中的元素必须是相同类型的。
  2. 每个数组都有一个shape(一个表示各维度大小的元组)和一个dtype(说明数组数据类型的对象)
1
data.shape

(2, 3)

1
data.dtype

dtype('float64')

创建ndarray对象

np.array函数

  • 最简单的办法是使用array函数,它接受一切序列型对象(字符串、列表、元组等,包括其他数组)
1
2
3
data1 = [6, 7.5, 8, 0, 1]
arr1 = np.array(data1)
arr1

array([ 6. , 7.5, 8. , 0. , 1. ])

1
data1

[6, 7.5, 8, 0, 1]

嵌套序列(比如等长的列表)将会被转化为一个多维数组

1
2
3
data2 = [[1,2,3,4],[5,6,7,8]]
arr2 = np.array(data2)
arr2
array([[1, 2, 3, 4],
       [5, 6, 7, 8]])
1
arr2.ndim

2

1
arr2.shape

(2, 4)

除非显示说明,np.array会尝试为新建的数组推断出一个较为合适的数据类型,保存在dtype对象中。

1
arr1.dtype

dtype('float64')

1
arr2.dtype

dtype('int32')

  • zerosones分别可以创建指定长度或形状的全部为0或1的数组。empty可以创建一个没有任何具体数值的数组。创建多维数组只需传入一个表示形状的元组即刻。
1
np.zeros((10,))

array([ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.])

1
np.zeros(10)

array([ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.])

1
np.zeros((3,2))
array([[ 0.,  0.],
       [ 0.,  0.],
       [ 0.,  0.]])
1
np.zeros((3,6))
array([[ 0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.]])
1
np.ones((3,2))
array([[ 1.,  1.],
       [ 1.,  1.],
       [ 1.,  1.]])
1
np.empty((2, 3, 2))
 array([[[ 0.,  0.],
         [ 0.,  0.],
         [ 0.,  0.]],

       [[ 0.,  0.],
        [ 0.,  0.],
        [ 0.,  0.]]])

np.empty并会总是返回全0数组,多数情况下,他返回的都是一些未初始化的垃圾值。

  • arrange函数(Python内置函数range的NumPy版)
1
np.arange(15)

array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14])

1
np.arange(0,15,2)

array([ 0, 2, 4, 6, 8, 10, 12, 14])

ndarray的数据类型

1
2
3
arr3 = np.array([1, 2, 3], dtype = np.float64)
arr4 = np.array([1, 2, 3], dtype = np.int32)
arr3.dtype

dtype('float64')

1
arr4.dtype

dtype('int32')

  • 可以通过ndarray的astype方法显式转换其dtype:
1
2
arr5 = np.array([1, 2, 3, 4, 5])
arr5.dtype

dtype('int32')

1
2
float_arr = arr5.astype(np.float)
float_arr.dtype

dtype('float64')

如果浮点数转换为整数,则小数部分会被截断。

1
2
arr6 = np.array([3.7, -1.2, -2.6, 0.5, 12.9, 10.1])
arr6

array([ 3.7, -1.2, -2.6, 0.5, 12.9, 10.1])

1
arr6.dtype

dtype('float64')

1
arr6.astype(np.int)

array([ 3, -1, -2, 0, 12, 10])

如果某字符串全是数字,也可以用astype将其转换为数值形式:

1
numeric_strings = np.array(['1.25', '-9.6', '42'], dtype=np.string_)
1
2
newarr = numeric_strings.astype(np.float)
newarr

array([ 1.25, -9.6 , 42. ])

1
numeric_strings

array([b'1.25', b'-9.6', b'42'], dtype='|S4')

1
newarr.dtype

dtype('float64')

  • 数组类型转化时,会返回新建数组,即原数组的一个copy,原数组并未发生改变。
  • 如果转化过程失败,会引发TypeError
  • 转换过程中可以不具体指定位数,譬如np.intnp.float等,NumPy会自动完成映射
1
2
calibers = np.array([.22, .270, .357, .380, .44, .50], dtype=np.float64)
int_array = np.arange(10)
  • dtype可以用于指定新数组与已有数组同数据类型
1
int_array.astype(calibers.dtype)

array([ 0., 1., 2., 3., 4., 5., 6., 7., 8., 9.])

切片和索引

基本的索引和切片

一维数组很简单,从表面上看,跟Python列表差不多。

1
2
arr8 = np.arange(10)
arr8

array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])

1
arr8[5]

5

1
arr8[5:8]

array([5, 6, 7])

  • 将一个标量值赋值给一个片段时,该值会自动传播(广播)到整个选区。
1
2
arr8[5:8] = 12
arr8

array([ 0, 1, 2, 3, 4, 12, 12, 12, 8, 9])

  • 跟列表的重要区别在于,数组切片是原始数组的视图,即数据不会被复制,视图上的任何修改都会直接反映到源数据上
1
2
3
arr8_slcie = arr8[5:8]
arr8_slcie[1] = 12345
arr8

array([ 0, 1, 2, 3, 4, 12, 12345, 12, 8, 9])

1
2
arr8_slcie[:] = 64
arr8

array([ 0, 1, 2, 3, 4, 64, 64, 64, 8, 9])

1
2
3
4
ls = list(range(10))
ls_slice = ls[5:8]
ls_slice[1] = 12
ls

[0, 1, 2, 3, 4, 5, 6, 7, 8, 9]

1
ls_slice

[5, 12, 7]

  • 如果需要得到ndarray切片的一份副本而非视图,则需进行显式地复制操作,例如arr8[5:8].copy()
  • 在二维数组中,各索引位置上的元素不是标量,而是一维数组
1
2
arr2d = np.array([[1,2,3],[4,5,6],[7,8,9]])
arr2d
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])
1
arr2d[2]

array([7, 8, 9])

1
arr2d[0][2]
3
1
arr2d[0,2]
3
  • 多维数组中,如果省略后面地索引,则返回对象是一个维度低一级地ndarray
1
2
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]]])
1
arr3d.shape
(2, 2, 3)
1
arr3d[0]
array([[1, 2, 3],
       [4, 5, 6]])
  • 标量值和数组都可以赋值给arr3d[0]
1
2
3
old_values = arr3d[0].copy()
arr3d[0] = 42
arr3d
array([[[42, 42, 42],
        [42, 42, 42]],

       [[ 7,  8,  9],
        [10, 11, 12]]])
1
2
arr3d[0] = old_values
arr3d
array([[[ 1,  2,  3],
        [ 4,  5,  6]],

       [[ 7,  8,  9],
        [10, 11, 12]]])
1
arr3d[1, 0]
array([7, 8, 9])

切片索引

  • ndarray地切片语法跟Python列表这样的一维对象差不多
  • 高维度对象花样更多,可以在一个或多个轴上进行切片,也可以跟整数索引混合使用。
1
2
arr8 = np.arange(10)
arr8[1:6]
array([1, 2, 3, 4, 5])
1
arr2d
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])
1
arr2d[1,2]
6
1
arr2d[:2]
array([[1, 2, 3],
       [4, 5, 6]])

可以看出,上述例子是沿着第0轴(即第一个轴)进行切片地。也就是说,切片是沿着一个轴向选取元素,因此可以一次传入多个切片,就像传入多个索引那样。

1
arr2d[:2,1:]
array([[2, 3],
       [5, 6]])

切片操作时,只能得到相同维数的数组视图。通过将整数索引和切片索引混合,可以得到低维度的切片。

1
arr2d[1,:2]
array([4, 5])
1
arr2d[2,:1]
array([7])
1
arr2d[2,0]
7
  • 只有冒号表示选取整个轴
  • 对切片进行赋值也会被扩散到整个选区
1
arr2d[:,0]
array([1, 4, 7])
1
2
arr2d[:2,1:] = 0
arr2d
array([[1, 0, 0],
       [4, 0, 0],
       [7, 8, 9]])

布尔型索引

1
2
names = np.array(['Bob', 'Joe', 'Will', 'Bob', 'Will', 'Joe', 'Joe'])
data = np.random.randn(7, 4) # 7行4列标准正态分布的随机数
1
cl = np.array(['A','B','A','C'])
1
names
array(['Bob', 'Joe', 'Will', 'Bob', 'Will', 'Joe', 'Joe'],
      dtype='<U4')
1
data
array([[ 0.0044272 , -0.52538554,  0.85946006,  2.39584017],
       [ 0.47111574,  2.19407251, -1.83018561,  1.42571829],
       [-0.25892768,  0.08458305,  0.6599042 ,  1.7552628 ],
       [-0.78157234, -2.04734931,  1.11835735,  1.59707815],
       [-1.21924958,  0.04668819, -1.09584675,  1.05554286],
       [ 0.31169982, -0.09411045, -0.13023334,  0.8459281 ],
       [-0.48763537, -1.05315113, -0.64942157,  0.73403295]])

假设每个名字对应data数组总的一行,如何选出对应"Bob"的所有行?跟算术运算一样,数组的比较运算(如==)也是矢量化的。

1
names == 'Bob'
array([ True, False, False,  True, False, False, False], dtype=bool)
  • 布尔型数组可以用于索引, 但布尔型数组的长度必须与被索引的轴长度一致。
  • 此外,还可以将布尔型数组跟切片、整数(或整数序列)混合使用
1
data[names=='Bob']
array([[ 0.0044272 , -0.52538554,  0.85946006,  2.39584017],
       [-0.78157234, -2.04734931,  1.11835735,  1.59707815]])
1
data[:,cl=='A']
array([[ 0.0044272 ,  0.85946006],
       [ 0.47111574, -1.83018561],
       [-0.25892768,  0.6599042 ],
       [-0.78157234,  1.11835735],
       [-1.21924958, -1.09584675],
       [ 0.31169982, -0.13023334],
       [-0.48763537, -0.64942157]])
1
data[names=='Joe', 2:]
array([[-1.83018561,  1.42571829],
       [-0.13023334,  0.8459281 ],
       [-0.64942157,  0.73403295]])
1
data[names=='Will', 3]
array([ 1.7552628 ,  1.05554286])

要选取除Bob以外的其他值,既可以使用不等于符号(!=),也可以通过波浪号(~)(以前版本是减号(-))对条件进行否定。

1
data[names != 'Bob']
array([[ 0.47111574,  2.19407251, -1.83018561,  1.42571829],
       [-0.25892768,  0.08458305,  0.6599042 ,  1.7552628 ],
       [-1.21924958,  0.04668819, -1.09584675,  1.05554286],
       [ 0.31169982, -0.09411045, -0.13023334,  0.8459281 ],
       [-0.48763537, -1.05315113, -0.64942157,  0.73403295]])
1
data[~(names == 'Bob')]
array([[ 0.47111574,  2.19407251, -1.83018561,  1.42571829],
       [-0.25892768,  0.08458305,  0.6599042 ,  1.7552628 ],
       [-1.21924958,  0.04668819, -1.09584675,  1.05554286],
       [ 0.31169982, -0.09411045, -0.13023334,  0.8459281 ],
       [-0.48763537, -1.05315113, -0.64942157,  0.73403295]])

选取三个名字中的两个则需要组合应用多个布尔条件,使用 &(和)、|(或)之类的布尔算术运算符。Python关键词and 和 or 在布尔型数组中无效。

1
2
mask = (names == 'Bob') | (names == 'Will')
mask
array([ True, False,  True,  True,  True, False, False], dtype=bool)
1
data[mask]
array([[ 0.0044272 , -0.52538554,  0.85946006,  2.39584017],
       [-0.25892768,  0.08458305,  0.6599042 ,  1.7552628 ],
       [-0.78157234, -2.04734931,  1.11835735,  1.59707815],
       [-1.21924958,  0.04668819, -1.09584675,  1.05554286]])
  • 布尔型索引选取数组中的数据,总是创建数据的副本,即使返回一模一样的数组
  • 可以通过布尔型数组设置值。譬如,将data中所有的负值设置为0:
1
2
data[data < 0] = 0
data
array([[ 0.0044272 ,  0.        ,  0.85946006,  2.39584017],
       [ 0.47111574,  2.19407251,  0.        ,  1.42571829],
       [ 0.        ,  0.08458305,  0.6599042 ,  1.7552628 ],
       [ 0.        ,  0.        ,  1.11835735,  1.59707815],
       [ 0.        ,  0.04668819,  0.        ,  1.05554286],
       [ 0.31169982,  0.        ,  0.        ,  0.8459281 ],
       [ 0.        ,  0.        ,  0.        ,  0.73403295]])
1
2
data[names != 'Joe'] = 7
data
array([[ 7.        ,  7.        ,  7.        ,  7.        ],
       [ 0.47111574,  2.19407251,  0.        ,  1.42571829],
       [ 7.        ,  7.        ,  7.        ,  7.        ],
       [ 7.        ,  7.        ,  7.        ,  7.        ],
       [ 7.        ,  7.        ,  7.        ,  7.        ],
       [ 0.31169982,  0.        ,  0.        ,  0.8459281 ],
       [ 0.        ,  0.        ,  0.        ,  0.73403295]])

花式索引(Fancy indexing)

NumPy术语,指利用整数数组进行索引。

1
2
3
arr9 = np.empty((8, 4))
for i in range(8):
arr9[i] = i
1
arr9
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.]])
  • 为了以特定顺序选取子集,只需传入一个用于指定顺序的整数列或ndarray即可
  • 使用负数将会从末尾开始选取行
1
arr9[1]
array([ 1.,  1.,  1.,  1.])
1
arr9[[4, 3, 0, 6]]
array([[ 4.,  4.,  4.,  4.],
       [ 3.,  3.,  3.,  3.],
       [ 0.,  0.,  0.,  0.],
       [ 6.,  6.,  6.,  6.]])
1
arr9[:,[2,1,0]]
array([[ 0.,  0.,  0.],
       [ 1.,  1.,  1.],
       [ 2.,  2.,  2.],
       [ 3.,  3.,  3.],
       [ 4.,  4.,  4.],
       [ 5.,  5.,  5.],
       [ 6.,  6.,  6.],
       [ 7.,  7.,  7.]])
1
arr9[[-3, -5, -7]]
array([[ 5.,  5.,  5.,  5.],
       [ 3.,  3.,  3.,  3.],
       [ 1.,  1.,  1.,  1.]])
  • 一次传入多个索引数组,将会返回一个一维数组,其中的元素对应各个索引元组
1
2
arr10 = np.arange(32).reshape((8, 4))
arr10
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]])
1
2
# 返回一维数组
arr10[[1,5,7,2],[0, 3,1,2]]
array([ 4, 23, 29, 10])
1
arr10[[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
arr10[np.ix_([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]])
  • 花式索引与切片不一样,总是将数据复制到新的数组中。

数组的重塑

reshape函数

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

转置和轴对换

  • 转置(transpose),是重塑的一种特殊形式,返回的是源数据的视图(不会进行任何复制操作)。
  • 此外,还有T属性,在进行矩阵计算时经常用到,譬如利用np.dot计算矩阵内积\(X^{T}X\)
1
arr14.transpose()
array([[ 0,  4,  8],
       [ 1,  5,  9],
       [ 2,  6, 10],
       [ 3,  7, 11]])
1
2
arr13 = np.arange(16).reshape((2, 2, 4))  # 三维数组,2x2x4,轴索引编号分别为:0,1,2
arr13
array([[[ 0,  1,  2,  3],
        [ 4,  5,  6,  7]],

       [[ 8,  9, 10, 11],
        [12, 13, 14, 15]]])
1
arr13[0,1,2]
6
1
arr13.transpose((1, 0 ,2))  # 此处的(1, 0, 2)为维度的索引
array([[[ 0,  1,  2,  3],
        [ 8,  9, 10, 11]],

       [[ 4,  5,  6,  7],
        [12, 13, 14, 15]]])

上面的变化中,实际上是把0轴和1轴对换了。譬如,元素6的索引位置为(0,1,2),把0轴和1轴对换后,将位于(1,0,2)位置,虽然变化后数组依然为2x2x4。再例如:

1
2
arr15 = np.arange(24).reshape((2,3,4))    # 2x3x4的三位数组
arr15
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]]])
1
arr15.transpose((1,2,0))   #变为3X4x2的三维数组
array([[[ 0, 12],
        [ 1, 13],
        [ 2, 14],
        [ 3, 15]],

       [[ 4, 16],
        [ 5, 17],
        [ 6, 18],
        [ 7, 19]],
    
       [[ 8, 20],
        [ 9, 21],
        [10, 22],
        [11, 23]]])
1
2
arr11 = np.arange(15).reshape((3, 5))
arr11
array([[ 0,  1,  2,  3,  4],
       [ 5,  6,  7,  8,  9],
       [10, 11, 12, 13, 14]])
1
arr11.T
array([[ 0,  5, 10],
       [ 1,  6, 11],
       [ 2,  7, 12],
       [ 3,  8, 13],
       [ 4,  9, 14]])
1
2
3
arr12 = np.random.randn(6, 3)
arr12
np.dot(arr12.T, arr12)
array([[ 0.2929244 ,  0.08163964,  0.92965335],
       [-1.4212667 ,  1.76679612,  0.64797326],
       [-0.42830942,  1.59407581,  1.13162298],
       [-0.34227908, -0.89853783, -1.69788152],
       [ 0.22750799,  0.19027131, -0.61655676],
       [ 0.76157891,  0.27457782, -0.71410471]])

array([[ 3.03817   , -2.6099803 , -1.2362788 ],
       [-2.6099803 ,  6.58827763,  4.23684643],
       [-1.2362788 ,  4.23684643,  6.33758469]])
  • 简单的转置可以使用T,其实质就是进行轴对换。此外,ndarray还有一个swapaxes方法,需接受一对轴编号。
1
arr13
array([[[ 0,  1,  2,  3],
        [ 4,  5,  6,  7]],

       [[ 8,  9, 10, 11],
        [12, 13, 14, 15]]])
1
arr13.swapaxes(1,2)     # 交换1轴和2轴,即由原来的(2,2,4)三维数组变换为(2,4,2)三维数组
array([[[ 0,  4],
        [ 1,  5],
        [ 2,  6],
        [ 3,  7]],

       [[ 8, 12],
        [ 9, 13],
        [10, 14],
        [11, 15]]])

Numpy运算

简单运算

  • NumPy数组可以不用写循环对数据执行批量计算,这通常叫作矢量化
  • 大小相等的数组之间的任何算术运算都会应用到元素级
  • 不同形状的数组之间的算数运算执行方式,称为广播( Broadcasting)
1
2
arr7 = np.array([[1., 2., 3.], [4., 5., 6.]])
arr7
array([[ 1.,  2.,  3.],
       [ 4.,  5.,  6.]])
1
arr7 * arr7
array([[  1.,   4.,   9.],
       [ 16.,  25.,  36.]])
1
arr7 - arr7
array([[ 0.,  0.,  0.],
       [ 0.,  0.,  0.]])
  • 数组与标量的的算术运算也会将标量值传播到各个元素
1
1 / arr7
array([[ 1.        ,  0.5       ,  0.33333333],
       [ 0.25      ,  0.2       ,  0.16666667]])

通用函数:快速的元素级数组函数

通用函数(即ufunc)是一种对ndarray中的数据执行元素级运算的函数。许多ufunc函数都是简单的元素级变体,譬如np.sqrtnp.exp

1
2
arr16 = np.arange(10)
np.sqrt(arr16)
array([ 0.        ,  1.        ,  1.41421356,  1.73205081,  2.        ,
        2.23606798,  2.44948974,  2.64575131,  2.82842712,  3.        ])
1
np.exp(arr16)
array([  1.00000000e+00,   2.71828183e+00,   7.38905610e+00,
         2.00855369e+01,   5.45981500e+01,   1.48413159e+02,
         4.03428793e+02,   1.09663316e+03,   2.98095799e+03,
         8.10308393e+03])
  • np.expnp.sqrt都是一元ufunc。另外一些,譬如np.addnp.maximum接受两个数组,因此也叫二元ufunc,其返回一个结果数组。
1
2
x = np.random.randn(8)
x
array([-0.73503915,  0.67857597, -0.93048223, -1.18871658, -0.58625309,
        0.48785485,  0.4478705 ,  2.07487731])
1
2
y = np.random.randn(8)
y
array([-1.32420268,  0.83166662, -1.78110982,  0.34907941,  1.27166819,
       -0.90842411, -0.18884914,  1.34234266])
1
np.maximum(x, y)
array([-0.73503915,  0.83166662, -0.93048223,  0.34907941,  1.27166819,
        0.48785485,  0.4478705 ,  2.07487731])
1
np.isnan(x)
array([False, False, False, False, False, False, False, False], dtype=bool)
1
np.add(x,y)
array([-2.05924184,  1.51024258, -2.71159205, -0.83963717,  0.6854151 ,
       -0.42056927,  0.25902137,  3.41721997])
1
x+y
array([-2.05924184,  1.51024258, -2.71159205, -0.83963717,  0.6854151 ,
       -0.42056927,  0.25902137,  3.41721997])
  • 虽然不常见,但有些ufunc可以返回多个数组。譬如,np.modf,它是Python内置函数divmod的矢量化版本,用于计算浮点数数组的小数和整数部分。
1
2
arr17 = np.random.randn(7) * 5
np.modf(arr17)
(array([ 0.15108012, -0.0060202 , -0.29412582, -0.63925294, -0.90022958,
         0.04633377,  0.22060716]), array([ 6., -5., -0., -6., -1.,  6.,  4.]))
1
arr17
array([ 6.15108012, -5.0060202 , -0.29412582, -6.63925294, -1.90022958,
        6.04633377,  4.22060716])
  1. 一元ufunc
函数 说明
abs,fabs 计算整数、浮点数或复数的绝对值
sqrt 计算个元素的平方根
square 计算各元素的平方
exp 指数
log, log10, log2, log1p 对数
sign 计算各元素的正负号
ceil ceiling 值
floor floor值
rint 四舍五入
modf 返回小数部分和整数部分两个数组
isnan 哪些元素不是数字的布尔型数组
isfinite,isinf 有穷,无穷的布尔型数组
cos,cosh, sin, sinh, tan, tanh 普通型和双曲型三角函数
  1. 二元ufunc
函数 说明
add 数组中对应元素相加
subtract 从第一个数组中减去第二个数组中的元素
multiply 对应元素相称
divide, floor_divide 除法或向下圆整除法
power 第一个数组元素的第二个数组元素次幂\(A^B\)
maximum,fmax 元素级最大值,fmax将忽略NaN

广播

  • 不同大小的数组之间的运算叫作广播(broadcasting)

对矩阵减去列平均值进行距平化处理

1
2
arr = np.random.rand(4,3)
arr.mean(0)

array([0.34866595, 0.33478123, 0.54635458])

1
2
demand = arr - arr.mean(0)
demand

array([[-0.12673372, -0.32884738, -0.3755222 ], [-0.14293272, -0.01459838, 0.12151164], [ 0.03661667, -0.09273143, 0.26623119], [ 0.23304977, 0.43617719, -0.01222063]])

广播的原则:

  • 如果两个数组的后缘维度(trailing dimension,即从末尾开始算起的维度)的轴长度相符或其中一方的长度为1,则认为他们是广播兼容的。
  • 广播会在缺失和(或)长度为1的维度上进行
1
arr

array([[0.22193223, 0.00593384, 0.17083238], [0.20573323, 0.32018285, 0.66786622], [0.38528262, 0.2420498 , 0.81258577], [0.58171572, 0.77095842, 0.53413395]])

1
2
row_means = arr.mean(1)
row_means

array([0.13289949, 0.39792743, 0.47997273, 0.62893603])

1
row_means.reshape((4,1))

array([[0.13289949], [0.39792743], [0.47997273], [0.62893603]])

1
arr - row_means.reshape((4,1))

array([[ 0.08903274, -0.12696564, 0.0379329 ], [-0.1921942 , -0.07774458, 0.26993879], [-0.09469011, -0.23792293, 0.33261304], [-0.04722031, 0.14202239, -0.09480208]])

更复杂的高维数组的广播请进一步阅读相关资料。

利用数组进行数据处理

NumPy可以将许多数据处理任务表述为简介的数组表达式,从而避免写循环,通常称之为(矢量化)。矢量化数组运算要比等价的纯Python方式快一两个数量级,甚至更多,尤其是各种数值计算。

1
points = np.arange(-5, 5, 0.01)    # 1000个间隔相等的点
1
xs, ys = np.meshgrid(points, points)    # 接受两个一维数组,产生两个二维数组
1
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]])
1
xs
array([[-5.  , -4.99, -4.98, ...,  4.97,  4.98,  4.99],
       [-5.  , -4.99, -4.98, ...,  4.97,  4.98,  4.99],
       [-5.  , -4.99, -4.98, ...,  4.97,  4.98,  4.99],
       ..., 
       [-5.  , -4.99, -4.98, ...,  4.97,  4.98,  4.99],
       [-5.  , -4.99, -4.98, ...,  4.97,  4.98,  4.99],
       [-5.  , -4.99, -4.98, ...,  4.97,  4.98,  4.99]])
1
2
3
# 到原点的距离
z = np.sqrt((xs**2 + ys**2))
z
array([[ 7.07106781,  7.06400028,  7.05693985, ...,  7.04988652,
         7.05693985,  7.06400028],
       [ 7.06400028,  7.05692568,  7.04985815, ...,  7.04279774,
         7.04985815,  7.05692568],
       [ 7.05693985,  7.04985815,  7.04278354, ...,  7.03571603,
         7.04278354,  7.04985815],
       ..., 
       [ 7.04988652,  7.04279774,  7.03571603, ...,  7.0286414 ,
         7.03571603,  7.04279774],
       [ 7.05693985,  7.04985815,  7.04278354, ...,  7.03571603,
         7.04278354,  7.04985815],
       [ 7.06400028,  7.05692568,  7.04985815, ...,  7.04279774,
         7.04985815,  7.05692568]])
1
2
3
4
5
%matplotlib inline
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")
<matplotlib.image.AxesImage at 0x6fb7910>

<matplotlib.colorbar.Colorbar at 0x718cc90>

<matplotlib.text.Text at 0x6f9c4f0>

将条件逻辑表述为数组运算

numpy.where函数是三元表达式x if condition els y的矢量化版本。

1
2
3
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])
1
list(zip(xarr,yarr,cond))
[(1.1000000000000001, 2.1000000000000001, True),
 (1.2, 2.2000000000000002, False),
 (1.3, 2.2999999999999998, True),
 (1.3999999999999999, 2.3999999999999999, True),
 (1.5, 2.5, False)]

假设根据cond的值选取xarr或yarr的值:当cond中的值为True时,选取xarr,否则从yarr选取。列表生成式写法如下:

1
2
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]
1
np.where(cond,xarr,yarr)
array([ 1.1,  2.2,  1.3,  1.4,  2.5])

上述写法存在几个问题:

  1. 对大数组的处理速度不是很快,因为都是由纯Python完成
  2. 无法用于多维数组。

若使用np.where则更简洁:

1
2
result = np.where(cond, xarr, 0)
result
array([ 1.1,  0. ,  1.3,  1.4,  0. ])
  • np.where的第二个和第三个参数不必是数组,都可以是标量。
  • np.where通常用于根据一个数组而产生一个新的数组。

假设有一个随机数组成的矩阵,将所有的正值替换为2,所有的负值替换为-2。则

1
2
arr18 = np.random.randn(4,4)
arr18
array([[ 1.25303209,  1.18849831, -0.61261399, -1.35428929],
       [-0.35775028, -0.68108611, -0.33610514,  1.73650395],
       [ 0.37363973,  0.62543044, -1.48153305,  1.40538291],
       [ 0.85732867, -0.5978973 , -0.66642207, -1.1582354 ]])
1
np.where(arr18 > 0, 2, -2)
array([[ 2,  2, -2, -2],
       [-2, -2, -2,  2],
       [ 2,  2, -2,  2],
       [ 2, -2, -2, -2]])
1
np.where(arr18 > 0, 2, arr18)
array([[ 2.        ,  2.        , -0.61261399, -1.35428929],
       [-0.35775028, -0.68108611, -0.33610514,  2.        ],
       [ 2.        ,  2.        , -1.48153305,  2.        ],
       [ 2.        , -0.5978973 , -0.66642207, -1.1582354 ]])
  • 传递给np.where的数组大小可以不相等,甚至是标量。
  • np.where可以嵌套,例如:np.where(cond1 & cond2, 0, (np.where(cond1, 1, np.where(cond2, 2, 3)))

数学和统计方法

  • 可以通过数组上的一组函数对整个数组或某个轴向的数据进行统计计算。
  • sum,mean以及标准差std等聚合计算既可以当作数组的实例方法调用,也可以当作顶级NumPy函数使用。
1
2
arr19 = np.random.randn(5, 4)
arr19
array([[ 1.467288  , -1.05983669,  0.05658104,  0.82467192],
       [-1.1483169 , -0.77066991,  0.23524319,  0.97062944],
       [-0.50571417, -0.39945756,  0.19360303,  0.43184963],
       [ 2.36430206, -0.85428373, -0.34905838,  2.8272828 ],
       [-0.44314677, -2.19239688,  2.03344435,  0.01128036]])
1
arr19.mean()
0.18466474273947711
1
np.mean(arr19)
0.18466474273947711
1
arr19.sum()
3.6932948547895421
  • meansum这类函数可以接受一个axis参数,用于计算该轴上的统计值。最终结果是一个少一维的数组。
1
arr19.mean(axis = 1)
array([ 0.32217607, -0.17827854, -0.06992977,  0.99706069, -0.14770473])
1
arr19.mean(axis=0)
array([ 0.34688244, -1.05532895,  0.43396265,  1.01314283])
1
arr19.sum(axis=0)
array([ 1.73441222, -5.27664477,  2.16981324,  5.06571416])
  • 其它如cumsumcumprod之类的方法则不聚合,而是产生一个由中间结果组成的数组。
1
2
3
arr20 = np.array([[0, 1, 2], [3, 4, 5], [6, 7, 8]])
arr20
arr20.cumsum(axis=0)
array([[0, 1, 2],
       [3, 4, 5],
       [6, 7, 8]])

>array([[ 0,  1,  2],
>[ 3,  5,  7],
>[ 9, 12, 15]], dtype=int32)
1
arr20.cumprod(axis = 1)
array([[  0,   0,   0],
       [  3,  12,  60],
       [  6,  42, 336]], dtype=int32)

基本数组统计方法表。

方法 说明
sum 对数组中全部或某轴向的元素求和。零长度数组的sum为0
mean 算术平均数
std、var 分别为标准差和方差
min、max 最小值和最大值
argmin、argmax 分别为最小、最大元素的索引
cumsum 所有元素的累计和
cumprod 所有元素的累计积

用于布尔型数组的方法

在前述方法中,布尔值会被强制转换为1(True)和0(False)。因此sum经常被用于对布尔型数组中的True值计数。

1
2
arr21 = np.random.randn(100)
(arr21 > 0).sum() # 正值的个数
49
1
(arr21> 0).sum()
49
  • 方法anyall对布尔型数组也非常有用:any用于测试数组中是否存在一个或多个True, 而all则检查数组中所有值是否都是True
  • 这两个方法也可以用于非布尔型数组,所有非零元素都会被当作True
1
2
bools = np.array([False, False, True, False])
bools.any()
True
1
bools.all()
False
1
2
data3 = np.array([1.1, 1.2, 0, 2, 3])
data3.any()
True
1
data3.all()
False

排序

跟Python内置的列表类型一样,NumPy数组也可以通过sort方法就地进行排序。

1
2
arr22 = np.random.randn(8)
arr22
array([ 1.37429684, -1.07242886,  1.70060711, -0.83261411,  1.72842932,  1.10081484,  0.47399465,  0.63373348])
1
2
arr22.sort()
arr22
array([-1.07242886, -0.83261411,  0.47399465,  0.63373348,  1.10081484,
        1.37429684,  1.70060711,  1.72842932])
  • 多维数组可以在任何一个轴上进行排序,只需将轴编号作为参数传入sort即可:
1
2
arr23 = np.random.randn(5,3)
arr23
array([[ 0.79301426, -0.29812058, -0.88582163],
       [ 1.61494158,  1.31626884,  0.3309185 ],
       [ 0.46263742, -0.58083798, -0.36151492],
       [ 0.42998483,  0.07636906, -1.00172838],
       [-1.37108813,  0.98236945,  0.2967868 ]])
1
2
sor= np.sort(arr23)
sor
array([[-0.88582163, -0.29812058,  0.79301426],
       [ 0.3309185 ,  1.31626884,  1.61494158],
       [-0.58083798, -0.36151492,  0.46263742],
       [-1.00172838,  0.07636906,  0.42998483],
       [-1.37108813,  0.2967868 ,  0.98236945]])
1
2
arr23.sort(axis=1)
arr23
array([[-0.88582163, -0.29812058,  0.79301426],
       [ 0.3309185 ,  1.31626884,  1.61494158],
       [-0.58083798, -0.36151492,  0.46263742],
       [-1.00172838,  0.07636906,  0.42998483],
       [-1.37108813,  0.2967868 ,  0.98236945]])
  • 顶级方法sort返回的是数组的已排序副本,而就地排序则会修改数组本身。
  • 计算数组分位数最简单的办法就是对数组进行排序,然后选取特定位置的值。
1
2
3
large_arr = np.random.randn(1000)
large_arr.sort()
large_arr[int(0.05 * len(large_arr))] # 5%分位数
-1.6694434719038056

唯一化以及其他的一些集合逻辑

NumPy提供一些针对一维ndarray的基本集合运算。譬如,np.unique,找出数组中的唯一值并返回已排序的结果。

1
names = np.array(['Bob', 'Joe', 'Will', 'Bob', 'Will', 'Joe', 'Joe'])
1
np.unique(names)
array(['Bob', 'Joe', 'Will'], dtype='<U4')
1
2
ints = np.array([3,3,3,2,2,1,1,4,4])
np.unique(ints)
array([1, 2, 3, 4])
1
sorted(set(names))
['Bob', 'Joe', 'Will']
  • 函数np.in1d用于测试一个数组中的值在另一个数组中的成员资格,返回一个布尔型数组。
1
2
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)

数组的集合运算函数表:

方法 说明
unique(x) 计算x中的唯一元素,并返回有序结果
intersect1d(x,y) 计算x和y的公共元素,并返回有序结果
union1d(x,y) 计算x和y的并集,并返回有序结果
in1d(x,y) 得到一个表示“x的元素是否包含于y”的布尔型数组
setdiff1d(x,y) 集合的差,即元素在x中且不在y中
setxor1d(x,y) 集合的对称差,即存在于一个数组中但不同时存在于两个数组中的元素

用于数组的文件输入输出

NumPy能够读写磁盘上的文本数据或二进制数据。后面讲到的pandas用于将表格型数据读取到内存中。

读取二进制文件

将数组以二进制格式保存到磁盘

np.savenp.load是读写磁盘数组的两个主要函数。默认情况下,数组是以未压缩的原始二进制格式保存在扩展名为 .npy 的文件中。

1
2
arr24 = np.arange(10)
np.save('arr24', arr24)
1
np.load('arr24.npy')
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
  • 通过np.savez可以将多个数组保存到一个压缩文件中,将数组以关键字参数的形式传入即可:
  • 加载 .npz 文件时,会得到一个类似字典的对象。
1
2
arr25 = np.random.randn(3, 4)
np.savez('arr_archived.npz', a = arr24, b = arr25)

读取二进制文件

1
2
arch = np.load('arr_archived.npz')
arch['b']
array([[ 0.22327739, -1.43571685, -1.51532615, -0.57348726],
       [-0.78856893,  1.93887747, -0.62795721, -0.89444234],
       [ 0.55478518, -1.09806106, -0.03176941, -0.49335069]])
1
arch['a']
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])

存取文本文件

Python也有文件的读写函数,但并不太易用。因此本课程主要介绍Pandas中的read_csvread_table函数,有时也用NumPy的np.loadtxt将数据加载到NumPy的数组中。

1
2
arr = np.loadtxt('arr_csv.csv', delimiter=',')
arr
array([[ 0.21588889,  0.86253552,  0.2659228 ],
       [ 0.38891274,  0.39542911,  0.52475854],
       [ 0.55477862,  0.1809702 ,  0.7752453 ],
       [ 0.10663222,  0.69513153,  0.06849945],
       [ 0.57416129,  0.61863014,  0.27521237],
       [ 0.63571405,  0.73961201,  0.00164165]])
1
np.savetxt(delimiter='\t',fname='newsave.csv', X=arr)

更多文件的读写,尤其是表格型数据,采用pandas方法。

Numpy常见库

线性代数库

线性代数(矩阵乘法、矩阵分解、行列式等)。NumPy提供了一个用于矩阵乘法的 dot 函数,既是一个数组的方法,也是numpy命名空间的一个函数。

1
2
3
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]])
1
y
array([[ 6., 23.],
       [-1.,  7.],
       [ 8.,  9.]])
1
x.dot(y)
array([[ 28.,  64.],
       [ 67., 181.]])
1
np.dot(x, y)

array([ 6., 15.])

1
np.dot(x, np.ones(3))

线性代数库linalg

  • numpy.linalg中有一组标准的矩阵分解运算,譬如求逆、行列式等。
1
2
3
X = np.random.randn(5,5)
mat = X.T.dot(X)
mat
array([[ 7.79231619,  3.1955974 , -0.68745836,  4.63211203, -0.0837485 ],
       [ 3.1955974 ,  5.06344351, -1.86563103,  2.04927385, -1.2547239 ],
       [-0.68745836, -1.86563103,  3.27810894,  0.13828478,  1.38250281],
       [ 4.63211203,  2.04927385,  0.13828478,  3.56053886, -1.04325734],
       [-0.0837485 , -1.2547239 ,  1.38250281, -1.04325734,  6.94635016]])
1
np.linalg.inv(mat)
array([[ 0.92312571,  0.01304268,  0.36260544, -1.29679853,> -0.25344566],
       [ 0.01304268,  0.3676009 ,  0.22858607, -0.24189231, -0.01526656],
       [ 0.36260544,  0.22858607,  0.61529209, -0.67960302, -0.17886571],
       [-1.29679853, -0.24189231, -0.67960302,  2.25503587,  0.41460951],
       [-0.25344566, -0.01526656, -0.17886571,  0.41460951,  0.23601542]])
1
mat.dot(np.linalg.inv(mat))
array([[ 1.00000000e+00, -1.33356867e-16, -4.47558657e-16,  3.33066907e-16, -2.15105711e-16],
       [ 1.11022302e-16,  1.00000000e+00,  1.11022302e-16,   -2.22044605e-16,  1.11022302e-16],
       [ 5.55111512e-17,  8.32667268e-17,  1.00000000e+00,    -6.66133815e-16, -5.55111512e-17],
       [ 1.66533454e-16, -2.42861287e-17, -3.05311332e-16,    1.00000000e+00, -1.11022302e-16],
       [-2.22044605e-16,  4.16333634e-17,  2.22044605e-16,      4.44089210e-16,  1.00000000e+00]])
1
2
from numpy.linalg import inv, qr
inv(mat)
array([[ 0.92312571,  0.01304268,  0.36260544, -1.29679853,  -0.25344566],
       [ 0.01304268,  0.3676009 ,  0.22858607, -0.24189231, -0.01526656],
       [ 0.36260544,  0.22858607,  0.61529209, -0.67960302, -0.17886571],
       [-1.29679853, -0.24189231, -0.67960302,  2.25503587,  0.41460951],
       [-0.25344566, -0.01526656, -0.17886571,  0.41460951,  0.23601542]])

常用的numpy.linalg函数

函数 说明
diag 以一维数组形式返回方阵的对角线元素;或将一维数组转换为方阵(非对角线元素为0)
dot 矩阵乘法
trace 计算对角线元素之和
det 计算矩阵行列式
eig 计算方阵的特征值和特征向量
inv 计算方阵的逆
solve 求解线性方程组\(Ax=b\)
lstsq 计算\(Ax=b\)的最小二乘解

随机数生成

Python的random模块

  • Python 的random 模块可以用于产生随机数
1
2
import random
random.random()

0.3176503054519141

1
random.uniform(1, 5)

2.4497935064174916

1
random.choice([1, 3, 5, 7, 9])

7

1
2
3
seq = [2, 4, 6, 8, 10]
random.shuffle(seq)
seq

[10, 6, 8, 4, 2]

1
2
samp = random.sample(seq, 3)
samp

[2, 8, 10]

Python的random函数

函数 说明
random [0,1)均匀分布随机数
randint 给定区间的整数
uniform 给定区间均匀分布随机数
shuffle 将给定序列随机打乱
choice 从非空序列中返回一个随机元素
sample 从总体序列或集合中抽取k个,无重复抽样

numpy.random

  • numpy.random模块对Python内置的random进行了补充,增加了一些高效生产多种概率分布的样本值的函数。
1
2
samples = np.random.normal(size=(4, 4))
samples
array([[ 0.17877605, -0.96452746,  0.90583362, -0.57433145],
       [ 0.90285923,  0.30519106, -0.85516783,  0.11017296],
       [ 0.02661199, -1.59621354,  0.75243553, -1.79864872],
       [ 0.45939344, -0.38791448,  0.62207872,  2.16328784]])
1
2
N = 1000000
%time test = [i**2 for i in range(N)]

CPU times: user 330 ms, sys: 8.93 ms, total: 339 ms Wall time: 344 ms

1
import random
1
%timeit sample1 = [random.random() for _ in range(N)]

168 ms ± 316 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

1
%timeit sample2 = np.random.random(N)

18.4 ms ± 2.37 ms per loop (mean ± std. dev. of 7 runs, 100 loops each)

部分numpy.random函数

函数 说明
seed 确定随机数生成器的种子
permutation 返回一个序列的随机排列或返回一个随机排列的范围
shuffle 对一个序列就地随机排列
rand 产生均匀分布的样本值
randint 从给定的上下限范围内随机生成整数
randn 标准正态分布抽样值
normal 正态分布样本值
chisquare 卡方分布样本值
uniform [0,1)均匀分布样本值
choice 从从一维数组中抽样

范例:随机漫步

从0开始,步长1和-1出现的概率相等,实现1000步的随机漫步。

1
2
3
4
5
# 步数
nSteps = 1000
draws = np.random.randint(0, 2, size=nSteps) # 1000个0或1的随机数
steps = np.where(draws > 0, 1, -1)
walk = steps.cumsum() # 累积求和,即当前位置
1
2
3
4
5
6
# 绘制图
from matplotlib import pyplot as plt

fig = plt.figure()
plt.plot(walk)
plt.show()