前言

如果没有numpy,Python内部只能用listarray来表示矩阵。假如用list来表示[1,2,3],由于list的元素可以是任何对象,因此list中所保存的是对象的指针,这样就需要有3个指针和三个整数对象,比较浪费内存和CPU计算时间。Python的arraylist不同,它直接保存数值,和C语言的一维数组比较类似,但是不支持多维,表达形式很简陋,写科学计算的算法很难受。numpy弥补了这些不足,核心贡献就是提供了ndarray这个存储单一数据类型的多维数组结构,实现上采用预编译好的C语言代码,性能上的表现十分不错。

ndarray内存结构

使用下面的代码可以生成本文讲解的ndarray案例:

from numpy import np

a = np.array([[0,1,2],[3,4,5],[6,7,8]], dtype=np.float32)

我们来看一下ndarray如何在内存中储存的:关于数组的描述信息保存在一个数据结构中,这个结构引用两个对象,一块用于保存数据的存储区域和一个用于描述元素类型的dtype对象。

数据存储区域保存着数组中所有元素的二进制数据,dtype对象则知道如何将元素的二进制数据转换为可用的值。数组的维数、大小等信息都保存在ndarray数组对象的数据结构中。

strides中保存的是当每个轴的下标增加1时,数据存储区中的指针所增加的字节数。例如图中的strides为12,4,即第0轴的下标增加1时,数据的地址增加12个字节:即a[1,0]的地址比a[0,0]的地址要高12个字节,正好是3个单精度浮点数的总字节数;第1轴下标增加1时,数据的地址增加4个字节,正好是单精度浮点数的字节数。

深入dtype

Python内置的基本数据类型,每个类别只有一种,比如只有一种整数int,只有一种浮点数float。对于不需要关心计算机如何存储数据的应用,只提供一种选择是非常方便的。但是对于科学计算而言,这是不够的,为了性能和精度往往需要根据具体场景控制更多的细节。为此,NumPy内置了24种基本类型,基本上可以和C语言的数据类型对应上,其中部分类型对应为Python内置的类型。下表列举了常用NumPy基本类型。

类型 注释
bool_ 兼容Python内置的bool类型
bool8 8位布尔
int_ 兼容Python内置的int类型
int8 8位整数
int16 16位整数
int32 32位整数
int64 64位整数
uint8 无符号8位整数
uint16 无符号16位整数
uint32 无符号32位整数
uint64 无符号64位整数
float_ 兼容Python内置的float类型
float16 16位浮点数
float32 32位浮点数
float64 64位浮点数
str_ 兼容Python内置的str类型

24个scalar types并不是dtype,但是可以作为参数传递给np.dtype()构造函数产生一个dtype对象,如np.dtype(np.int32)。在NumPy中所有需要dtype作为参数的函数都可以使用scalar types代替,会自动转化为对应的dtype类型。

Structured Array

由于numpy只支持单一数据类型,对于常见的表格型数据,我们需要通过numpy提供的Structrued Array机制自定义dtype

定义结构化数组有四种方式:1) string, 2) tuple, 3) list, or 4) dictionary。本文推荐使用后两种方式。

import numpy as np

# list方式:a list of tuples. Each tuple has 2 or 3 elements specifying: 1) The name of the field (‘’ is permitted), 2) the type of the field, and 3) the shape (optional)
persontype = [('name', np.str_), ('age', np.int16), ('weight', np.float32)]
  
# dict方式:需要指定的键值有names和formats
persontype = np.dtype({
        'names': ['name', 'age', 'weight'], 
        'formats': [np.str_, np.int16, np.float32]
    })
a = np.array([("Zhang", 32, 75.5), ("Wang", 24, 65.2)], dtype=persontype)

创建ndarray

共有三类基本方法:一是从Python内置的array-like数据结构转化得到;二是利用numpy提供的创建函数直接生成;三是使用genfromtxt()方法生成。

# -*- coding: utf-8 -*-
import numpy as np

# 从python的list转换
x = np.array([[1,2.0],[0,0],(1+1j,3.)])

np.zeros((2, 3))
np.arange(2, 3, 0.1) # start, end, step
np.linspace(1., 4., 6) # start, end, num
np.indices((3, 3)) # 返回一个array,元素0是行下标,元素1是列下标;行下标为一个3*3二维array,对应3*3矩阵的行下标;列下标为一个3*3二维array,对应3*3矩阵的列下标

ndtype=[('a',int), ('b', float), ('c', int)]
names = ["A", "B", "C"]
np.genfromtxt("file_name.txt", 
  delimiter=",",
  names=names,
  dtype=ndtype, 
  autostrip=True,
  comments="#",
  skip_header=3, 
  skip_footer=5,
  usecols=(0, -1))

# ndarray to list
a = np.array([[1, 2], [3, 4]])
a.tolist()

索引

ndarray索引的基本使用方法示例如下;

x = np.arange(10)
x.shape = (2,5) 
x[0][2]
x[:,0:5:2]

关于ndarray的索引方式,有以下几个重点需要记住:

  • 虽然x[0,2] = x[0][2],但是前者效率比后者高,因为后者在应用第一个索引后需要先创建一个temporary array,然后再应用第二个索引,最后找到目标值。
  • 分片操作不会引发copy操作,而是创建原ndarray的view;他们所指向的内存是同一片区域,无论是修改原ndarray还是修改view,都会同时改变二者的值。
  • index array和boolean index返回的是copy,不是view。

关于上面列举的分片操作不会引发copy操作,我们来进一步探讨一下。先看一下numpy的例子:

import numpy as np
a = np.arange(10)
a.shape = (2,5)
print "input a:", a

b = a[0, ]
print "copy b:", b

b[1] = -1
print "modify b:", b
print "final a:", a 
## input a: [[0 1 2 3 4]
##  [5 6 7 8 9]]
## copy b: [0 1 2 3 4]
## modify b: [ 0 -1  2  3  4]
## final a: [[ 0 -1  2  3  4]
##  [ 5  6  7  8  9]]

再来看一下R的例子:

a <- matrix(0:9, nrow=2, byrow=TRUE)
a
##      [,1] [,2] [,3] [,4] [,5]
## [1,]    0    1    2    3    4
## [2,]    5    6    7    8    9
b <- a[1,]
b[2] <- -1
b
## [1]  0 -1  2  3  4
a
##      [,1] [,2] [,3] [,4] [,5]
## [1,]    0    1    2    3    4
## [2,]    5    6    7    8    9

可以看到numpy和R在矩阵的分片操作有不同的设计理念:在R里分片操作会引起数据的复制,在numpy里不会。事实上,R的设计理念很多时候可以用一句话来概括:copy on modify,一旦对数据有修改就会引起内存上的复制操作,这个操作要花不少时间,因此经常会听到人们抱怨R费内存且速度慢。所以,我们可以看到numpy在处理这件事情上明显要用心很多,根据场景设计了不同的策略,不是简单地采用R的一刀切方式。当然,这也带来了一些学习成本,需要对numpy足够熟悉才能避免踩坑。R社区里对copy on modify的哲学也有诟病并在努力改变,比如同是data.frame操作库的data.table和dplyr,data.table性能比dplyr高很多,部分原因也是data.table规避了copy on modify的方式。

再来看一下Structured Array的索引效率。

import numpy as np
import time

persontype = [('name', np.str_), ('age', np.int16), ('weight', np.float32)]
a = np.array([("Zhang", 32, 75.5), ("Wang", 24, 65.2)], dtype=persontype)

start = time.clock()
a[1]
end = time.clock()
print "a[1] cost time: %f" % (end-start)

start = time.clock()
a['name']
end = time.clock()
print "a['name'] cost time: %f" % (end-start)

start = time.clock()
a[1]['name']
end = time.clock()
print "a[1]['name'] cost time: %f" % (end-start)

start = time.clock()
a['name'][1]
end = time.clock()
print "a['name'][1] cost time: %f" % (end-start)
## a[1] cost time: 0.000006
## a['name'] cost time: 0.000020
## a[1]['name'] cost time: 0.000037
## a['name'][1] cost time: 0.000002

从上面的结果,我们发现,获取相同的数据有多种操作,不同的操作性能差别很大。我做了一个推测,纯粹是瞎猜:numpy在建立结构化数组时,将整个结构体连续存储在一起,即按行存储,因此a[1]的速度最快;但是为了保证提取列的效率,对a['name']建立了索引,因此a['name']的效率也很高;但是这个索引只对整个a起作用,如果输入只有a的一部分,仍然需要遍历整个a,去提取出对应的数据,因此a[1]['name']a['name'][1]的效率差很多。

参考文献