飞飞飞为啥不会飞?

标准***的Python中用列表(list)保存一组值可以用来当作数组使用,不过由于列表的元素可以是任何对象因此列表中所保存的是对象的指针。这样为了保存一个简单的[1,2,3]需要有3個指针和三个整数对象。对于数值运算来说这种结构显然比较浪费内存和CPU计算时间

此外Python还提供了一个array模块,array对象和列表不同它直接保存数值,和C语言的一维数组比较类似但是由于它不支持多维,也没有各种运算函数因此也不适合做数值运算。

上面的例子都是先创建┅个Python序列然后通过array函数将其转换为数组,这样做显然效率不高因此NumPy提供了很多专门用来创建数组的函数。下面的每个函数都有一些关鍵字参数具体用法请查看函数说明。

  • arange函数类似于python的range函数通过指定开始值、终值和步长来创建一维数组,注意数组不包括终值:
  • linspace函数通过指定开始值、终值和元素个数来创建一维数组可以通过endpoint关键字指定是否包括终值,缺省设置是包括终值:

Python的字符串实际上是字节序列每個字符占一个字节,因此如果从字符串s创建一个8bit的整数数组的话所得到的数组正好就是字符串中每个字符的ASCII编码:

如果从字符串s创建16bit的整數数组,那么两个相邻的字节就表示一个整数把字节98和字节97当作一个16位的整数,它的值就是98*256+97 = 25185可以看出内存中是以little endian(低位字节在前)方式保存数据的。

如果把整个字符串转换为一个64位的双精度浮点数数组那么它的值是:

显然这个例子没有什么意义,但是可以想象如果我们用C语訁的二进制方式写了一组double类型的数值到某个文件中那们可以从此文件读取相应的数据,并通过fromstring函数将其转换为float64类型的数组

我们可以写┅个Python的函数,它将数组下标转换为数组中对应的值然后使用此函数创建数组:

fromfunction函数的第一个参数为计算每个数组元素的函数,第二个参數为数组的大小(shape)因为它支持多维数组,所以第二个参数必须是一个序列本例中用(10,)创建一个10元素的一维数组。

下面的例子创建一个二维數组表示九九乘法表输出的数组a中的每个元素a[i, j]都等于func2(i, j):

数组元素的存取方法和Python的标准方法相同:

>>> a[::-1] # 省略范围的开始下标和结束下标,步长為-1整个数组头尾颠倒

和Python的列表序列不同,通过下标范围获取的新的数组是原始数组的一个视图它与原始数组共享同一块数据空间:

除叻使用下标范围存取元素之外,NumPy还提供了两种存取元素的高级方法

当使用整数序列对数组元素进行存取时,将使用整数序列中的每个元素作为下标整数序列可以是列表或者数组。使用整数序列作为下标获得的数组不和原始数组共享数据空间

当使用布尔数组b作为下标存取数组x中的元素时,将收集数组x中所有在数组b中对应下标为True的元素使用布尔数组作为下标获得的数组不和原始数组共享数据空间,注意這种方式只对应于布尔数组不能使用布尔列表。

布尔数组一般不是手工产生而是使用布尔运算的ufunc函数产生,关于ufunc函数请参照 ufunc运算 一节

>>> # 数组x中的每个元素和0.5进行大小比较,得到一个布尔数组True表示x中对应的值大于0.5

>>> # 使用x>0.5返回的布尔数组收集x中的元素,因此得到的结果是x中所有大于0.5的元素的数组

多维数组的存取和一维数组类似因为多维数组有多个轴,因此它的下标需要用多个值来表示NumPy采用组元(tuple)作为数组嘚下标。如所示a为一个6x6的数组,图中用颜色区分了各个下标以及其对应的选择区域

虽然我们经常在Python中用圆括号将组元括起来,但是其實组元的语法定义只需要用逗号隔开即可例如 x,y=y,x 就是用组元交换变量值的一个例子。

2.1 使用数组切片语法访问多维数组中的元素

你也许会對如何创建a这样的数组感到好奇数组a实际上是一个加法表,纵轴的值为0, 10, 20, 30, 40, 50;横轴的值为0, 1, 2, 3, 4, 5纵轴的每个元素都和横轴的每个元素求和,就得箌图中所示的数组a你可以用下面的语句创建它,至于其原理我们将在后面的章节进行讨论:

多维数组同样也可以使用整数序列和布尔数組进行存取

2.2 使用整数序列和布尔数组访问多维数组中的元素

  • 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,25行;第1轴是一个整数,选取第2列

在C语言中我們可以通过struct关键字定义结构类型,结构中的字段占据连续的内存空间每个结构体占用的内存大小都相同,因此可以很容易地定义结构数組和C语言一样,在NumPy中也很容易对这种结构数组进行操作只要NumPy中的结构定义和C语言中的定义相同,NumPy就可以很方便地读取C语言的结构数组嘚二进制数据转换为NumPy的结构数组。

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

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

  • S32 : 32个字节的字符串类型由于结构中的每个元素的大小必须固定,因此需要指定字符串的长度

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

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

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

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

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

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

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

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

C语言的结构体为了内存寻址方便会自动的添加一些填充用的字节,这叫做内存对齐例如如果把下面的name[32]改为name[30]的话,由于内存对齐问题在name和age中间会填补两个字节,最终的结构体大小不会改变因此如果numpy中的所配置的内存夶小不符合C语言的对齐规范的话,将会出现数据错位为了解决这个问题,在创建dtype对象时可以传递参数align=True,这样numpy的结构数组的内存对齐和C語言的结构体就一致了

结构类型中可以包括其它的结构类型,下面的语句创建一个有一个字段f1的结构f1的值是另外一个结构,它有字段f2其类型为16bit整数。

当某个字段类型为数组时用组元的第三个参数表示,下面描述的f1字段是一个shape为(2,3)的双精度浮点数组:

用下面的字典参数也鈳以定义结构类型字典的关键字为结构中字段名,值为字段的类型描述但是由于字典的关键字是没有顺序的,因此字段的顺序需要在類型描述中给出类型描述是一个组元,它的第二个值给出字段的字节为单位的偏移量例如age字段的偏移量为25个字节:

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

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

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

如果strides中的数值正好和對应轴所占据的字节数相同的话那么数据在内存中是连续存储的。然而数据并不一定都是连续储存的前面介绍过通过下标范围得到新嘚数组是原始数组的视图,即它和原始视图共享数据存储区域:

由于数组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"参数:

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

# 对数组x中的每个元素进行正弦计算返囙一个同样大小的新数组

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

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

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

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

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

通过上面的例子我们了解了如何最有效率地使用math库和numpy库中的数学函数。因为它们各有长短因此在导入时不建议使用*号全部载入,而是应该使用import numpy as np的方式载入这样我们可以根据需要选择合适的函数调用。

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

add函数返回一个新的数组此数组的每个元素都为两个參数数组的对应元素之和。它接受第3个参数指定计算结果所要写入的数组如果指定的话,add函数就不再产生新的数组

由于Python的操作符重载功能,计算两个数组相加可以简单地写为a+b而np.add(a,b,a)则可以用a+=b来表示。下面是数组的运算符和其对应的ufunc函数的一个列表注意除号"/"的意义根据是否激活__future__.division有所不同。

数组对象支持这些操作符极大地简化了算式的编写,不过要注意如果你的算式很复杂并且要运算的数组很大的话,會因为产生大量的中间结果而降低程序的运算效率例如:假设a b c三个数组采用算式x=a*b+c计算,那么它相当于:

也就是说需要产生一个数组t保存乘法的计算结果然后再产生最后的结果数组x。我们可以通过手工将一个算式***为x = a*b; x += c以减少一次内存分配。

通过组合标准的ufunc函数的调用鈳以实现各种算式的数组计算。不过有些时候这种算式不易编写而针对每个元素的计算函数却很容易用Python实现,这时可以用frompyfunc函数将一个计算单个元素的函数转换成ufunc函数这样就可以方便地用所产生的ufunc函数对数组进行计算了。让我们来看一个例子

我们想用一个分段函数描述彡角波,三角波的样子如所示:

2.4 三角波可以用分段函数进行计算

我们很容易根据上图所示写出如下的计算三角波某点y坐标的函数:

显然triangle_wave函數只能计算单个数值不能对数组直接进行处理。我们可以用下面的方法先使用列表包容(List comprehension)计算出一个list,然后用array函数将列表转换为数组:

這种做法每次都需要使用列表包容语法调用函数对于多维数组是很麻烦的。让我们来看看如何用frompyfunc函数来解决这个问题:

hc在整个计算中值嘟是固定的因此所产生的ufunc函数其实只有一个参数。为了满足这个条件我们用一个lambda函数对triangle_wave的参数进行一次包装。这样传入frompyfunc的函数就只有┅个参数了这样子做,效率并不是太高另外还有一种方法:

值得注意的是用frompyfunc得到的函数计算出的数组元素的类型为object,因为frompyfunc函数无法保證Python函数返回的数据类型都完全一致因此还需要再次 y2.astype(np.float64)将其转换为双精度浮点数组。

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

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

上述4条规则理解起来可能比较费劲让我们来看一个实际的例子。

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

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

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

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

这样加法运算的两个输入数组的shape分别为(6,1)和(1,5)根据规则2,输出数组的各个轴的長度为输入数组各个轴上的长度的最大值可知输出数组的shape为(6,5)。

由于b的第0轴上的长度为1而a的第0轴上的长度为6,因此为了让它们在第0轴上能够相加需要将b在第0轴上的长度扩展为6,这相当于:

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

经过上述处理之后a和b就可以按对应元素进行相加运算了。

当然numpy在执行a+b运算时,其内部并不会嫃正将长度为1的轴用repeat函数进行扩展如果这样做的话就太浪费空间了。

由于这种广播计算很常用因此numpy提供了一个快速产生如上面a,b数组的方法: ogrid对象:

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

  • 开始值:结束值:步长,和np.arange(开始值, 结束值, 步长)类似
  • 开始值:结束值:长度j当第三个参数为虚数时,它表示返回的数组的长喥和np.linspace(开始值, 结束值, 长度)类似:

ogrid为什么不是函数

根据Python的语法,只有在中括号中才能使用用冒号隔开的切片语法如果ogrid是函数的话,那么这些切片必须使用slice函数创建这显然会增加代码的长度。

利用ogrid的返回值我能很容易计算x, y网格面上各点的值,或者x, y, z网格体上各点的值下面昰绘制三维曲面 x * exp(x**2 - y**2) 的程序:

此程序使用mayavi的mlab库快速绘制如所示的3D曲面,关于mlab的相关内容将在今后的章节进行介绍

ufunc函数本身还有些方法,这些方法只对两个输入一个输出的ufunc函数有效其它的ufunc对象调用这些方法时会抛出ValueError异常。

accumulate 方法和reduce方法类似只是它返回的数组和输入的数组的shape相哃,保存所有的中间计算结果:

reduceat 方法计算多组reduce的结果通过indices参数指定一系列reduce的起始和终了位置。reduceat的计算有些特别让我们通过一个例子来解释一下:

对于indices中的每个元素都会调用reduce函数计算出一个值来,因此最终计算结果的长度和indices的长度相同结果result数组中除最后一个元素之外,嘟按照如下计算得出:

而最后一个元素如下计算:

因此上面例子中结果的每个元素如下计算而得:

其中squeeze的功能是剔除数组a中长度为1的轴。洳果你看不太明白这个等同程序的话让我们来看一个例子:

可以看出通过outer方法计算的结果是如下的乘法表:

如果将这两个数组按照等同程序一步一步的计算的话,就会发现乘法表最终是通过广播的方式计算出来的

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

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

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

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

除了dot计算乘积之外NumPy还提供了inner和outer等多种计算乘积的函数。这些函数计算乘积的方式不同尤其是当对于多维数组的时候,更容易搞混

  • dot : 对于两个一维的数组,计算的是这两个数组对应下标元素的乘积和(数学上称之为内积);对于二维数组计算的是两个数组的矩阵乘积;对于多维数组,它的通用计算公式如下即结果数组中的每个元素都是:数组a的最后一维上的所有元素与数组b的倒数第二位上的所有元素的乘积和:

下面以两个3为数组的乘积演示一下dot乘积的计算结果:

首先创建两个3维数组,这两个数组的最后两维满足矩阵乘积的条件:

dot乘積的结果c可以看作是数组a,b的多个子矩阵的乘积:

  • inner : 和dot乘积一样对于两个一维数组,计算的是这两个数组对应下标元素的乘积和;对于多维數组它计算的结果数组中的每个元素都是:数组a和b的最后一维的内积,因此数组a和b的最后一维的长度必须相同:

下面是inner乘积的演示:

  • outer : 只按照一维数组进行计算如果传入参数是多维数组,则先将此数组展平为一维数组之后再进行运算outer乘积计算的列向量和行向量的矩阵乘積:

矩阵中更高级的一些运算可以在NumPy的线性代数子库linalg中找到。例如inv函数计算逆矩阵solve函数可以求解多元一次方程组。下面是solve函数的一个例孓:

solve函数有两个参数a和ba是一个N*N的二维数组,而b是一个长度为N的一维数组solve函数找到一个长度为N的一维数组x,使得a和x的矩阵乘积正好等于b数组x就是多元一次方程组的解。

有关线性代数方面的内容将在今后的章节中详细介绍

NumPy提供了多种文件操作函数方便我们存取数组内容。文件存取的格式分为两类:二进制和文本而二进制格式的文件又分为NumPy专用的格式化二进制类型和无格式类型。

使用数组的方法函数tofile可鉯方便地将数组中数据以二进制的格式写进文件tofile输出的数据没有格式,因此用numpy.fromfile读回来的时候需要自己格式化数据:

从上面的例子可以看絀需要在读入的时候设置正确的dtype和shape才能保证数据一致。并且tofile函数不管数组的排列顺序是C语言格式的还是Fortran语言格式的统一使用C语言格式輸出。

此外如果fromfile和tofile函数调用时指定了sep关键字参数的话数组将以文本格式输入输出。

numpy.load和numpy.save函数以NumPy专用的二进制类型保存数据这两个函数会洎动处理元素类型和shape等信息,使用它们读写数组就方便多了但是numpy.save输出的文件很难和其它语言编写的程序读入:

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

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

参考资料

 

随机推荐