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.
50,884 characters · 10 sections · 86 citation commands
Bayesian nonparametric graphical models for time-varying parameters VAR
\singlespacing
\singlespacing
Over the last decade, the availability of large datasets in economics and finance has allowed the introduction of high dimensional models. In particular, large datasets in macroeconomics help to improve the forecasts, while in finance some authors have investigated the use of large datasets to analyse financial crises, contagion effects and their impact on the real economy. In order to deal with high dimensional models, the introduction of Bayesian nonparametric techniques have become popular in different fields (such as statistics and machine learning), but only few attempts have been made in econometrics. In particular, Bayesian nonparametric approach allows to improve the estimation efficiency and the prediction accuracy in time series analysis.
Recently, time-varying parameter (TVP) models provide an interesting alternative to process multiple change points; for example, time-varying structural vector autoregressive (VAR) models have been used in Primiceri05Bayes_TVP_VAR for study monetary policy application; Dangl12predictive_TVP forecast equity returns by mean of TVP models; and in Belmonte14shrink_TVP the European inflation has been studied via a time-varying parameters model. As shown in Primiceri05Bayes_TVP_VAR, DelNegroPrimiceri15Bayes_TVP_VAR and Bitto19Shrinkage_Bayes_TVP_VAR, the advantage in capturing gradual changes is due to the flexibility of TVP models. We combine the ideas behind time-varying parameters models and Bayesian nonparametric techniques, thus allowing to model complex phenomena in a flexible and efficient manner. Moreover, we provide an innovative Bayesian nonparametric time-varying graphical framework for making inference in high-dimensional time series.
In this paper, we allow coefficients to be sparse, meaning that only a fraction of the time varying parameters have significant effects, but we retain flexibility in modelling non-zero coefficients, by including temporal dependence in the prior structure. In order to achieve these goals, we define a shrinkage prior on the VAR coefficients by means of a Bayesian nonparametric prior (BNP). This distribution is a spike-and-slab prior, where on the spike (parametric) distribution, we impose two different specifications: a Dirac spike and a “diffuse” spike. On the other hand, on the (non-parametric) slab distribution, we use a well know Bayesian nonparametric Lasso prior as in Billio19BNP_sparseVAR.
The prior previously described groups the time-varying parameter vector autoregressive (TVP-VAR) coefficients into clusters and shrinks the coefficients within a cluster toward common notation. Differently from Markov-switching approach Krolzig97MS and random walk processes Primiceri05Bayes_TVP_VAR,DelNegroPrimiceri15Bayes_TVP_VAR, we impose time variation on the distribution of the VAR coefficients. In the literature of time-varying coefficients, the VAR coefficients are represented as a direct dependence, in practice they can be represented as state-space models, where they are functions of the previous time. On the other hand, we introduce a different structure, the indirect dependence on the VAR coefficients. In this case, we have a dependence construct through the atoms of the Dirichlet process and not on a direct way. Thus, we include a Bayesian nonparametric dependent prior specification on the VAR coefficients and the covariance matrix by means of a time-series dependent Dirichlet process (tsDDP) as in NietoBarajas12TimeSeries_DDP.
Following Billio19BNP_sparseVAR, our hierarchical prior overcomes overparametrization and overfitting issues by clustering the VAR coefficients into groups and by shrinking the coefficients of each group toward a common location. This hierarchical prior allows to contemporaneously estimate the (potentially) sparse time-varying causal network structure and to cluster the corresponding coefficients. In our BNP-TVP-VAR model, time-varying coefficients allow to (i) estimate the temporal networks of contemporaneous and causal structures, (ii) identify different sources of time variation, from the size of shocks and/or the propagation mechanism, and (iii) accommodate for potential non-linearities.
We also contribute to the literature on financial and macroeconomic contagion (see Billioetal12GrangerNet; Bianchi19GraphicalSUR and Barigozzi19NETS_network_estimation) through the lens of Granger causality and graphs/networks representation. Our BNP prior is particularly suited for studying Granger causality from time series and in particular it allows to estimate the most significant time-varying dependence interrelationships between the variables of interest. As explained above, we can extract time-varying graphs by using the posterior random partition induced by the non-parametric (slab) distribution, which allows to cluster the edges into groups.
Since their introduction in macroeconomics (see Sims80VAR_macroeconomics), vector autoregessive (VAR) models have been extensively used in econometrics and time series statistics. Large VAR models have been used to analyse and forecast high-dimensional macroeconomic data (e.g., McCracken16FRED-MD_dataset) and financial panels (e.g., Barigozzi19NETS_network_estimation). Moreover, in recent years VAR models have been used for studying financial and macroeconomic contagion (e.g., Cogley05Bayes_TVP_VAR_macro, Stock07US_inflation_forecast, Diebold12GVAR_VolatilitySpillover and Bianchi19GraphicalSUR). Although, VAR models have been extensively used for assessing the impact and spread of external shocks (i.e., to perform impulse-response analysis), forecasting, estimating networks from Granger-causal relationships and to study systemic risk and financial contagion (e.g., Diebold09VolatilitySpillover, Billioetal12GrangerNet and Barigozzi19NETS_network_estimation).
Despite being a potentially very flexible statistical tools, the high number of parameters and the typical limited length of standard macroeconomic datasets make unrestricted inference daunting as the cross-sectional size increases. This has favoured the use of penalised regression and Bayesian methods for dealing with the problem of over-parametrisation. The general idea is to use informative priors to shrink the unrestricted model towards a more parsimonious setting, thereby reducing parameter uncertainty and improving forecast accuracy (see Karlsson13Forecasting_BVAR_survey, Koop10BVAR_survey for a survey).
In the Bayesian VAR (a.k.a. BVAR) literature, a plethora of different prior distributions have been proposed to perform sparse estimation (e.g., see Giannone14priorVAR). Starting from the well-known Minnesota prior (see DoanLitterman84BVAR_MinnesotaPrior_survey, Litterman86BVAR_MinnesotaPrior), which specifies an objective prior on the coefficient and covariance matrices of a VAR, several parametric approaches have been developed exploiting hierarchical structures and finite mixtures (e.g., Kalli14Bayes_sparse_TVP_VAR, Gefang14Bayes_DoubleElasticNet, Huber19AdaptiveShrinkage_VAR, Kastner18Sparse_largeVAR).
Among the recent contributions for dealing with large dimensional models, we distinguish two approaches: the first attempts to reduce the size of the data to handle or to process during each step of the inferential algorithm, while the second is concerned with the reduction of the size of the parameter space. Within the first class, we mention the Bayesian compressed VAR of KoopKorobilis18BayesianCompressedVAR, who tackled the dimensionality issue by using random projections to compress the data, and the Bayesian composite likelihood approach of Chan18CompLike_BayesVAR. On the other hand, Gefang19VariationalBayes_largeVAR and Koop18VariationalBayes_VAR adopted a variational Bayes approach for performing efficient approximate posterior inference in large parameter spaces. Also, Kastner18Sparse_largeVAR exploited factor models and hierarchical shrinkage priors for providing a parsimonious parametrisation of the covariance matrix which allows for equation-by-equation estimation. Additional contributions for estimating large VAR and VARMA models include Koop13Large_TVP_VAR, Korobilis16PanelVAR_priorShrinkage and Chan16Large_BayesVARMA.
In addition to the large cross-sectional dimensionality, also the temporal length of many economic and financial datasets is steadily increasing. Thus, the possible relations between different variables of interest can be described by static matrix of coefficients. This assumption can be elapsed by introducing a time-variation of the matrix of coefficients of the time series. In particular, the most common approach consists in specifying a process governing the evolution of the parameters of interest. According to the force driving this dynamics, we distinguish observation-driven and parameter-driven time-varying parameter (TVP) models. The first class is mainly represented by generalised autoregressive score models (GAS, see Creal13GAS), while the second one includes Markov switching (e.g., Hamilton89MS, Krolzig97MS), change point (e.g. Pesaran06ChangePoint_VAR) and random walk models (e.g. DelNegroPrimiceri15Bayes_TVP_VAR, Primiceri05Bayes_TVP_VAR). These processes are able to describe parameters whose evolution is subject to switching regimes, structural breaks or smooth changes, respectively.
In the Bayesian and frequentist literature, the use of parametric models has been widely studied by applying different shrinkage methods (such as the Least Absolute Shrinkage and Selection Operator, known as LASSO). In particular, important papers focus on sparse and efficient estimation in high-dimensional datasets. However, more recently increasing attention is being devoted to the issue of over-shrinkage and to the modelling of non-zero coefficients (e.g., Giannone18IllusionSparsity). Consequently, there is an increasing need for adequate statistical tools capable of flexibly model the dynamics described by a VAR process, allowing for sparsity without incurring into over-shrinking.
In this paper, we aim to contribute to the growing literature on the use of Bayesian nonparametrics in time series analysis. In particular, Bayesian nonparametric techniques are widely in statistics, machine learning and data analysis as powerful tools for flexible modelling of complex data structure. Only recently, Bayesian nonparametrics has increased popularity in econometrics and in economic time series modelling to capture observation clustering effects (see e.g. Bassetti14BetaPitmanYor, Kalli18BNP_VAR and Billio19BNP_sparseVAR).
Up to our knowledge, our paper is the first to provide sparse Bayesian nonparametric VAR model when the coefficients are time-varying and the proposed two-stage prior specification can be easily extended to other classes, such as the seemingly unrelated regression (SUR) models. We propose a novel Bayesian nonparametric prior structure, which provides a sparse estimation of the coefficient matrix of a VAR model. This representation allows to manage the flexibility of non-zero entries and most importantly, to manage the time-variation in the matrix of coefficients through the atoms of the Dirichlet process and not through a state-space representation.
Our approach substantially differs from the existing literature in two aspects as described below. First, we consider a spike-and-slab prior distribution for each entry of the coefficient matrix, where on the spike we have a parametric prior specification by mean of Dirac or diffuse prior. On the hand, the slab component has random nonparametric prior. Second, we impose prior dependence on the coefficients by specifying a Markov process for their random distribution. As a by-product of the estimation procedure, we are able to extract a time series of dependent Granger-causality graphs. This shows how the BNP-TVP-VAR contributes to the literature on the estimation of time-varying networks from economic and/or financial series.
The paper is organized as follows. Section (ref) introduces the modelling framework and presents the BNP-TVP prior structure, then Section (ref) presents posterior approximation and describes how to extract Granger-causal time varying graphs from time series. Finally, Section (ref) draws the conclusions.
Let $n$ be the number of units in a dataset and $\mathbf{y}_t = (y_{1,t},\dots,y_{n,t})$ a vector of $n$ variables available at time $t$. A time-varying parameters vector autoregressive model of order $p$ (TVP-VAR($p$)) is defined as
where $B_{t}$ is the $(n\times n)$ matrix of time-varying coefficients an $t=p,\ldots,T$ is the time period. We assume that the error terms $\boldsymbol{\epsilon}_t = (\epsilon_{1,t},\ldots, \epsilon_{n,t})'$ are i.i.d. for $t$ with Gaussian distribution $ \mathcal{N}(\mathbf{0},\boldsymbol{\Sigma})$. Eq. (ref) can be written in the more compact form as
where we define $\boldsymbol{\beta}_t = \operatorname{vec}(\mathbf{B}_t)$; $\mathbf{X}_t = (\mathbf{y}_{t-1}' \otimes \mathbf{I}_n)$; $\otimes$ is the Kronecker product and $\operatorname{vec}(\cdot)$ the column-wise vectorization operator that stacks the columns of a matrix into a column vector.
Let us consider the problem of defining a flexible prior for a time-varying parameter model and we define the following TVP-VAR($p$) with $p$ equal to $1$ as
In Eq. (ref), the number of parameters of the $n$-dimensional TVP-VAR($1$) model is $(T-1)n^2 + n(n+1)/2 = O(n^2)$, thus scales quadratically in $n$. In macroeconomic and financial applications, the number of variables of interest ranges from $n=3$ (small size model) to $n=20$ (large model) and even $n=100$ (huge model). This highlights the twofold need for shrinkage estimation methods and in particular, for the introduction of sparse estimation of the coefficient matrix. In fact, it is very hard both to provide a meaningful interpretation for a large VAR with full time-varying matrix $\mathbf{B}_t$ and to have an efficient and computationally feasible algorithm for an unrestricted estimation.
Motivated by this fact, we provide a prior distribution, which allows for sparse estimation in a time-varying parameter setting. For each coefficient of the matrix of parameters $B_t$, we introduce a mixture prior with independent location and scale parameters:
where $P(\beta_{j,t})$ is the probability distribution of the vector matrix of coefficients (for example, we can choose it as a Double Exponential or Laplace distribution). One of the most successful and widespread approach in the Bayesian literature consists in the use of (independent) spike-and-slab prior distributions (e.g., Mitchell88SpikeSlab_priors, George93VariableSelection, Smith96BNP_SpikeSlab_Dirac, George97SpikeSlabPrior) for each coefficient $\beta_{j,t}$. Based on it, we specify a spike-and-slab prior distribution for each $\beta_{j,t}$, with $j=1,\ldots,n^2$ and $t=2,\ldots,T$, of the form
where $R,Q$ correspond to the spike and slab distributions, respectively, and $\pi_t$ is the time-varying mixing probability (i.e., the prior probability of the spike component). In the literature, we have two commonly choices for $D$: a Dirac mass at $0$ such that $R(\beta_{j,t}) = \delta_{\lbrace 0 \rbrace}(\beta_{j,t})$; and a centered (in zero) Normal distribution $R(\beta_{j,t}) = \mathcal{N}(\beta_{j,t}|0,\tau_0)$. The Dirac spike is a degenerate distribution that allows for variable selection as a by-product of the estimation. Instead, the choice of a continuous, diffuse prior (like a Gaussian) allows for shrinkage of the coefficients and is computationally faster, but requires the post-processing specification of threshold for the sake of variable selection.
The standard choice for the slab component $Q$ is a heavy-tailed distribution belonging to the family of Generalised Hyperbolic distribution (e.g., Double Exponential, Cauchy, $t$-Student), since the aim of this component is to capture potentially large non-zero coefficients. In the case of Dirac spike, the prior for each coefficient, for $j=1,\ldots,n^2$ and $t=2,\ldots,T$, is given by
while for a diffuse (Gaussian) spike we have
As previously described, we study the performance of the two different specification of the spike component by comparing the performances of the two constructions in extracting time-varying Granger causality networks from time series data.
The literature on time-varying parameter (TVP) models is vast. Some common parametric specifications include the threshold AR (TAR, e.g., Tong80ThresholdAR), smooth transition AR (STAR, e.g.,Terasvirta94SmoothTransitionAR), along with their multivariate generalisations, Markov switching process (e.g., Hamilton89MS, Krolzig97MS), change point process (e.g. Pesaran06ChangePoint_VAR) and random walk process (e.g. DelNegroPrimiceri15Bayes_TVP_VAR, Primiceri05Bayes_TVP_VAR). The choice of the particular specifications has been motivated by the intent to capture a particular feature of the dynamic evolution of the coefficients, such as changing regimes, structural breaks or smooth variations.
Differently from the existing literature on TVP-VAR models, we model the temporal dependence of the autoregressive parameters $\beta_{j,t}$ via assuming that the underlying prior (random) distributions evolve according to a discrete-time Markov process. A standard approach in Bayesian nonparametrics involves the specification of a Dirichlet Process (a.k.a. DP, see Ferguson73BNP_DirichletProcess) or a Dirichlet Process mixture (a.k.a. DPM, see e.g., Lo84BNP_DirichletProcessMixture) prior for the distribution of the parameters of interest. The use of DP and related priors for a random probability measure $P$ allows for clustering of the variables $x_i \mathbin{\overset{iid}{\kern\z@\sim}} P$.
We proceed by introducing the prior temporal dependence between the random measures $P_2,\ldots,P_T$ via the time series Dirichlet Process (tsDDP) of NietoBarajas12TimeSeries_DDP. As stated in the paper, for time series models, it is convenient to use dependence on the weights and common location. In practice, we apply a common discretization over the sequence of random measures, while the assumption of common weights and dependent location will lead to a discretization over the probability scale. In opposite to Taddy10AR_DynamicSpatialPoisson, who was working with equally spaced time points, we accommodate for unequal time points. In our analysis a latent binomial process to induce the desired correlation has been used, differently from the stick-breaking random probability measures as in Taddy10AR_DynamicSpatialPoisson, which use a beta autoregression on the fractions of the stick-breaking constructions by mean of two sets of latent variables.
By exploiting the stick-breaking construction of Sethuraman94StcikBreaking_DP, the time series Dirichlet Process imposes a dependence for the random probability measures
where the locations are fixed $\boldsymbol{\theta}_{i,t} = \boldsymbol{\theta}_{i}$ and the weights $w_{i,t}$ vary over time. The dependence is described by a Markov process for each un-normalised stick-breaking weight $v_{i,t}$, with $i=1,\ldots,\infty$, via auxiliary variables $z_{i,t}$ (in the spirit of Pitt02AR_auxiliary_variables, Pitt05AR_auxiliary_variables_general), as follows
The hyper-parameter $m_{i,t}$ tunes the strength of the dependence between $P_t$ and $P_{t+1}$, such that $m_{i,t} =0$ implies $P_t \perp P_{t+1}$, while $m_{i,t} \to \infty$ implies $P_t = P_{t+1}$ with probability 1 (see NietoBarajas12TimeSeries_DDP). Note that this construction implies that at each time $t=2,\ldots,T$, the marginal distribution of each random measure is a Dirichlet Process, that is
with total mass parameter $\alpha$ and base measure $P_0$, such that the base measure defines the expectation and the mass parameter is interpreted as the precision parameter.
Eq. (ref) explains the joint distribution of $z_{i,t}$ and $v_{i,t}$ and it allows us to define the joint model for $(P_1,\ldots,P_T)$ as a $tsDDP(\alpha,P_0,\mathbf{m})$, where $\mathbf{m}$ is the sequence of the strength of dependence, $m_{i,t}$ for $i=1,2,\ldots$ and $t=1,\ldots,T$. In Figure (ref), we show a single draw of $(P_1,\ldots,P_T) \sim tsDDP(\alpha,P_0,\mathbf{m})$, with total mass parameter $\alpha$ equal to $10$; base measure $P_0$ as a normal distribution with zero mean and variance $\sqrt{2}$, i.e. $P_0 \sim \mathcal{N}(0,\sqrt{2})$ and total timing $T=6$.
In order to assess the dependence structure induced by the time series Dirichlet Process, we consider the correlation between two random probability measures $P_{t}$ and $P_{t+1}$. The following proposition is explaining this correlation:
We can summarize what we have described above in the following prior structure\footnote{We use the shape-scale parametrisation of the Gamma distribution (thus $\mathbb{E}[x] = ab$ and $\mathbb{V}[x] = ab^2$) and Inverse Gamma distribution, whose probability density functions are, respectively
The exponential distribution is obtained as a particular case when $a=1$, that is $\mathcal{E}xp(x|b) = \mathcal{G}a(x|1,b)$.}, for $j=1,\ldots,n^2$ and $t=2,\ldots,T$ as follows
where $R(\beta_{j,t})$ is either a Dirac mass at $0$ or a Normal distribution centered in zero and with variance $\tau_0$, which marginally has Inverse Gamma prior distribution $\tau_0 \sim \mathcal{IG}(\tau_0|a_0,b_0)$. If we marginalize over $\lambda_{j}$, we have a Double Exponential slab distribution for each entry of the coefficient matrix. Following the notation of Eq. (ref) we have $Q(\beta_{j,t}) = \mathcal{D}E(\beta_{j,t}|\mu_{j},\tau_{j})$, resulting in
For the covariance matrix $\Sigma$, we assume an Inverse Wishart prior distribution that is:
where $\nu$ and $\boldsymbol{\Psi}$ are the degrees of freedom and scale hyperparameters, respectively.
In summary, the observational model in Eq. (ref) together with the prior structure in Eqs. (ref) to (ref) lead to the BNP-TVP-VAR(1) model. Eq. (ref)--(ref) represent our hierarchical prior and in Figure (ref) we represent them through a Directed Acyclic Graph (DAG) for the Normal spike specification. The observable and non-observable random variables are indicated through shadow and empty circles, respectively. On the left side we have the prior for $\Sigma$, while on the right side, we have the hierarchial prior for $\beta_t$, with a description of the first and second stage of the hierarchy by means of the shrinking parameters $\mu$, $\tau$ and $\lambda$.
\tikzstyle{hyper}= [circle, fill=white, inner sep=1pt, dashed, minimum size=20pt, font=\fontsize{10}{10}\selectfont, draw=black] \tikzstyle{vertex} = [draw, black, circle, minimum size=25pt, inner sep=0pt] \tikzstyle{plate} = [draw, rectangle, rounded corners, minimum width=34pt, minimum height=34pt, fit=#1] \tikzstyle{dots} = [circle,text width=1.7em,minimum size=23pt,text centered, inner sep=0pt]
Sufficient conditions for stationarity of TVP autoregressive models are given for the univariate case (despite the proof is valid also in the multivariate setting) in Brandt86Stationary_TVP_AR, while Bourgerol92Stationary_TVP_AR_iid provides conditions for multivariate regressions where the coefficients are independent and identically distributed. The sufficient conditions given by Brandt86Stationary_TVP_AR is reported below.
Following NietoBarajas12TimeSeries_DDP, we assume $m_{i,t} = m$, for each $i=1,2,\ldots$ and $t=1,\ldots,T$. Higher values of $m$ strengthen the dependence between the un-normalised weights $v_{i,t}$, however when big $m$ may induce the prior to overcome the likelihood, especially when the sample size is small. For this reason they specify a Poisson prior distribution for $m$, truncated on $\lbrace 1,\ldots,5 \rbrace$. Given the complexity of our prior specification, we prefer to fix the value of $m=5$, which is sufficiently small to avoid overweighting of the prior\footnote{In our empirical application, the sample size is $T=248$, while NietoBarajas12TimeSeries_DDP have $T=8$.} and then check the robustness of the results to alternative values of $m$. We choose the following values for the hyper-parameters:
This choice amounts to assuming a uniform prior on each $\pi_t$ and a rather uninformative prior on the covariance matrix, $\Sigma$. The value of the concentration parameter $\alpha$ is set according to standard practice in Dirichlet Process literature. The hyper-parameters $(c,d,a_1,b_1)$ imply that for each new component of the Dirichlet Process the prior distribution of $\mu_j$ is centered at zero mean with medium-high variance, whereas the prior for $\tau_j$ has mean $2$. Instead, the values of $(a_0,b_0)$ imply that the prior variance of the (diffuse) spike distribution is $0.8$, reflecting that this component should account for coefficients $\beta_{j,t}$ not significantly different from zero.
Since the joint posterior distribution is not tractable and it is complex to be sample from, Bayesian estimator cannot be obtained analytically. In this paper, we rely on simulation based inference methods, and develop a Gibbs sampler algorithm for approximating the posterior distribution.
In order to deal with the finite mixture provided by the spike-and-slab prior and the infinite mixture given by the DPM, we exploited a data augmentation approach. For each $j=1,\ldots,n^2$ and $t=2,\ldots,T$, we introduce two sets of allocation variables $\gamma_{j,t},d_{j,t}$; a set of stick-breaking variables, $\mathbf{v}_{t} = \{v_{i,t}:i=1,2,\ldots\}$; a set of auxiliary variables $z_{i,t}$ (for $i=1,2,\ldots$) and a set of slice variables, $\mathbf{u}_{t} = \{u_{j,t}:j=1,\ldots,n^2\}$. The allocation variables, $\gamma_{j,t}$, selects the spike component $R(\cdot)$, when $\gamma_{j,t}$ is equal to zero and the slab component, when it is equal to one. The second allocation variable, $d_{j,t}$, selects the component of the Dirichlet mixture to which each single coefficient $\beta_{j,t}$ is allocated to. The sequence of stick-breaking variables defines the mixture weights, whereas the slice variable, $u_{j,t}$, allows us to deal with the infinite mixture components by identifying a finite number of stick-breaking variables to be sampled and an upper bound for the allocation variables $d_{j,t}$.
Finally, we obtain the following joint posterior distribution
where $\mathbf{U}=\{u_{j,t} : j= 1,\ldots, n^2; \mbox{ and } t = 2,\ldots,T\}$ and $\mathbf{V}=\{v_{i,t}: i =1,2,\ldots \mbox{ and } t=2,\ldots,T\}$ are the collections of slice variables and stick-breaking components, respectively; $\mathbf{Z} = \{z_{i,t}: i =1,2,\ldots \mbox{ and } t=2,\ldots,T\}$ and $\boldsymbol{\lambda} = \{\lambda_j: j = 1,\ldots,n^2\}$ are the auxiliary and latent variables, respectively; $\mathbf{D} = \{d_{j,t}: j= 1,\ldots, n^2; \mbox{ and } t = 2,\ldots,T\}$ and $\Gamma=\{\gamma_{j,t}: j= 1,\ldots, n^2; \mbox{ and } t = 2,\ldots,T \}$ are the allocation variables; $(\bm{\mu},\bm{\tau})= \{(\mu_{k},\tau_k): k =1,\ldots, k^{\ast}\}$ are the atoms, where $k$ ranges from $1$ to the number $k^*$ of allocated DP components; $\mathbf{B} = \{\bm{\beta}_t: t= 2,\ldots, T\}$ is the vector of VAR coefficients and $\bm{\pi} = \{\pi_t: t=2,\ldots,T \}$ are the specific probabilities of shrinking coefficients to zero.
We obtain random samples from the posterior distributions by Gibbs sampling. The Gibbs sampler is based on the algorithm of Hatjispyros11Dependent_DPM and on the slice sampler approach of Walker07SliceSampler_DPMixture and Kalli11SliceSampler_DPM for estimating the weights and locations of each random measure $P_t$. For improving the mixing of the MCMC, we introduced some Hamiltonian Monte Carlo (see Neal11HamiltonianMC) steps in spite of drawing from the full conditional posterior distribution. Hereafter, we show the iterative steps by using the conditional independence between variables, for $k=1,\ldots,k^{\ast}$, $j= 1,\ldots, n^2$, $i=1,2,\ldots$ and $t=2,\ldots,T$:
The detailed Gibbs sampler is described in (ref) and (ref).
Based on the Gibbs sampler previously described, we are able to extract time-varying Granger-causal graphs. In the literature, linkages and networks describing the relationships between variables of interest, such as macroeconomics and financial linkages (e.g. Billioetal12GrangerNet and Barigozzi19NETS_network_estimation) can be used to extract pairwise Granger causality. This approach is generating spurious causality effects and does not consider conditioning on variables of interest. The main problem relies on the high number of variables available relative to the number of data, thus it could lead to overparametrization and inefficiency in gauging the causal relationships. Our proposed prior can be used to extract the networks and pairwise Granger causality while reducing the overfitting and curse of dimensionality problems. Moreover, the introduction of our prior could lead to the extraction of edge-colored graphs, that allows us to identify stylized facts in financial or macroeconomics networks and to show the presence of communities, hubs and linkage heterogeneity.
From the MCMC output of the time-varying coefficient matrix $\mathbf{B}_t$, we are able to extract time-varying Granger-causal graphs. At each time $t=2,\ldots,T$, we use the posterior random partition induced by the nonparametric (slab) distribution to cluster the edges of the graph (i.e., the entries $\beta_{j,t}$, $j=1,\ldots,n^2$ of the vectorised coefficient matrix $\boldsymbol{\beta}_t$) into groups.
Formally, a graph $G$ is a pair $(V,E)$, where $V$ is a set of nodes and $E$ is a set of nodes pairs, named links or edges. The nodes are labeled and a link/edge is identified by the pair of nodes it connects, $(i,j)$. In particular, we have the existence of an edge if and only if the time-varying VAR coefficients of the variable $y_{i,t-1}$ in the equation of $y_{j,t}$ is not null. In our network analysis, we focus on the adjacency matrix constructed a posterior from the allocation variables and it allows to take both values between $0$ and $1$ if we apply a threshold, while if the values are allowed to vary between $0$ and $1$, we have a weighted graph. The purpose is to estimate the most significant time-varying dependence interrelationships (in terms of Granger-causality) between the $n$ variables of interest.
We proposed the BNP-TVP-VAR model for sparse, nonparametric inference in time-varying VAR models. The use of spike-and-slab priors with time-series dependent Dirichlet Process prior for the slab component allows to contemporaneously shrink the autoregressive coefficients and flexibly modelling time-varying non-zero entries. We applied the proposed methodology using two alternative spike distributions: a Dirac and a Normal distribution. The performance of the resulting models has been compared in terms of: (i) sparse estimation and variable selection, and (ii) clustering structure. Moreover, we showed how the BNP-TVP-VAR model can be used for extracting Granger-causal time-dependent graphs from multivariate time series.