?Numpy手寫機(jī)器學(xué)習(xí)算法，3萬(wàn)行代碼！

點(diǎn)擊下方卡片，關(guān)注“新機(jī)器視覺”公眾號(hào)
重磅干貨，第一時(shí)間送達(dá)

轉(zhuǎn)自：機(jī)器之心

NumPy 作為 Python 生態(tài)中最受歡迎的科學(xué)計(jì)算包，很多讀者已經(jīng)非常熟悉它了。

它為 Python 提供高效率的多維數(shù)組計(jì)算，并提供了一系列高等數(shù)學(xué)函數(shù)，我們可以快速搭建模型的整個(gè)計(jì)算流程。毫不負(fù)責(zé)任地說(shuō)，NumPy 就是現(xiàn)代深度學(xué)習(xí)框架的「爸爸」。

盡管目前使用 NumPy 寫模型已經(jīng)不是主流，但這種方式依然不失為是理解底層架構(gòu)和深度學(xué)習(xí)原理的好方法。最近，來(lái)自普林斯頓的一位博士后將 NumPy 實(shí)現(xiàn)的所有機(jī)器學(xué)習(xí)模型全部開源，超過(guò) 3 萬(wàn)行代碼、30 多個(gè)模型，并提供了相應(yīng)的論文和一些實(shí)現(xiàn)的測(cè)試效果。

項(xiàng)目地址：https://github.com/ddbourgin/numpy-ml

粗略估計(jì)，該項(xiàng)目大約有 30 個(gè)主要機(jī)器學(xué)習(xí)模型，此外還有 15 個(gè)用于預(yù)處理和計(jì)算的小工具，全部.py 文件數(shù)量有 62 個(gè)之多。平均每個(gè)模型的代碼行數(shù)在 500 行以上，在神經(jīng)網(wǎng)絡(luò)模型的 layer.py 文件中，代碼行數(shù)接近 4000。

這，應(yīng)該是目前用 NumPy 手寫機(jī)器學(xué)習(xí)模型的「最高境界」吧。
誰(shuí)用 NumPy 手推了一大波 ML 模型？
通過(guò)項(xiàng)目的代碼目錄，我們能發(fā)現(xiàn)，作者基本上把主流模型都實(shí)現(xiàn)了一遍，這個(gè)工作量簡(jiǎn)直驚為天人。作者 David Bourgin 是一位大神，于 2018 年獲得加州大學(xué)伯克利分校計(jì)算認(rèn)知科學(xué)博士學(xué)位，隨后在普林斯頓大學(xué)從事博士后研究。

盡管畢業(yè)不久，David 在頂級(jí)期刊與計(jì)算機(jī)會(huì)議上都發(fā)表了一些優(yōu)秀論文。在 ICML 2019 中，其關(guān)于認(rèn)知模型先驗(yàn)的研究就被接收為少有的 Oral 論文。


David Bourgin 就是用 NumPy 手寫 ML 模型、手推反向傳播的大神。這么多的工作量，當(dāng)然還是需要很多參考資源的，David 會(huì)理解這些資源或?qū)崿F(xiàn)，并以一種更易讀的方式寫出來(lái)。

他表示，從 autograd repo 學(xué)到了很多，但二者的不同之處在于，他顯式地進(jìn)行了所有梯度計(jì)算，以突出概念/數(shù)學(xué)的清晰性。當(dāng)然，這么做的缺點(diǎn)也很明顯，在每次需要微分一個(gè)新函數(shù)時(shí)，你都要寫出它的公式……

估計(jì) David Bourgin 在寫完這個(gè)項(xiàng)目后，機(jī)器學(xué)習(xí)基礎(chǔ)已經(jīng)極其牢固了。
項(xiàng)目總體介紹
這個(gè)項(xiàng)目最大的特點(diǎn)是作者把機(jī)器學(xué)習(xí)模型都用 NumPy 手寫了一遍，包括更顯式的梯度計(jì)算和反向傳播過(guò)程。可以說(shuō)它就是一個(gè)機(jī)器學(xué)習(xí)框架了，只不過(guò)代碼可讀性會(huì)強(qiáng)很多。

David Bourgin 表示他一直在慢慢寫或收集不同模型與模塊的純 NumPy 實(shí)現(xiàn)，它們跑起來(lái)可能沒那么快，但是模型的具體過(guò)程一定足夠直觀。每當(dāng)我們想了解模型 API 背后的實(shí)現(xiàn)，卻又不想看復(fù)雜的框架代碼，那么它可以作為快速的參考。

文章后面會(huì)具體介紹整個(gè)項(xiàng)目都有什么模型，這里先簡(jiǎn)要介紹它的整體結(jié)構(gòu)。如下所示為項(xiàng)目文件，不同的文件夾即不同種類的代碼集。


在每一個(gè)代碼集下，作者都會(huì)提供不同實(shí)現(xiàn)的參考資料，例如模型的效果示例圖、參考論文和參考鏈接等。如下所示，David 在實(shí)現(xiàn)神經(jīng)網(wǎng)絡(luò)層級(jí)的過(guò)程中，還提供了參考論文。


當(dāng)然如此龐大的代碼總會(huì)存在一些 Bug，作者也非常希望我們能一起完善這些實(shí)現(xiàn)。如果我們以前用純 NumPy 實(shí)現(xiàn)過(guò)某些好玩的模型，那也可以直接提交 PR 請(qǐng)求。因?yàn)閷?shí)現(xiàn)基本上都只依賴于 NumPy，那么環(huán)境配置就簡(jiǎn)單很多了，大家差不多都能跑得動(dòng)。
手寫 NumPy 全家福
作者在 GitHub 中提供了模型/模塊的實(shí)現(xiàn)列表，列表結(jié)構(gòu)基本就是代碼文件的結(jié)構(gòu)了。整體上，模型主要分為兩部分，即傳統(tǒng)機(jī)器學(xué)習(xí)模型與主流的深度學(xué)習(xí)模型。

其中淺層模型既有隱馬爾可夫模型和提升方法這樣的復(fù)雜模型，也包含了線性回歸或最近鄰等經(jīng)典方法。而深度模型則主要從各種模塊、層級(jí)、損失函數(shù)、最優(yōu)化器等角度搭建代碼架構(gòu)，從而能快速構(gòu)建各種神經(jīng)網(wǎng)絡(luò)。

除了模型外，整個(gè)項(xiàng)目還有一些輔助模塊，包括一堆預(yù)處理相關(guān)的組件和有用的小工具。

該 repo 的模型或代碼結(jié)構(gòu)如下所示：

1. 高斯混合模型
EM 訓(xùn)練
2. 隱馬爾可夫模型
維特比解碼
似然計(jì)算
通過(guò) Baum-Welch/forward-backward 算法進(jìn)行 MLE 參數(shù)估計(jì)
3. 隱狄利克雷分配模型（主題模型）
用變分 EM 進(jìn)行 MLE 參數(shù)估計(jì)的標(biāo)準(zhǔn)模型
用 MCMC 進(jìn)行 MAP 參數(shù)估計(jì)的平滑模型
4. 神經(jīng)網(wǎng)絡(luò)
4.1 層/層級(jí)運(yùn)算
Add
Flatten
Multiply
Softmax
全連接/Dense
稀疏進(jìn)化連接
LSTM
Elman 風(fēng)格的 RNN
最大+平均池化
點(diǎn)積注意力
受限玻爾茲曼機(jī) (w. CD-n training)
2D 轉(zhuǎn)置卷積 (w. padding 和 stride)
2D 卷積 (w. padding、dilation 和 stride)
1D 卷積 (w. padding、dilation、stride 和 causality)
4.2 模塊
雙向 LSTM
ResNet 風(fēng)格的殘差塊（恒等變換和卷積）
WaveNet 風(fēng)格的殘差塊（帶有擴(kuò)張因果卷積）
Transformer 風(fēng)格的多頭縮放點(diǎn)積注意力
4.3 正則化項(xiàng)
Dropout
歸一化
批歸一化（時(shí)間上和空間上）
層歸一化（時(shí)間上和空間上）
4.4 優(yōu)化器
SGD w/ 動(dòng)量
AdaGrad
RMSProp
Adam
4.5 學(xué)習(xí)率調(diào)度器
常數(shù)
指數(shù)
Noam/Transformer
Dlib 調(diào)度器
4.6 權(quán)重初始化器
Glorot/Xavier uniform 和 normal
He/Kaiming uniform 和 normal
標(biāo)準(zhǔn)和截?cái)嗾龖B(tài)分布初始化
4.7 損失
交叉熵
平方差
Bernoulli VAE 損失
帶有梯度懲罰的 Wasserstein 損失
4.8 激活函數(shù)
ReLU
Tanh
Affine
Sigmoid
Leaky ReLU
4.9 模型
Bernoulli 變分自編碼器
帶有梯度懲罰的 Wasserstein GAN
4.10 神經(jīng)網(wǎng)絡(luò)工具
col2im (MATLAB 端口)
im2col (MATLAB 端口)
conv1D
conv2D
deconv2D
minibatch
5. 基于樹的模型
決策樹 (CART)
[Bagging] 隨機(jī)森林
[Boosting] 梯度提升決策樹
6. 線性模型
嶺回歸
Logistic 回歸
最小二乘法
貝葉斯線性回歸 w/共軛先驗(yàn)
7.n 元序列模型
最大似然得分
Additive/Lidstone 平滑
簡(jiǎn)單 Good-Turing 平滑
8. 強(qiáng)化學(xué)習(xí)模型
使用交叉熵方法的智能體
首次訪問(wèn) on-policy 蒙特卡羅智能體
加權(quán)增量重要采樣蒙特卡羅智能體
Expected SARSA 智能體
TD-0 Q-learning 智能體
Dyna-Q / Dyna-Q+ 優(yōu)先掃描
9. 非參數(shù)模型
Nadaraya-Watson 核回歸
k 最近鄰分類與回歸
10. 預(yù)處理
離散傅立葉變換 (1D 信號(hào))
雙線性插值 (2D 信號(hào))
最近鄰插值 (1D 和 2D 信號(hào))
自相關(guān) (1D 信號(hào)）
信號(hào)窗口
文本分詞
特征哈希
特征標(biāo)準(zhǔn)化
One-hot 編碼/解碼
Huffman 編碼/解碼
詞頻逆文檔頻率編碼
11. 工具
相似度核
距離度量
優(yōu)先級(jí)隊(duì)列
Ball tree 數(shù)據(jù)結(jié)構(gòu)
項(xiàng)目示例
由于代碼量龐大，這里整理了一些示例。

例如，實(shí)現(xiàn)點(diǎn)積注意力機(jī)制：
class DotProductAttention(LayerBase):
    def __init__(self, scale=True, dropout_p=0, init="glorot_uniform", optimizer=None):
        super().__init__(optimizer)
        self.init = init
        self.scale = scale
        self.dropout_p = dropout_p
        self.optimizer = self.optimizer
        self._init_params()

    def _fwd(self, Q, K, V):
        scale = 1 / np.sqrt(Q.shape[-1]) if self.scale else 1
        scores = Q @ K.swapaxes(-2, -1) * scale  # attention scores
        weights = self.softmax.forward(scores)  # attention weights
        Y = weights @ V
        return Y, weights

    def _bwd(self, dy, q, k, v, weights):
        d_k = k.shape[-1]
        scale = 1 / np.sqrt(d_k) if self.scale else 1

        dV = weights.swapaxes(-2, -1) @ dy
        dWeights = dy @ v.swapaxes(-2, -1)
        dScores = self.softmax.backward(dWeights)
        dQ = dScores @ k * scale
        dK = dScores.swapaxes(-2, -1) @ q * scale
        return dQ, dK, dV

在以上代碼中，Q、K、V 三個(gè)向量輸入到「_fwd」函數(shù)中，用于計(jì)算每個(gè)向量的注意力分?jǐn)?shù)，并通過(guò) softmax 的方式得到權(quán)重。而「_bwd」函數(shù)則計(jì)算 V、注意力權(quán)重、注意力分?jǐn)?shù)、Q 和 K 的梯度，用于更新網(wǎng)絡(luò)權(quán)重。

在一些實(shí)現(xiàn)中，作者也進(jìn)行了測(cè)試，并給出了測(cè)試結(jié)果。如圖為隱狄利克雷（Latent Dirichlet allocation，LDA）實(shí)現(xiàn)進(jìn)行文本聚類的結(jié)果。左圖為詞語(yǔ)在特定主題中的分布熱力圖。右圖則為文檔在特定主題中的分布熱力圖。

圖注：隱狄利克雷分布實(shí)現(xiàn)的效果。