本文介绍如何使用 numpy,一个强大的 Python 科学计算库。

对于 Windows 系统的用户,推荐到加州大学尔湾分校的 Python whl 库下载站 这里下载与 Intel mkl 链接好的 numpy 版本;注意选择与用户计算机的 python & 操作系统相适应的 numpy 版本。链接了 mkl 与未链接 mkl 的 numpy 在性能上会有显著差异,你可以前往上述网址安装链接了 mkl 库的 numpy。请在安装其他依赖于 numpy 的包或库前,先安装好 numpy。

全文默认加载的 numpy 包:

import numpy as np

本文最后运行时的 numpy 版本:

print(np.__version__)
1.11.2

数组创建

numpy 最基础的操作就是创建数组。数组可以是高维(多于二维)的;在本文中,把二维的数组也叫做矩阵

用 np.array() 创建数组:

dt1 = [1, 2, 3, 4, 5, 6]
arr1 = np.array(dt1)
print(arr1)  # 一行六列的数组
[1 2 3 4 5 6]
dt2 = [[1, 2, 3], [4, 5, 6]]
arr2 = np.array(dt2)
print(arr2)  # 自动识别为两行三列
[[1 2 3]
 [4 5 6]]

用 arr.shape/ndim/size 来返回数组的尺寸:

print(arr2.shape, '\n',  # 返回结果是一个元组(tuple)
      arr2.ndim, '\n',   # 返回结果是数组的维度(即“轴”的数量,这里只有行和列)
      arr2.size)         # 返回结果是数组的总元素个数
(2, 3) 
 2 
 6

特殊矩阵

预分配矩阵:np.empty() / empty_like()

预分配矩阵只是初始化了矩阵尺寸,但是不保证元素值为0。如果想要生成全零阵,使用 np.zeros()。

np.empty([2, 3])
array([[ 0.,  0.,  0.],
       [ 0.,  0.,  0.]])
np.empty_like(arr2)  # 生成一个与 arr2 同尺寸的预分配矩阵
array([[0, 0, 0],
       [0, 0, 0]])

全0 / 全1 / 填充矩阵:np.zeros() / ones() / full()

全 0 与全 1 矩阵很简单:

np.ones([1, 4])
array([[ 1.,  1.,  1.,  1.]])

填充矩阵用指定的数值填充所有的元素:

np.full([2, 3], 1.2)
array([[ 1.2,  1.2,  1.2],
       [ 1.2,  1.2,  1.2]])

以上矩阵也有 like() 形式的命令。

np.zeros_like(arr2)
array([[0, 0, 0],
       [0, 0, 0]])

单位阵:np.identity() / eye()

严格的单位阵使用 np.identity() 命令产生:

np.identity(3)  # 严格正方形单位阵
array([[ 1.,  0.,  0.],
       [ 0.,  1.,  0.],
       [ 0.,  0.,  1.]])

如果是非严格的单位阵,使用 np.eye() 产生:

np.eye(2, 3)  # 不严格的、对角线为 1 的矩阵
array([[ 1.,  0.,  0.],
       [ 0.,  1.,  0.]])

也可以用参数 k 指定相对主对角线的偏移量。偏移方向以向右上为正,左下为负:

np.eye(3, 3, k=1)
array([[ 0.,  1.,  0.],
       [ 0.,  0.,  1.],
       [ 0.,  0.,  0.]])

对角阵:np.diagflat()

同样使用参数 k 指定偏移量。

np.diagflat([1, 2, 3], k=1)
array([[0, 1, 0, 0],
       [0, 0, 2, 0],
       [0, 0, 0, 3],
       [0, 0, 0, 0]])

下三角阵:np.tri()

下三角全为 1(包括对角线),其余为 0. 可以指定偏移量:

np.tri(3, 4, k=1)
array([[ 1.,  1.,  0.,  0.],
       [ 1.,  1.,  1.,  0.],
       [ 1.,  1.,  1.,  1.]])

要截取某个已知矩阵的上三角或者下三角,参考下文的内容。

范德蒙矩阵:np.vander()

范德蒙矩阵(Vandermonde matrix)的定义:

可以用参数 N 指定列数。

np.vander([1, 2, 3], N=4)
array([[ 1,  1,  1,  1],
       [ 8,  4,  2,  1],
       [27,  9,  3,  1]])

特殊数列

numpy 里有一系列可以快速生成数列的命令。

等差数列:np.arange()

生成的数列区间,默认是左闭右开 [start, stop)。

np.arange(start=1, stop=7, step=2)
array([1, 3, 5])

线性 / 对数空间:np.linspace() / logspace()

线性空间也是等差的;相比 arange() 命令只能指定步长,线性空间可以指定数组尺寸:

np.linspace(1, 2, num=5)
array([ 1.  ,  1.25,  1.5 ,  1.75,  2.  ])

也可以利用 endpoint 参数,产生类似 arange() 函数的效果:

np.linspace(1, 2, 5, endpoint=False)
array([ 1. ,  1.2,  1.4,  1.6,  1.8])

logspace() 是以从 $10^{start}$ 到 $10^{end}$ 的对数空间:

np.logspace(1, 2, num=4)
array([  10.        ,   21.5443469 ,   46.41588834,  100.        ])

从已有数据转换 / 复制

深度复制:np.copy()

arr3 = np.copy(arr2)
arr3 is arr2
False
arr3 = arr2.copy()
arr3 is arr2
False

类型转换:np.asarray() / asmatrix()

注意:这类函数不会生成一个复制;它们只生成一个同址引用

tmp = np.asarray(arr2) 
arr2[0, 0] = 100
tmp
array([[100,   2,   3],
       [  4,   5,   6]])
tmp is arr2
True

np.asmatrix() 函数非常灵活:

np.asmatrix("1, 2; 3, 4")
matrix([[1, 2],
        [3, 4]])

尺寸重塑:arr.reshape() / resize() / flatten() / ravel()

reshape() 生成一个重塑后的矩阵,但不改动原矩阵;resize() 则直接将原矩阵重塑。

arr2.reshape(3, 2)
array([[100,   2],
       [  3,   4],
       [  5,   6]])
arr2.resize(3, 2)
arr2
array([[100,   2],
       [  3,   4],
       [  5,   6]])

arr.flatten() 直接将矩阵转为一个向量。参数 order 指定了转换方式(此处只介绍两个参数):

  • ‘C’:默认。按行读取,行内从左到右。C-style.
  • ‘F’:按列读取,列内从上到下。Fortran-style.
arr2.flatten('F')
array([100,   3,   5,   2,   4,   6])

np.ravel() 的作用是一样的,不过会优先返回原值的引用(如果可能)。此外,ravel() 能够应用于非 ndarray 对象(比如一系列 ndarray 组成的列表),而 flatten() 不能。

转置:arr.T / arr.transpose()

对于二维数组而言,它们没有区别:

arr2.T
array([[100,   3,   5],
       [  2,   4,   6]])

从已有数据截取 / 拼接

切片

arr = np.linspace(1, 9, 9).reshape(3, 3)
arr[:, :2]
array([[ 1.,  2.],
       [ 4.,  5.],
       [ 7.,  8.]])
# 第三行
arr[2]
array([ 7.,  8.,  9.])

截取矩阵对角线:np.diag()

np.diag(arr)
array([ 1.,  5.,  9.])

截取矩阵上 / 下三角:np.triu() / np.tril()

np.triu() 函数用于截取矩阵的上三角(Upper Triangle)部分:

np.triu(arr)
array([[ 1.,  2.,  3.],
       [ 0.,  5.,  6.],
       [ 0.,  0.,  9.]])

可以设置偏移参数 k :

np.tril(arr, k=1)
array([[ 1.,  2.,  0.],
       [ 4.,  5.,  6.],
       [ 7.,  8.,  9.]])

拼接矩阵:np.vstack() / hstack()

竖向与横向两个方向:

a = np.arange(1, 5).reshape(2, 2)
b = np.arange(5, 9).reshape(2, 2)

np.vstack((a, b))  # 也可以写 np.vstack([a, b])
array([[1, 2],
       [3, 4],
       [5, 6],
       [7, 8]])
np.hstack([a, b])
array([[1, 2, 5, 6],
       [3, 4, 7, 8]])

拼接网格阵:np.meshgrid()

由两个列向量组成的二维坐标阵,常常用于绘图命令。

参数 sparse 为 True 表示尽可能地节省内存。

from matplotlib import pyplot as plt

x = np.arange(-5, 5, 0.1)
y = np.arange(-5, 5, 0.1)
xx, yy = np.meshgrid(x, y, sparse=True)
z = np.sin(xx**2 + yy**2) / (xx**2 + yy**2)

plt.contourf(x,y,z)
plt.show()

png

矩阵操作

矩阵基础运算

a = np.linspace(1, 6, 6).reshape(2, 3)
b = np.linspace(2, 7, 6).reshape(3, 2)
np.dot(a, b)  # 或者用 a.dot(b) 也可以
array([[ 28.,  34.],
       [ 64.,  79.]])
a + b.T
array([[  3.,   6.,   9.],
       [  7.,  10.,  13.]])

元素级运算

一元函数

直接从元素操作中移植的函数,比如幂指函数等等。一元函数调用形式是 np.function(arr) 简表如下:

一元函数 含义
abs/fabs 求模。对于非复数矩阵,fabs更快。
sqrt/square 平方根/平方
exp 指数
log/log10/log2/log1p 对数:自然对数/10/2为底,以及以2为底的log(1+x)
sign 正负号判断函数:返回 1, 0, 或者 -1
ceil/floor 向上取整/向下取整
rint 四舍五入到整数
isnan 返回关于非数值(np.nan)判断的布尔型数组
isfinite/isinf 判断非无限大值/无限大值
sin/cos/tan 三角函数
arcsin/arccos/arctan 反三角函数
sinh/cosh/tanh/arcsinh/arccosh/arctanh 以上三角函数的双曲形式

二元函数

以及矩阵的二元函数,调用形式是 np.function(x, y):

二元函数 含义
add/substract/multiply 对应元素相加/相减/相乘
divide/floor_divide 对应元素除法及向下整除法(弃余数)
power 计算 x(i,j)^y(i,j)
maximum/fmax 元素级的最大值。fmax 表示忽略 NaN
  注:在比较含有 NaN 的矩阵时可能出现问题,我尚不清楚 NumPy 做出了怎样的改变。
minimum/fmin 仿上
mod 取余
copysign 将y的符号传递给x中的对应元素

统计函数

一些与统计相关的函数,或者简单分析矩阵数据特征的函数:

统计函数 含义
average/mean 平均值。其中 average 还可以用 weight 参数指定权重
median 中位数
diff 一阶差分
cumsum/cumprod 累和/累积
sum/prod 求和 / 求积
std/var 标准差 / 方差

一些例子:

np.abs(np.array([3 + 4j, -3]))  # 求模
array([ 5.,  3.])
np.average(np.arange(6).reshape((3,2)), axis=1, weights=[1./4, 3./4])  # 加权平均值
array([ 0.75,  2.75,  4.75])
arr2 = np.array([7, 6])
arr3 = np.array([1, np.inf])

np.fmax(arr2, arr3)   # 求最大值
array([  7.,  inf])

线性代数

在 np.linalg 中有许多线性代数的函数。

逆矩阵:np.linalg.inv()

arr = np.array([[1, 2], [4, 3]])
np.linalg.inv(arr)
array([[-0.6,  0.4],
       [ 0.8, -0.2]])

迹:np.trace()

np.trace(arr)  # 矩阵的迹(主对角线元素之和)
4

线性方程组:np.linalg.solve() / lstsq()

前者是精确解,后者是最小二乘法的解。

Y = np.array([5, 10]).reshape(2, 1)
np.linalg.solve(arr, Y)  # 解线性方程组 arr * X = Y
array([[ 1.],
       [ 2.]])
np.linalg.lstsq(arr, Y)
(array([[ 1.],
        [ 2.]]),
 array([], dtype=float64),
 2,
 array([ 5.39834564,  0.92620968]))

QR 分解:np.linalg.qr()

返回值是元组。

np.linalg.qr(Y)
(array([[-0.4472136 ],
        [-0.89442719]]), array([[-11.18033989]]))

其他矩阵函数

一些矩阵函数并不是为矩阵运算准备的,可能只是为了筛选数据的需要,或者是一些其他的原因。

条件逻辑:np.where()

np.where(condition_arr, xarr, yarr) 表示:如果condition_arr 对应值为真,那么从 xarr 取值;否则从 yarr 取值。

xarr = np.array([[1, 2, 3],[-1, -2, 0], [4, 5, 6]])
np.where(xarr >= 0, xarr, np.nan)  # 把小于零的全部替换为 NaN
array([[  1.,   2.,   3.],
       [ nan,  nan,   0.],
       [  4.,   5.,   6.]])

布尔矩阵:arr.any() / all()

boolarr = np.array([False, True, False])
print(boolarr.any(),  # 是否存在 True
      boolarr.all())  # 是否全是 True
True False

唯一化:np.unique()

所有相同元素只保留一个,并排序,实质上返回了一个有序集合

xarr = np.array([[1, 2, 2, 5], [3, -1, 4, 4]])
yarr = np.array([[-1, 2, 2, 6], [-2, -1, 1, 3]])
np.unique(xarr)
array([-1,  1,  2,  3,  4,  5])

集合运算

将矩阵的同元素视作一个,进行类似集合的运算。

np.intersect1d(xarr, yarr),  # 交集
(array([-1,  1,  2,  3]),)
np.union1d(xarr, yarr),  # 并集
(array([-2, -1,  1,  2,  3,  4,  5,  6]),)
np.in1d(xarr, yarr),  # xarr 中的对应元素是否在 yarr 中
(array([ True,  True,  True, False,  True,  True, False, False], dtype=bool),)
np.setdiff1d(xarr, yarr),  # 在 xarr 中但不在 yarr 中的元素
(array([4, 5]),)
np.setxor1d(xarr, yarr)  # 仅存在于 xarr / yarr 之一的元素  
array([-2,  4,  5,  6])

保存与读取

# 保存
# np.save(filename_string, arr)  # 把单个数组保存到 filename_string.npy 中
# np.savez(zipname_string, a=arr1, b=arr2, ...)  # 把多个数组压缩保存到 zipname_string.npz 中

# 读取
# arr = np.load('filename_string.npy')  # 直接获取单个数组
# arrzip = np.load('zipname_string.npz')  # 将压缩包读取到字典 arrzip 中
# arr1 = arrzip['a']

# 保存 / 读取 txt
# arr = np.loadtxt('txtname.txt', delimeter=',')  # 保存 txt 用 savetxt