9.7 數(shù)組上的計(jì)算：廣播
本節(jié)是《Python 數(shù)據(jù)科學(xué)手冊(cè)》（Python Data Science Handbook）的摘錄。
譯者：飛龍
協(xié)議：CC BY-NC-SA 4.0
我們?cè)谏弦还?jié)中看到，NumPy 的通用函數(shù)如何用于向量化操作，從而消除緩慢的 Python 循環(huán)。向量化操作的另一種方法是使用 NumPy 的廣播功能。廣播只是一組規(guī)則，用于在不同大小的數(shù)組上應(yīng)用二元ufunc（例如，加法，減法，乘法等）。
廣播簡(jiǎn)介
回想一下，對(duì)于相同大小的數(shù)組，二元操作是逐元素執(zhí)行的：
import numpy as np
a = np.array([0, 1, 2])
b = np.array([5, 5, 5])
a + b
# array([5, 6, 7])
廣播允許在不同大小的數(shù)組上執(zhí)行這類二元操作 - 例如，我們可以輕松將數(shù)組和標(biāo)量相加（將其視為零維數(shù)組）：
a + 5
# array([5, 6, 7])
我們可以將此視為一個(gè)操作，將值5拉伸或復(fù)制為數(shù)組[5,5,5]，并將結(jié)果相加。
NumPy 廣播的優(yōu)勢(shì)在于，這種值的重復(fù)實(shí)際上并沒(méi)有發(fā)生，但是當(dāng)我們考慮廣播時(shí)，它是一種有用的心理模型。
我們可以類似地，將其擴(kuò)展到更高維度的數(shù)組。 將兩個(gè)二維數(shù)組相加時(shí)觀察結(jié)果：
M = np.ones((3, 3))
M
'''
array([[ 1., 1., 1.],
[ 1., 1., 1.],
[ 1., 1., 1.]])
'''
M + a
'''
array([[ 1., 2., 3.],
[ 1., 2., 3.],
[ 1., 2., 3.]])
'''
這里，一維數(shù)組a被拉伸，或者在第二維上廣播，來(lái)匹配M的形狀。
雖然這些示例相對(duì)容易理解，但更復(fù)雜的情況可能涉及兩個(gè)數(shù)組的廣播。請(qǐng)考慮以下示例：
a = np.arange(3)
b = np.arange(3)[:, np.newaxis]
print(a)
print(b)
'''
[0 1 2]
[[0]
[1]
[2]]
'''
a + b
'''
array([[0, 1, 2],
[1, 2, 3],
[2, 3, 4]])
'''
就像之前我們拉伸或廣播一個(gè)值來(lái)匹配另一個(gè)的形狀，這里我們拉伸a```和b``來(lái)匹配一個(gè)共同的形狀，結(jié)果是二維數(shù)組！
這些示例的幾何圖形為下圖（產(chǎn)生此圖的代碼可以在“附錄”中找到，并改編自 astroML 中發(fā)布的源碼，經(jīng)許可而使用）。
淺色方框代表廣播的值：同樣，這個(gè)額外的內(nèi)存實(shí)際上并沒(méi)有在操作過(guò)程中分配，但是在概念上想象它是有用的。
廣播規(guī)則
NumPy 中的廣播遵循一套嚴(yán)格的規(guī)則來(lái)確定兩個(gè)數(shù)組之間的交互：
- 規(guī)則 1：如果兩個(gè)數(shù)組的維數(shù)不同，則維數(shù)較少的數(shù)組的形狀，將在其左側(cè)填充。
- 規(guī)則 2：如果兩個(gè)數(shù)組的形狀在任何維度上都不匹配，則該維度中形狀等于 1 的數(shù)組將被拉伸來(lái)匹配其他形狀。
- 規(guī)則 3：如果在任何維度中，大小不一致且都不等于 1，則會(huì)引發(fā)錯(cuò)誤。
為了講清楚這些規(guī)則，讓我們?cè)敿?xì)考慮幾個(gè)例子。
廣播示例 1
讓我們看一下將二維數(shù)組和一維數(shù)組相加：
M = np.ones((2, 3))
a = np.arange(3)
讓我們考慮這兩個(gè)數(shù)組上的操作。數(shù)組的形狀是。
M.shape = (2, 3)a.shape = (3,)
我們?cè)谝?guī)則 1 中看到數(shù)組a的維數(shù)較少，所以我們?cè)谧筮吿畛渌?/p>
M.shape -> (2, 3)a.shape -> (1, 3)
根據(jù)規(guī)則 2，我們現(xiàn)在看到第一個(gè)維度不一致，因此我們將此維度拉伸來(lái)匹配：
M.shape -> (2, 3)a.shape -> (2, 3)
形狀匹配了，我們看到最終的形狀將是(2, 3)
M + a
'''
array([[ 1., 2., 3.],
[ 1., 2., 3.]])
'''
廣播示例 2
我們來(lái)看一個(gè)需要廣播兩個(gè)數(shù)組的例子：
a = np.arange(3).reshape((3, 1))
b = np.arange(3)
同樣，我們將首先寫出數(shù)組的形狀：
a.shape = (3, 1)b.shape = (3,)
規(guī)則 1 說(shuō)我們必須填充b的形狀：
a.shape -> (3, 1)b.shape -> (1, 3)
規(guī)則 2 告訴我們，我們更新這些中的每一個(gè)，來(lái)匹配另一個(gè)數(shù)組的相應(yīng)大小：
a.shape -> (3, 3)b.shape -> (3, 3)
因?yàn)榻Y(jié)果匹配，所以這些形狀是兼容的。我們?cè)谶@里可以看到：
a + b
'''
array([[0, 1, 2],
[1, 2, 3],
[2, 3, 4]])
'''
廣播示例 3
現(xiàn)在讓我們來(lái)看一個(gè)兩個(gè)數(shù)組不兼容的例子：
M = np.ones((3, 2))
a = np.arange(3)
這與第一個(gè)例子略有不同：矩陣M是轉(zhuǎn)置的。這對(duì)計(jì)算有何影響？數(shù)組的形狀是
M.shape = (3, 2)a.shape = (3,)
同樣，規(guī)則 1 告訴我們必須填充a的形狀：
M.shape -> (3, 2)a.shape -> (1, 3)
根據(jù)規(guī)則 2，a的第一個(gè)維度被拉伸來(lái)匹配M：
M.shape -> (3, 2)a.shape -> (3, 3)
現(xiàn)在我們到了規(guī)則 3 - 最終的形狀不匹配，所以這兩個(gè)數(shù)組是不兼容的，正如我們可以通過(guò)嘗試此操作來(lái)觀察：
M + a
'''
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<ipython-input-13-9e16e9f98da6> in <module>()
----> 1 M + a
ValueError: operands could not be broadcast together with shapes (3,2) (3,)
'''
注意這里潛在的混淆：你可以想象使a和M兼容，比如在右邊填充a的形狀，而不是在左邊。但這不是廣播規(guī)則的運(yùn)作方式！
在某些情況下，這種靈活性可能會(huì)有用，但這會(huì)導(dǎo)致潛在的二義性。如果在右側(cè)填充是你想要的，你可以通過(guò)數(shù)組的形狀調(diào)整，來(lái)明確地執(zhí)行此操作（我們將使用“NumPy 數(shù)組基礎(chǔ)”中介紹的np.newaxis關(guān)鍵字）：
a[:, np.newaxis].shape
# (3, 1)
M + a[:, np.newaxis]
'''
array([[ 1., 1.],
[ 2., 2.],
[ 3., 3.]])
'''
還要注意，雖然我們一直專注于+運(yùn)算符，但這些廣播規(guī)則適用于任何二元ufunc。
例如，這里是logaddexp(a, b)函數(shù)，它比原始方法更精確地計(jì)算log(exp(a) + exp(b))：
np.logaddexp(M, a[:, np.newaxis])
'''
array([[ 1.31326169, 1.31326169],
[ 1.69314718, 1.69314718],
[ 2.31326169, 2.31326169]])
'''
對(duì)于可用的通用函數(shù)的更多信息，請(qǐng)參閱“NumPy 數(shù)組上的計(jì)算：通用函數(shù)”。
實(shí)戰(zhàn)中的廣播
廣播操作是我們將在本書(shū)中看到的許多例子的核心。我們現(xiàn)在來(lái)看一些它們可能有用的簡(jiǎn)單示例。
數(shù)組中心化
在上一節(jié)中，我們看到ufunc允許 NumPy 用戶不再需要顯式編寫慢速 Python 循環(huán)。廣播擴(kuò)展了這種能力。一個(gè)常見(jiàn)的例子是數(shù)據(jù)數(shù)組的中心化。
想象一下，你有一組 10 個(gè)觀測(cè)值，每個(gè)觀測(cè)值由 3 個(gè)值組成。使用標(biāo)準(zhǔn)約定（參見(jiàn)“Scikit-Learn 中的數(shù)據(jù)表示”），我們將其存儲(chǔ)在10x3數(shù)組中：
X = np.random.random((10, 3))
我們可以使用第一維上的“均值”聚合，來(lái)計(jì)算每個(gè)特征的平均值：
Xmean = X.mean(0)
Xmean
# array([ 0.53514715, 0.66567217, 0.44385899])
現(xiàn)在我們可以通過(guò)減去均值（這是一個(gè)廣播操作）來(lái)中心化X數(shù)組：
X_centered = X - Xmean
要仔細(xì)檢查我們是否已正確完成此操作，我們可以檢查中心化的數(shù)組是否擁有接近零的均值：
X_centered.mean(0)
# array([ 2.22044605e-17, -7.77156117e-17, -1.66533454e-17])
在機(jī)器精度范圍內(nèi)，平均值現(xiàn)在為零。
繪制二維函數(shù)
廣播非常有用的一個(gè)地方是基于二維函數(shù)展示圖像。如果我們想要定義一個(gè)函數(shù)z = f(x, y)，廣播可用于在網(wǎng)格中計(jì)算函數(shù)：
# x 和 y 是從 0 到 5 的 50 步
x = np.linspace(0, 5, 50)
y = np.linspace(0, 5, 50)[:, np.newaxis]
z = np.sin(x) ** 10 + np.cos(10 + y * x) * np.cos(x)
我們將使用 Matplotlib 繪制這個(gè)二維數(shù)組（這些工具將在“密度和等高線圖”中完整討論）：
%matplotlib inline
import matplotlib.pyplot as plt
plt.imshow(z, origin='lower', extent=[0, 5, 0, 5],
cmap='viridis')
plt.colorbar();
結(jié)果是引人注目的二維函數(shù)的圖形。
