MCMC抽樣法（學習筆記）

前言：

后驗分布種給出了每個參數(shù)值的可能性镰吵，但很難看出完整的分布挺尿，也無法用分析解決册赛。于是我們需要MCMC法臀叙；
MCMC法有兩大要素：蒙特卡羅模型和馬爾可夫鏈

蒙特卡羅模型

蒙特卡羅模型：重復生成（大量）隨機數(shù)來估計固定參數(shù)（僅給出一個參數(shù)的近似值但不能直接計算）
- 舉例：估算不規(guī)則區(qū)域面積
  
  image.png
通過隨機生成20個點來估計圓面積（ $\approx \frac{15}{20}*100=75$ ）
應用：
1. 估算不規(guī)則區(qū)域面積
2. 模擬復雜過程（預測天氣笆豁、估計選舉結(jié)果......)

馬爾可夫鏈

歷史發(fā)展：

image.png
- 如圖郎汪，鐘形曲線被看作自然界中常見模式，但鐘形線甚至大多數(shù)規(guī)律所描述的一系列事件闯狱，彼此都是完全獨立的煞赢。因此，俄羅斯數(shù)學家和神學家帕維爾·涅克拉索夫（Pavel Nekrasov）認為：現(xiàn)實世界中的事物相互依存哄孤，不符合數(shù)學模型或分布照筑！
- 但是安德烈·馬爾科夫（Andrey Markov）試圖證明非獨立事件也可能符合數(shù)學模型。馬爾可夫鏈一個著名的應用案例是——從一部俄羅斯詩歌作品中輸出數(shù)千個“兩字母對”（two-character pairs）瘦陈。利用這些對凝危，計算每一字母的條件概率。也就是說晨逝，給定前面的字母或空格蛾默，就能預測下一個字母是A還是T，或是一個空格捉貌。利用這些概率支鸡，馬爾可夫可以模擬任意長的字符序列，這就是馬爾可夫鏈趁窃。當字符串一定長的時候牧挣，起始字母的選擇對后續(xù)字母的影響可以忽略不計即字母的分布遵循一種模式。也就是說醒陆，相互依賴的時間瀑构，如果它們受到固定概率的影響，也將遵循平均水平模式刨摩！
馬爾可夫鏈有概率相關的序列組成寺晌，每一事件來自一組結(jié)果，每個結(jié)果根據(jù)一組固定的概率確定下一個結(jié)果码邻。
一種重要特征——無記憶：預測下一個時間所需要的所有信息都能從當前狀態(tài)找到折剃，從歷史事件中找不到新的信息（類“蛇梯棋”）
- 數(shù)學定義：
$p(E_i|E_i,...,E_1)=p(E_i|E_{i-1})$
- 一個形象的比喻：假如每天的天氣是一個狀態(tài)的話，那今天是不是晴天只依賴于昨天的天氣而和前天的天氣沒有任何關系像屋！
- 非周期的馬爾可夫鏈具有收斂性質(zhì)：如果一個非周期的馬爾可夫鏈有狀態(tài)轉(zhuǎn)移矩陣P，且它的任何兩個狀態(tài)是連通的边篮，那么轉(zhuǎn)移n次后的狀態(tài)與初始狀態(tài)無關己莺！
  - 例^[1]：對于如下圖的股市模型：
    
    image.png
    
    共有三種狀態(tài)奏甫，定義狀態(tài)概率轉(zhuǎn)化矩陣為P，且P(i,j)=P(j|i)(從狀態(tài)i轉(zhuǎn)化到狀態(tài)j的概率)：
    
    $P= \begin{pmatrix} 0.9&0.075&0.025\\ 0.15&0.8&0.05\\ 0.25&0.25&0.5 \end{pmatrix}$
    我們會發(fā)現(xiàn)凌受，對于任意初始的股市概率阵子，假設初始概率分布：[0.3,0.4,0.3],代入轉(zhuǎn)移矩陣計算狀態(tài)2，3胜蛉，……挠进，n次狀態(tài)轉(zhuǎn)移（本例中約為6次）會發(fā)現(xiàn)股市概率分布穩(wěn)定在：[0.625,0.3125,0.0625]這組概率分布上（不論初始股市概率是什么）。
?
- 利用MCMC模型估計后驗分布形狀
  
  利用MCMC方法,有效從后驗分布種抽取樣本誊册，然后計算樣本中的平均值
  - 步驟：
    
    輸入選定的馬爾可夫鏈狀態(tài)轉(zhuǎn)移矩陣Q领突，設定狀態(tài)轉(zhuǎn)移次數(shù)閾值n₁ ，需要的樣本個數(shù)n₂ $\longrightarrow$ 從任意簡單概率分布取一個初始狀態(tài)值x₀（或自己給定） $\longrightarrow$ 重復X步驟n₁+n₂ -1次 $\longrightarrow$ 獲得需要的平穩(wěn)分布對應的樣本集( $x_{n1}, x_{n_1+1},...,x_{n_!+n_2+-1}$ ) $\longrightarrow$ 得到概率分布
    - X步驟：(t表示重復的次數(shù)):
      1. 從條件概率分布Q(x|x_t)中取得樣本x_*
      2. 均勻分布采樣u~uniform[0,1] （獲得一個隨機數(shù)u）
      3. 如果u<α（x_t ,x_*）(接受率案怯，α（x_t ,x~~）= $\pi$ (x~_)Q(x~, x_t)君旦；若為M-H分布，則α(i,j)=min{ $\frac{\pi (j)Q(j,i)}{\pi (i)Q(i,j)}$ ,1},提高了采樣效率)嘲碱，則接受轉(zhuǎn)移x_t $\rightarrow x_*$ 金砍，即x_t+1 =x~~
      4. 否則不接受轉(zhuǎn)移，即 $x_{t+1}=x_t$

利用MCMC模型估計后驗分布形狀

利用MCMC方法,有效從后驗分布種抽取樣本麦锯，然后計算樣本中的平均值

步驟：

輸入選定的馬爾可夫鏈狀態(tài)轉(zhuǎn)移矩陣Q恕稠，設定狀態(tài)轉(zhuǎn)移次數(shù)閾值n₁ ，需要的樣本個數(shù)n₂ $\longrightarrow$ 從任意簡單概率分布取一個初始狀態(tài)值x₀（或自己給定） $\longrightarrow$ 重復X步驟n₁+n₂ -1次
- X步驟：(t表示重復的次數(shù))
1. 從條件概率分布Q(x|x_t)中取得樣本x_*
2. 均勻分布采樣u~uniform[0,1] （獲得一個隨機數(shù)u）
3. 如果u<α（x_t ,x_*）(接受率扶欣，α（x_t ,x~~）= $\pi$ (x~_)Q(x~, x_t))谱俭，則接受轉(zhuǎn)移x_t $\rightarrow x_*$ ，即x_t+1 =x~~
4. 否則不接受轉(zhuǎn)移宵蛀，即 $x_{t+1}=x_t$
$\longrightarrow$ 獲得需要的平穩(wěn)分布對應的樣本集( $x_{n1}, x_{n_1+1},...,x_{n_!+n_2+-1}$ ) $\longrightarrow$ 得到概率分布
例：MCMC單參數(shù)運用——估計人類平均身高
1. 人類平均身高的后驗分布同時受到觀測曲線和先驗分布的影響昆著，假設人類平均身高的后驗分布如圖：
  
  image.png
2. 對單參數(shù)（x軸上的數(shù)值），MCMC法沿x軸開始隨機采樣术陶，在這里后驗概率影響Q凑懂，則：
  
  image.png
  
  （集中的點垂直疊加，便于觀察梧宫，事實上應該都是取在x軸上的）
3. 做直方圖接谨，獲取概率分布：
  
  image.png

REFERENCE

《告別數(shù)學公式，圖文解讀什么是馬爾可夫鏈蒙特卡羅方法（MCMC）》：https://zhuanlan.zhihu.com/p/32982140
《MCMC(三)MCMC采樣和M-H采樣》：https://www.cnblogs.com/pinard/p/6638955.html

來自wiki ?