
標(biāo)準(zhǔn)安裝的Python中用列表(list)保存一組值,可以用來當(dāng)作數(shù)組使用,不過由于列表的元素可以是任何對象,因此列表中所保存的是對象的指針。這樣為了保存一個(gè)簡單的[1,2,3],需要有3個(gè)指針和三個(gè)整數(shù)對象。對于數(shù)值運(yùn)算來說這種結(jié)構(gòu)顯然比較浪費(fèi)內(nèi)存和CPU計(jì)算時(shí)間。
此外Python還提供了一個(gè)array模塊,array對象和列表不同,它直接保存數(shù)值,和C語言的一維數(shù)組比較類似。但是由于它不支持多維,也沒有各種運(yùn)算函數(shù),因此也不適合做數(shù)值運(yùn)算。
NumPy的誕生彌補(bǔ)了這些不足,NumPy提供了兩種基本的對象:ndarray(N-dimensional array object)和 ufunc(universal function object)。ndarray(下文統(tǒng)一稱之為數(shù)組)是存儲單一數(shù)據(jù)類型的多維數(shù)組,而ufunc則是能夠?qū)?shù)組進(jìn)行處理的函數(shù)。

函數(shù)庫的導(dǎo)入
本書的示例程序假設(shè)用以下推薦的方式導(dǎo)入NumPy函數(shù)庫:
import numpy as np2.1.1 創(chuàng)建
首先需要?jiǎng)?chuàng)建數(shù)組才能對其進(jìn)行其它操作。
我們可以通過給array函數(shù)傳遞Python的序列對象創(chuàng)建數(shù)組,如果傳遞的是多層嵌套的序列,將創(chuàng)建多維數(shù)組(下例中的變量c):
>>> a = np.array([1, 2, 3, 4])>>> b = np.array((5, 6, 7, 8))>>> c = np.array([[1, 2, 3, 4],[4, 5, 6, 7], [7, 8, 9, 10]])>>> barray([5, 6, 7, 8])>>> carray([[1, 2, 3, 4], [4, 5, 6, 7], [7, 8, 9, 10]])>>> c.dtypedtype('int32')數(shù)組的大小可以通過其shape屬性獲得:
>>> a.shape(4,)>>> c.shape(3, 4)數(shù)組a的shape只有一個(gè)元素,因此它是一維數(shù)組。而數(shù)組c的shape有兩個(gè)元素,因此它是二維數(shù)組,其中第0軸的長度為3,第1軸的長度為4。還可以通過修改數(shù)組的shape屬性,在保持?jǐn)?shù)組元素個(gè)數(shù)不變的情況下,改變數(shù)組每個(gè)軸的長度。下面的例子將數(shù)組c的shape改為(4,3),注意從(3,4)改為(4,3)并不是對數(shù)組進(jìn)行轉(zhuǎn)置,而只是改變每個(gè)軸的大小,數(shù)組元素在內(nèi)存中的位置并沒有改變:
>>> c.shape = 4,3>>> carray([[ 1, 2, 3], [ 4, 4, 5], [ 6, 7, 7], [ 8, 9, 10]])當(dāng)某個(gè)軸的元素為-1時(shí),將根據(jù)數(shù)組元素的個(gè)數(shù)自動(dòng)計(jì)算此軸的長度,因此下面的程序?qū)?shù)組c的shape改為了(2,6):
>>> c.shape = 2,-1>>> carray([[ 1, 2, 3, 4, 4, 5], [ 6, 7, 7, 8, 9, 10]])使用數(shù)組的reshape方法,可以創(chuàng)建一個(gè)改變了尺寸的新數(shù)組,原數(shù)組的shape保持不變:
>>> d = a.reshape((2,2))>>> darray([[1, 2], [3, 4]])>>> aarray([1, 2, 3, 4])數(shù)組a和d其實(shí)共享數(shù)據(jù)存儲內(nèi)存區(qū)域,因此修改其中任意一個(gè)數(shù)組的元素都會同時(shí)修改另外一個(gè)數(shù)組的內(nèi)容:
>>> a[1] = 100 # 將數(shù)組a的第一個(gè)元素改為100>>> d # 注意數(shù)組d中的2也被改變了array([[ 1, 100], [ 3, 4]])數(shù)組的元素類型可以通過dtype屬性獲得。上面例子中的參數(shù)序列的元素都是整數(shù),因此所創(chuàng)建的數(shù)組的元素類型也是整數(shù),并且是32bit的長整型。可以通過dtype參數(shù)在創(chuàng)建時(shí)指定元素類型:
>>> np.array([[1, 2, 3, 4],[4, 5, 6, 7], [7, 8, 9, 10]], dtype=np.float)array([[ 1., 2., 3., 4.], [ 4., 5., 6., 7.], [ 7., 8., 9., 10.]])>>> np.array([[1, 2, 3, 4],[4, 5, 6, 7], [7, 8, 9, 10]], dtype=np.complex)array([[ 1.+0.j, 2.+0.j, 3.+0.j, 4.+0.j], [ 4.+0.j, 5.+0.j, 6.+0.j, 7.+0.j], [ 7.+0.j, 8.+0.j, 9.+0.j, 10.+0.j]])上面的例子都是先創(chuàng)建一個(gè)Python序列,然后通過array函數(shù)將其轉(zhuǎn)換為數(shù)組,這樣做顯然效率不高。因此NumPy提供了很多專門用來創(chuàng)建數(shù)組的函數(shù)。下面的每個(gè)函數(shù)都有一些關(guān)鍵字參數(shù),具體用法請查看函數(shù)說明。
arange函數(shù)類似于python的range函數(shù),通過指定開始值、終值和步長來創(chuàng)建一維數(shù)組,注意數(shù)組不包括終值:
>>> np.arange(0,1,0.1)array([ 0. , 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9])linspace函數(shù)通過指定開始值、終值和元素個(gè)數(shù)來創(chuàng)建一維數(shù)組,可以通過endpoint關(guān)鍵字指定是否包括終值,缺省設(shè)置是包括終值:
>>> np.linspace(0, 1, 12)array([ 0. , 0.09090909, 0.18181818, 0.27272727, 0.36363636, 0.45454545, 0.54545455, 0.63636364, 0.72727273, 0.81818182, 0.90909091, 1. ])logspace函數(shù)和linspace類似,不過它創(chuàng)建等比數(shù)列,下面的例子產(chǎn)生1(10^0)到100(10^2)、有20個(gè)元素的等比數(shù)列:
>>> np.logspace(0, 2, 20)array([ 1. , 1.27427499, 1.62377674, 2.06913808, 2.6366509 , 3.35981829, 4.2813324 , 5.45559478, 6.95192796, 8.8586679 , 11.28837892, 14.38449888, 18.32980711, 23.35721469, 29.76351442, 37.92690191, 48.32930239, 61.58482111, 78.47599704, 100. ])此外,使用frombuffer, fromstring, fromfile等函數(shù)可以從字節(jié)序列創(chuàng)建數(shù)組,下面以fromstring為例:
>>> s = "abcdefgh"Python的字符串實(shí)際上是字節(jié)序列,每個(gè)字符占一個(gè)字節(jié),因此如果從字符串s創(chuàng)建一個(gè)8bit的整數(shù)數(shù)組的話,所得到的數(shù)組正好就是字符串中每個(gè)字符的ASCII編碼:
>>> np.fromstring(s, dtype=np.int8)array([ 97, 98, 99, 100, 101, 102, 103, 104], dtype=int8)如果從字符串s創(chuàng)建16bit的整數(shù)數(shù)組,那么兩個(gè)相鄰的字節(jié)就表示一個(gè)整數(shù),把字節(jié)98和字節(jié)97當(dāng)作一個(gè)16位的整數(shù),它的值就是98*256+97 = 25185。可以看出內(nèi)存中是以little endian(低位字節(jié)在前)方式保存數(shù)據(jù)的。
>>> np.fromstring(s, dtype=np.int16)array([25185, 25699, 26213, 26727], dtype=int16)>>> 98*256+9725185如果把整個(gè)字符串轉(zhuǎn)換為一個(gè)64位的雙精度浮點(diǎn)數(shù)數(shù)組,那么它的值是:
>>> np.fromstring(s, dtype=np.float)array([ 8.54088322e+194])顯然這個(gè)例子沒有什么意義,但是可以想象如果我們用C語言的二進(jìn)制方式寫了一組double類型的數(shù)值到某個(gè)文件中,那們可以從此文件讀取相應(yīng)的數(shù)據(jù),并通過fromstring函數(shù)將其轉(zhuǎn)換為float64類型的數(shù)組。
我們可以寫一個(gè)Python的函數(shù),它將數(shù)組下標(biāo)轉(zhuǎn)換為數(shù)組中對應(yīng)的值,然后使用此函數(shù)創(chuàng)建數(shù)組:
>>> def func(i):... return i%4+1...>>> np.fromfunction(func, (10,))array([ 1., 2., 3., 4., 1., 2., 3., 4., 1., 2.])fromfunction函數(shù)的第一個(gè)參數(shù)為計(jì)算每個(gè)數(shù)組元素的函數(shù),第二個(gè)參數(shù)為數(shù)組的大小(shape),因?yàn)樗С侄嗑S數(shù)組,所以第二個(gè)參數(shù)必須是一個(gè)序列,本例中用(10,)創(chuàng)建一個(gè)10元素的一維數(shù)組。
下面的例子創(chuàng)建一個(gè)二維數(shù)組表示九九乘法表,輸出的數(shù)組a中的每個(gè)元素a[i, j]都等于func2(i, j):
>>> def func2(i, j):... return (i+1) * (j+1)...>>> a = np.fromfunction(func2, (9,9))>>> aarray([[ 1., 2., 3., 4., 5., 6., 7., 8., 9.], [ 2., 4., 6., 8., 10., 12., 14., 16., 18.], [ 3., 6., 9., 12., 15., 18., 21., 24., 27.], [ 4., 8., 12., 16., 20., 24., 28., 32., 36.], [ 5., 10., 15., 20., 25., 30., 35., 40., 45.], [ 6., 12., 18., 24., 30., 36., 42., 48., 54.], [ 7., 14., 21., 28., 35., 42., 49., 56., 63.], [ 8., 16., 24., 32., 40., 48., 56., 64., 72.], [ 9., 18., 27., 36., 45., 54., 63., 72., 81.]])2.1.2 存取元素
數(shù)組元素的存取方法和Python的標(biāo)準(zhǔn)方法相同:
>>> a = np.arange(10)>>> a[5] # 用整數(shù)作為下標(biāo)可以獲取數(shù)組中的某個(gè)元素5>>> a[3:5] # 用范圍作為下標(biāo)獲取數(shù)組的一個(gè)切片,包括a[3]不包括a[5]array([3, 4])>>> a[:5] # 省略開始下標(biāo),表示從a[0]開始array([0, 1, 2, 3, 4])>>> a[:-1] # 下標(biāo)可以使用負(fù)數(shù),表示從數(shù)組后往前數(shù)array([0, 1, 2, 3, 4, 5, 6, 7, 8])>>> a[2:4] = 100,101 # 下標(biāo)還可以用來修改元素的值>>> aarray([ 0, 1, 100, 101, 4, 5, 6, 7, 8, 9])>>> a[1:-1:2] # 范圍中的第三個(gè)參數(shù)表示步長,2表示隔一個(gè)元素取一個(gè)元素array([ 1, 101, 5, 7])>>> a[::-1] # 省略范圍的開始下標(biāo)和結(jié)束下標(biāo),步長為-1,整個(gè)數(shù)組頭尾顛倒array([ 9, 8, 7, 6, 5, 4, 101, 100, 1, 0])>>> a[5:1:-2] # 步長為負(fù)數(shù)時(shí),開始下標(biāo)必須大于結(jié)束下標(biāo)array([ 5, 101])和Python的列表序列不同,通過下標(biāo)范圍獲取的新的數(shù)組是原始數(shù)組的一個(gè)視圖。它與原始數(shù)組共享同一塊數(shù)據(jù)空間:
>>> b = a[3:7] # 通過下標(biāo)范圍產(chǎn)生一個(gè)新的數(shù)組b,b和a共享同一塊數(shù)據(jù)空間>>> barray([101, 4, 5, 6])>>> b[2] = -10 # 將b的第2個(gè)元素修改為-10>>> barray([101, 4, -10, 6])>>> a # a的第5個(gè)元素也被修改為10array([ 0, 1, 100, 101, 4, -10, 6, 7, 8, 9])除了使用下標(biāo)范圍存取元素之外,NumPy還提供了兩種存取元素的高級方法。
使用整數(shù)序列
當(dāng)使用整數(shù)序列對數(shù)組元素進(jìn)行存取時(shí),將使用整數(shù)序列中的每個(gè)元素作為下標(biāo),整數(shù)序列可以是列表或者數(shù)組。使用整數(shù)序列作為下標(biāo)獲得的數(shù)組不和原始數(shù)組共享數(shù)據(jù)空間。
>>> x = np.arange(10,1,-1)>>> xarray([10, 9, 8, 7, 6, 5, 4, 3, 2])>>> x[[3, 3, 1, 8]] # 獲取x中的下標(biāo)為3, 3, 1, 8的4個(gè)元素,組成一個(gè)新的數(shù)組array([7, 7, 9, 2])>>> b = x[np.array([3,3,-3,8])] #下標(biāo)可以是負(fù)數(shù)>>> b[2] = 100>>> barray([7, 7, 100, 2])>>> x # 由于b和x不共享數(shù)據(jù)空間,因此x中的值并沒有改變array([10, 9, 8, 7, 6, 5, 4, 3, 2])>>> x[[3,5,1]] = -1, -2, -3 # 整數(shù)序列下標(biāo)也可以用來修改元素的值>>> xarray([10, -3, 8, -1, 6, -2, 4, 3, 2])使用布爾數(shù)組
當(dāng)使用布爾數(shù)組b作為下標(biāo)存取數(shù)組x中的元素時(shí),將收集數(shù)組x中所有在數(shù)組b中對應(yīng)下標(biāo)為True的元素。使用布爾數(shù)組作為下標(biāo)獲得的數(shù)組不和原始數(shù)組共享數(shù)據(jù)空間,注意這種方式只對應(yīng)于布爾數(shù)組,不能使用布爾列表。
>>> x = np.arange(5,0,-1)>>> xarray([5, 4, 3, 2, 1])>>> x[np.array([True, False, True, False, False])]>>> # 布爾數(shù)組中下標(biāo)為0,2的元素為True,因此獲取x中下標(biāo)為0,2的元素array([5, 3])>>> x[[True, False, True, False, False]]>>> # 如果是布爾列表,則把True當(dāng)作1, False當(dāng)作0,按照整數(shù)序列方式獲取x中的元素array([4, 5, 4, 5, 5])>>> x[np.array([True, False, True, True])]>>> # 布爾數(shù)組的長度不夠時(shí),不夠的部分都當(dāng)作Falsearray([5, 3, 2])>>> x[np.array([True, False, True, True])] = -1, -2, -3>>> # 布爾數(shù)組下標(biāo)也可以用來修改元素>>> xarray([-1, 4, -2, -3, 1])布爾數(shù)組一般不是手工產(chǎn)生,而是使用布爾運(yùn)算的ufunc函數(shù)產(chǎn)生,關(guān)于ufunc函數(shù)請參照 ufunc運(yùn)算 一節(jié)。
>>> x = np.random.rand(10) # 產(chǎn)生一個(gè)長度為10,元素值為0-1的隨機(jī)數(shù)的數(shù)組>>> xarray([ 0.72223939, 0.921226 , 0.7770805 , 0.2055047 , 0.17567449, 0.95799412, 0.12015178, 0.7627083 , 0.43260184, 0.91379859])>>> x>0.5>>> # 數(shù)組x中的每個(gè)元素和0.5進(jìn)行大小比較,得到一個(gè)布爾數(shù)組,True表示x中對應(yīng)的值大于0.5array([ True, True, True, False, False, True, False, True, False, True], dtype=bool)>>> x[x>0.5]>>> # 使用x>0.5返回的布爾數(shù)組收集x中的元素,因此得到的結(jié)果是x中所有大于0.5的元素的數(shù)組array([ 0.72223939, 0.921226 , 0.7770805 , 0.95799412, 0.7627083 , 0.91379859])2.1.3 多維數(shù)組
多維數(shù)組的存取和一維數(shù)組類似,因?yàn)槎嗑S數(shù)組有多個(gè)軸,因此它的下標(biāo)需要用多個(gè)值來表示,NumPy采用組元(tuple)作為數(shù)組的下標(biāo)。如圖2.1所示,a為一個(gè)6x6的數(shù)組,圖中用顏色區(qū)分了各個(gè)下標(biāo)以及其對應(yīng)的選擇區(qū)域。
組元不需要圓括號
雖然我們經(jīng)常在Python中用圓括號將組元括起來,但是其實(shí)組元的語法定義只需要用逗號隔開即可,例如 x,y=y,x 就是用組元交換變量值的一個(gè)例子。
圖2.1 使用數(shù)組切片語法訪問多維數(shù)組中的元素
如何創(chuàng)建這個(gè)數(shù)組
你也許會對如何創(chuàng)建a這樣的數(shù)組感到好奇,數(shù)組a實(shí)際上是一個(gè)加法表,縱軸的值為0, 10, 20, 30, 40, 50;橫軸的值為0, 1, 2, 3, 4, 5。縱軸的每個(gè)元素都和橫軸的每個(gè)元素求和,就得到圖中所示的數(shù)組a。你可以用下面的語句創(chuàng)建它,至于其原理我們將在后面的章節(jié)進(jìn)行討論:
>>> 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]])多維數(shù)組同樣也可以使用整數(shù)序列和布爾數(shù)組進(jìn)行存取。
圖2.2 使用整數(shù)序列和布爾數(shù)組訪問多維數(shù)組中的元素
a[(0,1,2,3,4),(1,2,3,4,5)] : 用于存取數(shù)組的下標(biāo)和仍然是一個(gè)有兩個(gè)元素的組元,組元中的每個(gè)元素都是整數(shù)序列,分別對應(yīng)數(shù)組的第0軸和第1軸。從兩個(gè)序列的對應(yīng)位置取出兩個(gè)整數(shù)組成下標(biāo): a[0,1], a[1,2], ..., a[4,5]。a[3:, [0, 2, 5]] : 下標(biāo)中的第0軸是一個(gè)范圍,它選取第3行之后的所有行;第1軸是整數(shù)序列,它選取第0, 2, 5三列。a[mask, 2] : 下標(biāo)的第0軸是一個(gè)布爾數(shù)組,它選取第0,2,5行;第1軸是一個(gè)整數(shù),選取第2列。2.1.4 結(jié)構(gòu)數(shù)組
在C語言中我們可以通過struct關(guān)鍵字定義結(jié)構(gòu)類型,結(jié)構(gòu)中的字段占據(jù)連續(xù)的內(nèi)存空間,每個(gè)結(jié)構(gòu)體占用的內(nèi)存大小都相同,因此可以很容易地定義結(jié)構(gòu)數(shù)組。和C語言一樣,在NumPy中也很容易對這種結(jié)構(gòu)數(shù)組進(jìn)行操作。只要NumPy中的結(jié)構(gòu)定義和C語言中的定義相同,NumPy就可以很方便地讀取C語言的結(jié)構(gòu)數(shù)組的二進(jìn)制數(shù)據(jù),轉(zhuǎn)換為NumPy的結(jié)構(gòu)數(shù)組。
假設(shè)我們需要定義一個(gè)結(jié)構(gòu)數(shù)組,它的每個(gè)元素都有name, age和weight字段。在NumPy中可以如下定義:
import numpy as nppersontype = np.dtype({ 'names':['name', 'age', 'weight'], 'formats':['S32','i', 'f']})a = np.array([("Zhang",32,75.5),("Wang",24,65.2)], dtype=persontype)我們先創(chuàng)建一個(gè)dtype對象persontype,通過其字典參數(shù)描述結(jié)構(gòu)類型的各個(gè)字段。字典有兩個(gè)關(guān)鍵字:names,formats。每個(gè)關(guān)鍵字對應(yīng)的值都是一個(gè)列表。names定義結(jié)構(gòu)中的每個(gè)字段名,而formats則定義每個(gè)字段的類型:
S32 : 32個(gè)字節(jié)的字符串類型,由于結(jié)構(gòu)中的每個(gè)元素的大小必須固定,因此需要指定字符串的長度i : 32bit的整數(shù)類型,相當(dāng)于np.int32f : 32bit的單精度浮點(diǎn)數(shù)類型,相當(dāng)于np.float32然后我們調(diào)用array函數(shù)創(chuàng)建數(shù)組,通過關(guān)鍵字參數(shù) dtype=persontype, 指定所創(chuàng)建的數(shù)組的元素類型為結(jié)構(gòu)persontype。運(yùn)行上面程序之后,我們可以在IPython中執(zhí)行如下的語句查看數(shù)組a的元素類型
>>> a.dtypedtype([('name', '|S32'), ('age', '<i4'), ('weight', '<f4')])這里我們看到了另外一種描述結(jié)構(gòu)類型的方法: 一個(gè)包含多個(gè)組元的列表,其中形如 (字段名, 類型描述) 的組元描述了結(jié)構(gòu)中的每個(gè)字段。類型描述前面為我們添加了 '|', '<' 等字符,這些字符用來描述字段值的字節(jié)順序:
| : 忽視字節(jié)順序< : 低位字節(jié)在前> : 高位字節(jié)在前結(jié)構(gòu)數(shù)組的存取方式和一般數(shù)組相同,通過下標(biāo)能夠取得其中的元素,注意元素的值看上去像是組元,實(shí)際上它是一個(gè)結(jié)構(gòu):
>>> a[0]('Zhang', 32, 75.5)>>> a[0].dtypedtype([('name', '|S32'), ('age', '<i4'), ('weight', '<f4')])a[0]是一個(gè)結(jié)構(gòu)元素,它和數(shù)組a共享內(nèi)存數(shù)據(jù),因此可以通過修改它的字段,改變原始數(shù)組中的對應(yīng)字段:
>>> c = a[1]>>> c["name"] = "Li">>> a[1]["name"]"Li"結(jié)構(gòu)像字典一樣可以通過字符串下標(biāo)獲取其對應(yīng)的字段值:
>>> a[0]["name"]'Zhang'我們不但可以獲得結(jié)構(gòu)元素的某個(gè)字段,還可以直接獲得結(jié)構(gòu)數(shù)組的字段,它返回的是原始數(shù)組的視圖,因此可以通過修改b[0]改變a[0]["age"]:
>>> b=a[:]["age"] # 或者a["age"]>>> barray([32, 24])>>> b[0] = 40>>> a[0]["age"]40通過調(diào)用a.tostring或者a.tofile方法,可以直接輸出數(shù)組a的二進(jìn)制形式:
>>> a.tofile("test.bin")利用下面的C語言程序可以將test.bin文件中的數(shù)據(jù)讀取出來。
內(nèi)存對齊
C語言的結(jié)構(gòu)體為了內(nèi)存尋址方便,會自動(dòng)的添加一些填充用的字節(jié),這叫做內(nèi)存對齊。例如如果把下面的name[32]改為name[30]的話,由于內(nèi)存對齊問題,在name和age中間會填補(bǔ)兩個(gè)字節(jié),最終的結(jié)構(gòu)體大小不會改變。因此如果numpy中的所配置的內(nèi)存大小不符合C語言的對齊規(guī)范的話,將會出現(xiàn)數(shù)據(jù)錯(cuò)位。為了解決這個(gè)問題,在創(chuàng)建dtype對象時(shí),可以傳遞參數(shù)align=True,這樣numpy的結(jié)構(gòu)數(shù)組的內(nèi)存對齊和C語言的結(jié)構(gòu)體就一致了。
#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();}結(jié)構(gòu)類型中可以包括其它的結(jié)構(gòu)類型,下面的語句創(chuàng)建一個(gè)有一個(gè)字段f1的結(jié)構(gòu),f1的值是另外一個(gè)結(jié)構(gòu),它有字段f2,其類型為16bit整數(shù)。
>>> np.dtype([('f1', [('f2', np.int16)])])dtype([('f1', [('f2', '<i2')])])當(dāng)某個(gè)字段類型為數(shù)組時(shí),用組元的第三個(gè)參數(shù)表示,下面描述的f1字段是一個(gè)shape為(2,3)的雙精度浮點(diǎn)數(shù)組:
>>> np.dtype([('f0', 'i4'), ('f1', 'f8', (2, 3))])dtype([('f0', '<i4'), ('f1', '<f8', (2, 3))])用下面的字典參數(shù)也可以定義結(jié)構(gòu)類型,字典的關(guān)鍵字為結(jié)構(gòu)中字段名,值為字段的類型描述,但是由于字典的關(guān)鍵字是沒有順序的,因此字段的順序需要在類型描述中給出,類型描述是一個(gè)組元,它的第二個(gè)值給出字段的字節(jié)為單位的偏移量,例如age字段的偏移量為25個(gè)字節(jié):
>>> np.dtype({'surname':('S25',0),'age':(np.uint8,25)})dtype([('surname', '|S25'), ('age', '|u1')])2.1.5 內(nèi)存結(jié)構(gòu)
下面讓我們來看看ndarray數(shù)組對象是如何在內(nèi)存中儲存的。如圖2.3所示,關(guān)于數(shù)組的描述信息保存在一個(gè)數(shù)據(jù)結(jié)構(gòu)中,這個(gè)結(jié)構(gòu)引用兩個(gè)對象:一塊用于保存數(shù)據(jù)的存儲區(qū)域和一個(gè)用于描述元素類型的dtype對象。
圖2.3 ndarray數(shù)組對象在內(nèi)存中的儲存方式
數(shù)據(jù)存儲區(qū)域保存著數(shù)組中所有元素的二進(jìn)制數(shù)據(jù),dtype對象則知道如何將元素的二進(jìn)制數(shù)據(jù)轉(zhuǎn)換為可用的值。數(shù)組的維數(shù)、大小等信息都保存在ndarray數(shù)組對象的數(shù)據(jù)結(jié)構(gòu)中。圖中顯示的是如下數(shù)組的內(nèi)存結(jié)構(gòu):
>>> a = np.array([[0,1,2],[3,4,5],[6,7,8]], dtype=np.float32)strides中保存的是當(dāng)每個(gè)軸的下標(biāo)增加1時(shí),數(shù)據(jù)存儲區(qū)中的指針?biāo)黾拥淖止?jié)數(shù)。例如圖中的strides為12,4,即第0軸的下標(biāo)增加1時(shí),數(shù)據(jù)的地址增加12個(gè)字節(jié):即a[1,0]的地址比a[0,0]的地址要高12個(gè)字節(jié),正好是3個(gè)單精度浮點(diǎn)數(shù)的總字節(jié)數(shù);第1軸下標(biāo)增加1時(shí),數(shù)據(jù)的地址增加4個(gè)字節(jié),正好是單精度浮點(diǎn)數(shù)的字節(jié)數(shù)。
如果strides中的數(shù)值正好和對應(yīng)軸所占據(jù)的字節(jié)數(shù)相同的話,那么數(shù)據(jù)在內(nèi)存中是連續(xù)存儲的。然而數(shù)據(jù)并不一定都是連續(xù)儲存的,前面介紹過通過下標(biāo)范圍得到新的數(shù)組是原始數(shù)組的視圖,即它和原始視圖共享數(shù)據(jù)存儲區(qū)域:
>>> b = a[::2,::2]>>> barray([[ 0., 2.], [ 6., 8.]], dtype=float32)>>> b.strides(24, 8)由于數(shù)組b和數(shù)組a共享數(shù)據(jù)存儲區(qū),而b中的第0軸和第1軸都是數(shù)組a中隔一個(gè)元素取一個(gè),因此數(shù)組b的strides變成了24,8,正好都是數(shù)組a的兩倍。 對照前面的圖很容易看出數(shù)據(jù)0和2的地址相差8個(gè)字節(jié),而0和6的地址相差24個(gè)字節(jié)。
元素在數(shù)據(jù)存儲區(qū)中的排列格式有兩種:C語言格式和Fortan語言格式。在C語言中,多維數(shù)組的第0軸是最上位的,即第0軸的下標(biāo)增加1時(shí),元素的地址增加的字節(jié)數(shù)最多;而Fortan語言的多維數(shù)組的第0軸是最下位的,即第0軸的下標(biāo)增加1時(shí),地址只增加一個(gè)元素的字節(jié)數(shù)。在NumPy中,元素在內(nèi)存中的排列缺省是以C語言格式存儲的,如果你希望改為Fortan格式的話,只需要給數(shù)組傳遞order="F"參數(shù):
>>> c = np.array([[0,1,2],[3,4,5],[6,7,8]], dtype=np.float32, order="F")>>> c.strides(4, 12)2.2 ufunc運(yùn)算
ufunc是universal function的縮寫,它是一種能對數(shù)組的每個(gè)元素進(jìn)行操作的函數(shù)。NumPy內(nèi)置的許多ufunc函數(shù)都是在C語言級別實(shí)現(xiàn)的,因此它們的計(jì)算速度非常快。讓我們來看一個(gè)例子:
>>> x = np.linspace(0, 2*np.pi, 10)# 對數(shù)組x中的每個(gè)元素進(jìn)行正弦計(jì)算,返回一個(gè)同樣大小的新數(shù)組>>> y = np.sin(x)>>> yarray([ 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產(chǎn)生一個(gè)從0到2*PI的等距離的10個(gè)數(shù),然后將其傳遞給sin函數(shù),由于np.sin是一個(gè)ufunc函數(shù),因此它對x中的每個(gè)元素求正弦值,然后將結(jié)果返回,并且賦值給y。計(jì)算之后x中的值并沒有改變,而是新創(chuàng)建了一個(gè)數(shù)組保存結(jié)果。如果我們希望將sin函數(shù)所計(jì)算的結(jié)果直接覆蓋到數(shù)組x上去的話,可以將要被覆蓋的數(shù)組作為第二個(gè)參數(shù)傳遞給ufunc函數(shù)。例如::
>>> t = np.sin(x,x)>>> xarray([ 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)Truesin函數(shù)的第二個(gè)參數(shù)也是x,那么它所做的事情就是對x中的每給值求正弦值,并且把結(jié)果保存到x中的對應(yīng)的位置中。此時(shí)函數(shù)的返回值仍然是整個(gè)計(jì)算的結(jié)果,只不過它就是x,因此兩個(gè)變量的id是相同的(變量t和變量x指向同一塊內(nèi)存區(qū)域)。
我用下面這個(gè)小程序,比較了一下numpy.math和Python標(biāo)準(zhǔn)庫的math.sin的計(jì)算速度::
import timeimport mathimport numpy as npx = [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() - startx = [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在我的電腦上計(jì)算100萬次正弦值,numpy.sin比math.sin快10倍多。這得利于numpy.sin在C語言級別的循環(huán)計(jì)算。numpy.sin同樣也支持對單個(gè)數(shù)值求正弦,例如:numpy.sin(0.5)。不過值得注意的是,對單個(gè)數(shù)的計(jì)算math.sin則比numpy.sin快得多了,讓我們看下面這個(gè)測試程序:
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的計(jì)算速度只有math.sin的1/5。這是因?yàn)閚umpy.sin為了同時(shí)支持?jǐn)?shù)組和單個(gè)值的計(jì)算,其C語言的內(nèi)部實(shí)現(xiàn)要比math.sin復(fù)雜很多,如果我們同樣在Python級別進(jìn)行循環(huán)的話,就會看出其中的差別了。此外,numpy.sin返回的數(shù)的類型和math.sin返回的類型有所不同,math.sin返回的是Python的標(biāo)準(zhǔn)float類型,而numpy.sin則返回一個(gè)numpy.float64類型:
>>> type(math.sin(0.5))<type 'float'>>>> type(np.sin(0.5))<type 'numpy.float64'>通過上面的例子我們了解了如何最有效率地使用math庫和numpy庫中的數(shù)學(xué)函數(shù)。因?yàn)樗鼈兏饔虚L短,因此在導(dǎo)入時(shí)不建議使用*號全部載入,而是應(yīng)該使用import numpy as np的方式載入,這樣我們可以根據(jù)需要選擇合適的函數(shù)調(diào)用。
NumPy中有眾多的ufunc函數(shù)為我們提供各式各樣的計(jì)算。除了sin這種單輸入函數(shù)之外,還有許多多個(gè)輸入的函數(shù),add函數(shù)就是一個(gè)最常用的例子。先來看一個(gè)例子:
>>> a = np.arange(0,4)>>> aarray([0, 1, 2, 3])>>> b = np.arange(1,5)>>> barray([1, 2, 3, 4])>>> np.add(a,b)array([1, 3, 5, 7])>>> np.add(a,b,a)array([1, 3, 5, 7])>>> aarray([1, 3, 5, 7])add函數(shù)返回一個(gè)新的數(shù)組,此數(shù)組的每個(gè)元素都為兩個(gè)參數(shù)數(shù)組的對應(yīng)元素之和。它接受第3個(gè)參數(shù)指定計(jì)算結(jié)果所要寫入的數(shù)組,如果指定的話,add函數(shù)就不再產(chǎn)生新的數(shù)組。
由于Python的操作符重載功能,計(jì)算兩個(gè)數(shù)組相加可以簡單地寫為a+b,而np.add(a,b,a)則可以用a+=b來表示。下面是數(shù)組的運(yùn)算符和其對應(yīng)的ufunc函數(shù)的一個(gè)列表,注意除號"/"的意義根據(jù)是否激活__future__.division有所不同。
| y = x1 + x2: | add(x1, x2 [, y]) |
| y = x1 - x2: | subtract(x1, x2 [, y]) |
| y = x1 * x2: | multiply (x1, x2 [, y]) |
| y = x1 / x2: | divide (x1, x2 [, y]), 如果兩個(gè)數(shù)組的元素為整數(shù),那么用整數(shù)除法 |
| y = x1 / x2: | true divide (x1, x2 [, y]), 總是返回精確的商 |
| y = x1 // x2: | floor divide (x1, x2 [, y]), 總是對返回值取整 |
| y = -x: | negative(x [,y]) |
| y = x1**x2: | power(x1, x2 [, y]) |
| y = x1 % x2: | remainder(x1, x2 [, y]), mod(x1, x2, [, y]) |
數(shù)組對象支持這些操作符,極大地簡化了算式的編寫,不過要注意如果你的算式很復(fù)雜,并且要運(yùn)算的數(shù)組很大的話,會因?yàn)楫a(chǎn)生大量的中間結(jié)果而降低程序的運(yùn)算效率。例如:假設(shè)a b c三個(gè)數(shù)組采用算式x=a*b+c計(jì)算,那么它相當(dāng)于:
t = a * bx = t + cdel t也就是說需要產(chǎn)生一個(gè)數(shù)組t保存乘法的計(jì)算結(jié)果,然后再產(chǎn)生最后的結(jié)果數(shù)組x。我們可以通過手工將一個(gè)算式分解為x = a*b; x += c,以減少一次內(nèi)存分配。
通過組合標(biāo)準(zhǔn)的ufunc函數(shù)的調(diào)用,可以實(shí)現(xiàn)各種算式的數(shù)組計(jì)算。不過有些時(shí)候這種算式不易編寫,而針對每個(gè)元素的計(jì)算函數(shù)卻很容易用Python實(shí)現(xiàn),這時(shí)可以用frompyfunc函數(shù)將一個(gè)計(jì)算單個(gè)元素的函數(shù)轉(zhuǎn)換成ufunc函數(shù)。這樣就可以方便地用所產(chǎn)生的ufunc函數(shù)對數(shù)組進(jìn)行計(jì)算了。讓我們來看一個(gè)例子。
我們想用一個(gè)分段函數(shù)描述三角波,三角波的樣子如圖2.4所示:
圖2.4 三角波可以用分段函數(shù)進(jìn)行計(jì)算
我們很容易根據(jù)上圖所示寫出如下的計(jì)算三角波某點(diǎn)y坐標(biāo)的函數(shù):
def triangle_wave(x, c, c0, hc): x = x - int(x) # 三角波的周期為1,因此只取x坐標(biāo)的小數(shù)部分進(jìn)行計(jì)算 if x >= c: r = 0.0 elif x < c0: r = x / c0 * hc else: r = (c-x) / (c-c0) * hc return r顯然triangle_wave函數(shù)只能計(jì)算單個(gè)數(shù)值,不能對數(shù)組直接進(jìn)行處理。我們可以用下面的方法先使用列表包容(List comprehension),計(jì)算出一個(gè)list,然后用array函數(shù)將列表轉(zhuǎn)換為數(shù)組:
x = np.linspace(0, 2, 1000)y = np.array([triangle_wave(t, 0.6, 0.4, 1.0) for t in x])這種做法每次都需要使用列表包容語法調(diào)用函數(shù),對于多維數(shù)組是很麻煩的。讓我們來看看如何用frompyfunc函數(shù)來解決這個(gè)問題:
triangle_ufunc = np.frompyfunc( lambda x: triangle_wave(x, 0.6, 0.4, 1.0), 1, 1)y2 = triangle_ufunc(x)frompyfunc的調(diào)用格式為frompyfunc(func, nin, nout),其中func是計(jì)算單個(gè)元素的函數(shù),nin是此函數(shù)的輸入?yún)?shù)的個(gè)數(shù),nout是此函數(shù)的返回值的個(gè)數(shù)。雖然triangle_wave函數(shù)有4個(gè)參數(shù),但是由于后三個(gè)c, c0, hc在整個(gè)計(jì)算中值都是固定的,因此所產(chǎn)生的ufunc函數(shù)其實(shí)只有一個(gè)參數(shù)。為了滿足這個(gè)條件,我們用一個(gè)lambda函數(shù)對triangle_wave的參數(shù)進(jìn)行一次包裝。這樣傳入frompyfunc的函數(shù)就只有一個(gè)參數(shù)了。這樣子做,效率并不是太高,另外還有一種方法:
def triangle_func(c, c0, hc): def trifunc(x): x = x - int(x) # 三角波的周期為1,因此只取x坐標(biāo)的小數(shù)部分進(jìn)行計(jì)算 if x >= c: r = 0.0 elif x < c0: r = x / c0 * hc else: r = (c-x) / (c-c0) * hc return r # 用trifunc函數(shù)創(chuàng)建一個(gè)ufunc函數(shù),可以直接對數(shù)組進(jìn)行計(jì)算, 不過通過此函數(shù) # 計(jì)算得到的是一個(gè)Object數(shù)組,需要進(jìn)行類型轉(zhuǎn)換 return np.frompyfunc(trifunc, 1, 1)y2 = triangle_func(0.6, 0.4, 1.0)(x)我們通過函數(shù)triangle_func包裝三角波的三個(gè)參數(shù),在其內(nèi)部定義一個(gè)計(jì)算三角波的函數(shù)trifunc,trifunc函數(shù)在調(diào)用時(shí)會采用triangle_func的參數(shù)進(jìn)行計(jì)算。最后triangle_func返回用frompyfunc轉(zhuǎn)換結(jié)果。
值得注意的是用frompyfunc得到的函數(shù)計(jì)算出的數(shù)組元素的類型為object,因?yàn)閒rompyfunc函數(shù)無法保證Python函數(shù)返回的數(shù)據(jù)類型都完全一致。因此還需要再次 y2.astype(np.float64)將其轉(zhuǎn)換為雙精度浮點(diǎn)數(shù)組。
2.2.1 廣播
當(dāng)我們使用ufunc函數(shù)對兩個(gè)數(shù)組進(jìn)行計(jì)算時(shí),ufunc函數(shù)會對這兩個(gè)數(shù)組的對應(yīng)元素進(jìn)行計(jì)算,因此它要求這兩個(gè)數(shù)組有相同的大小(shape相同)。如果兩個(gè)數(shù)組的shape不同的話,會進(jìn)行如下的廣播(broadcasting)處理:
讓所有輸入數(shù)組都向其中shape最長的數(shù)組看齊,shape中不足的部分都通過在前面加1補(bǔ)齊輸出數(shù)組的shape是輸入數(shù)組shape的各個(gè)軸上的最大值如果輸入數(shù)組的某個(gè)軸和輸出數(shù)組的對應(yīng)軸的長度相同或者其長度為1時(shí),這個(gè)數(shù)組能夠用來計(jì)算,否則出錯(cuò)當(dāng)輸入數(shù)組的某個(gè)軸的長度為1時(shí),沿著此軸運(yùn)算時(shí)都用此軸上的第一組值上述4條規(guī)則理解起來可能比較費(fèi)勁,讓我們來看一個(gè)實(shí)際的例子。
先創(chuàng)建一個(gè)二維數(shù)組a,其shape為(6,1):
>>> a = np.arange(0, 60, 10).reshape(-1, 1)>>> aarray([[ 0], [10], [20], [30], [40], [50]])>>> a.shape(6, 1)再創(chuàng)建一維數(shù)組b,其shape為(5,):
>>> b = np.arange(0, 5)>>> barray([0, 1, 2, 3, 4])>>> b.shape(5,)計(jì)算a和b的和,得到一個(gè)加法表,它相當(dāng)于計(jì)算a,b中所有元素組的和,得到一個(gè)shape為(6,5)的數(shù)組:
>>> c = a + b>>> carray([[ 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屬性)不同,根據(jù)規(guī)則1,需要讓b的shape向a對齊,于是將b的shape前面加1,補(bǔ)齊為(1,5)。相當(dāng)于做了如下計(jì)算:
>>> b.shape=1,5>>> barray([[0, 1, 2, 3, 4]])這樣加法運(yùn)算的兩個(gè)輸入數(shù)組的shape分別為(6,1)和(1,5),根據(jù)規(guī)則2,輸出數(shù)組的各個(gè)軸的長度為輸入數(shù)組各個(gè)軸上的長度的最大值,可知輸出數(shù)組的shape為(6,5)。
由于b的第0軸上的長度為1,而a的第0軸上的長度為6,因此為了讓它們在第0軸上能夠相加,需要將b在第0軸上的長度擴(kuò)展為6,這相當(dāng)于:
>>> b = b.repeat(6,axis=0)>>> barray([[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軸上的長度擴(kuò)展為5,這相當(dāng)于:
>>> a = a.repeat(5, axis=1)>>> aarray([[ 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]])經(jīng)過上述處理之后,a和b就可以按對應(yīng)元素進(jìn)行相加運(yùn)算了。
當(dāng)然,numpy在執(zhí)行a+b運(yùn)算時(shí),其內(nèi)部并不會真正將長度為1的軸用repeat函數(shù)進(jìn)行擴(kuò)展,如果這樣做的話就太浪費(fèi)空間了。
由于這種廣播計(jì)算很常用,因此numpy提供了一個(gè)快速產(chǎn)生如上面a,b數(shù)組的方法: ogrid對象:
>>> x,y = np.ogrid[0:5,0:5]>>> xarray([[0], [1], [2], [3], [4]])>>> yarray([[0, 1, 2, 3, 4]])ogrid是一個(gè)很有趣的對象,它像一個(gè)多維數(shù)組一樣,用切片組元作為下標(biāo)進(jìn)行存取,返回的是一組可以用來廣播計(jì)算的數(shù)組。其切片下標(biāo)有兩種形式:
開始值:結(jié)束值:步長,和np.arange(開始值, 結(jié)束值, 步長)類似
開始值:結(jié)束值:長度j,當(dāng)?shù)谌齻€(gè)參數(shù)為虛數(shù)時(shí),它表示返回的數(shù)組的長度,和np.linspace(開始值, 結(jié)束值, 長度)類似:
>>> x, y = np.ogrid[0:1:4j, 0:1:3j]>>> xarray([[ 0. ], [ 0.33333333], [ 0.66666667], [ 1. ]])>>> yarray([[ 0. , 0.5, 1. ]])ogrid為什么不是函數(shù)
根據(jù)Python的語法,只有在中括號中才能使用用冒號隔開的切片語法,如果ogrid是函數(shù)的話,那么這些切片必須使用slice函數(shù)創(chuàng)建,這顯然會增加代碼的長度。
利用ogrid的返回值,我能很容易計(jì)算x, y網(wǎng)格面上各點(diǎn)的值,或者x, y, z網(wǎng)格體上各點(diǎn)的值。下面是繪制三維曲面 x * exp(x**2 - y**2) 的程序:
import numpy as npfrom enthought.mayavi import mlabx, 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庫快速繪制如圖2.5所示的3D曲面,關(guān)于mlab的相關(guān)內(nèi)容將在今后的章節(jié)進(jìn)行介紹。
圖2.5 使用ogrid創(chuàng)建的三維曲面
2.2.2 ufunc的方法
ufunc函數(shù)本身還有些方法,這些方法只對兩個(gè)輸入一個(gè)輸出的ufunc函數(shù)有效,其它的ufunc對象調(diào)用這些方法時(shí)會拋出ValueError異常。
reduce 方法和Python的reduce函數(shù)類似,它沿著axis軸對array進(jìn)行操作,相當(dāng)于將<op>運(yùn)算符插入到沿axis軸的所有子數(shù)組或者元素當(dāng)中。
<op>.reduce (array=, axis=0, dtype=None)例如:
>>> np.add.reduce([1,2,3]) # 1 + 2 + 36>>> np.add.reduce([[1,2,3],[4,5,6]], axis=1) # 1,4 + 2,5 + 3,6array([ 6, 15])accumulate 方法和reduce方法類似,只是它返回的數(shù)組和輸入的數(shù)組的shape相同,保存所有的中間計(jì)算結(jié)果:
>>> np.add.accumulate([1,2,3])array([1, 3, 6])>>> np.add.accumulate([[1,2,3],[4,5,6]], axis=1)array([[ 1, 3, 6], [ 4, 9, 15]])reduceat 方法計(jì)算多組reduce的結(jié)果,通過indices參數(shù)指定一系列reduce的起始和終了位置。reduceat的計(jì)算有些特別,讓我們通過一個(gè)例子來解釋一下:
>>> a = np.array([1,2,3,4])>>> result = np.add.reduceat(a,indices=[0,1,0,2,0,3,0])>>> resultarray([ 1, 2, 3, 3, 6, 4, 10])對于indices中的每個(gè)元素都會調(diào)用reduce函數(shù)計(jì)算出一個(gè)值來,因此最終計(jì)算結(jié)果的長度和indices的長度相同。結(jié)果result數(shù)組中除最后一個(gè)元素之外,都按照如下計(jì)算得出:
if indices[i] < indices[i+1]: result[i] = np.reduce(a[indices[i]:indices[i+1]])else: result[i] = a[indices[i]而最后一個(gè)元素如下計(jì)算:
np.reduce(a[indices[-1]:])因此上面例子中,結(jié)果的每個(gè)元素如下計(jì)算而得:
1 : a[0] = 12 : a[1] = 23 : a[0] + a[1] = 1 + 23 : a[2] = 36 : a[0] + a[1] + a[2] = 1 + 2 + 3 = 64 : a[3] = 410: a[0] + a[1] + a[2] + a[4] = 1+2+3+4 = 10可以看出result[::2]和a相等,而result[1::2]和np.add.accumulate(a)相等。
outer 方法,<op>.outer(a,b)方法的計(jì)算等同于如下程序:
>>> a.shape += (1,)*b.ndim>>> <op>(a,b)>>> a = a.squeeze()其中squeeze的功能是剔除數(shù)組a中長度為1的軸。如果你看不太明白這個(gè)等同程序的話,讓我們來看一個(gè)例子:
>>> np.multiply.outer([1,2,3,4,5],[2,3,4])array([[ 2, 3, 4], [ 4, 6, 8], [ 6, 9, 12], [ 8, 12, 16], [10, 15, 20]])可以看出通過outer方法計(jì)算的結(jié)果是如下的乘法表:
# 2, 3, 4# 1# 2# 3# 4# 5如果將這兩個(gè)數(shù)組按照等同程序一步一步的計(jì)算的話,就會發(fā)現(xiàn)乘法表最終是通過廣播的方式計(jì)算出來的。
2.3 矩陣運(yùn)算
NumPy和Matlab不一樣,對于多維數(shù)組的運(yùn)算,缺省情況下并不使用矩陣運(yùn)算,如果你希望對數(shù)組進(jìn)行矩陣運(yùn)算的話,可以調(diào)用相應(yīng)的函數(shù)。
matrix對象
numpy庫提供了matrix類,使用matrix類創(chuàng)建的是矩陣對象,它們的加減乘除運(yùn)算缺省采用矩陣方式計(jì)算,因此用法和matlab十分類似。但是由于NumPy中同時(shí)存在ndarray和matrix對象,因此用戶很容易將兩者弄混。這有違Python的“顯式優(yōu)于隱式”的原則,因此并不推薦在較復(fù)雜的程序中使用matrix。下面是使用matrix的一個(gè)例子:
>>> a = np.matrix([[1,2,3],[5,5,6],[7,9,9]])>>> a*a**-1matrix([[ 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]])因?yàn)閍是用matrix創(chuàng)建的矩陣對象,因此乘法和冪運(yùn)算符都變成了矩陣運(yùn)算,于是上面計(jì)算的是矩陣a和其逆矩陣的乘積,結(jié)果是一個(gè)單位矩陣。
矩陣的乘積可以使用dot函數(shù)進(jìn)行計(jì)算。對于二維數(shù)組,它計(jì)算的是矩陣乘積,對于一維數(shù)組,它計(jì)算的是其點(diǎn)積。當(dāng)需要將一維數(shù)組當(dāng)作列矢量或者行矢量進(jìn)行矩陣運(yùn)算時(shí),推薦先使用reshape函數(shù)將一維數(shù)組轉(zhuǎn)換為二維數(shù)組:
>>> a = array([1, 2, 3])>>> a.reshape((-1,1))array([[1], [2], [3]])>>> a.reshape((1,-1))array([[1, 2, 3]])除了dot計(jì)算乘積之外,NumPy還提供了inner和outer等多種計(jì)算乘積的函數(shù)。這些函數(shù)計(jì)算乘積的方式不同,尤其是當(dāng)對于多維數(shù)組的時(shí)候,更容易搞混。
dot : 對于兩個(gè)一維的數(shù)組,計(jì)算的是這兩個(gè)數(shù)組對應(yīng)下標(biāo)元素的乘積和(數(shù)學(xué)上稱之為內(nèi)積);對于二維數(shù)組,計(jì)算的是兩個(gè)數(shù)組的矩陣乘積;對于多維數(shù)組,它的通用計(jì)算公式如下,即結(jié)果數(shù)組中的每個(gè)元素都是:數(shù)組a的最后一維上的所有元素與數(shù)組b的倒數(shù)第二位上的所有元素的乘積和:
dot(a, b)[i,j,k,m] = sum(a[i,j,:] * b[k,:,m])下面以兩個(gè)3為數(shù)組的乘積演示一下dot乘積的計(jì)算結(jié)果:
首先創(chuàng)建兩個(gè)3維數(shù)組,這兩個(gè)數(shù)組的最后兩維滿足矩陣乘積的條件:
>>> a = np.arange(12).reshape(2,3,2)>>> b = np.arange(12,24).reshape(2,2,3)>>> c = np.dot(a,b)dot乘積的結(jié)果c可以看作是數(shù)組a,b的多個(gè)子矩陣的乘積:
>>> np.alltrue( c[0,:,0,:] == np.dot(a[0],b[0]) )True>>> np.alltrue( c[1,:,0,:] == np.dot(a[1],b[0]) )True>>> np.alltrue( c[0,:,1,:] == np.dot(a[0],b[1]) )True>>> np.alltrue( c[1,:,1,:] == np.dot(a[1],b[1]) )Trueinner : 和dot乘積一樣,對于兩個(gè)一維數(shù)組,計(jì)算的是這兩個(gè)數(shù)組對應(yīng)下標(biāo)元素的乘積和;對于多維數(shù)組,它計(jì)算的結(jié)果數(shù)組中的每個(gè)元素都是:數(shù)組a和b的最后一維的內(nèi)積,因此數(shù)組a和b的最后一維的長度必須相同:
inner(a, b)[i,j,k,m] = sum(a[i,j,:]*b[k,m,:])下面是inner乘積的演示:
>>> a = np.arange(12).reshape(2,3,2)>>> b = np.arange(12,24).reshape(2,3,2)>>> c = np.inner(a,b)>>> c.shape(2, 3, 2, 3)>>> c[0,0,0,0] == np.inner(a[0,0],b[0,0])True>>> c[0,1,1,0] == np.inner(a[0,1],b[1,0])True>>> c[1,2,1,2] == np.inner(a[1,2],b[1,2])Trueouter : 只按照一維數(shù)組進(jìn)行計(jì)算,如果傳入?yún)?shù)是多維數(shù)組,則先將此數(shù)組展平為一維數(shù)組之后再進(jìn)行運(yùn)算。outer乘積計(jì)算的列向量和行向量的矩陣乘積:
>>> np.outer([1,2,3],[4,5,6,7])array([[ 4, 5, 6, 7], [ 8, 10, 12, 14], [12, 15, 18, 21]])矩陣中更高級的一些運(yùn)算可以在NumPy的線性代數(shù)子庫linalg中找到。例如inv函數(shù)計(jì)算逆矩陣,solve函數(shù)可以求解多元一次方程組。下面是solve函數(shù)的一個(gè)例子:
>>> a = np.random.rand(10,10)>>> b = np.random.rand(10)>>> x = np.linalg.solve(a,b)>>> np.sum(np.abs(np.dot(a,x) - b))3.1433189384699745e-15solve函數(shù)有兩個(gè)參數(shù)a和b。a是一個(gè)N*N的二維數(shù)組,而b是一個(gè)長度為N的一維數(shù)組,solve函數(shù)找到一個(gè)長度為N的一維數(shù)組x,使得a和x的矩陣乘積正好等于b,數(shù)組x就是多元一次方程組的解。
有關(guān)線性代數(shù)方面的內(nèi)容將在今后的章節(jié)中詳細(xì)介紹。
2.4 文件存取
NumPy提供了多種文件操作函數(shù)方便我們存取數(shù)組內(nèi)容。文件存取的格式分為兩類:二進(jìn)制和文本。而二進(jìn)制格式的文件又分為NumPy專用的格式化二進(jìn)制類型和無格式類型。
使用數(shù)組的方法函數(shù)tofile可以方便地將數(shù)組中數(shù)據(jù)以二進(jìn)制的格式寫進(jìn)文件。tofile輸出的數(shù)據(jù)沒有格式,因此用numpy.fromfile讀回來的時(shí)候需要自己格式化數(shù)據(jù):
>>> a = np.arange(0,12)>>> a.shape = 3,4>>> aarray([[ 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類型讀入數(shù)據(jù)>>> b # 讀入的數(shù)據(jù)是錯(cuò)誤的array([ 2.12199579e-314, 6.36598737e-314, 1.06099790e-313, 1.48539705e-313, 1.90979621e-313, 2.33419537e-313])>>> a.dtype # 查看a的dtypedtype('int32')>>> b = np.fromfile("a.bin", dtype=np.int32) # 按照int32類型讀入數(shù)據(jù)>>> b # 數(shù)據(jù)是一維的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]])從上面的例子可以看出,需要在讀入的時(shí)候設(shè)置正確的dtype和shape才能保證數(shù)據(jù)一致。并且tofile函數(shù)不管數(shù)組的排列順序是C語言格式的還是Fortran語言格式的,統(tǒng)一使用C語言格式輸出。
此外如果fromfile和tofile函數(shù)調(diào)用時(shí)指定了sep關(guān)鍵字參數(shù)的話,數(shù)組將以文本格式輸入輸出。
numpy.load和numpy.save函數(shù)以NumPy專用的二進(jìn)制類型保存數(shù)據(jù),這兩個(gè)函數(shù)會自動(dòng)處理元素類型和shape等信息,使用它們讀寫數(shù)組就方便多了,但是numpy.save輸出的文件很難和其它語言編寫的程序讀入:
>>> np.save("a.npy", a)>>> c = np.load( "a.npy" )>>> carray([[ 0, 1, 2, 3], [ 4, 5, 6, 7], [ 8, 9, 10, 11]])如果你想將多個(gè)數(shù)組保存到一個(gè)文件中的話,可以使用numpy.savez函數(shù)。savez函數(shù)的第一個(gè)參數(shù)是文件名,其后的參數(shù)都是需要保存的數(shù)組,也可以使用關(guān)鍵字參數(shù)為數(shù)組起一個(gè)名字,非關(guān)鍵字參數(shù)傳遞的數(shù)組會自動(dòng)起名為arr_0, arr_1, ...。savez函數(shù)輸出的是一個(gè)壓縮文件(擴(kuò)展名為npz),其中每個(gè)文件都是一個(gè)save函數(shù)保存的npy文件,文件名對應(yīng)于數(shù)組名。load函數(shù)自動(dòng)識別npz文件,并且返回一個(gè)類似于字典的對象,可以通過數(shù)組名作為關(guān)鍵字獲取數(shù)組的內(nèi)容:
>>> 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"] # 數(shù)組aarray([[1, 2, 3], [4, 5, 6]])>>> r["arr_1"] # 數(shù)組barray([ 0. , 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9])>>> r["sin_array"] # 數(shù)組carray([ 0. , 0.09983342, 0.19866933, 0.29552021, 0.38941834, 0.47942554, 0.56464247, 0.64421769, 0.71735609, 0.78332691])如果你用解壓軟件打開result.npz文件的話,會發(fā)現(xiàn)其中有三個(gè)文件:arr_0.npy, arr_1.npy, sin_array.npy,其中分別保存著數(shù)組a, b, c的內(nèi)容。
使用numpy.savetxt和numpy.loadtxt可以讀寫1維和2維的數(shù)組:
>>> a = np.arange(0,12,0.5).reshape(4,-1)>>> np.savetxt("a.txt", a) # 缺省按照'%.18e'格式保存數(shù)據(jù),以空格分隔>>> 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=",") #改為保存為整數(shù),以逗號分隔>>> np.loadtxt("a.txt",delimiter=",") # 讀入的時(shí)候也需要指定逗號分隔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.]])文件名和文件對象
本節(jié)介紹所舉的例子都是傳遞的文件名,也可以傳遞已經(jīng)打開的文件對象,例如對于load和save函數(shù)來說,如果使用文件對象的話,可以將多個(gè)數(shù)組儲存到一個(gè)npy文件中:
>>> a = np.arange(8)>>> b = np.add.accumulate(a)>>> c = a + b>>> f = file("result.npy", "wb")>>> np.save(f, a) # 順序?qū),b,c保存進(jìn)文件對象f>>> np.save(f, b)>>> np.save(f, c)>>> f.close()>>> f = file("result.npy", "rb")>>> np.load(f) # 順序從文件對象f中讀取內(nèi)容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])
|
新聞熱點(diǎn)
疑難解答
圖片精選