R語言coda貝葉斯MCMC Metropolis-Hastings采樣鏈分析和收斂診斷可視化
原文鏈接：http://tecdat.cn/?p=27228?
原文出處：拓端數據部落公眾號
?
作為先決條件，我們將使用幾行代碼，在代碼中，我們創(chuàng)建了一些測試數據，其中因變量 y 線性依賴于自變量 x（預測變量）；定義線性模型擬合數據的可能性和先驗；并實現一個簡單的 Metropolis-Hastings MCMC 從該模型的后驗分布中采樣。
x = (-(sleze-1)/2):((sple-1)/2)y =? treA * x + tuB + rnorm(n=sapeize,mean=0,sd=tuSd)
所以，讓我們運行 MCMC：
stavalue = c(4,2,8)cn = rmtrisMCC(avae, 10000)
由 coda 促成的鏈的一些簡單總結
好吧，coda 是一個 R 包，它提供了許多用于繪制和分析后驗樣本的標準函數。為了使這些功能起作用，您需要將輸出作為“mcmc”或“mcmc.list”類的對象，我們將在后面討論。
擁有一個 coda 對象的好處是我們通常想要用鏈做的很多事情都已經實現了，所以例如我們可以簡單地 summary() 和 plot() 輸出
summary(chn)plot(cn)
它提供了一些關于控制臺的有用信息和一個大致如下所示的圖：
圖：?一個 coda 對象的 plot() 函數的結果
對 plot() 函數的結果：每一行對應一個參數，因此每個參數有兩個圖。左邊的圖稱為軌跡圖——它顯示了參數在鏈運行時所取的值。右圖通常稱為邊際密度圖。基本上，它是軌跡圖中值的（平滑的）直方圖，即參數值在鏈中的分布。
邊際密度隱藏了相關性
邊際密度是參數取值與所有其他“邊緣化”參數的平均值，即其他參數根據其后驗概率具有任何值。通常，邊際密度被視為貝葉斯分析的主要輸出（例如，通過報告它們的均值和標準差），但我強烈建議不要進一步分析這種做法。原因是邊際密度“隱藏”了參數之間的相關性，如果存在相關性，參數的不確定性在邊際中似乎要大得多。
Plot(data.frame(can))
在我們的例子中，應該沒有大的相關性，因為我以這種方式設置了示例
x = (-(samee-1)/2):((smeie-1)/2) + 20
再次運行 MCMC 并檢查相關性應該會給你一個完全不同的畫面。
圖：?不平衡 x 值擬合的邊際密度（對角線）、配對密度（下圖）和相關系數（上圖）
您可以看到第一個和第二個參數（斜率和截距）之間的強相關性，并且您還可以看到每個參數 X2 的邊際不確定性增加了。
請注意，我們在這里只檢查了配對相關性，可能仍然有更高階的交互不會出現在這樣的分析中，所以你可能仍然遺漏了一些東西。
收斂診斷
現在，到收斂：一個 MCMC 從后驗分布創(chuàng)建一個樣本，我們通常想知道這個樣本是否足夠接近后驗以用于分析。有幾種標準方法可以檢查這一點，但我建議使用 Gelman-Rubin 診斷。
cha2 = runmeooi_MCMC(arvue, 10000)cominchns = mcmc.list(cai, ain2)plot(coinchns)
結果圖應該是這樣的
圖：?結果
diag 為您提供每個參數的比例縮減因子。因子 1 意味著方差和鏈內方差相等，較大的值意味著鏈之間仍然存在顯著差異。
改善收斂/混合
那么，如果還沒有收斂怎么辦？當然，你總是可以讓 MCMC 運行更長時間，但另一個選擇是讓它收斂得更快?？赡軙l(fā)生兩件事：
與我們從中抽樣的分布相比，您的提議函數很窄——接受率高，但我們沒有得到任何結果，混合不好
與我們從中抽樣的分布相比，您的提議函數太寬了——接受率低，大部分時間我們都呆在原地
最受歡迎的見解
1.使用R語言進行METROPLIS-IN-GIBBS采樣和MCMC運行
2.R語言中的Stan概率編程MCMC采樣的貝葉斯模型
3.R語言實現MCMC中的Metropolis–Hastings算法與吉布斯采樣
4.R語言BUGS JAGS貝葉斯分析 馬爾科夫鏈蒙特卡洛方法（MCMC）采樣
5.R語言中的block Gibbs吉布斯采樣貝葉斯多元線性回歸
6.R語言Gibbs抽樣的貝葉斯簡單線性回歸仿真分析
7.R語言用Rcpp加速Metropolis-Hastings抽樣估計貝葉斯邏輯回歸模型的參數
8.R語言使用Metropolis- Hasting抽樣算法進行邏輯回歸
9.R語言中基于混合數據抽樣(MIDAS)回歸的HAR-RV模型預測GDP增長
