Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
37,660 characters · 5 sections · 28 citation commands
Approximation Properties of Variational Bayes for Vector Autoregressions
Variational Bayes (VB) is an approximate method of Bayesian inference originally developed in machine learning jordan1999introduction,attias2000variational, finding its way into statistics Ormerod2010,Blei2017 and econometrics koop2018variational,gunawan2017fast. The merit of VB is that it is fast and scalable and therefore useful in big data settings where MCMC might be difficult or infeasible (e.g. forecasting with large vector autoregressions). A drawback of the VB is that its approximation properties are still not well-understood. \\ There are a number of papers on properties of VB including but not limited to wang2005inadequacy, wang2006convergence, hall2011asymptotic, bickel2013asymptotic, you2014variational, pati2017statistical and wang2018frequentist. These papers focus mainly on asymptotic properties of VB. Few analytical results, if any, are available on VB's finite sample approximation error. Here we study the approximation error of VB for a linear Gaussian multiple equation regression. We focus on vector autoregressions (VAR) since large VARs are of enormous interest in macro-econometric forecasting koop2013forecasting,Giannone2015 and VB is likely to be very useful in this context koop2018variational although our results are valid for any linear Gaussian multi-equation regression.\\ We derive the divergence between the mean and the mode of VB posterior and those of exact posterior, the magnitude of underestimation in posterior variance and the magnitude of the KL divergence. We show that the approximation error decreases with an increase in sample size and strength of a prior, but increases with the size of the model. We also derive the divergence in the predictive densities and demonstrate that VB estimates predictive densities well. This is important given the enormous interest in forecasting with large VARs for which MCMC is difficult or time consuming.\\ The paper is organized as follows: Section (ref) provides an overview of VB while Section (ref) derives the mean-field VB posterior for both conjugate and independent prior VARs. In Section (ref), we provide the main results on approximation errors and Section (ref) concludes the paper.
Let $\mathbf{y}$ be a vector of observations, $p(\mathbf{y}|\mathbf{\theta })$ the likelihood for a model and $p(\mathbf{\theta})$ a prior where the objective is to learn about the parameter vector $\mathbf{\theta }$ given the observations and the prior. Bayes theorem can be used to derive the posterior kernel easily but learning about $\mathbf{\theta}$ (its posterior density and moments) requires integration of the posterior kernel. For complex models, integration is often done through MCMC but even MCMC can sometimes become intractable (e.g. in big data settings). In VB the idea is to approximate the posterior $p(\mathbf{\theta }|\mathbf{y})$ with a density $q(\mathbf{\theta })$ that is more tractable. The optimal $q(\mathbf{\theta })$ is obtained by minimizing the KL divergence \footnote{For a more detailed review of VB see Ormerod2010 or Blei2017.} \[KL(q||p)=\int_{\mathbf{\theta }\in \,\Theta }{q(\mathbf{\theta })\ln [{q(\mathbf{\theta })}/{p(\mathbf{\theta }|\mathbf{y})}\;]d\mathbf{\theta }}\] over all $q(\mathbf{\theta })$s that belong to a family of densities $\Lambda .$ A commonly used family includes all $q(\mathbf{\theta })$s that can be written as $q(\mathbf{\theta})=\prod\limits_{k=1}^{M}{{{q}_{k}}(}{{\mathbf{\theta }}_{k}})$ where $\{{{\mathbf{\theta }}_{1}},....,{{\mathbf{\theta }}_{M}}\}$ is some partition of $\mathbf{\theta }$ and the ${{q}_{k}}$s are probability density functions. This factorized form corresponds to an approximation framework known as the mean field; thus, this version of VB is often referred to as mean-field variational Bayes (MFVB). \\ Using calculus of variations, it has been shown that the optimal ${{q}_{k}}$s are obtainable through the following iterative procedure
where $p(\mathbf{y},\mathbf{\theta })=p(\mathbf{y}|\mathbf{\theta })p(\mathbf{\theta })$, ${{E}_{{{\theta }_{-k}}}}=\int{\log p(\mathbf{y},\mathbf{\theta }){{\prod\nolimits_{j\ne k}{{{q}_{j}}({{\mathbf{\theta }}_{j}})}}^{old}}}$ and ${{q}_{1}}{{({{\mathbf{\theta }}_{1}})}^{old}},....,{{q}_{M}}{{({{\mathbf{\theta }}_{M}})}^{old}}$ are current estimates for $q$. The iteration is stopped when the increase in $\ln \underline{ML}$ (defined below) becomes negligible. In practice the partition of $\mathbf{\theta }$ is often chosen in a way such that $q_{j}^{{}}({{\mathbf{\theta }}_{j}})$s are known density functions that depend on the moments of the remaining parameters. As a result, it is often sufficient to iterate over the appropriate moments extracted from the distributions in (2.1).\\ It can be shown that
where $ML$ denotes the marginal likelihood. If we define
then given that $KL(q||p)$ is positive, $\ln \underline{ML}$ provides a lower bound to the log-marginal-likelihood. Note that minimization of $KL(q||p)$ is equivalent to maximization of $\ln \underline{ML}$ therefore $\ln \underline{ML}$ (which can be computed using VB) is often used to monitor the convergence of (2.1).
Since the seminal paper of sims1980macroeconomics, vector autoregressions have become an essential tool in macroeconomic analysis and forecasting. They are extensions of univariate autoregressive models to a vector of variables \footnote{For a detailed review of Bayesian VARs see e.g. Karlsson2013 or koop2010bayesian.}. A VAR(d) can be written as
where $t=1,...,T$ indexes time, $d$ is the lag order, ${{\mathbf{y}}_{t}}$ is an $M\times 1$ vector of observations at time $t$, ${\mathbf{\gamma }}_{0}$ is an $M\times 1$ vector of constants, the ${{\mathbf{\Gamma }}_{i}}$s are $M\times M$ matrices of autoregressive coefficients and ${\mathbf{\varepsilon}}_{t}$ is an $M\times 1$ error vector with mean $\mathbf{0}$ and covariance matrix $\mathbf{\Sigma}$. Define $p=Md+1$ (the number of parameters in each equation), a matrix ${{\mathbf{Y}}_{T\times M}}=({{\mathbf{y}}_{1}},...,{{\mathbf{y}}_{T}}{)}'$, a $1\times p$ vector ${{\mathbf{x}}_{t}}=(1,{{\mathbf{{y}'}}_{T-1}},...,{{\mathbf{{y}'}}_{T-d}})$, ${{\mathbf{X}}_{T\times p}}=({{\mathbf{{x}'}}_{1}},...,{{\mathbf{{x}'}}_{T}}{)}'$, a matrix of coefficients ${{\mathbf{\Gamma }}_{p\,\,\times M}}=({{\mathbf{\gamma }}_{0}},{{\mathbf{\Gamma }}_{1}},...,{{\mathbf{\Gamma }}_{d}}{)}'$ and a matrix ${{\mathbf{E}}_{T\times M}}=({{\mathbf{\varepsilon }}_{1}},...,{{\mathbf{\varepsilon }}_{T}}{)}'$, then equation (3.1) can be written as:
where $\mathcal{MN}$ denotes a matricvariate normal distribution (see Appendix-(ref) for a definition). We first consider a natural-conjugate prior where the exact and the VB approximate posterior and the error in the MFVB approximation can be derived analytically.\\ Normal-Wishart Conjugate Prior: For a VAR model, the natural-conjugate prior is given in the following normal-Wishart form:
where $\mathbf{\underline{\Gamma }}$ is a $p\times M$ matrix of prior means, $\mathbf{\underline{V}}$ and $\mathbf{\underline{S}}$ are $p\times p$ and $M\times M$ positive-definite matrices, $\mathcal{W}$ and $\underline{\upsilon }$ denote the Wishart distribution and its degrees of freedom. In this case, it is well-known that the exact posterior is a matricvariate normal-Wishart distribution Karlsson2013:
Other quantities of interest are marginal posteriors, posterior moments, marginal likelihoods and predictive densities that are available in analytical forms and are given in Table-(ref).\\ Now let’s consider MFVB inference where the VB approximate posterior density is factorized into $q(\mathbf{\Gamma },{{\mathbf{\Sigma }}^{-1}})=q(\mathbf{\Gamma })q({{\mathbf{\Sigma }}^{-1}})$. Using (2.1) it can be shown (see the Appendix-(ref)) that the VB posterior is
without a need for an iterative process !. Related quantities of interest such as marginal posteriors, posterior moments, marginal likelihoods and predictive densities are given in Table-(ref).\\ Normal-Wishart Independent Prior: A conjugate prior has the advantage that one can derive the approximation errors analytically. The potential drawback is that the results might be crucially dependent on the conjugacy assumption. Therefore, we also consider a VAR with an independent prior. In this case we cannot derive the approximation error analytically, but we use a real data example to see the extent to which the results differ from the conjugate VAR’s results.\\ To discuss VAR with an independent prior we appeal to another way of writing a VAR model koop2010bayesian:\\
$\And \ {{\mathbf{Z}}_{t}}=\left(
\right)$ where if $\mathbf{z'}_{mt}=\mathbf{x}_t$ for $m=1,...,M$, we obtain the same equation as (3.1). In a more compact form:
\, where $\mathcal{N}$ denotes a multivariate normal distribution. For the independent case, we assume that the parameter vector $\mathbf{\beta }$ and the precision matrix $\mathbf{\Sigma }_{{}}^{-1}$ are a priori independent:
It is well-known that under this prior, the posterior can be obtained through a Gibbs sampler Karlsson2013 defined by
Marginal likelihoods and predictive densities are not available in analytical form in this case, but can be obtained through simulation Karlsson2013.\\ \\ Now let's assume $q(\mathbf{\beta },{{\mathbf{\Sigma }}^{-1}})=q(\mathbf{\beta})q({{\mathbf{\Sigma }}^{-1}})$, then MFVB posterior can be obtained using (2.1) as
where $\{{{\mathbf{\bar{\beta}}}_{q}},{{\mathbf{\bar{V}}}_{q}},\mathbf{\bar{S}}_{q}^{-1}\}$ are obtained by using a starting value for $\overline{\mathbf{\Sigma }_{q}^{-1}}$ and iterating the following until convergence
$\ln \underline{ML}={{E}_{q}}\ln p(\mathbf{y,\theta })-{{E}_{q}}\ln q(\mathbf{\theta })$, which is used for monitoring of the convergence, can be obtained using (2.3), (A.1), and the entropy of normal and Wishart distributions, as
where ${\Gamma }_{M}$ is the multivariate gamma function, $\left| \ .\ \right|$ is the determinant and $tr\{.\}$ denotes the trace of a matrix.\\ While calculating the exact predictive density requires simulation, using arguments similar to those used in Appendix-(ref), the MFVB predictive density can be obtained as sum of a normal and a t-distribution with the following mean and variance
Table-(ref) compares various features of MFVB and of the exact posterior for the conjugate prior case. As this table shows, the posterior means for both $\mathbf{\Gamma }$ and ${{\mathbf{\Sigma }}^{-1}}$ are equal but MFVB underestimates variances of the parameters, confirming the folklore. The variance of ${{\mathbf{\Gamma }}_{q}}$ is smaller than the variance of $\mathbf{\Gamma }$ by a factor of $1-\frac{M+1}{T+\underline{\upsilon }}$ which is close to one unless the number of equations is large (see Table-(ref)). Variances of elements of $\mathbf{\Sigma }_{q}^{-1}$ are smaller by a factor of $1-({p+1)}/{(T+\underline{\upsilon })}\;$ indicating small underestimation when $p$ (the number of parameters in each equation) is small compared to $T+\underline{\upsilon }$. However, it could be substantial for large VARs (see Table-(ref)). The mode of the approximate posterior is equal to the mode of the exact posterior for $\mathbf{\Gamma }$ but for $\mathbf{\Sigma }_{q}^{-1}$ it is smaller by a factor of $1-{p}/{(T+\underline{\upsilon }+p)}\;$. The KL divergence which is a measure of closeness of the approximate posterior to the exact posterior can be derived as \[KL=\ln ML-\ln \underline{ML}=-\dfrac{Mp}{2}(\ln 2+1)+\dfrac{M}{2}\ln \left( \dfrac{{{{\bar{\upsilon }}}_{q}}^{{{{\bar{\upsilon }}}_{q}}}}{{{{\bar{\upsilon }}}_{{}}}^{{{{\bar{\upsilon }}}_{{}}}}} \right)-\ln \left( \dfrac{{{\Gamma }_{M}}({{{{\bar{\upsilon }}}_{q}}}/{2}\;)}{{{\Gamma }_{M}}\left( {{{{\bar{\upsilon }}}_{{}}}}/{2}\; \right)} \right)\] and using Stirling’s approximation to gamma functions as \[KL\sim \dfrac{1}{2}\sum\limits_{j=1}^{M}{\ln \left( \dfrac{{{\left( \bar{\upsilon }-j+1 \right)}^{(\bar{\upsilon }-j)}}{{{\bar{\upsilon }}}_{q}}^{{{{\bar{\upsilon }}}_{q}}}}{{{\left( {{{\bar{\upsilon }}}_{q}}-j+1 \right)}^{({{{\bar{\upsilon }}}_{q}}-j)}}{{{\bar{\upsilon }}}_{{}}}^{{\bar{\upsilon }}}} \right)}\] Note that KL increases as the size of the model (i.e. the numbers of parameters in each equation) increases and decreases with an increase in the number of observations $T$, and with the strength of the prior for ${{\mathbf{\Sigma }}^{-1}}$ (i.e. $\underline{\upsilon }$). The KL does not depend on the data and the other hyper-parameters !\\ Note that
where $i,j=1,..,M$ and $h=1,...,Mp$. The first line shows that with an increase in the number of observations, variances of parameters based on MFVB approaches variances of parameters from the exact Bayes. The second line shows the same results with an increase in strength of the prior for $\Sigma^{-1}$. \\Also note that the predictive densities based on the exact and MFVB posteriors share the same mean and an almost equal variance. The size of the model does not seem to affect predictive densities (see Tables (ref) and (ref)).\\ To study the approximation error for the VAR with an independent prior, we estimate a small (M=3) and a medium VAR (with M=7) using the same data set as in Giannone2015 assuming four lags (d = 4). The data set is based on the benchmark US data set of Stock2008 and consists of 200 observations on US macro quarterly variables from the first quarter of 1959 to the last quarter of 2008. The seven variables are real GDP, federal funds rate, GDP deflator, real investment, real consumption, hours worked and real compensation per hour. Hyper-parameters (i.e. $\mathbf{\underline{\Gamma }},\mathbf{\underline{S}},\mathbf{\underline{V}},\underline{\upsilon },\mathbf{\underline{\beta }}$) have been chosen based on the so called Minnesota prior as described in e.g. Karlsson2013. \\
Table-(ref) reports $\ln ML$, $\ln \underline{ML}$ and KL divergence for both the conjugate and independent prior models. As it can be seen, the KL divergences for conjugate and independent VARs, although not equal, are of the same order of magnitude and they both increase with the size of the model. Our other results (not reported here) show that KL divergence decreases with increasing $\underline{\upsilon }$ and $T$. \\ Since the number of variables in VAR models is large, we cannot present comparisons for each of the parameters here, but Table-(ref) gives means and variances of the diagonal elements of the precision matrix for the small VAR. The results suggest that MFVB's mean of ${{\mathbf{\Sigma }}^{-1}}$ is very well estimated but its variance is underestimated and the amount of underestimation is of a similar order to that from the the conjugate case (see Table-(ref) for the medium VAR’s estimation results).
Table-(ref) gives the ratios of means and variances from VB to those from the exact (or MCMC) predictive densities for the small VAR [see Table-(ref) for medium VAR]. For the conjugate case, the means are equal but the variances are smaller for VB, although the difference is negligible. A similar statement can be made for the independent prior case although underestimation seems slightly higher (perhaps more a result of simulation errors).
Tables (ref) to (ref) suggest VB's approximation errors derived for the conjugate prior case also holds for the independent prior case although not in an exact manner.
We derive analytical expressions for the MFVB's approximation error for a VAR with a conjugate prior. We also show that errors of a similar order hold for a VAR with independent priors based on a real data example. Our results confirm the claim that MFVB underestimates variances of the parameters and shows that the magnitude of underestimation (and also the KL divergence) change directly with the size of the model but inversely with the number of observations and the strength of the prior for variance. MFVB, at least for standard VARs, is as good as MCMC for estimating posterior means, and for forecasting koop2018variational. MFVB also provides an ideal candidate for reciprocal importance sampling hajargasht2018accurate. Estimation of variances can also be improved giordano2017covariances.