資料來(lái)源：https://github.com/BrambleXu/pydata-notebook
在這篇附錄中，我會(huì)深入NumPy庫(kù)的數(shù)組計(jì)算。這會(huì)包括ndarray更內(nèi)部的細(xì)節(jié)，和更高級(jí)的數(shù)組操作和算法。
本章包括了一些雜亂的章節(jié)，不需要仔細(xì)研究。
A.1 ndarray對(duì)象的內(nèi)部機(jī)理
NumPy的ndarray提供了一種將同質(zhì)數(shù)據(jù)塊（可以是連續(xù)或跨越）解釋為多維數(shù)組對(duì)象的方式。正如你之前所看到的那樣，數(shù)據(jù)類(lèi)型（dtype）決定了數(shù)據(jù)的解釋方式，比如浮點(diǎn)數(shù)、整數(shù)、布爾值等。
ndarray如此強(qiáng)大的部分原因是所有數(shù)組對(duì)象都是數(shù)據(jù)塊的一個(gè)跨度視圖（strided view）。你可能想知道數(shù)組視圖arr[::2,::-1]不復(fù)制任何數(shù)據(jù)的原因是什么。簡(jiǎn)單地說(shuō)，ndarray不只是一塊內(nèi)存和一個(gè)dtype，它還有跨度信息，這使得數(shù)組能以各種步幅（step size）在內(nèi)存中移動(dòng)。更準(zhǔn)確地講，ndarray內(nèi)部由以下內(nèi)容組成：
- 一個(gè)指向數(shù)據(jù)（內(nèi)存或內(nèi)存映射文件中的一塊數(shù)據(jù)）的指針。
- 數(shù)據(jù)類(lèi)型或dtype，描述在數(shù)組中的固定大小值的格子。
- 一個(gè)表示數(shù)組形狀（shape）的元組。
- 一個(gè)跨度元組（stride），其中的整數(shù)指的是為了前進(jìn)到當(dāng)前維度下一個(gè)元素需要“跨過(guò)”的字節(jié)數(shù)。
圖A-1簡(jiǎn)單地說(shuō)明了ndarray的內(nèi)部結(jié)構(gòu)。
例如，一個(gè)10×5的數(shù)組，其形狀為(10,5)：
In [10]: np.ones((10, 5)).shape
Out[10]: (10, 5)
一個(gè)典型的（C順序，稍后將詳細(xì)講解）3×4×5的float64（8個(gè)字節(jié)）數(shù)組，其跨度為(160,40,8) —— 知道跨度是非常有用的，通常，跨度在一個(gè)軸上越大，沿這個(gè)軸進(jìn)行計(jì)算的開(kāi)銷(xiāo)就越大：
In [11]: np.ones((3, 4, 5), dtype=np.float64).strides
Out[11]: (160, 40, 8)
雖然NumPy用戶(hù)很少會(huì)對(duì)數(shù)組的跨度信息感興趣，但它們卻是構(gòu)建非復(fù)制式數(shù)組視圖的重要因素?？缍壬踔量梢允秦?fù)數(shù)，這樣會(huì)使數(shù)組在內(nèi)存中后向移動(dòng)，比如在切片obj[::-1]或obj[:,::-1]中就是這樣的。
NumPy數(shù)據(jù)類(lèi)型體系
你可能偶爾需要檢查數(shù)組中所包含的是否是整數(shù)、浮點(diǎn)數(shù)、字符串或Python對(duì)象。因?yàn)楦↑c(diǎn)數(shù)的種類(lèi)很多（從float16到float128），判斷dtype是否屬于某個(gè)大類(lèi)的工作非常繁瑣。幸運(yùn)的是，dtype都有一個(gè)超類(lèi)（比如np.integer和np.floating），它們可以跟np.issubdtype函數(shù)結(jié)合使用：
In [12]: ints = np.ones(10, dtype=np.uint16)
In [13]: floats = np.ones(10, dtype=np.float32)
In [14]: np.issubdtype(ints.dtype, np.integer)
Out[14]: True
In [15]: np.issubdtype(floats.dtype, np.floating)
Out[15]: True
調(diào)用dtype的mro方法即可查看其所有的父類(lèi)：
In [16]: np.float64.mro()
Out[16]:
[numpy.float64,
numpy.floating,
numpy.inexact,
numpy.number,
numpy.generic,
float,
object]
然后得到：
In [17]: np.issubdtype(ints.dtype, np.number)
Out[17]: True
大部分NumPy用戶(hù)完全不需要了解這些知識(shí)，但是這些知識(shí)偶爾還是能派上用場(chǎng)的。圖A-2說(shuō)明了dtype體系以及父子類(lèi)關(guān)系。
A.2 高級(jí)數(shù)組操作
除花式索引、切片、布爾條件取子集等操作之外，數(shù)組的操作方式還有很多。雖然pandas中的高級(jí)函數(shù)可以處理數(shù)據(jù)分析工作中的許多重型任務(wù)，但有時(shí)你還是需要編寫(xiě)一些在現(xiàn)有庫(kù)中找不到的數(shù)據(jù)算法。
數(shù)組重塑
多數(shù)情況下，你可以無(wú)需復(fù)制任何數(shù)據(jù)，就將數(shù)組從一個(gè)形狀轉(zhuǎn)換為另一個(gè)形狀。只需向數(shù)組的實(shí)例方法reshape傳入一個(gè)表示新形狀的元組即可實(shí)現(xiàn)該目的。例如，假設(shè)有一個(gè)一維數(shù)組，我們希望將其重新排列為一個(gè)矩陣（結(jié)果見(jiàn)圖A-3）：
In [18]: arr = np.arange(8)
In [19]: arr
Out[19]: array([0, 1, 2, 3, 4, 5, 6, 7])
In [20]: arr.reshape((4, 2))
Out[20]:
array([[0, 1],
[2, 3],
[4, 5],
[6, 7]])
多維數(shù)組也能被重塑：
In [21]: arr.reshape((4, 2)).reshape((2, 4))
Out[21]:
array([[0, 1, 2, 3],
[4, 5, 6, 7]])
作為參數(shù)的形狀的其中一維可以是－1，它表示該維度的大小由數(shù)據(jù)本身推斷而來(lái)：
In [22]: arr = np.arange(15)
In [23]: arr.reshape((5, -1))
Out[23]:
array([[ 0, 1, 2],
[ 3, 4, 5],
[ 6, 7, 8],
[ 9, 10, 11],
[12, 13, 14]])
與reshape將一維數(shù)組轉(zhuǎn)換為多維數(shù)組的運(yùn)算過(guò)程相反的運(yùn)算通常稱(chēng)為扁平化（flattening）或散開(kāi)（raveling）：
In [27]: arr = np.arange(15).reshape((5, 3))
In [28]: arr
Out[28]:
array([[ 0, 1, 2],
[ 3, 4, 5],
[ 6, 7, 8],
[ 9, 10, 11],
[12, 13, 14]])
In [29]: arr.ravel()
Out[29]: array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14])
如果結(jié)果中的值與原始數(shù)組相同，ravel不會(huì)產(chǎn)生源數(shù)據(jù)的副本。flatten方法的行為類(lèi)似于ravel，只不過(guò)它總是返回?cái)?shù)據(jù)的副本：
In [30]: arr.flatten()
Out[30]: array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14])
數(shù)組可以被重塑或散開(kāi)為別的順序。這對(duì)NumPy新手來(lái)說(shuō)是一個(gè)比較微妙的問(wèn)題，所以在下一小節(jié)中我們將專(zhuān)門(mén)講解這個(gè)問(wèn)題。
C和Fortran順序
NumPy允許你更為靈活地控制數(shù)據(jù)在內(nèi)存中的布局。默認(rèn)情況下，NumPy數(shù)組是按行優(yōu)先順序創(chuàng)建的。在空間方面，這就意味著，對(duì)于一個(gè)二維數(shù)組，每行中的數(shù)據(jù)項(xiàng)是被存放在相鄰內(nèi)存位置上的。另一種順序是列優(yōu)先順序，它意味著每列中的數(shù)據(jù)項(xiàng)是被存放在相鄰內(nèi)存位置上的。
由于一些歷史原因，行和列優(yōu)先順序又分別稱(chēng)為C和Fortran順序。在FORTRAN 77中，矩陣全都是列優(yōu)先的。
像reshape和reval這樣的函數(shù)，都可以接受一個(gè)表示數(shù)組數(shù)據(jù)存放順序的order參數(shù)。一般可以是'C'或'F'（還有'A'和'K'等不常用的選項(xiàng)，具體請(qǐng)參考NumPy的文檔）。圖A-3對(duì)此進(jìn)行了說(shuō)明：
In [31]: arr = np.arange(12).reshape((3, 4))
In [32]: arr
Out[32]:
array([[ 0, 1, 2, 3],
[ 4, 5, 6, 7],
[ 8, 9, 10, 11]])
In [33]: arr.ravel()
Out[33]: array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11])
In [34]: arr.ravel('F')
Out[34]: array([ 0, 4, 8, 1, 5, 9, 2, 6, 10, 3, 7, 11])
二維或更高維數(shù)組的重塑過(guò)程比較令人費(fèi)解（見(jiàn)圖A-3）。C和Fortran順序的關(guān)鍵區(qū)別就是維度的行進(jìn)順序：
- C/行優(yōu)先順序：先經(jīng)過(guò)更高的維度（例如，軸1會(huì)先于軸0被處理）。
- Fortran/列優(yōu)先順序：后經(jīng)過(guò)更高的維度（例如，軸0會(huì)先于軸1被處理）。
數(shù)組的合并和拆分
numpy.concatenate可以按指定軸將一個(gè)由數(shù)組組成的序列（如元組、列表等）連接到一起：
In [35]: arr1 = np.array([[1, 2, 3], [4, 5, 6]])
In [36]: arr2 = np.array([[7, 8, 9], [10, 11, 12]])
In [37]: np.concatenate([arr1, arr2], axis=0)
Out[37]:
array([[ 1, 2, 3],
[ 4, 5, 6],
[ 7, 8, 9],
[10, 11, 12]])
In [38]: np.concatenate([arr1, arr2], axis=1)
Out[38]:
array([[ 1, 2, 3, 7, 8, 9],
[ 4, 5, 6, 10, 11, 12]])
對(duì)于常見(jiàn)的連接操作，NumPy提供了一些比較方便的方法（如vstack和hstack）。因此，上面的運(yùn)算還可以表達(dá)為：
In [39]: np.vstack((arr1, arr2))
Out[39]:
array([[ 1, 2, 3],
[ 4, 5, 6],
[ 7, 8, 9],
[10, 11, 12]])
In [40]: np.hstack((arr1, arr2))
Out[40]:
array([[ 1, 2, 3, 7, 8, 9],
[ 4, 5, 6, 10, 11, 12]])
與此相反，split用于將一個(gè)數(shù)組沿指定軸拆分為多個(gè)數(shù)組：
In [41]: arr = np.random.randn(5, 2)
In [42]: arr
Out[42]:
array([[-0.2047, 0.4789],
[-0.5194, -0.5557],
[ 1.9658, 1.3934],
[ 0.0929, 0.2817],
[ 0.769 , 1.2464]])
In [43]: first, second, third = np.split(arr, [1, 3])
In [44]: first
Out[44]: array([[-0.2047, 0.4789]])
In [45]: second
Out[45]:
array([[-0.5194, -0.5557],
[ 1.9658, 1.3934]])
In [46]: third
Out[46]:
array([[ 0.0929, 0.2817],
[ 0.769 , 1.2464]])
傳入到np.split的值[1,3]指示在哪個(gè)索引處分割數(shù)組。
表A-1中列出了所有關(guān)于數(shù)組連接和拆分的函數(shù)，其中有些是專(zhuān)門(mén)為了方便常見(jiàn)的連接運(yùn)算而提供的。
堆疊輔助類(lèi)：r_和c_
NumPy命名空間中有兩個(gè)特殊的對(duì)象——r_和c_，它們可以使數(shù)組的堆疊操作更為簡(jiǎn)潔：
In [47]: arr = np.arange(6)
In [48]: arr1 = arr.reshape((3, 2))
In [49]: arr2 = np.random.randn(3, 2)
In [50]: np.r_[arr1, arr2]
Out[50]:
array([[ 0. , 1. ],
[ 2. , 3. ],
[ 4. , 5. ],
[ 1.0072, -1.2962],
[ 0.275 , 0.2289],
[ 1.3529, 0.8864]])
In [51]: np.c_[np.r_[arr1, arr2], arr]
Out[51]:
array([[ 0. , 1. , 0. ],
[ 2. , 3. , 1. ],
[ 4. , 5. , 2. ],
[ 1.0072, -1.2962, 3. ],
[ 0.275 , 0.2289, 4. ],
[ 1.3529, 0.8864, 5. ]])
它還可以將切片轉(zhuǎn)換成數(shù)組：
In [52]: np.c_[1:6, -10:-5]
Out[52]:
array([[ 1, -10],
[ 2, -9],
[ 3, -8],
[ 4, -7],
[ 5, -6]])
r_和c_的具體功能請(qǐng)參考其文檔。
元素的重復(fù)操作：tile和repeat
對(duì)數(shù)組進(jìn)行重復(fù)以產(chǎn)生更大數(shù)組的工具主要是repeat和tile這兩個(gè)函數(shù)。repeat會(huì)將數(shù)組中的各個(gè)元素重復(fù)一定次數(shù)，從而產(chǎn)生一個(gè)更大的數(shù)組：
In [53]: arr = np.arange(3)
In [54]: arr
Out[54]: array([0, 1, 2])
In [55]: arr.repeat(3)
Out[55]: array([0, 0, 0, 1, 1, 1, 2, 2, 2])
筆記：跟其他流行的數(shù)組編程語(yǔ)言（如MATLAB）不同，NumPy中很少需要對(duì)數(shù)組進(jìn)行重復(fù)（replicate）。這主要是因?yàn)閺V播（broadcasting，我們將在下一節(jié)中講解該技術(shù)）能更好地滿足該需求。
默認(rèn)情況下，如果傳入的是一個(gè)整數(shù)，則各元素就都會(huì)重復(fù)那么多次。如果傳入的是一組整數(shù)，則各元素就可以重復(fù)不同的次數(shù)：
In [56]: arr.repeat([2, 3, 4])
Out[56]: array([0, 0, 1, 1, 1, 2, 2, 2, 2])
對(duì)于多維數(shù)組，還可以讓它們的元素沿指定軸重復(fù)：
In [57]: arr = np.random.randn(2, 2)
In [58]: arr
Out[58]:
array([[-2.0016, -0.3718],
[ 1.669 , -0.4386]])
In [59]: arr.repeat(2, axis=0)
Out[59]:
array([[-2.0016, -0.3718],
[-2.0016, -0.3718],
[ 1.669 , -0.4386],
[ 1.669 , -0.4386]])
注意，如果沒(méi)有設(shè)置軸向，則數(shù)組會(huì)被扁平化，這可能不會(huì)是你想要的結(jié)果。同樣，在對(duì)多維進(jìn)行重復(fù)時(shí)，也可以傳入一組整數(shù)，這樣就會(huì)使各切片重復(fù)不同的次數(shù)：
In [60]: arr.repeat([2, 3], axis=0)
Out[60]:
array([[-2.0016, -0.3718],
[-2.0016, -0.3718],
[ 1.669 , -0.4386],
[ 1.669 , -0.4386],
[ 1.669 , -0.4386]])
In [61]: arr.repeat([2, 3], axis=1)
Out[61]:
array([[-2.0016, -2.0016, -0.3718, -0.3718, -0.3718],
[ 1.669 , 1.669 , -0.4386, -0.4386, -0.4386]])
tile的功能是沿指定軸向堆疊數(shù)組的副本。你可以形象地將其想象成“鋪瓷磚”：
In [62]: arr
Out[62]:
array([[-2.0016, -0.3718],
[ 1.669 , -0.4386]])
In [63]: np.tile(arr, 2)
Out[63]:
array([[-2.0016, -0.3718, -2.0016, -0.3718],
[ 1.669 , -0.4386, 1.669 , -0.4386]])
第二個(gè)參數(shù)是瓷磚的數(shù)量。對(duì)于標(biāo)量，瓷磚是水平鋪設(shè)的，而不是垂直鋪設(shè)。它可以是一個(gè)表示“鋪設(shè)”布局的元組：
In [64]: arr
Out[64]:
array([[-2.0016, -0.3718],
[ 1.669 , -0.4386]])
In [65]: np.tile(arr, (2, 1))
Out[65]:
array([[-2.0016, -0.3718],
[ 1.669 , -0.4386],
[-2.0016, -0.3718],
[ 1.669 , -0.4386]])
In [66]: np.tile(arr, (3, 2))
Out[66]:
array([[-2.0016, -0.3718, -2.0016, -0.3718],
[ 1.669 , -0.4386, 1.669 , -0.4386],
[-2.0016, -0.3718, -2.0016, -0.3718],
[ 1.669 , -0.4386, 1.669 , -0.4386],
[-2.0016, -0.3718, -2.0016, -0.3718],
[ 1.669 , -0.4386, 1.669 , -0.4386]])
花式索引的等價(jià)函數(shù)：take和put
在第4章中我們講過(guò)，獲取和設(shè)置數(shù)組子集的一個(gè)辦法是通過(guò)整數(shù)數(shù)組使用花式索引：
In [67]: arr = np.arange(10) * 100
In [68]: inds = [7, 1, 2, 6]
In [69]: arr[inds]
Out[69]: array([700, 100, 200, 600])
ndarray還有其它方法用于獲取單個(gè)軸向上的選區(qū)：
In [70]: arr.take(inds)
Out[70]: array([700, 100, 200, 600])
In [71]: arr.put(inds, 42)
In [72]: arr
Out[72]: array([ 0, 42, 42, 300, 400, 500, 42, 42,800, 900])
In [73]: arr.put(inds, [40, 41, 42, 43])
In [74]: arr
Out[74]: array([ 0, 41, 42, 300, 400, 500, 43, 40, 800, 900])
要在其它軸上使用take，只需傳入axis關(guān)鍵字即可：
In [75]: inds = [2, 0, 2, 1]
In [76]: arr = np.random.randn(2, 4)
In [77]: arr
Out[77]:
array([[-0.5397, 0.477 , 3.2489, -1.0212],
[-0.5771, 0.1241, 0.3026, 0.5238]])
In [78]: arr.take(inds, axis=1)
Out[78]:
array([[ 3.2489, -0.5397, 3.2489, 0.477 ],
[ 0.3026, -0.5771, 0.3026, 0.1241]])
put不接受axis參數(shù)，它只會(huì)在數(shù)組的扁平化版本（一維，C順序）上進(jìn)行索引。因此，在需要用其他軸向的索引設(shè)置元素時(shí)，最好還是使用花式索引。
A.3 廣播
廣播（broadcasting）指的是不同形狀的數(shù)組之間的算術(shù)運(yùn)算的執(zhí)行方式。它是一種非常強(qiáng)大的功能，但也容易令人誤解，即使是經(jīng)驗(yàn)豐富的老手也是如此。將標(biāo)量值跟數(shù)組合并時(shí)就會(huì)發(fā)生最簡(jiǎn)單的廣播：
In [79]: arr = np.arange(5)
In [80]: arr
Out[80]: array([0, 1, 2, 3, 4])
In [81]: arr * 4
Out[81]: array([ 0, 4, 8, 12, 16])
這里我們說(shuō)：在這個(gè)乘法運(yùn)算中，標(biāo)量值4被廣播到了其他所有的元素上。
看一個(gè)例子，我們可以通過(guò)減去列平均值的方式對(duì)數(shù)組的每一列進(jìn)行距平化處理。這個(gè)問(wèn)題解決起來(lái)非常簡(jiǎn)單：
In [82]: arr = np.random.randn(4, 3)
In [83]: arr.mean(0)
Out[83]: array([-0.3928, -0.3824, -0.8768])
In [84]: demeaned = arr - arr.mean(0)
In [85]: demeaned
Out[85]:
array([[ 0.3937, 1.7263, 0.1633],
[-0.4384, -1.9878, -0.9839],
[-0.468 , 0.9426, -0.3891],
[ 0.5126, -0.6811, 1.2097]])
In [86]: demeaned.mean(0)
Out[86]: array([-0., 0., -0.])
圖A-4形象地展示了該過(guò)程。用廣播的方式對(duì)行進(jìn)行距平化處理會(huì)稍微麻煩一些。幸運(yùn)的是，只要遵循一定的規(guī)則，低維度的值是可以被廣播到數(shù)組的任意維度的（比如對(duì)二維數(shù)組各列減去行平均值）。
于是就得到了：
雖然我是一名經(jīng)驗(yàn)豐富的NumPy老手，但經(jīng)常還是得停下來(lái)畫(huà)張圖并想想廣播的原則。再來(lái)看一下最后那個(gè)例子，假設(shè)你希望對(duì)各行減去那個(gè)平均值。由于arr.mean(0)的長(zhǎng)度為3，所以它可以在0軸向上進(jìn)行廣播：因?yàn)閍rr的后緣維度是3，所以它們是兼容的。根據(jù)該原則，要在1軸向上做減法（即各行減去行平均值），較小的那個(gè)數(shù)組的形狀必須是(4,1)：
In [87]: arr
Out[87]:
array([[ 0.0009, 1.3438, -0.7135],
[-0.8312, -2.3702, -1.8608],
[-0.8608, 0.5601, -1.2659],
[ 0.1198, -1.0635, 0.3329]])
In [88]: row_means = arr.mean(1)
In [89]: row_means.shape
Out[89]: (4,)
In [90]: row_means.reshape((4, 1))
Out[90]:
array([[ 0.2104],
[-1.6874],
[-0.5222],
[-0.2036]])
In [91]: demeaned = arr - row_means.reshape((4, 1))
In [92]: demeaned.mean(1)
Out[92]: array([ 0., -0., 0., 0.])
圖A-5說(shuō)明了該運(yùn)算的過(guò)程。
圖A-6展示了另外一種情況，這次是在一個(gè)三維數(shù)組上沿0軸向加上一個(gè)二維數(shù)組。
沿其它軸向廣播
高維度數(shù)組的廣播似乎更難以理解，而實(shí)際上它也是遵循廣播原則的。如果不然，你就會(huì)得到下面這樣一個(gè)錯(cuò)誤：
In [93]: arr - arr.mean(1)
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<ipython-input-93-7b87b85a20b2> in <module>()
----> 1 arr - arr.mean(1)
ValueError: operands could not be broadcast together with shapes (4,3) (4,)
人們經(jīng)常需要通過(guò)算術(shù)運(yùn)算過(guò)程將較低維度的數(shù)組在除0軸以外的其他軸向上廣播。根據(jù)廣播的原則，較小數(shù)組的“廣播維”必須為1。在上面那個(gè)行距平化的例子中，這就意味著要將行平均值的形狀變成(4,1)而不是(4,)：
In [94]: arr - arr.mean(1).reshape((4, 1))
Out[94]:
array([[-0.2095, 1.1334, -0.9239],
[ 0.8562, -0.6828, -0.1734],
[-0.3386, 1.0823, -0.7438],
[ 0.3234, -0.8599, 0.5365]])
對(duì)于三維的情況，在三維中的任何一維上廣播其實(shí)也就是將數(shù)據(jù)重塑為兼容的形狀而已。圖A-7說(shuō)明了要在三維數(shù)組各維度上廣播的形狀需求。
于是就有了一個(gè)非常普遍的問(wèn)題（尤其是在通用算法中），即專(zhuān)門(mén)為了廣播而添加一個(gè)長(zhǎng)度為1的新軸。雖然reshape是一個(gè)辦法，但插入軸需要構(gòu)造一個(gè)表示新形狀的元組。這是一個(gè)很郁悶的過(guò)程。因此，NumPy數(shù)組提供了一種通過(guò)索引機(jī)制插入軸的特殊語(yǔ)法。下面這段代碼通過(guò)特殊的np.newaxis屬性以及“全”切片來(lái)插入新軸：
In [95]: arr = np.zeros((4, 4))
In [96]: arr_3d = arr[:, np.newaxis, :]
In [97]: arr_3d.shape
Out[97]: (4, 1, 4)
In [98]: arr_1d = np.random.normal(size=3)
In [99]: arr_1d[:, np.newaxis]
Out[99]:
array([[-2.3594],
[-0.1995],
[-1.542 ]])
In [100]: arr_1d[np.newaxis, :]
Out[100]: array([[-2.3594, -0.1995, -1.542 ]])
因此，如果我們有一個(gè)三維數(shù)組，并希望對(duì)軸2進(jìn)行距平化，那么只需要編寫(xiě)下面這樣的代碼就可以了：
In [101]: arr = np.random.randn(3, 4, 5)
In [102]: depth_means = arr.mean(2)
In [103]: depth_means
Out[103]:
array([[-0.4735, 0.3971, -0.0228, 0.2001],
[-0.3521, -0.281 , -0.071 , -0.1586],
[ 0.6245, 0.6047, 0.4396, -0.2846]])
In [104]: depth_means.shape
Out[104]: (3, 4)
In [105]: demeaned = arr - depth_means[:, :, np.newaxis]
In [106]: demeaned.mean(2)
Out[106]:
array([[ 0., 0., -0., -0.],
[ 0., 0., -0., 0.],
[ 0., 0., -0., -0.]])
有些讀者可能會(huì)想，在對(duì)指定軸進(jìn)行距平化時(shí)，有沒(méi)有一種既通用又不犧牲性能的方法呢？實(shí)際上是有的，但需要一些索引方面的技巧：
def demean_axis(arr, axis=0):
means = arr.mean(axis)
# This generalizes things like [:, :, np.newaxis] to N dimensions
indexer = [slice(None)] * arr.ndim
indexer[axis] = np.newaxis
return arr - means[indexer]
通過(guò)廣播設(shè)置數(shù)組的值
算術(shù)運(yùn)算所遵循的廣播原則同樣也適用于通過(guò)索引機(jī)制設(shè)置數(shù)組值的操作。對(duì)于最簡(jiǎn)單的情況，我們可以這樣做：
In [107]: arr = np.zeros((4, 3))
In [108]: arr[:] = 5
In [109]: arr
Out[109]:
array([[ 5., 5., 5.],
[ 5., 5., 5.],
[ 5., 5., 5.],
[ 5., 5., 5.]])
但是，假設(shè)我們想要用一個(gè)一維數(shù)組來(lái)設(shè)置目標(biāo)數(shù)組的各列，只要保證形狀兼容就可以了：
In [110]: col = np.array([1.28, -0.42, 0.44, 1.6])
In [111]: arr[:] = col[:, np.newaxis]
In [112]: arr
Out[112]:
array([[ 1.28, 1.28, 1.28],
[-0.42, -0.42, -0.42],
[ 0.44, 0.44, 0.44],
[ 1.6 , 1.6 , 1.6 ]])
In [113]: arr[:2] = [[-1.37], [0.509]]
In [114]: arr
Out[114]:
array([[-1.37 , -1.37 , -1.37 ],
[ 0.509, 0.509, 0.509],
[ 0.44 , 0.44 , 0.44 ],
[ 1.6 , 1.6 , 1.6 ]])
A.4 ufunc高級(jí)應(yīng)用
雖然許多NumPy用戶(hù)只會(huì)用到通用函數(shù)所提供的快速的元素級(jí)運(yùn)算，但通用函數(shù)實(shí)際上還有一些高級(jí)用法能使我們丟開(kāi)循環(huán)而編寫(xiě)出更為簡(jiǎn)潔的代碼。
ufunc實(shí)例方法
NumPy的各個(gè)二元ufunc都有一些用于執(zhí)行特定矢量化運(yùn)算的特殊方法。表A-2匯總了這些方法，下面我將通過(guò)幾個(gè)具體的例子對(duì)它們進(jìn)行說(shuō)明。
reduce接受一個(gè)數(shù)組參數(shù)，并通過(guò)一系列的二元運(yùn)算對(duì)其值進(jìn)行聚合（可指明軸向）。例如，我們可以用np.add.reduce對(duì)數(shù)組中各個(gè)元素進(jìn)行求和：
In [115]: arr = np.arange(10)
In [116]: np.add.reduce(arr)
Out[116]: 45
In [117]: arr.sum()
Out[117]: 45
起始值取決于ufunc（對(duì)于add的情況，就是0）。如果設(shè)置了軸號(hào)，約簡(jiǎn)運(yùn)算就會(huì)沿該軸向執(zhí)行。這就使你能用一種比較簡(jiǎn)潔的方式得到某些問(wèn)題的答案。在下面這個(gè)例子中，我們用np.logical_and檢查數(shù)組各行中的值是否是有序的：
In [118]: np.random.seed(12346) # for reproducibility
In [119]: arr = np.random.randn(5, 5)
In [120]: arr[::2].sort(1) # sort a few rows
In [121]: arr[:, :-1] < arr[:, 1:]
Out[121]:
array([[ True, True, True, True],
[False, True, False, False],
[ True, True, True, True],
[ True, False, True, True],
[ True, True, True, True]], dtype=bool)
In [122]: np.logical_and.reduce(arr[:, :-1] < arr[:, 1:], axis=1)
Out[122]: array([ True, False, True, False, True], dtype=bool)
注意，logical_and.reduce跟all方法是等價(jià)的。
ccumulate跟reduce的關(guān)系就像cumsum跟sum的關(guān)系那樣。它產(chǎn)生一個(gè)跟原數(shù)組大小相同的中間“累計(jì)”值數(shù)組：
In [123]: arr = np.arange(15).reshape((3, 5))
In [124]: np.add.accumulate(arr, axis=1)
Out[124]:
array([[ 0, 1, 3, 6, 10],
[ 5, 11, 18, 26, 35],
[10, 21, 33, 46, 60]])
outer用于計(jì)算兩個(gè)數(shù)組的叉積：
In [125]: arr = np.arange(3).repeat([1, 2, 2])
In [126]: arr
Out[126]: array([0, 1, 1, 2, 2])
In [127]: np.multiply.outer(arr, np.arange(5))
Out[127]:
array([[0, 0, 0, 0, 0],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 2, 4, 6, 8],
[0, 2, 4, 6, 8]])
outer輸出結(jié)果的維度是兩個(gè)輸入數(shù)據(jù)的維度之和：
In [128]: x, y = np.random.randn(3, 4), np.random.randn(5)
In [129]: result = np.subtract.outer(x, y)
In [130]: result.shape
Out[130]: (3, 4, 5)
最后一個(gè)方法reduceat用于計(jì)算“局部約簡(jiǎn)”，其實(shí)就是一個(gè)對(duì)數(shù)據(jù)各切片進(jìn)行聚合的groupby運(yùn)算。它接受一組用于指示如何對(duì)值進(jìn)行拆分和聚合的“面元邊界”：
In [131]: arr = np.arange(10)
In [132]: np.add.reduceat(arr, [0, 5, 8])
Out[132]: array([10, 18, 17])
最終結(jié)果是在arr[0:5]、arr[5:8]以及arr[8:]上執(zhí)行的約簡(jiǎn)。跟其他方法一樣，這里也可以傳入一個(gè)axis參數(shù)：
In [133]: arr = np.multiply.outer(np.arange(4), np.arange(5))
In [134]: arr
Out[134]:
array([[ 0, 0, 0, 0, 0],
[ 0, 1, 2, 3, 4],
[ 0, 2, 4, 6, 8],
[ 0, 3, 6, 9, 12]])
In [135]: np.add.reduceat(arr, [0, 2, 4], axis=1)
Out[135]:
array([[ 0, 0, 0],
[ 1, 5, 4],
[ 2, 10, 8],
[ 3, 15, 12]])
表A-2總結(jié)了部分的ufunc方法。
編寫(xiě)新的ufunc
有多種方法可以讓你編寫(xiě)自己的NumPy ufuncs。最常見(jiàn)的是使用NumPy C API，但它超越了本書(shū)的范圍。在本節(jié)，我們講純粹的Python ufunc。
numpy.frompyfunc接受一個(gè)Python函數(shù)以及兩個(gè)分別表示輸入輸出參數(shù)數(shù)量的參數(shù)。例如，下面是一個(gè)能夠?qū)崿F(xiàn)元素級(jí)加法的簡(jiǎn)單函數(shù)：
In [136]: def add_elements(x, y):
.....: return x + y
In [137]: add_them = np.frompyfunc(add_elements, 2, 1)
In [138]: add_them(np.arange(8), np.arange(8))
Out[138]: array([0, 2, 4, 6, 8, 10, 12, 14], dtype=object)
用frompyfunc創(chuàng)建的函數(shù)總是返回Python對(duì)象數(shù)組，這一點(diǎn)很不方便。幸運(yùn)的是，還有另一個(gè)辦法，即numpy.vectorize。雖然沒(méi)有frompyfunc那么強(qiáng)大，但可以讓你指定輸出類(lèi)型：
In [139]: add_them = np.vectorize(add_elements, otypes=[np.float64])
In [140]: add_them(np.arange(8), np.arange(8))
Out[140]: array([ 0., 2., 4., 6., 8., 10., 12., 14.])
雖然這兩個(gè)函數(shù)提供了一種創(chuàng)建ufunc型函數(shù)的手段，但它們非常慢，因?yàn)樗鼈冊(cè)谟?jì)算每個(gè)元素時(shí)都要執(zhí)行一次Python函數(shù)調(diào)用，這就會(huì)比NumPy自帶的基于C的ufunc慢很多：
In [141]: arr = np.random.randn(10000)
In [142]: %timeit add_them(arr, arr)
4.12 ms +- 182 us per loop (mean +- std. dev. of 7 runs, 100 loops each)
In [143]: %timeit np.add(arr, arr)
6.89 us +- 504 ns per loop (mean +- std. dev. of 7 runs, 100000 loops each)
本章的后面，我會(huì)介紹使用Numba（http://numba.pydata.org/），創(chuàng)建快速Python ufuncs。
A.5 結(jié)構(gòu)化和記錄式數(shù)組
你可能已經(jīng)注意到了，到目前為止我們所討論的ndarray都是一種同質(zhì)數(shù)據(jù)容器，也就是說(shuō)，在它所表示的內(nèi)存塊中，各元素占用的字節(jié)數(shù)相同（具體根據(jù)dtype而定）。從表面上看，它似乎不能用于表示異質(zhì)或表格型的數(shù)據(jù)。結(jié)構(gòu)化數(shù)組是一種特殊的ndarray，其中的各個(gè)元素可以被看做C語(yǔ)言中的結(jié)構(gòu)體（struct，這就是“結(jié)構(gòu)化”的由來(lái)）或SQL表中帶有多個(gè)命名字段的行：
In [144]: dtype = [('x', np.float64), ('y', np.int32)]
In [145]: sarr = np.array([(1.5, 6), (np.pi, -2)], dtype=dtype)
In [146]: sarr
Out[146]:
array([( 1.5 , 6), ( 3.1416, -2)],
dtype=[('x', '<f8'), ('y', '<i4')])
定義結(jié)構(gòu)化dtype（請(qǐng)參考NumPy的在線文檔）的方式有很多。最典型的辦法是元組列表，各元組的格式為(field_name,field_data_type)。這樣，數(shù)組的元素就成了元組式的對(duì)象，該對(duì)象中各個(gè)元素可以像字典那樣進(jìn)行訪問(wèn)：
In [147]: sarr[0]
Out[147]: ( 1.5, 6)
In [148]: sarr[0]['y']
Out[148]: 6
字段名保存在dtype.names屬性中。在訪問(wèn)結(jié)構(gòu)化數(shù)組的某個(gè)字段時(shí)，返回的是該數(shù)據(jù)的視圖，所以不會(huì)發(fā)生數(shù)據(jù)復(fù)制：
In [149]: sarr['x']
Out[149]: array([ 1.5 , 3.1416])
嵌套dtype和多維字段
在定義結(jié)構(gòu)化dtype時(shí)，你可以再設(shè)置一個(gè)形狀（可以是一個(gè)整數(shù)，也可以是一個(gè)元組）：
In [150]: dtype = [('x', np.int64, 3), ('y', np.int32)]
In [151]: arr = np.zeros(4, dtype=dtype)
In [152]: arr
Out[152]:
array([([0, 0, 0], 0), ([0, 0, 0], 0), ([0, 0, 0], 0), ([0, 0, 0], 0)],
dtype=[('x', '<i8', (3,)), ('y', '<i4')])
在這種情況下，各個(gè)記錄的x字段所表示的是一個(gè)長(zhǎng)度為3的數(shù)組：
In [153]: arr[0]['x']
Out[153]: array([0, 0, 0])
這樣，訪問(wèn)arr['x']即可得到一個(gè)二維數(shù)組，而不是前面那個(gè)例子中的一維數(shù)組：
In [154]: arr['x']
Out[154]:
array([[0, 0, 0],
[0, 0, 0],
[0, 0, 0],
[0, 0, 0]])
這就使你能用單個(gè)數(shù)組的內(nèi)存塊存放復(fù)雜的嵌套結(jié)構(gòu)。你還可以嵌套dtype，作出更復(fù)雜的結(jié)構(gòu)。下面是一個(gè)簡(jiǎn)單的例子：
In [155]: dtype = [('x', [('a', 'f8'), ('b', 'f4')]), ('y', np.int32)]
In [156]: data = np.array([((1, 2), 5), ((3, 4), 6)], dtype=dtype)
In [157]: data['x']
Out[157]:
array([( 1., 2.), ( 3., 4.)],
dtype=[('a', '<f8'), ('b', '<f4')])
In [158]: data['y']
Out[158]: array([5, 6], dtype=int32)
In [159]: data['x']['a']
Out[159]: array([ 1., 3.])
pandas的DataFrame并不直接支持該功能，但它的分層索引機(jī)制跟這個(gè)差不多。
為什么要用結(jié)構(gòu)化數(shù)組
跟pandas的DataFrame相比，NumPy的結(jié)構(gòu)化數(shù)組是一種相對(duì)較低級(jí)的工具。它可以將單個(gè)內(nèi)存塊解釋為帶有任意復(fù)雜嵌套列的表格型結(jié)構(gòu)。由于數(shù)組中的每個(gè)元素在內(nèi)存中都被表示為固定的字節(jié)數(shù)，所以結(jié)構(gòu)化數(shù)組能夠提供非?？焖俑咝У拇疟P(pán)數(shù)據(jù)讀寫(xiě)（包括內(nèi)存映像）、網(wǎng)絡(luò)傳輸?shù)裙δ堋?/p>
結(jié)構(gòu)化數(shù)組的另一個(gè)常見(jiàn)用法是，將數(shù)據(jù)文件寫(xiě)成定長(zhǎng)記錄字節(jié)流，這是C和C++代碼中常見(jiàn)的數(shù)據(jù)序列化手段（業(yè)界許多歷史系統(tǒng)中都能找得到）。只要知道文件的格式（記錄的大小、元素的順序、字節(jié)數(shù)以及數(shù)據(jù)類(lèi)型等），就可以用np.fromfile將數(shù)據(jù)讀入內(nèi)存。這種用法超出了本書(shū)的范圍，知道這點(diǎn)就可以了。
A.6 更多有關(guān)排序的話題
跟Python內(nèi)置的列表一樣，ndarray的sort實(shí)例方法也是就地排序。也就是說(shuō)，數(shù)組內(nèi)容的重新排列是不會(huì)產(chǎn)生新數(shù)組的：
In [160]: arr = np.random.randn(6)
In [161]: arr.sort()
In [162]: arr
Out[162]: array([-1.082 , 0.3759, 0.8014, 1.1397, 1.2888, 1.8413])
在對(duì)數(shù)組進(jìn)行就地排序時(shí)要注意一點(diǎn)，如果目標(biāo)數(shù)組只是一個(gè)視圖，則原始數(shù)組將會(huì)被修改：
In [163]: arr = np.random.randn(3, 5)
In [164]: arr
Out[164]:
array([[-0.3318, -1.4711, 0.8705, -0.0847, -1.1329],
[-1.0111, -0.3436, 2.1714, 0.1234, -0.0189],
[ 0.1773, 0.7424, 0.8548, 1.038 , -0.329 ]])
In [165]: arr[:, 0].sort() # Sort first column values in-place
In [166]: arr
Out[166]:
array([[-1.0111, -1.4711, 0.8705, -0.0847, -1.1329],
[-0.3318, -0.3436, 2.1714, 0.1234, -0.0189],
[ 0.1773, 0.7424, 0.8548, 1.038 , -0.329 ]])
相反，numpy.sort會(huì)為原數(shù)組創(chuàng)建一個(gè)已排序副本。另外，它所接受的參數(shù)（比如kind）跟ndarray.sort一樣：
In [167]: arr = np.random.randn(5)
In [168]: arr
Out[168]: array([-1.1181, -0.2415, -2.0051, 0.7379, -1.0614])
In [169]: np.sort(arr)
Out[169]: array([-2.0051, -1.1181, -1.0614, -0.2415, 0.7379])
In [170]: arr
Out[170]: array([-1.1181, -0.2415, -2.0051, 0.7379, -1.0614])
這兩個(gè)排序方法都可以接受一個(gè)axis參數(shù)，以便沿指定軸向?qū)Ω鲏K數(shù)據(jù)進(jìn)行單獨(dú)排序：
In [171]: arr = np.random.randn(3, 5)
In [172]: arr
Out[172]:
array([[ 0.5955, -0.2682, 1.3389, -0.1872, 0.9111],
[-0.3215, 1.0054, -0.5168, 1.1925, -0.1989],
[ 0.3969, -1.7638, 0.6071, -0.2222, -0.2171]])
In [173]: arr.sort(axis=1)
In [174]: arr
Out[174]:
array([[-0.2682, -0.1872, 0.5955, 0.9111, 1.3389],
[-0.5168, -0.3215, -0.1989, 1.0054, 1.1925],
[-1.7638, -0.2222, -0.2171, 0.3969, 0.6071]])
你可能注意到了，這兩個(gè)排序方法都不可以被設(shè)置為降序。其實(shí)這也無(wú)所謂，因?yàn)閿?shù)組切片會(huì)產(chǎn)生視圖（也就是說(shuō)，不會(huì)產(chǎn)生副本，也不需要任何其他的計(jì)算工作）。許多Python用戶(hù)都很熟悉一個(gè)有關(guān)列表的小技巧：values[::-1]可以返回一個(gè)反序的列表。對(duì)ndarray也是如此：
In [175]: arr[:, ::-1]
Out[175]:
array([[ 1.3389, 0.9111, 0.5955, -0.1872, -0.2682],
[ 1.1925, 1.0054, -0.1989, -0.3215, -0.5168],
[ 0.6071, 0.3969, -0.2171, -0.2222, -1.7638]])
間接排序：argsort和lexsort
在數(shù)據(jù)分析工作中，常常需要根據(jù)一個(gè)或多個(gè)鍵對(duì)數(shù)據(jù)集進(jìn)行排序。例如，一個(gè)有關(guān)學(xué)生信息的數(shù)據(jù)表可能需要以姓和名進(jìn)行排序（先姓后名）。這就是間接排序的一個(gè)例子，如果你閱讀過(guò)有關(guān)pandas的章節(jié)，那就已經(jīng)見(jiàn)過(guò)不少高級(jí)例子了。給定一個(gè)或多個(gè)鍵，你就可以得到一個(gè)由整數(shù)組成的索引數(shù)組（我親切地稱(chēng)之為索引器），其中的索引值說(shuō)明了數(shù)據(jù)在新順序下的位置。argsort和numpy.lexsort就是實(shí)現(xiàn)該功能的兩個(gè)主要方法。下面是一個(gè)簡(jiǎn)單的例子：
In [176]: values = np.array([5, 0, 1, 3, 2])
In [177]: indexer = values.argsort()
In [178]: indexer
Out[178]: array([1, 2, 4, 3, 0])
In [179]: values[indexer]
Out[179]: array([0, 1, 2, 3, 5])
一個(gè)更復(fù)雜的例子，下面這段代碼根據(jù)數(shù)組的第一行對(duì)其進(jìn)行排序：
In [180]: arr = np.random.randn(3, 5)
In [181]: arr[0] = values
In [182]: arr
Out[182]:
array([[ 5. , 0. , 1. , 3. , 2. ],
[-0.3636, -0.1378, 2.1777, -0.4728, 0.8356],
[-0.2089, 0.2316, 0.728 , -1.3918, 1.9956]])
In [183]: arr[:, arr[0].argsort()]
Out[183]:
array([[ 0. , 1. , 2. , 3. , 5. ],
[-0.1378, 2.1777, 0.8356, -0.4728, -0.3636],
[ 0.2316, 0.728 , 1.9956, -1.3918, -0.2089]])
lexsort跟argsort差不多，只不過(guò)它可以一次性對(duì)多個(gè)鍵數(shù)組執(zhí)行間接排序（字典序）。假設(shè)我們想對(duì)一些以姓和名標(biāo)識(shí)的數(shù)據(jù)進(jìn)行排序：
In [184]: first_name = np.array(['Bob', 'Jane', 'Steve', 'Bill', 'Barbara'])
In [185]: last_name = np.array(['Jones', 'Arnold', 'Arnold', 'Jones', 'Walters'])
In [186]: sorter = np.lexsort((first_name, last_name))
In [187]: sorter
Out[187]: array([1, 2, 3, 0, 4])
In [188]: zip(last_name[sorter], first_name[sorter])
Out[188]: <zip at 0x7fa203eda1c8>
剛開(kāi)始使用lexsort的時(shí)候可能會(huì)比較容易頭暈，這是因?yàn)殒I的應(yīng)用順序是從最后一個(gè)傳入的算起的。不難看出，last_name是先于first_name被應(yīng)用的。
筆記：Series和DataFrame的sort_index以及Series的order方法就是通過(guò)這些函數(shù)的變體（它們還必須考慮缺失值）實(shí)現(xiàn)的。
其他排序算法
穩(wěn)定的（stable）排序算法會(huì)保持等價(jià)元素的相對(duì)位置。對(duì)于相對(duì)位置具有實(shí)際意義的那些間接排序而言，這一點(diǎn)非常重要：
In [189]: values = np.array(['2:first', '2:second', '1:first', '1:second',
.....: '1:third'])
In [190]: key = np.array([2, 2, 1, 1, 1])
In [191]: indexer = key.argsort(kind='mergesort')
In [192]: indexer
Out[192]: array([2, 3, 4, 0, 1])
In [193]: values.take(indexer)
Out[193]:
array(['1:first', '1:second', '1:third', '2:first', '2:second'],
dtype='<U8')
mergesort（合并排序）是唯一的穩(wěn)定排序，它保證有O(n log n)的性能（空間復(fù)雜度），但是其平均性能比默認(rèn)的quicksort（快速排序）要差。表A-3列出了可用的排序算法及其相關(guān)的性能指標(biāo)。大部分用戶(hù)完全不需要知道這些東西，但了解一下總是好的。
部分排序數(shù)組
排序的目的之一可能是確定數(shù)組中最大或最小的元素。NumPy有兩個(gè)優(yōu)化方法，numpy.partition和np.argpartition，可以在第k個(gè)最小元素劃分的數(shù)組：
In [194]: np.random.seed(12345)
In [195]: arr = np.random.randn(20)
In [196]: arr
Out[196]:
array([-0.2047, 0.4789, -0.5194, -0.5557, 1.9658, 1.3934, 0.0929,
0.2817, 0.769 , 1.2464, 1.0072, -1.2962, 0.275 , 0.2289,
1.3529, 0.8864, -2.0016, -0.3718, 1.669 , -0.4386])
In [197]: np.partition(arr, 3)
Out[197]:
array([-2.0016, -1.2962, -0.5557, -0.5194, -0.3718, -0.4386, -0.2047,
0.2817, 0.769 , 0.4789, 1.0072, 0.0929, 0.275 , 0.2289,
1.3529, 0.8864, 1.3934, 1.9658, 1.669 , 1.2464])
當(dāng)你調(diào)用partition(arr, 3)，結(jié)果中的頭三個(gè)元素是最小的三個(gè)，沒(méi)有特定的順序。numpy.argpartition與numpy.argsort相似，會(huì)返回索引，重排數(shù)據(jù)為等價(jià)的順序：
In [198]: indices = np.argpartition(arr, 3)
In [199]: indices
Out[199]:
array([16, 11, 3, 2, 17, 19, 0, 7, 8, 1, 10, 6, 12, 13, 14, 15, 5,
4, 18, 9])
In [200]: arr.take(indices)
Out[200]:
array([-2.0016, -1.2962, -0.5557, -0.5194, -0.3718, -0.4386, -0.2047,
0.2817, 0.769 , 0.4789, 1.0072, 0.0929, 0.275 , 0.2289,
1.3529, 0.8864, 1.3934, 1.9658, 1.669 , 1.2464])
numpy.searchsorted：在有序數(shù)組中查找元素
searchsorted是一個(gè)在有序數(shù)組上執(zhí)行二分查找的數(shù)組方法，只要將值插入到它返回的那個(gè)位置就能維持?jǐn)?shù)組的有序性：
In [201]: arr = np.array([0, 1, 7, 12, 15])
In [202]: arr.searchsorted(9)
Out[202]: 3
你可以傳入一組值就能得到一組索引：
In [203]: arr.searchsorted([0, 8, 11, 16])
Out[203]: array([0, 3, 3, 5])
從上面的結(jié)果中可以看出，對(duì)于元素0，searchsorted會(huì)返回0。這是因?yàn)槠淠J(rèn)行為是返回相等值組的左側(cè)索引：
In [204]: arr = np.array([0, 0, 0, 1, 1, 1, 1])
In [205]: arr.searchsorted([0, 1])
Out[205]: array([0, 3])
In [206]: arr.searchsorted([0, 1], side='right')
Out[206]: array([3, 7])
再來(lái)看searchsorted的另一個(gè)用法，假設(shè)我們有一個(gè)數(shù)據(jù)數(shù)組（其中的值在0到10000之間），還有一個(gè)表示“面元邊界”的數(shù)組，我們希望用它將數(shù)據(jù)數(shù)組拆分開(kāi)：
In [207]: data = np.floor(np.random.uniform(0, 10000, size=50))
In [208]: bins = np.array([0, 100, 1000, 5000, 10000])
In [209]: data
Out[209]:
array([ 9940., 6768., 7908., 1709., 268., 8003., 9037., 246.,
4917., 5262., 5963., 519., 8950., 7282., 8183., 5002.,
8101., 959., 2189., 2587., 4681., 4593., 7095., 1780.,
5314., 1677., 7688., 9281., 6094., 1501., 4896., 3773.,
8486., 9110., 3838., 3154., 5683., 1878., 1258., 6875.,
7996., 5735., 9732., 6340., 8884., 4954., 3516., 7142.,
5039., 2256.])
然后，為了得到各數(shù)據(jù)點(diǎn)所屬區(qū)間的編號(hào)（其中1表示面元[0,100)），我們可以直接使用searchsorted：
In [210]: labels = bins.searchsorted(data)
In [211]: labels
Out[211]:
array([4, 4, 4, 3, 2, 4, 4, 2, 3, 4, 4, 2, 4, 4, 4, 4, 4, 2, 3, 3, 3, 3, 4,
3, 4, 3, 4, 4, 4, 3, 3, 3, 4, 4, 3, 3, 4, 3, 3, 4, 4, 4, 4, 4, 4, 3,
3, 4, 4, 3])
通過(guò)pandas的groupby使用該結(jié)果即可非常輕松地對(duì)原數(shù)據(jù)集進(jìn)行拆分：
In [212]: pd.Series(data).groupby(labels).mean()
Out[212]:
2 498.000000
3 3064.277778
4 7389.035714
dtype: float64
A.7 用Numba編寫(xiě)快速NumPy函數(shù)
Numba是一個(gè)開(kāi)源項(xiàng)目，它可以利用CPUs、GPUs或其它硬件為類(lèi)似NumPy的數(shù)據(jù)創(chuàng)建快速函數(shù)。它使用了LLVM項(xiàng)目（http://llvm.org/），將Python代碼轉(zhuǎn)換為機(jī)器代碼。
為了介紹Numba，來(lái)考慮一個(gè)純粹的Python函數(shù)，它使用for循環(huán)計(jì)算表達(dá)式(x - y).mean()：
import numpy as np
def mean_distance(x, y):
nx = len(x)
result = 0.0
count = 0
for i in range(nx):
result += x[i] - y[i]
count += 1
return result / count
這個(gè)函數(shù)很慢：
In [209]: x = np.random.randn(10000000)
In [210]: y = np.random.randn(10000000)
In [211]: %timeit mean_distance(x, y)
1 loop, best of 3: 2 s per loop
In [212]: %timeit (x - y).mean()
100 loops, best of 3: 14.7 ms per loop
NumPy的版本要比它快過(guò)100倍。我們可以轉(zhuǎn)換這個(gè)函數(shù)為編譯的Numba函數(shù)，使用numba.jit函數(shù)：
In [213]: import numba as nb
In [214]: numba_mean_distance = nb.jit(mean_distance)
也可以寫(xiě)成裝飾器：
@nb.jit
def mean_distance(x, y):
nx = len(x)
result = 0.0
count = 0
for i in range(nx):
result += x[i] - y[i]
count += 1
return result / count
它要比矢量化的NumPy快：
In [215]: %timeit numba_mean_distance(x, y)
100 loops, best of 3: 10.3 ms per loop
Numba不能編譯Python代碼，但它支持純Python寫(xiě)的一個(gè)部分，可以編寫(xiě)數(shù)值算法。
Numba是一個(gè)深厚的庫(kù)，支持多種硬件、編譯模式和用戶(hù)插件。它還可以編譯NumPy Python API的一部分，而不用for循環(huán)。Numba也可以識(shí)別可以便以為機(jī)器編碼的結(jié)構(gòu)體，但是若調(diào)用CPython API，它就不知道如何編譯。Numba的jit函數(shù)有一個(gè)選項(xiàng)，nopython=True，它限制了可以被轉(zhuǎn)換為Python代碼的代碼，這些代碼可以編譯為L(zhǎng)LVM，但沒(méi)有任何Python C API調(diào)用。jit(nopython=True)有一個(gè)簡(jiǎn)短的別名numba.njit。
前面的例子，我們還可以這樣寫(xiě)：
from numba import float64, njit
@njit(float64(float64[:], float64[:]))
def mean_distance(x, y):
return (x - y).mean()
我建議你學(xué)習(xí)Numba的線上文檔（http://numba.pydata.org/）。下一節(jié)介紹一個(gè)創(chuàng)建自定義Numpy ufunc對(duì)象的例子。
用Numba創(chuàng)建自定義numpy.ufunc對(duì)象
numba.vectorize創(chuàng)建了一個(gè)編譯的NumPy ufunc，它與內(nèi)置的ufunc很像?？紤]一個(gè)numpy.add的Python例子：
from numba import vectorize
@vectorize
def nb_add(x, y):
return x + y
現(xiàn)在有：
In [13]: x = np.arange(10)
In [14]: nb_add(x, x)
Out[14]: array([ 0., 2., 4., 6., 8., 10., 12., 14., 16., 18.])
In [15]: nb_add.accumulate(x, 0)
Out[15]: array([ 0., 1., 3., 6., 10., 15., 21., 28., 36., 45.])
A.8 高級(jí)數(shù)組輸入輸出
我在第4章中講過(guò)，np.save和np.load可用于讀寫(xiě)磁盤(pán)上以二進(jìn)制格式存儲(chǔ)的數(shù)組。其實(shí)還有一些工具可用于更為復(fù)雜的場(chǎng)景。尤其是內(nèi)存映像（memory map），它使你能處理在內(nèi)存中放不下的數(shù)據(jù)集。
內(nèi)存映像文件
內(nèi)存映像文件是一種將磁盤(pán)上的非常大的二進(jìn)制數(shù)據(jù)文件當(dāng)做內(nèi)存中的數(shù)組進(jìn)行處理的方式。NumPy實(shí)現(xiàn)了一個(gè)類(lèi)似于ndarray的memmap對(duì)象，它允許將大文件分成小段進(jìn)行讀寫(xiě)，而不是一次性將整個(gè)數(shù)組讀入內(nèi)存。另外，memmap也擁有跟普通數(shù)組一樣的方法，因此，基本上只要是能用于ndarray的算法就也能用于memmap。
要?jiǎng)?chuàng)建一個(gè)內(nèi)存映像，可以使用函數(shù)np.memmap并傳入一個(gè)文件路徑、數(shù)據(jù)類(lèi)型、形狀以及文件模式：
In [214]: mmap = np.memmap('mymmap', dtype='float64', mode='w+',
.....: shape=(10000, 10000))
In [215]: mmap
Out[215]:
memmap([[ 0., 0., 0., ..., 0., 0., 0.],
[ 0., 0., 0., ..., 0., 0., 0.],
[ 0., 0., 0., ..., 0., 0., 0.],
...,
[ 0., 0., 0., ..., 0., 0., 0.],
[ 0., 0., 0., ..., 0., 0., 0.],
[ 0., 0., 0., ..., 0., 0., 0.]])
對(duì)memmap切片將會(huì)返回磁盤(pán)上的數(shù)據(jù)的視圖：
In [216]: section = mmap[:5]
如果將數(shù)據(jù)賦值給這些視圖：數(shù)據(jù)會(huì)先被緩存在內(nèi)存中（就像是Python的文件對(duì)象），調(diào)用flush即可將其寫(xiě)入磁盤(pán)：
In [217]: section[:] = np.random.randn(5, 10000)
In [218]: mmap.flush()
In [219]: mmap
Out[219]:
memmap([[ 0.7584, -0.6605, 0.8626, ..., 0.6046, -0.6212, 2.0542],
[-1.2113, -1.0375, 0.7093, ..., -1.4117, -0.1719, -0.8957],
[-0.1419, -0.3375, 0.4329, ..., 1.2914, -0.752 , -0.44 ],
...,
[ 0. , 0. , 0. , ..., 0. , 0. , 0. ],
[ 0. , 0. , 0. , ..., 0. , 0. , 0. ],
[ 0. , 0. , 0. , ..., 0. , 0. , 0. ]])
In [220]: del mmap
只要某個(gè)內(nèi)存映像超出了作用域，它就會(huì)被垃圾回收器回收，之前對(duì)其所做的任何修改都會(huì)被寫(xiě)入磁盤(pán)。當(dāng)打開(kāi)一個(gè)已經(jīng)存在的內(nèi)存映像時(shí)，仍然需要指明數(shù)據(jù)類(lèi)型和形狀，因?yàn)榇疟P(pán)上的那個(gè)文件只是一塊二進(jìn)制數(shù)據(jù)而已，沒(méi)有任何元數(shù)據(jù)：
In [221]: mmap = np.memmap('mymmap', dtype='float64', shape=(10000, 10000))
In [222]: mmap
Out[222]:
memmap([[ 0.7584, -0.6605, 0.8626, ..., 0.6046, -0.6212, 2.0542],
[-1.2113, -1.0375, 0.7093, ..., -1.4117, -0.1719, -0.8957],
[-0.1419, -0.3375, 0.4329, ..., 1.2914, -0.752 , -0.44 ],
...,
[ 0. , 0. , 0. , ..., 0. , 0. , 0. ],
[ 0. , 0. , 0. , ..., 0. , 0. , 0. ],
[ 0. , 0. , 0. , ..., 0. , 0. , 0. ]])
內(nèi)存映像可以使用前面介紹的結(jié)構(gòu)化或嵌套dtype。
HDF5及其他數(shù)組存儲(chǔ)方式
PyTables和h5py這兩個(gè)Python項(xiàng)目可以將NumPy的數(shù)組數(shù)據(jù)存儲(chǔ)為高效且可壓縮的HDF5格式（HDF意思是“層次化數(shù)據(jù)格式”）。你可以安全地將好幾百GB甚至TB的數(shù)據(jù)存儲(chǔ)為HDF5格式。要學(xué)習(xí)Python使用HDF5，請(qǐng)參考pandas線上文檔。
A.9 性能建議
使用NumPy的代碼的性能一般都很不錯(cuò)，因?yàn)閿?shù)組運(yùn)算一般都比純Python循環(huán)快得多。下面大致列出了一些需要注意的事項(xiàng)：
- 將Python循環(huán)和條件邏輯轉(zhuǎn)換為數(shù)組運(yùn)算和布爾數(shù)組運(yùn)算。
- 盡量使用廣播。
- 避免復(fù)制數(shù)據(jù)，盡量使用數(shù)組視圖（即切片）。
- 利用ufunc及其各種方法。
如果單用NumPy無(wú)論如何都達(dá)不到所需的性能指標(biāo)，就可以考慮一下用C、Fortran或Cython（等下會(huì)稍微介紹一下）來(lái)編寫(xiě)代碼。我自己在工作中經(jīng)常會(huì)用到Cython（http://cython.org），因?yàn)樗挥没ㄙM(fèi)我太多精力就能得到C語(yǔ)言那樣的性能。
連續(xù)內(nèi)存的重要性
雖然這個(gè)話題有點(diǎn)超出本書(shū)的范圍，但還是要提一下，因?yàn)樵谀承?yīng)用場(chǎng)景中，數(shù)組的內(nèi)存布局可以對(duì)計(jì)算速度造成極大的影響。這是因?yàn)樾阅懿顒e在一定程度上跟CPU的高速緩存（cache）體系有關(guān)。運(yùn)算過(guò)程中訪問(wèn)連續(xù)內(nèi)存塊（例如，對(duì)以C順序存儲(chǔ)的數(shù)組的行求和）一般是最快的，因?yàn)閮?nèi)存子系統(tǒng)會(huì)將適當(dāng)?shù)膬?nèi)存塊緩存到超高速的L1或L2CPU Cache中。此外，NumPy的C語(yǔ)言基礎(chǔ)代碼（某些）對(duì)連續(xù)存儲(chǔ)的情況進(jìn)行了優(yōu)化處理，這樣就能避免一些跨越式的內(nèi)存訪問(wèn)。
一個(gè)數(shù)組的內(nèi)存布局是連續(xù)的，就是說(shuō)元素是以它們?cè)跀?shù)組中出現(xiàn)的順序（即Fortran型（列優(yōu)先）或C型（行優(yōu)先））存儲(chǔ)在內(nèi)存中的。默認(rèn)情況下，NumPy數(shù)組是以C型連續(xù)的方式創(chuàng)建的。列優(yōu)先的數(shù)組（比如C型連續(xù)數(shù)組的轉(zhuǎn)置）也被稱(chēng)為Fortran型連續(xù)。通過(guò)ndarray的flags屬性即可查看這些信息：
In [225]: arr_c = np.ones((1000, 1000), order='C')
In [226]: arr_f = np.ones((1000, 1000), order='F')
In [227]: arr_c.flags
Out[227]:
C_CONTIGUOUS : True
F_CONTIGUOUS : False
OWNDATA : True
WRITEABLE : True
ALIGNED : True
UPDATEIFCOPY : False
In [228]: arr_f.flags
Out[228]:
C_CONTIGUOUS : False
F_CONTIGUOUS : True
OWNDATA : True
WRITEABLE : True
ALIGNED : True
UPDATEIFCOPY : False
In [229]: arr_f.flags.f_contiguous
Out[229]: True
在這個(gè)例子中，對(duì)兩個(gè)數(shù)組的行進(jìn)行求和計(jì)算，理論上說(shuō)，arr_c會(huì)比arr_f快，因?yàn)閍rr_c的行在內(nèi)存中是連續(xù)的。我們可以在IPython中用%timeit來(lái)確認(rèn)一下：
In [230]: %timeit arr_c.sum(1)
784 us +- 10.4 us per loop (mean +- std. dev. of 7 runs, 1000 loops each)
In [231]: %timeit arr_f.sum(1)
934 us +- 29 us per loop (mean +- std. dev. of 7 runs, 1000 loops each)
如果想從NumPy中提升性能，這里就應(yīng)該是下手的地方。如果數(shù)組的內(nèi)存順序不符合你的要求，使用copy并傳入'C'或'F'即可解決該問(wèn)題：
In [232]: arr_f.copy('C').flags
Out[232]:
C_CONTIGUOUS : True
F_CONTIGUOUS : False
OWNDATA : True
WRITEABLE : True
ALIGNED : True
UPDATEIFCOPY : False
注意，在構(gòu)造數(shù)組的視圖時(shí)，其結(jié)果不一定是連續(xù)的：
In [233]: arr_c[:50].flags.contiguous
Out[233]: True
In [234]: arr_c[:, :50].flags
Out[234]:
C_CONTIGUOUS : False
F_CONTIGUOUS : False
OWNDATA : False
WRITEABLE : True
ALIGNED : True
UPDATEIFCOPY : False
