原文連結：http://tecdat.cn/?p=23236

原文出處：拓端資料部落公眾號

什麼是頻率學派？

在頻率學派中，觀察樣本是隨機的，而引數是固定的、未知的數量。

概率被解釋為一個隨機過程的許多觀測的預期頻率。

有一種想法是 "真實的"，例如，在預測魚的生活環境時，鹽度和溫度之間的相互作用有一個迴歸係數？

什麼是貝葉斯學派？

在貝葉斯方法中，概率被解釋為對信念的主觀衡量。

所有的變數--因變數、引數和假設都是隨機變數。我們用資料來確定一個估計的確定性（可信度）。

這種鹽度X溫度的相互作用反映的不是絕對的，而是我們對魚的生活環境所瞭解的東西（本質上是草率的）。

目標

頻率學派

保證正確的誤差概率，同時考慮到抽樣、樣本大小和模型。

• 缺點：需要對置信區間、第一類和第二類錯誤進行復雜的解釋。

• 優點：更具有內在的 "客觀性 "和邏輯上的一致性。

貝葉斯學派

分析更多的資訊能在多大程度上提高我們對一個系統的認識。

• 缺點：這都是關於信仰的問題! ...有重大影響。

• 優點: 更直觀的解釋和實施，例如，這是這個假設的概率，這是這個引數等於這個值的概率。可能更接近於人類自然地解釋世界的方式。

實際應用中：為什麼用貝葉斯

• 具有有限資料的複雜模型，例如層次模型，其中

• 實際的先驗知識非常少

貝葉斯法則：

一些典型的貝葉斯速記法。

注意:

• 貝葉斯的最大問題在於確定先驗分佈。先驗應該是什麼？它有什麼影響？

目標:

計算引數的後驗分佈：π（θ|X）。

點估計是後驗的平均值。

一個可信的區間是

你可以把它解釋為一個引數在這個區間內的概率 。

計算

皮埃爾-西蒙-拉普拉斯（1749-1827）（見：Sharon Bertsch McGrayne: The Theory That Would Not Die)

• 有些問題是可分析的，例如二項式似然-貝塔先驗。
  • 如果你有幾個引數，而且是奇數分佈，你可以用數值乘以/整合先驗和似然（又稱網格近似）。
    • 但如果你有很多引數，這是不可能完成的操作
• 儘管該理論可以追溯到1700年，甚至它對推理的解釋也可以追溯到19世紀初，但它一直難以更廣泛地實施，直到馬爾科夫鏈蒙特卡洛技術的發展。

MCMC

MCMC的思想是對引數值θi進行 "抽樣"。

回顧一下，馬爾科夫鏈是一個隨機過程，它只取決於它的前一個狀態，而且（如果是遍歷的），會生成一個平穩的分佈。

技巧 "是找到漸進地接近正確分佈的抽樣規則（MCMC演算法）。

有幾種這樣的（相關）演算法。

• Metropolis-Hastings抽樣
• Gibbs 抽樣
• No U-Turn Sampling (NUTS)
• Reversible Jump

一個不斷髮展的文獻和工作體系!

Metropolis-Hastings 演算法

1. 開始:
2. 跳到一個新的候選位置:
3. 計算後驗:
4. 如果
5. 如果
6. 轉到第2步

Metropolis-Hastings: 硬幣例子

你丟擲了5個正面。你對θ的最初 "猜測 "是

MCMC:

2. p.old <- prior *likelihood
3. while(length(thetas) <= n){
4. theta.new <- theta + rnorm(1,0,0.05)
5. p.new <- prior *likelihood
6. if(p.new > p.old | runif(1) < p.new/p.old){
7. theta <- theta.new
8. p.old <- p.new
9. }

畫圖:

1. hist(thetas[-(1:100)] )
2. curve(6*x^5 )

取樣鏈：調整、細化、多鏈

• 那個 "朝向 "平穩的初始過渡被稱為 "預燒期"，必須加以修整。
  • 怎麼做？用眼睛看
• 取樣過程（顯然）是自相關的。
  • 如何做？通常是用眼看，用acf()作為指導。
    

• 為了保證你收斂到正確的分佈，你通常會從不同的位置獲得多條鏈（例如4條）。

• 有效樣本量

MCMC 診斷法

R軟體包幫助分析MCMC鏈。一個例子是線性迴歸的貝葉斯擬合（α,β,σ

plot(line)

預燒部分:

plot(line[[1]], start=10)

MCMC診斷法

檢視後驗分佈（同時評估收斂性）。

density(line)

引數之間的關聯性，以及鏈內的自相關關係

1. levelplot(line[[2]])
2. acfplot(line)

統計摘要

執行MCMC的工具（在R內部）

邏輯Logistic迴歸：嬰兒出生體重低

logitmcmc(low~age+as.factor(race)+smoke )

plot(mcmc)

MCMC與GLM邏輯迴歸的比較

MCMC與GLM邏輯迴歸的比較

對於這個應用，沒有很好的理由使用貝葉斯建模，除非--你是 "貝葉斯主義者"。 你有關於迴歸係數的真正先驗資訊（這基本上是不太可能的）。

一個主要的缺點是 先驗分佈棘手的調整引數。

但是，MCMC可以擬合的一些更復雜的模型（例如，層次的logit MCMChlogit）。

Metropolis-Hastings

Metropolis-Hastings很好，很簡單，很普遍。但是對迴圈次數很敏感。而且可能太慢，因為它最終會拒絕大量的迴圈。

Gibbs 取樣

在Gibbs吉布斯抽樣中，你不是用適當的概率接受/拒絕，而是用適當的條件概率在引數空間中行進。 並從該分佈中抽取一次。

然後你從新的條件分佈中抽取下一個引數。

比Metropolis-Hastings快得多。有效樣本量要高得多!

BUGS（OpenBUGS，WinBUGS）是使用吉布斯取樣器的貝葉斯推理。

JAGS是 "吉布斯取樣器"

其他取樣器

漢密爾頓蒙特卡洛（HMC）--是一種梯度的Metropolis-Hastings，因此速度更快，對引數之間的關聯性更好。

No-U Turn Sampler（NUTS）--由於不需要固定的長度，它的速度更快。這是STAN使用的方法（見http://arxiv.org/pdf/1111.4246v1.pdf）。

(Hoffman and Gelman 2011)

其他工具

你可能想建立你自己的模型，使用貝葉斯MC進行擬合，而不是依賴現有的模型。為此，有幾個工具可以選擇。

• BUGS / WinBUGS / OpenBUGS (Bayesian inference Using Gibbs Sampling) - 貝葉斯抽樣工具的鼻祖（自1989年起）。WinBUGS是專有的。OpenBUGS的支援率很低。

• JAGS（Just Another Gibbs Sampler）接受一個用類似於R語言的語法編寫的模型字串，並使用吉布斯抽樣從這個模型中編譯和生成MCMC樣本。可以在R中使用rjags包。

• Stan（以Stanislaw Ulam命名）是一個類似於JAGS的相當新的程式--速度更快，更強大，發展迅速。從偽R/C語法生成C++程式碼。安裝：http://mc-stan.org/rstan.html**

• Laplace’s Demon所有的貝葉斯工具都在R中： http://www.bayesian-inference.com/software

STAN

要用STAN擬合一個模型，步驟是:

1. 為模型生成一個STAN語法虛擬碼（在JAGS和BUGS中相同
2. 執行一個R命令，用C++語言編譯該模型
3. 使用生成的函式來擬合你的資料

STAN示例--線性迴歸

STAN程式碼是R（例如，具有分佈函式）和C（即你必須宣告你的變數）之間的一種混合。每個模型定義都有三個塊。

1.資料塊:

1. int n; //
2. vector[n] y; // Y 向量

這指定了你要輸入的原始資料。在本例中，只有Y和X，它們都是長度為n的（數字）向量，是一個不能小於0的整數。

2. 引數塊

  real beta1;  // slope

這些列出了你要估計的引數：截距、斜率和方差。

3. 模型塊

1. sigma ~ inv_gamma(0.001, 0.001);
3. yhat[i] <- beta0 + beta1 * (x[i] - mean(x));}
4. y ~ normal(yhat, sigma);

注意:

• 你可以向量化，但迴圈也同樣快
• 有許多分佈（和 "平均值 "等函式）可用

請經常參閱手冊！https://github.com/stan-dev/stan/releases/download/v2.9.0/stan-reference-2.9.0.pdf

2. 在R中編譯模型

你把你的模型儲存在一個單獨的檔案中， 然後用stan_model()命令編譯這個模型。

這個命令是把你描述的模型，用C++編碼和編譯一個NUTS取樣器。相信我，自己編寫C++程式碼是一件非常非常痛苦的事情（如果沒有很多經驗的話），而且它保證比R中的同等程式碼快得多。

注意：這一步可能會很慢。

3. 在R中執行該模型

這裡的關鍵函式是sampling()。還要注意的是，為了給你的模型提供資料，它必須是列表的形式

模擬一些資料。

2. X <- runif(100,0,20)
3. Y <- rnorm(100, beta0+beta1*X, sigma)

進行取樣!

sampling(stan, Data)

這裡有大量的輸出，因為它計算了

print(fit, digits = 2)

MCMC診斷法

為了應用coda系列的診斷工具，你需要從STAN擬合物件中提取鏈，並將其重新建立為mcmc.list。

1. extract(stan.fit
2. alply(chains, 2, mcmc)

最受歡迎的見解

1.matlab使用貝葉斯優化的深度學習

2.matlab貝葉斯隱馬爾可夫hmm模型實現

3.R語言Gibbs抽樣的貝葉斯簡單線性迴歸模擬

4.R語言中的block Gibbs吉布斯取樣貝葉斯多元線性迴歸

5.R語言中的Stan概率程式設計MCMC取樣的貝葉斯模型

6.Python用PyMC3實現貝葉斯線性迴歸模型

7.R語言使用貝葉斯 層次模型進行空間資料分析

8.R語言隨機搜尋變數選擇SSVS估計貝葉斯向量自迴歸（BVAR）模型

9.matlab貝葉斯隱馬爾可夫hmm模型實現