原文鏈接：http://tecdat.cn/?p=9390

介紹

向量自回歸（VAR）模型的一般缺點(diǎn)是，估計(jì)系數(shù)的數(shù)量與滯后的數(shù)量成比例地增加。因此，隨著滯后次數(shù)的增加，每個(gè)參數(shù)可用的信息較少。在貝葉斯VAR文獻(xiàn)中，減輕這種所謂的維數(shù)詛咒的一種方法是隨機(jī)搜索變量選擇（SSVS），由George等人提出（2008）。SSVS的基本思想是將通常使用的先驗(yàn)方差分配給應(yīng)包含在模型中的參數(shù)，將不相關(guān)參數(shù)的先驗(yàn)方差接近零。這樣，通常就可以估算出相關(guān)參數(shù)，并且無關(guān)變量的后驗(yàn)值接近于零，因此它們對(duì)預(yù)測(cè)和沖激響應(yīng)沒有顯著影響。這是通過在模型之前添加層次結(jié)構(gòu)來實(shí)現(xiàn)的，其中在采樣算法的每個(gè)步驟中評(píng)估變量的相關(guān)性。

這篇文章介紹了使用SSVS估計(jì)貝葉斯向量自回歸（BVAR）模型。它使用Lütkepohl（2007）的數(shù)據(jù)集E1，其中包含有關(guān)1960Q1至1982Q4德國(guó)固定投資，可支配收入和消費(fèi)支出的數(shù)據(jù)。加載數(shù)據(jù)并生成數(shù)據(jù)：

1. 


2. # 加載和轉(zhuǎn)換數(shù)據(jù)

3. e1 <- diff(log(e1))

4. 


5. # 生成VAR

6. data <- gen_var(e1, p = 4, deterministic = "const")

7. 


8. # 獲取數(shù)據(jù)矩陣

9. y <- data$Y[, 1:71]

10. x <- data$Z[, 1:71]

估算值

根據(jù)George等人所述的半自動(dòng)方法來設(shè)置參數(shù)的先驗(yàn)方差（2008）。對(duì)于所有變量，先驗(yàn)包含概率設(shè)置為0.5。誤差方差-協(xié)方差矩陣的先驗(yàn)信息不足。

1. # 重置隨機(jī)數(shù)提高可重復(fù)性

2. set.seed(1234567)

3. 


4. t <- ncol(y) # 觀察數(shù)

5. k <- nrow(y) # 內(nèi)生變量數(shù)

6. m <- k * nrow(x) # 估計(jì)系數(shù)數(shù)

7. 


8. # 系數(shù)先驗(yàn)

9. a_mu_prior <- matrix(0, m) # 先驗(yàn)均值的向量

10. 


11. # SSVS先驗(yàn)（半自動(dòng)方法）

12. ols <- tcrossprod(y, x) %*% solve(tcrossprod(x)) # OLS估計(jì)

13. sigma_ols <- tcrossprod(y - ols %*% x) / (t - nrow(x)) # OLS誤差協(xié)方差矩陣

14. cov_ols <- kronecker(solve(tcrossprod(x)), sigma_ols)

15. se_ols <- matrix(sqrt(diag(cov_ols))) # OLS標(biāo)準(zhǔn)誤

16. 


17. 


18. # 先驗(yàn)參數(shù)

19. prob_prior <- matrix(0.5, m)

20. 


21. # ?方差-協(xié)方差矩陣

22. u_sigma_df_prior <- 0 # 方差-協(xié)方差矩陣

23. u_sigma_scale_prior <- diag(0, k) # 先驗(yàn)協(xié)方差矩陣

24. u_sigma_df_post <- t + u_sigma_df_prior # 后驗(yàn)自由度

初始參數(shù)值設(shè)置為零，這意味著在Gibbs采樣器的第一步中應(yīng)相對(duì)自由地估算所有參數(shù)。

可以直接將SSVS添加到VAR模型的標(biāo)準(zhǔn)Gibbs采樣器算法中。在此示例中，常數(shù)項(xiàng)從SSVS中排除，這可以通過指定來實(shí)現(xiàn)include = 1:36。具有SSVS的Gibbs采樣器的輸出可以用通常的方式進(jìn)一步分析。因此，可以通過計(jì)算參數(shù)的繪制方式獲得點(diǎn)估計(jì)：

1. 


2. ## ? ? ? ? ?invest income ? cons

3. ## invest.1 -0.102 ?0.011 -0.002

4. ## income.1 ?0.044 -0.031 ?0.168

5. ## cons.1 ? ?0.074 ?0.140 -0.287

6. ## invest.2 -0.013 ?0.002 ?0.004

7. ## income.2 ?0.015 ?0.004 ?0.315

8. ## cons.2 ? ?0.027 -0.001 ?0.006

9. ## invest.3 ?0.033 ?0.000 ?0.000

10. ## income.3 -0.008 ?0.021 ?0.013

11. ## cons.3 ? -0.043 ?0.007 ?0.019

12. ## invest.4 ?0.250 ?0.001 -0.005

13. ## income.4 -0.064 -0.010 ?0.025

14. ## cons.4 ? -0.023 ?0.001 ?0.000

15. ## const ? ? 0.014 ?0.017 ?0.014

還可以通過計(jì)算變量的均值來獲得每個(gè)變量的后驗(yàn)概率。從下面的輸出中可以看出，在VAR（4）模型中似乎只有幾個(gè)變量是相關(guān)的。常數(shù)項(xiàng)的概率為100％，因?yàn)樗鼈円褟腟SVS中排除。

1. 


2. ## ? ? ? ? ?invest income cons

3. ## invest.1 ? 0.43 ? 0.23 0.10

4. ## income.1 ? 0.10 ? 0.18 0.67

5. ## cons.1 ? ? 0.11 ? 0.40 0.77

6. ## invest.2 ? 0.11 ? 0.09 0.14

7. ## income.2 ? 0.08 ? 0.07 0.98

8. ## cons.2 ? ? 0.07 ? 0.06 0.08

9. ## invest.3 ? 0.19 ? 0.07 0.06

10. ## income.3 ? 0.06 ? 0.13 0.10

11. ## cons.3 ? ? 0.09 ? 0.07 0.12

12. ## invest.4 ? 0.78 ? 0.09 0.16

13. ## income.4 ? 0.13 ? 0.09 0.18

14. ## cons.4 ? ? 0.09 ? 0.07 0.06

15. ## const ? ? ?1.00 ? 1.00 1.00

給定這些值，研究人員可以按照常規(guī)方式進(jìn)行操作，并根據(jù)Gibbs采樣器的輸出獲得預(yù)測(cè)和脈沖響應(yīng)。這種方法的優(yōu)勢(shì)在于它不僅考慮了參數(shù)不確定性，而且還考慮了模型不確定性。這可以通過系數(shù)的直方圖來說明，該直方圖描述了收入的第一個(gè)滯后項(xiàng)與消費(fèi)當(dāng)前值之間的關(guān)系。

hist(draws_a[6,],

?

通過兩個(gè)峰描述模型不確定性，并通過右峰在它們周圍的分布來描述參數(shù)不確定性。

但是，如果研究人員不希望使用模型，變量的相關(guān)性可能會(huì)從采樣算法的一個(gè)步驟更改為另一個(gè)步驟，那么另一種方法將是僅使用高概率的模型。這可以通過進(jìn)一步的模擬來完成，在該模擬中，對(duì)于不相關(guān)的變量使用非常嚴(yán)格的先驗(yàn)，而對(duì)于相關(guān)參數(shù)則使用沒有信息的先驗(yàn)。

后方抽取的均值類似于Lütkepohl（2007，5.2.10節(jié)）中的OLS估計(jì)值：

1. 


2. ## ? ? ? ? ?invest income ? cons

3. ## invest.1 -0.219 ?0.001 -0.001

4. ## income.1 ?0.000 ?0.000 ?0.262

5. ## cons.1 ? ?0.000 ?0.238 -0.334

6. ## invest.2 ?0.000 ?0.000 ?0.001

7. ## income.2 ?0.000 ?0.000 ?0.329

8. ## cons.2 ? ?0.000 ?0.000 ?0.000

9. ## invest.3 ?0.000 ?0.000 ?0.000

10. ## income.3 ?0.000 ?0.000 ?0.000

11. ## cons.3 ? ?0.000 ?0.000 ?0.000

12. ## invest.4 ?0.328 ?0.000 -0.001

13. ## income.4 ?0.000 ?0.000 ?0.000

14. ## cons.4 ? ?0.000 ?0.000 ?0.000

15. ## const ? ? 0.015 ?0.015 ?0.014

評(píng)價(jià)

bvar功能可用于將Gibbs采樣器的相關(guān)輸出收集到標(biāo)準(zhǔn)化對(duì)象中，例如predict獲得預(yù)測(cè)或irf進(jìn)行脈沖響應(yīng)分析。

?

hin(bvar_est, thin = 5)

預(yù)測(cè)

可以使用函數(shù)獲得置信區(qū)間的預(yù)測(cè)predict。

plot(bvar_pred)

?

脈沖響應(yīng)分析

1. 


2. plot(OIR

?

最受歡迎的見解

1.matlab使用貝葉斯優(yōu)化的深度學(xué)習(xí)

2.matlab貝葉斯隱馬爾可夫hmm模型實(shí)現(xiàn)

3.R語言Gibbs抽樣的貝葉斯簡(jiǎn)單線性回歸仿真

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

5.R語言中的Stan概率編程MCMC采樣的貝葉斯模型

6.Python用PyMC3實(shí)現(xiàn)貝葉斯線性回歸模型

7.R語言使用貝葉斯 層次模型進(jìn)行空間數(shù)據(jù)分析

8.R語言隨機(jī)搜索變量選擇SSVS估計(jì)貝葉斯向量自回歸（BVAR）模型

9.matlab貝葉斯隱馬爾可夫hmm模型實(shí)現(xiàn)

?