Metropolis-Hastings 演算法和 Gibbs sampling 演算法
我們有一個概率分布(pi)
，需要生成滿足這個分布的樣本。 如果樣本維度很低，只有一兩維，我們可以用反切法、拒絕採樣和重要性採樣等方法。 但是樣本維度很高，成千上百，這些方法就不適用了。這時我們就要使用一些高檔的演算法，比如下面要介紹的 Metropolis-Hasting 演算法和 Gibbs sampling 演算法。
Metropolis-Hasting 演算法和 Gibbs sampling 演算法是馬爾科夫鏈蒙特卡洛（Markov Chain Mento Carlo, MCMC）方法。我們先介紹 MCMC 方法。
1. 馬爾科夫蒙特卡洛方法
MCMC 方法是用蒙特卡洛方法去體現馬爾科夫鏈的方法。馬爾科夫鏈是狀態空間的轉換關係，下一個狀態只和當前的狀態有關。比如下圖就是一個馬爾科夫鏈的示意圖。
圖中轉移關係可以用一個概率轉換矩陣 p 表示，
egin
P = egin
0 &1 &0 \
0&0.1 &0.9 \
0.6 &0.4 &0
end
如果當前狀態分布為(u(x) = (0.3,0.2,0.5))
, 那麼下一個矩陣的狀態就是( u(x)p ), 再下一個就是(u(x)p^2),... 最後會收斂到一個平穩分布(pi)。這個平穩分布(pi)只和概率轉移矩陣 p 有關，而和初始狀態分布 u 是什麼沒有關係。
如何判斷一個馬爾科夫鏈是否能收斂到平穩分布，以及如何判斷一個狀態分布是不是一個馬爾科夫鏈的平穩分布呢？我們有下面定理。
細緻平衡條件: 已知各態歷經的的馬爾科夫鏈有概率轉移矩陣(p)
，以及已知狀態分布(pi)。如果對於任意兩個狀態 i 和 j，下面公式成立，則馬爾科夫鏈能夠收斂到(pi)。
pi(i) p(j|i) = pi(j) p(i|j) onumber
這裡的各態歷經是指任意兩個狀態之間可以通過有限步到達。
怎麼證明細緻平衡條件呢？我也不知道啊。
MCMC 方法的基本原理是利用細緻平衡條件構建一個概率轉移矩陣，使得目標概率就是概率轉移矩陣的平穩分布。怎麼構建這樣一個概率轉移矩陣呢？ Metropolis-Hasting 和 Gibbs sampling 演算法本質上就是構建概率轉移矩陣的不同方法。
2. Metropolis-Hastings 演算法
Metropolis-Hastings 演算法先提出一個可能不符合條件的概率轉移矩陣 q, 然後再進行調整。比如我們提出的 q 是均勻概率，即從任意狀態到任意狀態的概率是相等的。顯然在絕大部分情況下，q 的穩定概率不是目標概率(pi)
，即不滿足細緻平衡條件。
pi(i) q(j|i) e pi(j)q(i|j) onumber
如何讓這個不等式轉變成等式呢？根據對稱性，我們容易得到下面的等式。
label
pi(i) q(j|i) pi(j)q(i|j) = pi(j)q(i|j) pi(i)q(j|i)
這時整個概率轉移矩陣滿足細緻平衡條件。從 i 狀態轉到 j 狀態的概率是(q(j|i) pi(j)q(i|j))
，實現這個轉移概率的方式是 i 狀態以 q(j|i) 概率跳轉到 j 狀態，然後以(pi(j)q(i|j))接受跳轉 (拒絕跳轉就退回 i 狀態)。這樣整個 Metropolis-Hasting 演算法的框架就建立起來了。
這個原始的 Metropoli-Hasting 演算法的有一個小問題。 跳轉接受概率(pi(j)q(i|j))
和(pi(i)q(j|i))的值很小，演算法進行過程充斥著跳轉拒絕。為了改進這點，Metropoli-Hasting 演算法的方法是公式兩邊同時乘以一個係數，使得(pi(j)q(i|j))和(pi(i)q(j|i))中大的一項 scale 到 1，得到下面的公式。
pi(i) q(j|i) frac{pi(j)q(i|j)}{pi(i)q(j|i)} &=& pi(j)q(i|j) ;; when ;; pi(i)q(j|i) > pi(j)q(i|j) onumber \
oronumber \
pi(i) q(j|i) &=& pi(j)q(i|j)frac{pi(i)q(j|i)}{pi(j)q(i|j)} ;; when ;; pi(i)q(j|i) le pi(j)q(i|j) onumber \
這個公式可以進一步簡化為公式 2
pi(i) q(j|i) a(j|i) &=& pi(j)q(i|j) a(i|j) onumber \
a(j|i) &=& min{frac{pi(j)q(i|j)}{pi(i)q(j|i)},1} onumber \
a(i|j) &=& min{frac{pi(i)q(j|i)}{pi(j)q(i|j)},1}
根據上面的推導，我們容易得到 Metropolis-Hasting 演算法的流程。
3. Gibbs sampling 演算法
Gibbs sampling 演算法是 Metropolis-Hasting 演算法的一個特例。很雞賊的一個特例。m 維的一個樣本跳轉到另一個樣本的過程，可以拆解為 m 個子過程，每一個子過程對應一個維度。這時概率轉移矩陣可以看成是 m 個子概率轉移矩陣之積，即(p = prod_^ p_k )
其中(p_k)
表示第 k 維的變化概率。在(p_k)中，兩個狀態之間只有 k 維不同，其跳轉概率如下所示；不然為 0。
p_k(pmb_{dashv k, k=v_2}|pmb_{dashv k, k=v_1}) = frac{pi(pmb_{dashv k, k=v_2})}{sum_pi(pmb_{dashv k, k=v})} onumber
其中(pmb_{dashv k, k=v_2})
表示樣本第 k 維數據為(v_2)，其他沒有變化。這時候我們發現如下公式
&& pi(pmb_{dashv k, k=v_1}) p(pmb_{dashv k, k=v_2}|pmb_{dashv k, k=v_1}) onumber \
&=& pi(pmb_{dashv k, k=v_1}) frac{pi(pmb_{dashv k, k=v_2})}{sum_pi(pmb_{dashv k, k=v})} onumber \
&=& pi(pmb_{dashv k, k=v_2}) frac{pi(pmb_{dashv k, k=v_1})}{sum_pi(pmb_{dashv k, k=v})} onumber \
&=& pi(pmb_{dashv k, k=v_2}) p(pmb_{dashv k, k=v_1}|pmb_{dashv k, k=v_2}) onumber \
即(p_k)
和(pi)滿足細緻平衡條件的等式。那麼(p_k)就是我們要構建的概率轉移矩陣嘛？答案是否定的。因為完整的細緻平衡條件需要各態歷經，而一個狀態在概率轉移矩陣(p_k)下永遠不能到達另一個第 k-1 維數據不同的狀態。我們最終構建的概率轉移矩陣
p = prod_^ p_k
我們很容易證明(p)
依然滿足細緻平衡條件中的等式，同時還滿足各態歷經。根據這些推理，我們得到 Gibbs sampling 的演算法過程。
3. 總結
獨家解密：探月「網紅」兔背後的推手海報設計靈感：簡約獨特的圖形圖案排版 by Quim Marin量子通信衛星上天，VR遊戲時代的黎明也快到了？廣告人挖腦溝想創意，有人成就經典，有人只是挖溝5G關鍵技術之波束成型
TAG:推酷 |
※Bridge to the digital world——AR演算法技術分享
※資深演算法工程師眼中的深度學習：Ian Goodfellow 和Yoshua Bengio的「Deep Learning」讀書分享
※標籤傳播演算法（Label Propagation）及 Python 實現
※谷歌投資「演算法商店」 Algorithmia ，要打造AI版的Play Store？
※TensorFlow Agents日前開源，輕鬆在TensorFlow中構建並行強化學習演算法
※Facebook開源Zstandard新型壓縮演算法代替Zlib 簡單使用
※解析Monte-Carlo演算法
※《Alteration》：這部 VR 電影用上了 Facebook 的風格轉移演算法
※神奇的HyperLogLog演算法
※MetaTrader 5 演算法交易
※Science專訪谷歌Magenta負責人：AI創作焦點是機器學習演算法
※SlopOne推薦演算法（附Python源碼）
※機器學習演算法實踐：決策樹 (Decision Tree)
※DeepMind提出Rainbow：整合DQN演算法中的六種變體
※Nutella生產商Ferrero用一種演算法來設計包裝
※Facebook研製新工具，Google演算法則可能推送假新聞
※常見的排序演算法總結（JavaScript）
※「解密阿老師」從 AlphaGo 到Master，最大優勢是通用演算法
※遞歸演算法中一種不可忽略的技術（Memoization）