手寫 30 個主流機器學習演算法，程式碼超 3 萬行，全都開源了！

NumPy 作為 Python 生態中最受歡迎的科學計算包，很多讀者已經非常熟悉它了。它為 Python 提供高效率的多維陣列計算，並提供了一系列高等數學函式，我們可以快速搭建模型的整個計算流程。毫不負責任地說，NumPy 就是現代深度學習框架的「爸爸」。

儘管目前使用 NumPy 寫模型已經不是主流，但這種方式依然不失為是理解底層架構和深度學習原理的好方法。最近，來自普林斯頓的一位博士後將 NumPy 實現的所有機器學習模型全部開源，並提供了相應的論文和一些實現的測試效果。

• 專案地址：https://github.com/ddbourgin/numpy-ml

根據機器之心的粗略估計，該專案大約有 30 個主要機器學習模型，此外還有 15 個用於預處理和計算的小工具，全部.py 檔案數量有 62 個之多。平均每個模型的程式碼行數在 500 行以上，在神經網路模型的 layer.py 檔案中，程式碼行數接近 4000。

這，應該是目前用 NumPy 手寫機器學習模型的「最高境界」吧。

誰用 NumPy 手推了一大波 ML 模型

通過專案的程式碼目錄，我們能發現，作者基本上把主流模型都實現了一遍，這個工作量簡直驚為天人。我們發現作者 David Bourgin 也是一位大神，他於 2018 年獲得加州大學伯克利分校計算認知科學博士學位，隨後在普林斯頓大學從事博士後研究。

儘管畢業不久，David 在頂級期刊與計算機會議上都發表了一些優秀論文。在最近結束的 ICML 2019 中，其關於認知模型先驗的研究就被接收為少有的 Oral 論文。

David Bourgin 小哥哥就是用 NumPy 手寫 ML 模型、手推反向傳播的大神。這麼多的工作量，當然還是需要很多參考資源的，David 會理解這些資源或實現，並以一種更易讀的方式寫出來。

正如 reddit 讀者所質疑的：在 autograd repo 中已經有很多這樣的例子，為什麼你還要做這個專案？

作者表示，他的確從 autograd repo 學到了很多，但二者的不同之處在於，他顯式地進行了所有梯度計算，以突出概念/數學的清晰性。當然，這麼做的缺點也很明顯，在每次需要微分一個新函式時，你都要寫出它的公式……

估計 David Bourgin 小哥哥在寫完這個專案後，機器學習基礎已經極其牢固了。最後，David 表示下一步會新增文件和示例，以方便大家使用。

專案總體介紹

這個專案最大的特點是作者把機器學習模型都用 NumPy 手寫了一遍，包括更顯式的梯度計算和反向傳播過程。可以說它就是一個機器學習框架了，只不過程式碼可讀性會強很多。

David Bourgin 表示他一直在慢慢寫或收集不同模型與模組的純 NumPy 實現，它們跑起來可能沒那麼快，但是模型的具體過程一定足夠直觀。每當我們想了解模型 API 背後的實現，卻又不想看複雜的框架程式碼，那麼它可以作為快速的參考。

文章後面會具體介紹整個專案都有什麼模型，這裡先簡要介紹它的整體結構。如下所示為專案檔案，不同的資料夾即不同種類的程式碼集。

在每一個程式碼集下，作者都會提供不同實現的參考資料，例如模型的效果示例圖、參考論文和參考連結等。如下所示，David 在實現神經網路層級的過程中，還提供了參考論文。

當然如此龐大的程式碼總會存在一些 Bug，作者也非常希望我們能一起完善這些實現。如果我們以前用純 NumPy 實現過某些好玩的模型，那也可以直接提交 PR 請求。因為實現基本上都只依賴於 NumPy，那麼環境配置就簡單很多了，大家差不多都能跑得動。

手寫 NumPy 全家福

作者在 GitHub 中提供了模型/模組的實現列表，列表結構基本就是程式碼檔案的結構了。整體上，模型主要分為兩部分，即傳統機器學習模型與主流的深度學習模型。

其中淺層模型既有隱馬爾可夫模型和提升方法這樣的複雜模型，也包含了線性迴歸或最近鄰等經典方法。而深度模型則主要從各種模組、層級、損失函式、最優化器等角度搭建程式碼架構，從而能快速構建各種神經網路。

除了模型外，整個專案還有一些輔助模組，包括一堆預處理相關的元件和有用的小工具。

該 repo 的模型或程式碼結構如下所示：

1. 高斯混合模型

• EM 訓練

2. 隱馬爾可夫模型

• 維特比解碼

• 似然計算

• 通過 Baum-Welch/forward-backward 演算法進行 MLE 引數估計

3. 隱狄利克雷分配模型（主題模型）

• 用變分 EM 進行 MLE 引數估計的標準模型

• 用 MCMC 進行 MAP 引數估計的平滑模型

4. 神經網路

4.1 層/層級運算

• Add

• Flatten

• Multiply

• Softmax

• 全連線/Dense

• 稀疏進化連線

• LSTM

• Elman 風格的 RNN

• 最大+平均池化

• 點積注意力

• 受限玻爾茲曼機 (w. CD-n training)

• 2D 轉置卷積 (w. padding 和 stride)

• 2D 卷積 (w. padding、dilation 和 stride)

• 1D 卷積 (w. padding、dilation、stride 和 causality)

4.2 模組

• 雙向 LSTM

• ResNet 風格的殘差塊（恆等變換和卷積）

• WaveNet 風格的殘差塊（帶有擴張因果卷積）

• Transformer 風格的多頭縮放點積注意力

4.3 正則化項

• Dropout

• 歸一化

• 批歸一化（時間上和空間上）

• 層歸一化（時間上和空間上）

4.4 優化器

• SGD w/ 動量

• AdaGrad

• RMSProp

• Adam

4.5 學習率排程器

• 常數

• 指數

• Noam/Transformer

• Dlib 排程器

4.6 權重初始化器

• Glorot/Xavier uniform 和 normal

• He/Kaiming uniform 和 normal

• 標準和截斷正態分佈初始化

4.7 損失

• 交叉熵

• 平方差

• Bernoulli VAE 損失

• 帶有梯度懲罰的 Wasserstein 損失

4.8 啟用函式

• ReLU

• Tanh

• Affine

• Sigmoid

• Leaky ReLU

4.9 模型

• Bernoulli 變分自編碼器

• 帶有梯度懲罰的 Wasserstein GAN

4.10 神經網路工具

• col2im (MATLAB 埠)

• im2col (MATLAB 埠)

• conv1D

• conv2D

• deconv2D

• minibatch

5. 基於樹的模型

• 決策樹 (CART)

• [Bagging] 隨機森林

• [Boosting] 梯度提升決策樹

6. 線性模型

• 嶺迴歸

• Logistic 迴歸

• 最小二乘法

• 貝葉斯線性迴歸 w/共軛先驗

7.n 元序列模型

• 最大似然得分

• Additive/Lidstone 平滑

• 簡單 Good-Turing 平滑

8. 強化學習模型

• 使用交叉熵方法的智慧體

• 首次訪問 on-policy 蒙特卡羅智慧體

• 加權增量重要取樣蒙特卡羅智慧體

• Expected SARSA 智慧體

• TD-0 Q-learning 智慧體

• Dyna-Q / Dyna-Q+ 優先掃描

9. 非引數模型

• Nadaraya-Watson 核迴歸

• k 最近鄰分類與迴歸

10. 預處理

• 離散傅立葉變換 (1D 訊號)

• 雙線性插值 (2D 訊號)

• 最近鄰插值 (1D 和 2D 訊號)

• 自相關 (1D 訊號）

• 訊號視窗

• 文字分詞

• 特徵雜湊

• 特徵標準化

• One-hot 編碼/解碼

• Huffman 編碼/解碼

• 詞頻逆文件頻率編碼

11. 工具

• 相似度核

• 距離度量

• 優先順序佇列

• Ball tree 資料結構

專案示例

由於程式碼量龐大，機器之心在這裡整理了一些示例。

例如，實現點積注意力機制：

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 三個向量輸入到「_fwd」函式中，用於計算每個向量的注意力分數，並通過 softmax 的方式得到權重。而「_bwd」函式則計算 V、注意力權重、注意力分數、Q 和 K 的梯度，用於更新網路權重。

在一些實現中，作者也進行了測試，並給出了測試結果。如圖為隱狄利克雷（Latent Dirichlet allocation，LDA）實現進行文字聚類的結果。左圖為詞語在特定主題中的分佈熱力圖。右圖則為文件在特定主題中的分佈熱力圖。

圖注：隱狄利克雷分佈實現的效果。

- END -

如果看到這裡，說明你喜歡這篇文章，請轉發、點贊。掃描下方二維碼或者微信搜尋「perfect_iscas」，新增好友後即可獲得10套程式設計師全棧課程+1000套PPT和簡歷模板，向我私聊「進群」二字即可進入高質量交流群。

↓掃描二維碼進群↓

喜歡文章，點個在看