凡是搞计量经济的，都关注这个号了
稿件：econometrics666@sina.cn
所有计量经济圈方法论丛的code程序, 宏微观数据库和各种软件都放在社群里.欢迎到计量经济圈社群交流访问.
关于VAR方法，我们引荐了1.R软件中的时间序列分析程序包纵览，2.时间序列分析的各种程序, 38页集结整理成文档，3.时间序列数据分析的思维导图一览, 金融经济学者必备工具，4.送书: 应用时间序列分析(经典)，5.为啥时间序列模型比较难学？时间序列的正名路，6.时间序列中的协整检验和VECM，以及回归后的系列估计操作，7.时间序列模型分解，季节调整分析基础，8.空间和时间的计量，关注二位国人，9.TVP-VAR时变参数VAR系列文献和估计程序，10.向量自回归VAR模型操作指南针,为微观面板VAR铺基石，11.VAR宏观计量模型演进与发展，无方向确认推断更好，12.应用VAR模型时的15个注意点，总结得相当地道，13.面板数据单位根检验软件操作和解读全在这里，14.动态面板回归和软件操作，单位根和协整检验（Dynamic Panel Data），15.面板向量自回归PVAR是什么? 数据, 程序和解读一步到位，16.ARDL, ARIMA, VAR, (G)ARCH时间数据模型讲解及软件操作，17.动态因子模型是什么, 又怎么去实现？18.SVAR模型的起源、识别、估计与应用, 系统讲述，19.平滑转移自回归模型(STAR)应用与在R软件的操作，20.Copula函数，21.GVAR, 全局VAR模型是什么？该如何用软件实现, 有哪些研究文献和最新进展！
关于下方文字内容，作者：韦忠吉，清华大学经管学院，通信邮箱：weizhj.18@sem.tsinghua.edu.cn
Nikolas Kuschnig, Lukas Vashold. BVAR: Bayesian Vector Autoregressions with Hierarchical Prior Selection in RVector autoregression (VAR) models are widely used models for multivariate time series analysis, but often suffer from their dense parameterization. Bayesian methods are commonly employed as a remedy by imposing shrinkage on the model coefficients via informative priors, thereby reducing parameter uncertainty. The subjective choice of the informativeness of these priors is often criticized and can be alleviated via hierarchical modeling. This paper introduces BVAR, an R package dedicated to the estimation of Bayesian VAR models in a hierarchical fashion. It incorporates functionalities that permit addressing a wide range of research problems while retaining an easy-to-use and transparent interface. It features the most commonly used priors in the context of multivariate time series analysis as well as an extensive set of standard methods for analysis. Further functionalities include a framework for defining custom dummy-observation priors, the computation of impulse response functions, forecast error variance decompositions and forecasts.
VAR想必大家都知道，那么BVAR是什么呢？BVAR，简单来说，其实就是用贝叶斯统计的方法实现了VAR回归的功能。Nikolas Kuschnig和Lukas Vashold（2019）这篇工作论文，BVAR: Bayesian Vector Autoregressions with Hierarchical Prior Selection in R，将为我们答疑解惑，并教你如何在R中使用它。
VAR（Vector AutoRegression，向量自回归）被广泛运用于多个学科的多变量的时间序列分析中，例如Crespo Cuaresma, Feldkircher, 和Huber （2016）等。通过Sims（1980）的推广，VAR已经成了宏观实证研究的主要工具(Kilian and Lütkepohl，2017)。但是，大量的参数和宏观经济数据有限的时间可用性，可能导致过度参数化（over-parameterization）问题（Koop和Korobilis 2010）。因此，我们引入了BVAR（贝叶斯向量自回归）模型。信息性先验（informative priors）为模型添加了额外的结构，并将其缩小为简洁的、经过验证的基准。这样得到的模型，能降低参数不确定性并显著增强的样本外预测性能(Koop 2013)。这些先验的选择和他们的信息构成了一个挑战，仍然是讨论和批评的支点。Giannone，Lenza和Primiceri(2015)从层次建模思想出发，提出了基于数据的且有理论基础的先验信息设置方法。它们减轻了设定超参数（hyperparameters）的主观性，并在常见的分析中展示了出色的性能。BVAR是实现这些分层贝叶斯VAR模型的第一个R包(hierarchical Bayesian VAR)，它为评估和分析提供了一个完整且易于使用的工具包。下面让我们一起走近BVAR吧！我们首先介绍贝叶斯统计软件的发展现状。不过，在介绍贝叶斯统计软件的现状之前，先来简单介绍一下贝叶斯统计。英国学者托马斯·贝叶斯在《论有关机遇问题的求解》中提出一种归纳推理的理论，后被一些统计学者发展为一种系统的统计推断方法，称为贝叶斯方法。采用这种方法作统计推断所得的全部结果，构成贝叶斯统计的内容。贝叶斯统计的原理，主要依赖先验分布和后验分布两个概念。先验分布是总体分布参数θ的一个概率分布。贝叶斯学派的根本观点，是认为在关于θ的任何统计推断问题中，除了使用样本X所提供的信息外，还必须对θ规定一个先验分布，它是在进行推断时不可或缺的一个要素。贝叶斯学派把先验分布解释为在抽样前就有的关于θ的先验信息的概率表述，先验分布不必有客观的依据，它可以部分地或完全地基于主观信念。例如，某甲怀疑自己患有一种疾病A，在就诊时医生对他测了诸如体温、血压等指标，其结果构成样本X。引进参数θ：有病时，θ=1；无病时，θ=0。X的分布取决于θ是0还是1，因而知道了X有助于推断θ是否为1。按传统（频率）学派的观点，医生诊断时，只使用X提供的信息；而按贝叶斯学派观点，则认为只有在规定了一个介于0与1之间的数p作为事件{θ=1}的先验概率时，才能对甲是否有病（即θ是否为1）进行推断。p这个数刻画了本问题的先验分布，且可解释为疾病A的发病率。先验分布的规定对推断结果有影响，如在此例中，若疾病A的发病率很小，医生将倾向于只有在样本X显示出很强的证据时，才诊断甲有病。在这里先验分布的使用看来是合理的，但贝叶斯学派并不是基于 “p是发病率”这样一个解释而使用它的，事实上即使对本病的发病率毫无所知，也必须规定这样一个p，否则问题就无法求解。根据样本 X 的分布Pθ及θ的先验分布π（θ），用概率论中求条件概率分布的方法，可算出在已知X=x的条件下，θ的条件分布 π（θ|x）。因为这个分布是在抽样以后才得到的，故称为后验分布。贝叶斯学派认为：这个分布综合了样本X及先验分布π（θ）所提供的有关的信息。抽样的全部目的，就在于完成由先验分布到后验分布的转换。如上例，设p=P（θ=1)=0.001，而π（θ=1|x)=0.86，则贝叶斯学派解释为：在某甲的指标量出之前，他患病的可能性定为0.001，而在得到X后，认识发生了变化：其患病的可能性提高为0.86，这一点的实现既与X有关，也离不开先验分布。计算后验分布的公式本质上就是概率论中著名的贝叶斯公式：
这公式正是上面提到的贝叶斯1763年的文章的一个重要内容。贝叶斯推断方法的关键在于所作出的任何推断都必须也只须根据后验分布π（θ│X），而不能再涉及X的样本分布Pθ。贝叶斯统计软件是一个相当活跃的领域，日益增长的算力，促进了马尔科夫链蒙特卡洛（Markov chain Monte Carlo，MCMC)方法的发展。已有软件建立在Gibbs抽样——Metropolis-Hastings算法的一种变体之上，包括了BUGS（Lunn, Spiegelhalter, Thomas, 和Best，2009等)和JAGS（Plummer，2003）。此外还有Stan（Guo, Li和Riddell，2017等)，一种统计模型的命令式概率编程语言，采用无U型采样器（No-U-Turn Sampler，汉密尔顿蒙特卡洛方法的一种变体，见Hoffman and Gelman 2014）来提升鲁棒性和效率。这些包提供了灵活可扩展的贝叶斯推断工具，并且能够在几个平台之间跨平台交互（例如R和Python）。其他软件包提供了类似或更具体的实现，包括Python的pymc3 (Salvatier, Wiecki和Fonnesbeck 2016)， MCMCglmm (Hadfield 2010)， greta (Golding et al. 2018)和R的bvarsv (Krueger 2015)。集成嵌套拉普拉斯逼近（Integrated nested Laplace approximation，见Rue, Martino, and Chopin 2009)是一种具有显著计算优势的近似贝叶斯推断的替代方法，广泛应用于空间建模中(见Blangiardo, Cameletti, Baio和 Rue 2013)，在R- INLA (Rue, Martino, Lindgren et al. 2015)中有一个R实现。尽管有各种各样的贝叶斯推理软件，但没有专门用于贝叶斯VAR建模的综合选项。vars (Pfaff 2008)包代表了频率学多变量时间序列分析领域的基石，提供了一套完整的var相关功能。然而，尽管贝叶斯VAR模型很流行，但目前还没有相同的模型。(没办法，谁让频率学派才是主流呢)前面提到的包涵盖了各种可能的应用程序，并为建模和后续分析提供了强大而高效的工具。然而，贝叶斯变量的工作通常是通过MATLAB和R脚本完成的。因此，代码被循环利用并频繁地拼凑在一起，而没有任何中央存储库或版本控制促进重现性(Ram 2013)，而许多其他脚本只使用过一次。在本文中，我们提出了一个R软件包BVAR，它协调了更完整的贝叶斯推理工具包的优点和特别脚本的特殊性。该包实现了Giannoneet等(2015)提出的层级模型（hierarchical modeling）的方法，选择常用先验的信息量。当前版本实现了一系列常见分析的功能，并允许使用现有的框架，例如用于MCMC诊断的coda命令 (Plummer, Best, Cowles和Vines 2006)。絮絮叨叨这么多，下面我们要正式介绍今天的主角——BVAR了！我们会先简单地描述它背后的计量原理，然后再介绍大家在实证中更感兴趣的实操环节，把BVAR这个包安利给大家，并给一个实例说说BVAR到底怎么用。在这个部分，我们将介绍一下贝叶斯背景下的VAR模型，并根据Giannone等人（2015勾勒出实现分层先验选择的方法。对于VAR模型的贝叶斯估计的细节解释，我们将参照Koop和Korobilis(2010)以及Kilian和Lutkepohl(2017)。VAR模型是单变量自回归（AR）模型的一般化，常常是检验经济冲击的保留工具。它们基于模型中所有变量滞后值之间动态行为的概念。一个有限p阶的VAR模型通常写作VAR（p），并可以表示成：
这里面yt和a0分别是表示内生变量和常数项的M1向量，而Ap是表示系数的MM矩阵，εt是表示外生冲击的M1向量。需要顾及的系数一共有M+pM^2个，因此随着所包含变量和/或滞后的数量急剧上升。这叫“维数灾难”，对频率估计来说尤其如此。在贝叶斯方法中，设置先验信念（prior beliefs）被施加于模型参数中，绕过了这个诅咒，并允许更大的模型，且显著提高样本外预测精度估计，相关论述可见Bańbura和Giannone Reichlin（2010），Koop（2013）。恰当地告知先前的信念是至关重要的，而且很自然地会接受批评。Giannone等(2015)提出以基于数据的方式设置先验参数，即，方法是将它们视为附加参数。他们通过对该模式的边际似然值(Marginal Likelihood，ML)进行积分来实现，并用它作为探索参数空间的判断准则。在他们的论文中，作者从理论上证明了这种方法，并预测了脉冲响应函数估计的准确性。通过实证例子，Giannone等人(2015)也证明了预测准确性与标准模型、绩效以及因素模型。他们的方法被广泛使用(见Alquist, Kilian，和Vigfusson 2013;Miranda-Agrippino和Rey 2015;Altavilla, Giannone和Lenza 2014)。他们考虑了常用的高斯-反-维斯哈特族（Gaussian-inverse-Wishartf）的先验分布:
这里面b，Ω，ψ和d是超参数γ的低维向量的函数。由于（2）（3）式的共轭性（conjugacy），模型的ML可以作为一个函数，在封闭形式（closed form）下有效计算(见Giannone等 2015, p. 439)。在他们的论文中，作者考虑了三个具体先验——明尼苏达先验(Litterman)先验，它被用作基线；系数和先验和单一单位根先验。 明尼苏达先验(Litterman 1980)本质上强加了个体变量遵循随机漫步过程的假设。这个简约的标准通常在宏观经济时间序列方面表现良好(Kilian和Lutkepohl（2017）,p.356)，并且经常被用作评估准确性的基准。它由下面的方法所刻画:
其关键参数为λ，控制着先验的松紧度，即，衡量先验和数据的相对重要性。当λ趋近0时，先验是准确施加的，而当λ趋向无穷时，后验估计将趋近于普通最小二乘（OLS）估计。控制滞后阶数增加的方差衰减，可控制远距离观测的惩罚。最后，ψ控制着变量上先验滞后项的标准差。 明尼苏达先验的改进通常作为附加先验来实现，试图“降低VAR模型对初始观测的估计条件所隐含的确定性成分的重要性”(Giannone 等，2015, p.440)。首先是系数总和(sum-of-coefficients,SOC) 先验(Doan, Litterman，和Sims 1984)，强加了一个概念，即无变化预测在时间序列的开始是最优的。它是通过在数据矩阵上添加人工模型观测值的泰尔混合估计（Theil mixed estimation）实现的，其构造方法如下:
这里面是一个关于每个变量的前p个观察值的平均值的M*1向量。核心变量μ控制着方差，因此控制了先验的松紧度，即，μ趋向无穷时，先验变得无信息性。μ趋向0时，模型以同样多的单位根作为变量，并且没有协整。这促进了单一单位根(SUR)先验(Sims 1993;Sims和Zha1998)，它在数据中建立了协整关系。先验将变量推向无条件平均值或至少存在一个单位根。与之对应的虚拟观测值为:
设置这些先验的参数已被广泛讨论，并已有一些启发式方法(例如，Litterman 1980;Doan等，1984;Bańbura等，2010)。Giannone等（2015）写道，从贝叶斯视角来看，选择参数的过程在概念上是和模型中其他任何参数的推断是一样的。他们表明，由于数据的边际似然性，给定先验参数在具有共轭先验的封闭形式的VAR模型中可用，有可能把模型看成分层的(Giannone等 2015, p. 437)。 下面我们要来正式介绍BVAR的用法了！我们将用R包BVAR来实现。 BVAR包将Giannone等（2015）的分层方法在R中实现了（R Core Team，2019），并给了使用者一种容易上手且灵活的使用分层贝叶斯VAR的方法。它可以主要用于现代、(宏观)经济、多元时间序列分析。它对先验选择的分层方法可以防止有问题的参数选择。由于它的易操作性和灵活性，在Bańbura等 (2010) 或Koop (2013)看来，它很适合用来做经济分析，并可能用于像是模型的一致性检验。它也可以作为贝叶斯范式在多变量经济时间序列建模中的引入。 这个包除了基本款的和mvtnorm (Genz，Bretz，Miwa,Mi,Leisch.Scheipl和Hothorn 2019)之外没有其他依赖关系，因此可以跨平台使用，甚至可以在最小的安装上使用。它是在原装 R中实现的，以提高透明度，并降低贡献和/或适应的门槛。无论如何，一个函数方法的包的结构方便了未来的计算密集的步骤端口，例如，C。完整的文档、访问多种设置的helper函数和用于分析的标准方法使软件包易于操作，而不牺牲灵活性。 BVAR提供了广泛的定制选项，涉及到要使用的先验、它们的参数和它们的层次处理。对作为基线的明尼苏达先验来说，它的所有参数都是可调的，可以分层处理。包括SOC和/或SUR 先验的选项也很容易获得。此外，灵活的实现允许用户指定自己的自定义先验，只要它们是通过混合估计实现的，就像之前的虚拟观测对象的先验一样。超参数适当的起始值是用optim()命令来获取的（R Core Team 2019），用的是有限内存BFGS的准牛顿方法。MCMC，特别是Metropolis-Hastings (MH)算法，提供了进一步的选择。自然，删去和保存的数量是可调的，可以进行削减。为了正确地探索后验分布，先验参数的提议范围（proposal range）是至关重要的，因此可以单独设置。借助变量目标（variable target）和调整率（adjustment rates）,适当的接受率（acceptance rate）可以通过在老化阶段（burn-in phase）启用自动提议率（automatic proposal rate）调整来实现。 贝叶斯VAR模型的主要应用是利用脉冲响应函数(IRF)对(宏观)经济系统进行结构分析。这些函数作为冲击经济系统的一个表示，并用于分析模型变量的反应。这些冲击的准确传播是非常有趣的，但适当的识别是必要的，以便以一种有意义的方式解释它们。BVAR目前具有两种最常见的识别方案——即短期零限制（short-term zero restrictions）和符号限制（sign restrictions）。前者也称为递归标识，是通过Cholesky方差-协方差矩阵∑的分解实现的（见Kilian和Lütkepohl 2017，第8章）。这种方法的计算成本低，实现了精确的识别，而不需要太多的基于理论的关于可变行为的假设。只有变量的顺序是关键的，因为某些变量的同时期反应是有限的。符号限制(参见Kilian和Lutkepohl 2017，第13章)是另一种流行的识别方法，可以打包使用，遵循Rubio-Ramirez、Waggoner和Zha(2010)的方法。这个识别方案使得一些关于变量在一定冲击后行为的假定成为必要。随着模型维数的增加，这可能很快变得具有挑战性。此外，通过符号限制进行识别的代价是增加参数的不确定性，并导致IRF的精度损失。利用VAR模型可以解决的另一个问题是，冲击后哪些变量驱动某个变量的路径。为了帮助分析这一点，可以用预测误差方差分解(forecast error variance decompositions，FEVD)来实现。它们允许对决定经济系统行为的过程进行更详细的结构分析。贝叶斯VAR在预测练习中也表现得很好。它们已被证明优于许多其他方法(例如，Carriero、Kapetanios和Marcellino 2009;Koop 2013)，不需要像结构模型（structural models）等模型一样推断特殊的限制参数。此时，BVAR可用于进行无条件预测，其预测结果可与因子模型的预测结果相媲美(Giannone等，2015)。条件预测或情景分析，即，其中假设一个或多个内生变量的未来路径是已知的，将在未来实现。 使用R包估计贝叶斯VAR模型可以通过各种辅助函数和参数轻松定制。选择默认值是为了提供一个合理的起点，并允许逐步采用包。分析对R用户来说很好用——除了通用函数的实现，还可以选择绘制轨迹、密度、剩余、预测和脉冲响应，包括summary()、predict()、irf()命令等。最终和中间输出将以通常的格式和特征，用print()的方式提供，以提供透明的研究过程。因此，现有的框架可以用于进一步的分析——例如，coda (Plummer et al. 2006)用于检查收敛性，ggplot2 (Hadley 2016)用于绘图。 除了上面介绍的特性外，BVAR还包括FRED-QD数据集(McCracken和Ng 2016)，该数据集在修改后的开放数据共享属性许可(ODC-BY 1.0)下获得许可。它构成了描述战后时期美国经济的宏观经济变量的最大数据库之一，非常适合进行多元时间序列分析。该数据以季度为基础提供了248个宏观经济指标，最早可追溯到1959年第一季度，其中234个指标包含在一揽子数据中。它会定期更新，目前包含在套装中的版本将持续到2018年第四季度。FREDQD适用于广泛的经济现象研究，并经常用于开发新模型(如Huber，Koop和Onorante 2019)的基准测试。 说了这么多，你可能还是不知道BVAR怎么用，那么请看下面的例子吧。在这个例子中，我们将用上面提到过的FRED-QD数据集的一部分。当然，之后你也可以用它的其他部分。工作流程是这样的：(1)数据准备，(2)模型先验配置等，(3)模型估计，最后(4)对输出和绘制结果进行评估。 最主要的函数bvar()要求数据可以被表示成一个没有任何缺失的数值矩阵。我们从FRED-QD中提取了6个时间序列：以2012年为基准，十亿美元为单位的真实GDP，工业生产指数，以千人为单位的非农业就业人数，制造业平均每周工作时间，城市消费者CPI以及以百分数为单位的有效联邦基金利率。由于我们希望得到稳定的AR过程，因此我们对除联邦基金利率以外的所有变量进行了转换，以实现平稳性。对于GDP和CPI，我们使用年对数差分法，对于工业生产，使用非法对数差分法，而对于平均每周小时数，使用季度对数差分法，代码如下。图1提供了转换后的时间序列的概述。R> df <- fred_qd[, c("GDPC1", "INDPRO", "PAYEMS", + "CES0600000007", "CPIAUCSL", "FEDFUNDS")]R> for (i in c("GDPC1", "CPIAUCSL")) + df[5:nrow(df), i] <- diff(log(df[, i]), lag = 4) * 100 R> for(i in c("INDPRO", "PAYEMS", "CES0600000007")) + df[2:nrow(df), i] <- diff(log(df[, i]), lag = 1) * 100R> df <- df[5:nrow(df), ]
自然，这个选择过程也可以扩展到FRED-QD以及其他数据集的其他变量中。也可以参见Kilian和Lutkepohl(2017)的第2章和第19章，了解更多关于变量转换的信息 在准备了数据之后，我们需要指定先验，并调整模型的其他配置。注意，与预先设置相关的函数的前缀是bv_。因此，用户可以快速、简单地访问选项及其文档，从而方便了它们的使用，并降低了新用户的进入门槛。分析的方法和功能，与习惯用语密切相关。 如前所述，明尼苏达先验被用作基线，并且它的关键参数的拟和被包含在每个默认的层次建模练习中。我们使用bv_minnesota()来指定参数λ和α更精确的先验分布，以及提议值的上界和下界。按照Giannone等（2015），分布时伽马密度的，并且边界被用于抛弃难以置信的和不合理的数值，它们是从用于MH步骤的高斯提议分布（Gaussian proposal distribution）中得到的。参数var指定了模型中常数项的先验方差，通常被设置为相当分散的，即以较大的值。代码如下：+ lambda = bv_lambda(mode = 0.2, sd = 0.4, min = 0.0001, max = 5),+ alpha = bv_alpha(mode = 2, sd = 0.25, min = 1, max = 3),有了SOC和SUR先验，我们也把两个提前构建好的虚拟观测先验（pre-constructed dummy-observation）包括进来。关键参数的先验也被假定为是伽马分布的，并且这个过程跟λ和α一样。定制的虚拟观测先验是通过bv_dummy()类似地创建的，只需要一个额外的函数来构建观察。R> soc <- bv_soc(mode = 1, sd = 1, min = 1e-04, max = 50)R> sur <- bv_sur(mode = 1, sd = 1, min = 1e-04, max = 50)一旦指定了先验，我们就将它们提供给bv_priors()，对任何伪观察先验使用省略号参数(…)。通过参数hyper，我们选择哪一个先验参数应该被分层处理。它的缺省设置“auto”包括了λ，以及所有可选的虚拟观测先验的关键参数，这和提供特征向量c（“lambda”，“soc”，“sur”）是等价的。不被包含在分层步骤中的先验参数被当成固定的，并设置成和它们的mode参数相等。R> priors <- bv_priors(hyper = "auto", mn = mn, soc = soc, sur = sur)我们可以通过函数bv_irf()和bv_fcast()来调整IRF的计算或预测。考虑的范围可以针对这两种情况进行调整。FEVD是在IRF的一个额外步骤中计算的，可以通过fevd进行切换。 冲击的识别是通过identification命令实现的。若要跳过任一项的计算并加速估计，可将参数设为NULL;两者都可以事后计算。R> irfs <- bv_irf(horizon = 12, fevd = TRUE, identification = TRUE) 最终，我们调整MH步骤，以达到合适的接受率，并因此准确地找出模型参数的后验分布。这是通过bv_metropolis()实现的，它需要一个初始参数scale_hess，这是一种数值标量或向量，用于缩放逆Hessian，它是用于绘制分层处理参数的提议（proposal）。这也可以用设置adjust_acc=TRUE来补充，以在老化（burn-in）期间启动自动缩放调整。自动调整由acc_change百分比迭代完成，直到达到acc_lower和acc_upper之间的通过率。代码如下：R> mh <- bv_metropolis(scale_hess = 0.005, adjust_acc = TRUE, + acc_lower = 0.25,+ acc_upper = 0.35, acc_change = 0.02)这许多可供选择的设置让用户能够定制它们的模型和所有他们个人需要的组成部分。这对于解决一系列广泛的不同研究问题是必要的。但是，也可以使用更简单和更快的方法——默认设置应该能够满足广泛的应用程序。这使得用户能够(1)专注于他们模型的关键部分，(2)轻松地使用BVAR并逐步调整他们的模型。 BVAR的核心是bvar()函数。在准备数据和可选调整的各种设置后，可进行估计。除了设置对象，我们提供了延迟的数量，包括在我们的模型和一些关于MCMC迭代的选项。在n_save中，我们定义了迭代的总次数，在n_burn中我们设置了要丢弃的初始迭代的数量，通过n_thin我们给出要存储的抽取分数的分母。更进一步，verbose=TRUE, 提示打印中间结果，并在MCMC步骤中启用进度条。代码和结果如下：R> run <- bvar(df, lags = 5, n_draw = 25000, n_burn = 10000, n_thin = 1, + priors = priors, mh = mh, fcast = fcasts, irf = irfs, verbose = TRUE)Optimisation concluded. Posterior marginal likelihood: -1123.907 Parameters: lambda = 0.52; soc = 0.9; sur = 0.69 |====================================================| 100%Finished after 2.58 mins. 函数的返回值是一个类的对象bvar——一个有数个输出的命名了的清单。这总是包括了最开始感兴趣的参数，即，VAR的系数的后验提取，方差-协方差矩阵的后验提取，按层次结构处理的超参数的后验提取。其他内容包括每次提取的边际似然值，优化过的参数的起始值（通过optim得到），预先设置和自动设置，以及最开始要的bvar函数。各种各样的元信息也包括在内，例如，可接受的绘制数量，变量名和计算花费的时间。只有在实际计算的情况下才添加一些输出，即IRF、FEVD和预测结果。BVAR提供了用于bvar类型及衍生物的print(),plot()和summary()工具。print()工具提供了一些元信息，通过optim()对先验参数和优化结果进行了详细的分层处理。summary()工具在vars（Ptaff 2008）中提供了模仿了它的对手，并包括了关于估计模型的对数似然性、VAR模型的系数及其残差的方差-协方差矩阵的信息。这些也可以用logLik(),coef()和vcov()来实现。更多的标准工具例如fitted(),density()和residuals()也是可以用的。IRF、FEVD和预测结果的评估是相似的，后面将讨论。输入这行代码，就会出现下面这些结果：
收敛性对于MCMC算法的稳定性至关重要，因此需要特别关注。除了适当的MH步骤的接受率，跟踪和密度图的参数通常被用来评估收敛。plot()工具通过缺省提供了ML的画图和分层处理参数。参数边界以虚线灰线绘制。画图功能也提供了type参数，有轨迹，密度或者都画两个选项，还有参数的子集绘制超参数(vars)或绘图系数值作为代替(vars_response和vars_impulse)。可以看下面的代码，图4给出了相应的图象。图3和图4显示了关键超参数的收敛性。由图可见，没有明显的异常值是可识别的，根据明尼苏达先验，有轻微的收缩，这是因为绝大多数λ的概率在0.45和0.60之间。但是，可能需要使用额外的收敛诊断。这是通过一个coda (Plummer等，2006)通用as.mcmc()函数的方法实现的。代码如下：R> plot(run, type = "full", vars = "lambda", mfrow = c(2, 1))
Figure 3所有分层处理参数和ML的轨迹图和密度图
Figure 4明尼苏达先验法的关键参数——拟合曲线的轨迹图和密度图通过解释脉冲响应的结构分析是使用BVAR的一种便利且直接的方式。通用函数irf()可用于从bvar对象检索或计算irf。产生的东西具有绘图、打印和汇总自身的方法。plot()方法可以通过名称或位置将图形子集为特定的脉冲和/或响应。下面可以看到一个示例，代码如下，相应的输出如图5所示。R> plot(irf(run), vars_impulse = c("GDPC1", "FEDFUNDS"), + vars_response = c(1:5))
Figure 5 GDP增长、工业生产、非农就业、平均每周工作时间的脉冲响应，通货膨胀、综合经济冲击(左图)和货币政策冲击(右图)。灰色线表示第16和第84可信集。预测误差方差分解是结构分析的另一个重要工具。FEVD绘图在bvar对象中作为数组提供。汇总它们的一种方法是计算所有保存的绘图和时间段的中值，如下面的通用函数fevd()所示。
为了帮助进行预测工作，BVAR提供了一种predict()工具以及相关的绘图、打印和汇总方法。通过在predict()里面用bvar()或者先算完bvar()再用predict()就行。两种情况下，设置都是通过bv_fcast()来实现。在下面的例子中，我们把预测添加到我们之前的bvar对象中，把predict()的输出分配给run[["fcast"]]。然后，我们使用plot()方法，通过向vars提供它们的名称(或者它们的位置)来可视化两个包含的变量的预测。代码如下，相关的输出见图6。R> run[["fcast"]] <- predict(run, horizon = 8)R> plot(predict(run), vars = c("GDPC1", "FEDFUNDS"),+ orientation = "vertical")predict()和irf()都会提取已有的预测或IRF，或者基于bvar中已有的预测结果来计算。这样，它们通过使用不同的设置，而无需重新估计整个模型来计算结果。此外，参数conf_bands允许调整可信区间周围的中值预测/ IRF，这在其他方法中有使用到。下面的例子给出了增加水平和调整可信间隔后的IRF的事后计算。然后绘制结果，如图7所示。代码如下：R> plot(irf(run, conf_bands = 0.05, horizon = 20, fevd = TRUE),+ vars_impulse = c("GDPC1", "FEDFUNDS"), vars_response = c(1:5))
Figure 6对GDP增长和联邦基金利率的无条件预测。灰色线表示第16和第84可信集
Figure 7前后计算的GDP增长、工业生产、非农就业、平均每周工作时间和CPI通胀对20个季度的总需求冲击(左面板)和货币政策冲击(右面板)的脉冲响应。灰色线表示第5和第95可信集 本文介绍了R的贝叶斯向量自回归的基本估计方法。它提供了一种灵活、结构化和透明的方法来评估多元时间序列分析领域的广泛研究问题。通过对先验的规范，减少对主观选择的需要，它抵消了对贝叶斯方法的一些主要批评。通过提供几种方法和功能，该包允许快速评估各种模型输出。通过一个应用示例，我们说明了这个包的用法，并解释了它的实现和配置。 BVAR比大多数通用的贝叶斯统计软件有更低的进入成本，并包括了VAR模型后续的一系列分析。它对一系列问题是有用的，可以应用于更广泛的研究问题，而不像类似的包集中在贝叶斯VAR模型。R中的实现使包易于使用和扩展。下面这些短链接文章属于合集，可以收藏起来阅读，不然以后都找不到了。
计量经济圈组织了一个计量社群，有如下特征：热情互助最多、前沿趋势最多、社科资料最多、社科数据最多、科研牛人最多、海外名校最多。因此，建议积极进取和有强烈研习激情的中青年学者到社群交流探讨，始终坚信优秀是通过感染优秀而互相成就彼此的。