Science

numpy--快速科学计算

NumPy提供了两种基本的对象：

  1. ndarray（N-dimensional array object）
  2. ufunc（universal function object）。
    • ndarray(下文统一称之为数组)是存储单一数据类型的多维数组
    • ufunc则是能够对数组进行处理的函数。

ndarray

组成：
  1. 实际的数据
  2. 描述这些数据的元数据
存取元素的高级方法
  1. 使用整数序列
    当使用整数序列对数组元素进行存取时，将使用整数序列中的每个元素作为下标，整数序列可以是列 表或者数组。使用整数序列作为下标获得的数组不和原始数组共享数据空间。
  2. 使用布尔数组
    当使用布尔数组b作为下标存取数组x中的元素时，将收集数组x中所有在数组b中对应下标为True的 元素。使用布尔数组作为下标获得的数组不和原始数组共享数据空间，注意这种方式只对应于布尔数 组，不能使用布尔列表
1
2
3
4
5
6
7
8
9
10
11
12
13
>>> x = np.arange(10,1,-1) 
>>> x 
array([10, 9, 8, 7, 6, 5, 4, 3, 2]) 
>>> x[[3, 3, 1, 8]] # 获取x中的下标为3, 3, 1, 8的4个元素，组成一个新的数组 array([7, 7, 9, 2]) 
>>> b = x[np.array([3,3,-3,8])] #下标可以是负数 
>>> b[2] = 100 
>>> b 
array([7, 7, 100, 2]) 
>>> x # 由于b和x不共享数据空间，因此x中的值并没有改变 
array([10, 9, 8, 7, 6, 5, 4, 3, 2]) 
>>> x[[3,5,1]] = -1, -2, -3 # 整数序列下标也可以用来修改元素的值 
>>> x 
array([10, -3, 8, -1, 6, -2, 4, 3, 2])
1
2
3
4
5
6
7
8
9
10
11
12
13
>>> x = np.arange(5,0,-1) 
>>> x array([5, 4, 3, 2, 1]) 
>>> x[np.array([True, False, True, False, False])] 
>>> # 布尔数组中下标为0，2的元素为True，因此获取x中下标为0,2的元素 
array([5, 3]) 
>>> x[[True, False, True, False, False]] 
>>> # 如果是布尔列表，则把True当作1, False当作0，按照整数序列方式获取x中的元素 
array([4, 5, 4, 5, 5]) 
>>> x[np.array([True, False, True, True])] 
>>> # 布尔数组的长度不够时，不够的部分都当作False array([5, 3, 2]) 
>>> x[np.array([True, False, True, True])] = -1, -2, -3 
>>> # 布尔数组下标也可以用来修改元素 
>>> x array([-1, 4, -2, -3, 1])
1
2
3
4
5
6
7
8
>>> x = np.random.rand(10) # 产生一个长度为10，元素值为0-1的随机数的数组 
>>> x array([ 0.72223939, 0.921226 , 0.7770805 , 0.2055047 , 0.17567449, 0.95799412, 0.12015178, 0.7627083 , 0.43260184, 0.91379859]) 
>>> x>0.5 
>>> # 数组x中的每个元素和0.5进行大小比较，得到一个布尔数组，True表示x中对应的值大于0.5 
array([ True, True, True, False, False, True, False, True, False, True], dtype=bool) 
>>> x[x>0.5] 
>>> # 使用x>0.5返回的布尔数组收集x中的元素，因此得到的结果是x中所有大于0.5的元素的数组 
array([ 0.72223939, 0.921226 , 0.7770805 , 0.95799412, 0.7627083 , 0.91379859])

多维数组

数组a实际上是一个加法表，纵轴的值为0, 10, 20, 30, 40, 50；横轴的值为0, 1, 2, 3, 4, 5。纵轴的每个元素都和横轴的每个元素求和，就得到图中 所示的数组a。你可以用下面的语句创建它

1
2
>>> np.arange(0, 60, 10).reshape(-1, 1) + np.arange(0, 6) 
array([[ 0, 1, 2, 3, 4, 5], [10, 11, 12, 13, 14, 15], [20, 21, 22, 23, 24, 25], [30, 31, 32, 33, 34, 35], [40, 41, 42, 43, 44, 45], [50, 51, 52, 53, 54, 55]])

  • a[(0,1,2,3,4),(1,2,3,4,5)] : 用于存取数组的下标和仍然是一个有两个元素的组元，组元中的每个 元素都是整数序列，分别对应数组的第0轴和第1轴。从两个序列的对应位置取出两个整数组成下标： a[0,1], a[1,2], …, a[4,5]。
  • a[3:, [0, 2, 5]] : 下标中的第0轴是一个范围，它选取第3行之后的所有行；第1轴是整数序列，它 选取第0, 2, 5三列。
  • a[mask, 2] : 下标的第0轴是一个布尔数组，它选取第0，2，5行；第1轴是一个整数，选取第2 列

结构数组

假设我们需要定义一个结构数组，它的每个元素都有name, age和weight字段。在NumPy中可以如下 定义

1
2
3
import numpy as np 
persontype = np.dtype({ 'names':['name', 'age', 'weight'], 'formats':['S32','i', 'f']}) 
a = np.array([("Zhang",32,75.5),("Wang",24,65.2)], dtype=persontype)

我们先创建一个dtype对象persontype，通过其字典参数描述结构类型的各个字段。字典有两个关 键字：names，formats。每个关键字对应的值都是一个列表。names定义结构中的每个字段名，而 formats则定义每个字段的类型：

  • S32 : 32个字节的字符串类型，由于结构中的每个元素的大小必须固定，因此需要指定字符串的 长度
  • i : 32bit的整数类型，相当于np.int32
  • f : 32bit的单精度浮点数类型，相当于np.float32

然后我们调用array函数创建数组，通过关键字参数 dtype=persontype， 指定所创建的数组的元素类 型为结构persontype。运行上面程序之后，我们可以在IPython中执行如下的语句查看数组a的元素类型

1
2
>>> a.dtype 
dtype([('name', '|S32'), ('age', '<i4'), ('weight', '<f4')])

这里我们看到了另外一种描述结构类型的方法： 一个包含多个组元的列表，其中形如 (字段名, 类型描 述) 的组元描述了结构中的每个字段。类型描述前面为我们添加了 |, < 等字符，这些字符用来描述字 段值的字节顺序：

  • | : 忽视字节顺序
  • < : 低位字节在前
  • < : 低位字节在前

结构数组的存取方式和一般数组相同，通过下标能够取得其中的元素，注意元素的值看上去像是组 元，实际上它是一个结构：

1
2
3
>>> a[0] ('Zhang', 32, 75.5) 
>>> a[0].dtype 
dtype([('name', '|S32'), ('age', '<i4'), ('weight', '<f4')])

a[0]是一个结构元素，它和数组a共享内存数据，因此可以通过修改它的字段，改变原始数组中的对应字段：

1
2
>>> c = a[1] >>> c["name"] = "Li" 
>>> a[1]["name"] "Li"

结构像字典一样可以通过字符串下标获取其对应的字段值：

1
2
>>> a[0]["name"] 
'Zhang'

我们不但可以获得结构元素的某个字段，还可以直接获得结构数组的字段，它返回的是原始数组的视 图，因此可以通过修改b[0]改变a[0][‘’age’’]：

1
2
3
4
5
>>> b=a[:]["age"] # 或者a["age"] 
>>> b array([32, 24]) 
>>> b[0] = 40 
>>> a[0]["age"] 
40

通过调用a.tostring或者a.tofile方法，可以直接输出数组a的二进制形式

1
>>> a.tofile("test.bin")

利用下面的C语言程序可以将test.bin文件中的数据读取出来

  • 内存对齐 C语言的结构体为了内存寻址方便，会自动的添加一些填充用的字节，这叫做内存对齐。例如 如果把下面的name[32]改为name[30]的话，由于内存对齐问题，在name和age中间会填补两 个字节，最终的结构体大小不会改变。因此如果numpy中的所配置的内存大小不符合C语言的 对齐规范的话，将会出现数据错位。为了解决这个问题，在创建dtype对象时，可以传递参数 align=True，这样numpy的结构数组的内存对齐和C语言的结构体就一致了。
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
#include <stdio.h>
struct person 
{ char name[32]; 
  int age;
 float weight; };
struct person p[2];
void main () {
  FILE *fp;
  int i;
  fp=fopen("test.bin","rb");
  fread(p, sizeof(struct person), 2, fp);
  fclose(fp);
  for(i=0;i<2;i++)
    printf("%s %d %f\n", p[i].name, p[i].age, p[i].weight);
  getchar();
}

内存结构

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

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

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

strides中保存的是当每个轴的下标增加1时，数据存储区中的指针所增加的字节数。例如图中的 strides为12,4，即第0轴的下标增加1时，数据的地址增加12个字节：即a[1,0]的地址比a[0,0]的地址要 高12个字节，正好是3个单精度浮点数的总字节数；第1轴下标增加1时，数据的地址增加4个字节，正 好是单精度浮点数的字节数。
如果strides中的数值正好和对应轴所占据的字节数相同的话，那么数据在内存中是连续存储的。然而 数据并不一定都是连续储存的，前面介绍过通过下标范围得到新的数组是原始数组的视图，即它和原 始视图共享数据存储区域：

1
2
3
4
5
>>> b = a[::2,::2] 
>>> b
array([[ 0., 2.], [ 6., 8.]], dtype=float32) 
>>> b.strides 
(24, 8)

由于数组b和数组a共享数据存储区，而b中的第0轴和第1轴都是数组a中隔一个元素取一个，因此数组 b的strides变成了24,8，正好都是数组a的两倍。 对照前面的图很容易看出数据0和2的地址相差8个字 节，而0和6的地址相差24个字节。

元素在数据存储区中的排列格式有两种：C语言格式和Fortan语言格式。在C语言中，多维数组的第 0轴是最上位的，即第0轴的下标增加1时，元素的地址增加的字节数最多；而Fortan语言的多维数 组的第0轴是最下位的，即第0轴的下标增加1时，地址只增加一个元素的字节数。在NumPy中，元 素在内存中的排列缺省是以C语言格式存储的，如果你希望改为Fortan格式的话，只需要给数组传递 order=’’F’’参数：

1
2
3
>>> c = np.array([[0,1,2],[3,4,5],[6,7,8]], dtype=np.float32, order="F") 
>>> c.strides 
(4, 12)

ufunc运算

ufunc是universal function的缩写，它是一种能对数组的每个元素进行操作的函数。NumPy内置的许 多ufunc函数都是在C语言级别实现的，因此它们的计算速度非常快。

1
2
3
4
5
6
>>> x = np.linspace(0, 2*np.pi, 10) # 对数组x中的每个元素进行正弦计算，返回一个同样大小的新数组 
>>> y = np.sin(x) 
>>> y array([ 0.00000000e+00, 6.42787610e-01, 9.84807753e-01,
             8.66025404e-01, 3.42020143e-01, -3.42020143e-01, 
             -8.66025404e-01, -9.84807753e-01, -6.42787610e-01,
             -2.44921271e-16])

先用linspace产生一个从0到2*PI的等距离的10个数，然后将其传递给sin函数，由于np.sin是一个 ufunc函数，因此它对x中的每个元素求正弦值，然后将结果返回，并且赋值给y。计算之后x中的值并 没有改变，而是新创建了一个数组保存结果。如果我们希望将sin函数所计算的结果直接覆盖到数组x上 去的话，可以将要被覆盖的数组作为第二个参数传递给ufunc函数

1
2
3
4
5
6
7
8
>>> t = np.sin(x,x) 
>>> x 
array([ 0.00000000e+00, 6.42787610e-01, 9.84807753e-01,
        8.66025404e-01, 3.42020143e-01, -3.42020143e-01, 
       -8.66025404e-01, -9.84807753e-01, -6.42787610e-01, 
       -2.44921271e-16]) 
>>> id(t) == id(x) 
True

sin函数的第二个参数也是x，那么它所做的事情就是对x中的每给值求正弦值，并且把结果保存到x中的 对应的位置中。此时函数的返回值仍然是整个计算的结果，只不过它就是x，因此两个变量的id是相同 的(变量t和变量x指向同一块内存区域)。

我用下面这个小程序，比较了一下numpy.math和Python标准库的math.sin的计算速度

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
import time 
import math 
import numpy as np
x = [i * 0.001 for i in xrange(1000000)] 
start = time.clock() 
for i, t in enumerate(x): 
    x[i] = math.sin(t) 
print "math.sin:", time.clock() - start
    
x = [i * 0.001 for i in xrange(1000000)] 
x = np.array(x) 
start = time.clock() 
np.sin(x,x) 
print "numpy.sin:", time.clock() - start
# 输出 
# math.sin: 1.15426932753 
# numpy.sin: 0.0882399858083

在我的电脑上计算100万次正弦值，numpy.sin比math.sin快10倍多。这得利于numpy.sin在C语言级 别的循环计算。numpy.sin同样也支持对单个数值求正弦，例如：numpy.sin(0.5)。不过值得注意的 是，对单个数的计算math.sin则比numpy.sin快得多了，让我们看下面这个测试程序

1
2
3
4
5
6
7
x = [i * 0.001 for i in xrange(1000000)] 
start = time.clock() 
for i, t in enumerate(x): 
    x[i] = np.sin(t) 
print "numpy.sin loop:", time.clock() - start
# 输出 
# numpy.sin loop: 5.72166965355

请注意numpy.sin的计算速度只有math.sin的1/5。这是因为numpy.sin为了同时支持数组和单个值的 计算，其C语言的内部实现要比math.sin复杂很多，如果我们同样在Python级别进行循环的话，就会 看出其中的差别了。此外，numpy.sin返回的数的类型和math.sin返回的类型有所不同，math.sin返 回的是Python的标准float类型，而numpy.sin则返回一个numpy.float64类型：

1
2
3
4
>>> type(math.sin(0.5)) 
<type 'float'> 
>>> type(np.sin(0.5)) 
<type 'numpy.float64'>

NumPy中有众多的ufunc函数为我们提供各式各样的计算。除了sin这种单输入函数之外，还有许多多 个输入的函数，add函数就是一个最常用的例子。

1
2
3
4
5
6
7
8
9
10
11
>>> a = np.arange(0,4) 
>>> a array([0, 1, 2, 3]) 
>>> b = np.arange(1,5) 
>>> b 
array([1, 2, 3, 4]) 
>>> np.add(a,b) 
array([1, 3, 5, 7]) 
>>> np.add(a,b,a) 
array([1, 3, 5, 7]) 
>>> a 
array([1, 3, 5, 7])

add函数返回一个新的数组，此数组的每个元素都为两个参数数组的对应元素之和。它接受第3个参数 指定计算结果所要写入的数组，如果指定的话，add函数就不再产生新的数组。
由于Python的操作符重载功能，计算两个数组相加可以简单地写为a+b，而np.add(a,b,a)则可以用a +=b来表示。

广播

当我们使用ufunc函数对两个数组进行计算时，ufunc函数会对这两个数组的对应元素进行计算，因 此它要求这两个数组有相同的大小(shape相同)。如果两个数组的shape不同的话，会进行如下的广播 (broadcasting)处理:

  • 让所有输入数组都向其中shape最长的数组看齐，shape中不足的部分都通过在前面加1补齐
  • 输出数组的shape是输入数组shape的各个轴上的最大值
  • 如果输入数组的某个轴和输出数组的对应轴的长度相同或者其长度为1时，这个数组能够用来计 算，否则出错
  • 当输入数组的某个轴的长度为1时，沿着此轴运算时都用此轴上的第一组值

Example:

先创建一个二维数组a，其shape为(6,1):

1
2
3
4
>>> a = np.arange(0, 60, 10).reshape(-1, 1) 
>>> a 
array([[ 0], [10], [20], [30], [40], [50]]) 
>>> a.shape (6, 1)

再创建一维数组b，其shape为(5,)：

1
2
3
4
5
>>> b = np.arange(0, 5) 
>>> b 
array([0, 1, 2, 3, 4]) 
>>> b.shape 
(5,)

计算a和b的和，得到一个加法表，它相当于计算a,b中所有元素组的和，得到一个shape为(6,5)的数 组：

1
2
3
4
5
6
7
8
9
10
>>> c = a + b 
>>> c 
array([[ 0, 1, 2, 3, 4], 
       [10, 11, 12, 13, 14], 
       [20, 21, 22, 23, 24], 
       [30, 31, 32, 33, 34],
       [40, 41, 42, 43, 44], 
       [50, 51, 52, 53, 54]]) 
>>> c.shape 
(6, 5)

由于a和b的shape长度(也就是ndim属性)不同，根据规则1，需要让b的shape向a对齐，于是将b的 shape前面加1，补齐为(1,5)。相当于做了如下计算：

1
2
3
>>> b.shape=1,5 
>>> b 
array([[0, 1, 2, 3, 4]])

这样加法运算的两个输入数组的shape分别为(6,1)和(1,5)，根据规则2，输出数组的各个轴的长度为输 入数组各个轴上的长度的最大值，可知输出数组的shape为(6,5)。
由于b的第0轴上的长度为1，而a的第0轴上的长度为6，因此为了让它们在第0轴上能够相加，需要将b 在第0轴上的长度扩展为6，这相当于：

1
2
3
4
5
6
7
8
>>> b = b.repeat(6,axis=0) 
>>> b 
array([[0, 1, 2, 3, 4], 
       [0, 1, 2, 3, 4], 
       [0, 1, 2, 3, 4], 
       [0, 1, 2, 3, 4], 
       [0, 1, 2, 3, 4], 
       [0, 1, 2, 3, 4]])

由于a的第1轴的长度为1，而b的第一轴长度为5，因此为了让它们在第1轴上能够相加，需要将a在第1 轴上的长度扩展为5，这相当于：

1
2
3
4
5
6
7
8
>>> a = a.repeat(5, axis=1) 
>>> a 
array([[ 0, 0, 0, 0, 0], 
       [10, 10, 10, 10, 10], 
       [20, 20, 20, 20, 20], 
       [30, 30, 30, 30, 30], 
       [40, 40, 40, 40, 40], 
       [50, 50, 50, 50, 50]])

经过上述处理之后，a和b就可以按对应元素进行相加运算了。
当然，numpy在执行a+b运算时，其内部并不会真正将长度为1的轴用repeat函数进行扩展，如果这 样做的话就太浪费空间了。
由于这种广播计算很常用，因此numpy提供了一个快速产生如上面a,b数组的方法： ogrid对象：

1
2
3
4
5
6
7
8
9
>>> x,y = np.ogrid[0:5,0:5] 
>>> x 
array([[0], 
       [1], 
       [2], 
       [3], 
       [4]]) 
>>> y 
array([[0, 1, 2, 3, 4]])

ogrid是一个很有趣的对象，它像一个多维数组一样，用切片组元作为下标进行存取，返回的是一组可 以用来广播计算的数组。其切片下标有两种形式：

  • 开始值:结束值:步长，和np.arange(开始值, 结束值, 步长)类似
  • 开始值:结束值:长度j，当第三个参数为虚数时，它表示返回的数组的长度，和np.linspace(开始 值, 结束值, 长度)类似：
    1
    2
    3
    4
    5
    6
    7
    8
       >>> x, y = np.ogrid[0:1:4j, 0:1:3j] 
       >>> x 
       array([[ 0. ], 
    [ 0.33333333], 
    [ 0.66666667], 
    [ 1. ]]) 
    >>> y
    array([[ 0. , 0.5, 1. ]])
1
2
3
4
5
6
7
8
9
10
11
12
* ogrid为什么不是函数 根据Python的语法，只有在中括号中才能使用用冒号隔开的切片语法，如果ogrid是函数的话， 那么这些切片必须使用slice函数创建，这显然会增加代码的长度。
利用ogrid的返回值，我能很容易计算x, y网格面上各点的值，或者x, y, z网格体上各点的值。下面是绘 制三维曲面 x * exp(x**2 - y**2) 的程序：
​```python
import numpy as np 
from enthought.mayavi import mlab
x, y = np.ogrid[-2:2:20j, -2:2:20j] 
z = x * np.exp( - x**2 - y**2)
pl = mlab.surf(x, y, z, warp_scale="auto")
mlab.axes(xlabel='x', ylabel='y', zlabel='z') 
mlab.outline(pl)

此程序使用mayavi的mlab库快速绘制如图所示的3D曲面

矩阵运算

NumPy和Matlab不一样，对于多维数组的运算，缺省情况下并不使用矩阵运算，如果你希望对数组进 行矩阵运算的话，可以调用相应的函数。

matrix对象

numpy库提供了matrix类，使用matrix类创建的是矩阵对象，它们的加减乘除运算缺省采用矩阵 方式计算，因此用法和matlab十分类似。但是由于NumPy中同时存在ndarray和matrix对象，因 此用户很容易将两者弄混。这有违Python的“显式优于隐式”的原则，因此并不推荐在较复杂的 程序中使用matrix。下面是使用matrix的一个例子：

1
2
3
4
5
>>> a = np.matrix([[1,2,3],[5,5,6],[7,9,9]]) 
>>> a*a**-1 
matrix([[ 1.00000000e+00, 1.66533454e-16, -8.32667268e-17], 
        [ -2.77555756e-16, 1.00000000e+00, -2.77555756e-17], 
        [ 1.66533454e-16, 5.55111512e-17, 1.00000000e+00]])

因为a是用matrix创建的矩阵对象，因此乘法和幂运算符都变成了矩阵运算，于是上面计算的是矩 阵a和其逆矩阵的乘积，结果是一个单位矩阵。

矩阵的乘积可以使用dot函数进行计算。对于二维数组，它计算的是矩阵乘积，对于一维数组，它计算 的是其点积。当需要将一维数组当作列矢量或者行矢量进行矩阵运算时，推荐先使用reshape函数将一维数组转换为二维数组：

1
2
3
4
5
>>> a = array([1, 2, 3]) 
>>> a.reshape((-1,1)) 
array([[1], [2], [3]]) 
>>> a.reshape((1,-1)) 
array([[1, 2, 3]])

文件存取

NumPy提供了多种文件操作函数方便我们存取数组内容。文件存取的格式分为两类：二进制和文本。 而二进制格式的文件又分为NumPy专用的格式化二进制类型和无格式类型。
使用数组的方法函数tofile可以方便地将数组中数据以二进制的格式写进文件。tofile输出的数据没有格 式，因此用numpy.fromfile读回来的时候需要自己格式化数据

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
>>> a = np.arange(0,12) 
>>> a.shape = 3,4 
>>> a 
array([[ 0, 1, 2, 3], 
       [ 4, 5, 6, 7], 
       [ 8, 9, 10, 11]]) 
>>> a.tofile("a.bin") 
>>> b = np.fromfile("a.bin", dtype=np.float) # 按照float类型读入数据 
>>> b # 读入的数据是错误的 
array([ 2.12199579e-314, 6.36598737e-314, 1.06099790e-313, 
       1.48539705e-313, 1.90979621e-313, 2.33419537e-313]) 
>>> a.dtype # 查看a的dtype dtype('int32') 
>>> b = np.fromfile("a.bin", dtype=np.int32) # 按照int32类型读入数据
>>> b # 数据是一维的 
array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11]) 
>>> b.shape = 3, 4 # 按照a的shape修改b的shape 
>>> b # 这次终于正确了 
array([[ 0, 1, 2, 3], [ 4, 5, 6, 7], [ 8, 9, 10, 11]])

从上面的例子可以看出，需要在读入的时候设置正确的dtype和shape才能保证数据一致。并且tofile 函数不管数组的排列顺序是C语言格式的还是Fortran语言格式的，统一使用C语言格式输出。
此外如果fromfile和tofile函数调用时指定了sep关键字参数的话，数组将以文本格式输入输出。
numpy.load和numpy.save函数以NumPy专用的二进制类型保存数据，这两个函数会自动处理元素类 型和shape等信息，使用它们读写数组就方便多了，但是numpy.save输出的文件很难和其它语言编写 的程序读入：

1
2
3
4
5
6
>>> np.save("a.npy", a) 
>>> c = np.load( "a.npy" ) 
>>> c 
array([[ 0, 1, 2, 3], 
       [ 4, 5, 6, 7], 
       [ 8, 9, 10, 11]])

如果你想将多个数组保存到一个文件中的话，可以使用numpy.savez函数。savez函数的第一个参数是 文件名，其后的参数都是需要保存的数组，也可以使用关键字参数为数组起一个名字，非关键字参数 传递的数组会自动起名为arr_0, arr_1, …。savez函数输出的是一个压缩文件(扩展名为npz)，其中每个 文件都是一个save函数保存的npy文件，文件名对应于数组名。load函数自动识别npz文件，并且返回 一个类似于字典的对象，可以通过数组名作为关键字获取数组的内容：

1
2
3
4
5
6
7
8
9
10
11
12
13
>>> a = np.array([[1,2,3],[4,5,6]]) 
>>> b = np.arange(0, 1.0, 0.1) 
>>> c = np.sin(b) 
>>> np.savez("result.npz", a, b, sin_array = c) 
>>> r = np.load("result.npz") 
>>> r["arr_0"] # 数组a 
array([[1, 2, 3], 
       [4, 5, 6]]) 
>>> r["arr_1"] # 数组b 
array([ 0. , 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9]) 
>>> r["sin_array"] # 数组c 
array([ 0. , 0.09983342, 0.19866933, 0.29552021, 0.38941834, 
       0.47942554, 0.56464247, 0.64421769, 0.71735609, 0.78332691])

如果你用解压软件打开result.npz文件的话，会发现其中有三个文件：arr_0.npy， arr_1.npy， sin_array.npy，其中分别保存着数组a, b, c的内容

使用numpy.savetxt和numpy.loadtxt可以读写1维和2维的数组：

1
2
3
4
5
6
7
8
9
10
11
12
13
>>> a = np.arange(0,12,0.5).reshape(4,-1) 
>>> np.savetxt("a.txt", a) # 缺省按照'%.18e'格式保存数据，以空格分隔 
>>> np.loadtxt("a.txt") 
array([[ 0. , 0.5, 1. , 1.5, 2. , 2.5], 
       [ 3. , 3.5, 4. , 4.5, 5. , 5.5], 
       [ 6. , 6.5, 7. , 7.5, 8. , 8.5], 
       [ 9. , 9.5, 10. , 10.5, 11. , 11.5]]) 
>>> np.savetxt("a.txt", a, fmt="%d", delimiter=",") #改为保存为整数，以逗号分隔 
>>> np.loadtxt("a.txt",delimiter=",") # 读入的时候也需要指定逗号分隔 
array([[ 0., 0., 1., 1., 2., 2.], 
       [ 3., 3., 4., 4., 5., 5.], 
       [ 6., 6., 7., 7., 8., 8.], 
       [ 9., 9., 10., 10., 11., 11.]])

文件名和文件对象 本节介绍所举的例子都是传递的文件名，也可以传递已经打开的文件对象，例如对于load和save 函数来说，如果使用文件对象的话，可以将多个数组储存到一个npy文件中：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
>>> a = np.arange(8) 
>>> b = np.add.accumulate(a) 
>>> c = a + b 
>>> f = file("result.npy", "wb") 
>>> np.save(f, a) # 顺序将a,b,c保存进文件对象f 
>>> np.save(f, b) 
>>> np.save(f, c) 
>>> f.close() 
>>> f = file("result.npy", "rb") 
>>> np.load(f) # 顺序从文件对象f中读取内容 
array([0, 1, 2, 3, 4, 5, 6, 7]) 
>>> np.load(f) 
array([ 0, 1, 3, 6, 10, 15, 21, 28]) 
>>> np.load(f) 
array([ 0, 2, 5, 9, 14, 20, 27, 35])

一些函数和属性

linspace

1
numpy.linspace(start, stop, num=50, endpoint=True, retstep=False, dtype=None)

Example

1
2
3
4
5
6
>>> np.linspace(2.0, 3.0, num=5)
    array([ 2.  ,  2.25,  2.5 ,  2.75,  3.  ])
>>> np.linspace(2.0, 3.0, num=5, endpoint=False)
    array([ 2. ,  2.2,  2.4,  2.6,  2.8])
>>> np.linspace(2.0, 3.0, num=5, retstep=True)
    (array([ 2.  ,  2.25,  2.5 ,  2.75,  3.  ]), 0.25)

meshgrid

1
numpy.meshgrid(*xi, **kwargs)

Exampl

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
>>> nx, ny = (3, 2)
>>> x = np.linspace(0, 1, nx)
>>> y = np.linspace(0, 1, ny)
>>> xv, yv = meshgrid(x, y)
>>> xv
array([[ 0. ,  0.5,  1. ],
       [ 0. ,  0.5,  1. ]])
>>> yv
array([[ 0.,  0.,  0.],
       [ 1.,  1.,  1.]])
>>> xv, yv = meshgrid(x, y, sparse=True)  # make sparse output arrays
>>> xv
array([[ 0. ,  0.5,  1. ]])
>>> yv
array([[ 0.],
       [ 1.]])

r_

将切片拼接到一起至第一轴

Example

1
2
3
4
>>> np.r_[np.array([1,2,3]), 0, 0, np.array([4,5,6])]
array([1, 2, 3, 0, 0, 4, 5, 6])
>>> np.r_[-1:1:6j, [0]*3, 5, 6]
array([-1. , -0.6, -0.2,  0.2,  0.6,  1. ,  0. ,  0. ,  0. ,  5. ,  6. ])

c_

将切片拼接到一起至第二轴

Example

1
2
>>> np.c_[np.array([[1,2,3]]), 0, 0, np.array([[4,5,6]])]
array([[1, 2, 3, 0, 0, 4, 5, 6]])
Compartir