numpy入门实战教程(进阶篇)

Abstract:numpy入门实战教程进阶篇,附代码。

1.广播法则(Broadcast Rule)

广播法则能使通用函数有意义地处理不具有相同形状的输入

  • 第一法则:若所有的输入数组维度不都相同,一个“1”将被重复地添加在维度较小的数组上直至所有数据都拥有相同的维度
  • 第二法则:确定长度为1的数组沿着特殊方向表现地好像它有沿着那个方向最大形状的大小。对数组来说,沿着那个维度的数组元素的值理应相同

应用广播法则后,所有数组的大小必须匹配

2.花哨的索引和索引技巧

numpy比普通的python序列提供更多的索引功能。数组的索引功能:索引整数、切片、被整数数组和布尔数组索引。

2.1 通过数组索引(Indexed by Array)

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
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
# indexed via array
a = arange(12)**2
print a
i = array([1,1,3,8,5])
print i
print a[i] # # the elements of a at the positions i

j = array([[3,4],[9,7]])
print a[j]

# 当被索引数组a是多维的时,每一个唯一的索引数列指向a的第一维5。
# 以下示例通过将图片标签用调色版转换成色彩图像展示了这种行为
palette = array([
[0,0,0],
[255,0,0],
[0,255,0],
[0,0,255],
[255,0,255],
[255,255,255]])
image = array([[0,1,2,0],
[0,3,4,0]])
print palette[image]

# 我们也可以给出不不止一维的索引,每一维的索引数组必须有相同的形状
a = arange(12).reshape(3,4)
print a
i = array([
[0,1],
[1,2]])
j = array([[2,1],[3,3]])
print a[i,j]
print a[i,2]
print a[:,j]
print a[i,:]
l = [i,j]
print a[l] # 与 a[i,j] 相等

# 搜索时间序列最大值6
time = linspace(20, 145, 5) # time scale
data = sin(arange(20)).reshape(5,4) # 4 time-dependent series
print time
print data

ind = data.argmax(axis=0) # index of the maxima for each series
print ind

time_max = time[ind]
data_max = data[ind, range(data.shape[1])]
print time_max
print data_max
print all(data_max == data.max(axis=0)) # True
a = arange(5)
print a
# 当一个索引列表包含重复时,赋值被多次完成,保留最后的值
a = arange(5)
a[[0,0,2]] = [1,2,3]
print a

2.2 通过布尔数组索引(Indexed by Boolean Function)

当我们使用整数数组索引数组时,我们提供一个索引列表去选择。通过布尔数组索引法使我们可以显式地选择数组中我们想要和不想要的元素

我们能想到的使用布尔数组的索引最自然方式是使用和源数组一样形状的布尔数组

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
# indexed via boolean array
a = arange(12).reshape(3,4)
b = a > 4
print b
print a[b]
a[b] = 0 # All elements of 'a' higher than 4 become 0(取反)
print a

# 通过布尔来索引的方法更近似于整数索引;对数组的每个维度我们给一个一维布尔数组来选择我们想要的切片
a = arange(12).reshape(3,4)
b1 = array([False,True,True])
b2 = array([True,False,True,False])
print a[b1,:]
print a[b1]
print a[:,b2]
print a[b1,b2]
# 注意一维数组的长度必须和你想要切片的维度或轴的长度一致,
# 在之前的例子中,b1是一个秩为1长度为三的数组(a的行数),b2(长度为4)与a的第二秩(列)相一致

2.3 ix()函数

ix()函数可为了获得多元组的结果而用来结合不同向量。如,若你想要用所有向量a、b和c元素组成的 三元组来计算a+b*c:

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
# ix_() function
a = array([2,3,4,5])
b = array([8,5,4])
c = array([5,4,6,8,3])
ax,bx,cx = ix_(a,b,c)
print ax,bx,cx
print ax.shape,bx.shape,cx.shape

result = ax+bx*cx
print result
print result[3,2,4]

result1 = ax*bx+cx
print result1
print result1[3,2,4]
print result[3,2,4]
print a[3] + b[2] + c[4]

# second method
def ufunc_reduce(ufct, *vectors):
vs = ix_(*vectors)
r = ufct.identity
for v in vs:
r = ufct(r,v)
return r
print ufunc_reduce(add,a,b,c)
# 这个reduce与ufunc.reduce(比如说add.reduce)相比的优势在于它利用了广播法则,
# 避免了创建一个输出大小乘以向量个数的参数数组。

2.4 字符串索引(Indexed by String)

参见Structured arrays

3.线性代数(Linear Algebra)

3.1 简单数组运算

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
# coding:utf-8
from numpy import *
from numpy.linalg import *

# Linear Algebra
a = array([[1.0,2.0],[3.0,4.0]])
print (a)
print a.transpose()
print inv(a) # 计算矩阵的(乘法)逆
u = eye(2) # unit 2X2 matrix; "eye" represents "I",指元数组
print u
# j = array([[0.0, -1.0],[1.0, 0.0]])
# dot(i,j) # matrix product; dot():两个数组的点积
print trace(u) # trace():沿数组的对角线返回总和
y = array([[5.],[7.]])
print solve(a, y) # solve():求解线性矩阵方程或线性标量方程组

3.2 矩阵类

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
# coding:utf-8
from numpy import *
from numpy.linalg import *

# matrix
A = matrix('1.0 2.0; 3.0 4.0')
print A
print type(A)
print A.T # transpose
X = matrix('5.0 7.0')
Y = X.T
print Y
matrix([[5.],[7.]])
print (A * Y)
print (A.I) # inverse
print solve(A,Y)

3.3 索引:比较矩阵和二维数组

注意numpy中数组和矩阵有重要区别。numpy提供了2个基本对象:(其他对象都是建构在它们之上的)

  • 一个N维数组对象

  • 一个通用函数对象

矩阵:继承自numpy数组对象的二维数组对象。对数组和矩阵,索引都必须包含合适的一个或多个这些组合:整数标量、省略号、整数列表;布尔值,整数,布尔值构成的元组,和一个一维整数或布尔值数组。矩阵可被用作矩阵的索引,但通常需要数组、列表或其他形式来完成这个任务。

  • 索引从0开始

  • 使用矩阵的行和列表示一个二维数组或矩阵。沿0轴的方向被穿过的称作行,沿1轴的方向被穿过的是列。

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
32
33
34
35
36
37
38
39
40
41
42
# coding:utf-8
from numpy import *
from numpy.linalg import *

# index:compare matrix with 2d array
A = arange(12)
print A
A.shape = (3,4)
M = mat(A.copy()) # mat():将输入解释为矩阵
print (type(A)," ",type(M))
print (A)
print (M)
print (A[:])
print (A[:].shape)
print (M[:])
print (M[:].shape)
print (A[:,1])
print (A[:,1].shape)
# 注意最后两个结果的不同。对二维数组使用一个冒号产生一个一维数组,然而矩阵产生了一个二维矩阵。
# 例如,一个M[2,:]切片产生了一个形状为(1,4)的矩阵,相比之下,一个数组的切片总是产生一个最低可能维度11的数组。
# 例如,如果C是一个三维数组,C[...,1]产生一个二维的数组而C[1,:,1]产生一个一维数组。
# 从这时开始,如果相应的矩阵切片结果是相同的话,我们将只展示数组切片的结果。

# 假如我们想要一个数组的第一列和第三列,一种方法是使用列表切片:
print A[:,[1,3]]
# 稍微复杂点的方法是使用 take() 方法 method:
print A[:,].take([1,3],axis=1) # take():从轴沿一个数组中取元素
# 如果我们想跳过第一行,我们可以这样:
print A[1:,].take([1,3],axis=1)
print A[ix_((1,2),(1,3))]

# 保留第一行大于1的列。一种方法是创建布尔索引:
print A[0,:] > 1
print A[:,A[0,:]] > 1
print M[0,:] > 1
print M[:,M.A[0,:]>1]

# 如果我们想要在矩阵两个方向有条件地切片,我们必须稍微调整策略,代之以:
print A[A[:,0] > 2, A[0,:] > 1]
print M[M.A[:,0]>2,M.A[0,:]>1]
print A[ix_(A[:,0] > 2, A[0,:] > 1)]
print M[ix_(M.A[:,0] > 2, M.A[0,:] > 1)]

4.技巧和Tipps

4.1 “自动”改变形状

更改数组的维度,你可以省略一个尺寸,它将被自动推导出来。

1
2
3
4
5
# change the shape automatically
a = arange(30)
a.shape = 2,-1,3 # -1 means "whatever is needed"
print a.shape
print a

4.2 向量组合(Vector Stacking)

我们如何用两个相同尺寸的行向量列表构建一个二维数组?在 MATLAB 中这非常简单:如果 xy 是两个相同长度的向量,你仅仅需要做 m=[x;y]。在 Numpy 中这个过程通过函数 column_stackdstackhstackvstack 来完成,取决于你想要在那个维度上组合。例如:

1
2
3
4
5
6
7
8
9
10
11
12
# vector stacking
x = arange(0,10,2)
print x

y = arange(5)
print y

m = vstack([x,y]) # 按行顺序堆叠数组
print m

xy = hstack([x,y]) # 按列顺序组合数组
print xy

4.3 直方图(histogram)

Numpyhistogram函数应用到一个数组返回一对变量:直方图数组和箱式向量。注意:matplotlib也有一个用来建立直方图的函数(叫作hist,正如matlab中一样)与Numpy中的不同。
主要的差别是pylab.hist自动绘制直方图,而numpy.histogram仅仅产生数据

1
2
3
4
5
6
7
8
9
10
11
12
13
# histogram
# Build a vector of 10000 normal deviates with variance 0.5^2 and mean 2
mu, sigma = 2, 0.5
v = numpy.random.normal(mu,sigma,10000)
# Plot a normalized histogram with 50 bins
pylab.hist(v, bins=50, normed=1)
pylab.title('Matplotlib Version')# matplotlib version (plot)
pylab.show()
# Compute the histogram with numpy and then plot it
(n, bins) = numpy.histogram(v, bins=50, normed=True) # NumPy version (no plot)
pylab.plot(.5*(bins[1:]+bins[:-1]), n)
pylab.title('Numpy Version')
pylab.show()

学习资料

numpy v1.11手册

numpy参考手册

numpy快速入门教程

numpy通用索引:所有函数、类、术语

numpy全目录

kesci: Numpy快速上手指南 —- 进阶篇

Reference

【Python数据分析】Numpy的详细教程

Scipy: Quickstart tutorial

Thanks!