拓端tecdat：R語言RStan MCMC：NUTS取樣演算法用LASSO 構建貝葉斯線性迴歸模型分析職業聲望資料
原文連結：http://tecdat.cn/?p=24456
原文出處：拓端資料部落公眾號
如果你正在進行統計分析：想要加一些先驗資訊，最終你想要的是預測。所以你決定使用貝葉斯。
但是，你沒有共軛先驗。你可能會花費很長時間編寫 Metropolis-Hastings 程式碼，優化接受率和提議分佈，或者你可以使用 RStan。
Hamiltonian Monte Carlo（HMC）
HMC 是一種為 MH 演算法生成提議分佈的方法，該提議分佈被接受的概率很高。具體演算法過程請檢視參考文獻。
打個比方：
給粒子一些動量。
它在滑冰場周圍滑行，大部分時間都在密度高的地方。
拍攝這條軌跡的快照為後驗分佈提供了一個建議樣本。
然後我們使用 Metropolis-Hastings 進行校正。
NUTS取樣器（No-U-turn Sampler）
HMC，像RWMH一樣，需要對步驟的數量和大小進行一些調整。
No-U-Turn Sampler "或NUTs（Hoffman和Gelman（2014）），對這些進行了自適應的優化。
NUTS建立了一組可能的候選點，並在軌跡開始自相矛盾時立即停止。
Stan 的優點
可以產生高維度的提議，這些提議被接受的概率很高，而不需要花時間進行調整。
有內建的診斷程式來分析MCMC的輸出。
在C++中構建，所以執行迅速，輸出到R。
示例
如何使用 LASSO 構建貝葉斯線性迴歸模型。
構建 Stan 模型
資料：n、p、Y、X 先驗引數，超引數
引數：
模型：高斯似然、拉普拉斯和伽瑪先驗。
輸出：後驗樣本，後驗預測樣本。
資料
- int<lwer=0> n;
- vectr[n] y;
- rel<loer=0> a;
引數
- vetor[p+1] beta;
- real<lowr=0> siga;
轉換後的引數（可選）
- vectr[n] liped;
- lnpred = X*bea;
模型
- bta ~ dolexneial(0,w);
- siga ~ gama(a,b);
或沒有向量化，
- for(i in 1:n){
- y[i]~noral(X[i,]*beta,siga);
- }
生成的數量（可選）
- vecor[n] yprict;
- for(i in 1:n){
- prdit[i] = nrmlrng(lnprd[i],siga);
對後驗樣本的每一個元素都要評估一次這個程式碼。
職業聲望資料集
這裡我們使用職業聲望資料集，它有以下變數
教育：職業在職者的平均教育程度，年。
收入：在職者的平均收入，元。
女性：在職者中女性的百分比。
威望：Pineo-Porter的職業聲望得分，來自一項社會調查。
普查：人口普查的職業程式碼。
型別：職業的型別
bc: 藍領
prof: 專業、管理和技術
wc: 白領
在R中執行
- library(rstan)
- stan(file="byLASO",iter=50000)
在3.5秒內執行25000次預熱和25000次取樣。
第一次編譯c++程式碼，所以可能需要更長的時間。
繪製後驗分佈圖
- par(mrow=c(1,2))
- plot(denty(prs$bea)
預測分佈
- plot(density)
鏈診斷
splas[[1][1:5,]
鏈診斷
trac("beta" )鏈診斷
pa(pars="beta")更多鏈診斷
Stan 還可以從鏈中提取各種其他診斷，如置信區間、有效樣本量和馬爾可夫鏈平方誤差。
鏈的值與各種鏈屬性、對數似然、接受率和步長之間的比較圖。
Stan 出錯
stan使用的步驟太大。
可以通過手動增加期望的平均接受度來解決。
adapt_delta，高於其預設的0.8
stan(cntl = list(datta = 0.99, mxrh = 15))
這會減慢你的鏈的速度，但可能會產生更好的樣本。
自制函式
Stan 也相容自制函式。
如果你的先驗或似然函式不標準，則很有用。
- model {
- beta ~ doubexp(0,w);
- for(i in 1:n){
- logprb(‐0.5*fs(1‐(exp(normalog(
- siga))/yde));
- }
- }
結論
不要浪費時間編碼和調整 RWMH.
Stan 執行得更快，會自動調整，並且應該會產生較好的樣本。
參考文獻
Alder, Berni J, and T E Wainwright. 1959. “Studies in Molecular Dynamics. I. General Method.” The Journal of Chemical Physics 31 (2). AIP: 459–66.
Hoffman, Matthew D, and Andrew Gelman. 2014. “The No-U-Turn Sampler: Adaptively Setting Path Lengths in Hamiltonian Monte Carlo.” Journal of Machine Learning Research 15 (1): 1593–1623.
最受歡迎的見解
4.R語言中的block Gibbs吉布斯取樣貝葉斯多元線性迴歸
8.R語言隨機搜尋變數選擇SSVS估計貝葉斯向量自迴歸（BVAR）模型
▍關注我們 【大資料部落】第三方資料服務提供商,提供全面的統計分析與資料探勘諮詢服務,為客戶定製個性化的資料解決方案與行業報告等。 ▍諮詢連結：http://y0.cn/teradat ▍聯絡郵箱：[email protected]