Convergence Diagnostic in Bayesian VAR

[Jules89]

Hi Tom,

Is there any convergence diagnostic for the estimation of the Bayesian VAR models?
In general if you for example use a independent Minnesota-Wishart-Prior for a VAR and apply a Gibbs-Sampler you would discard a number of draws as burn-in-draws and then use the remaining ones for monte carlo integration. If you would take enough draws the Gibbs kernel would converge to its target distribution.
I would like to know whether there is a way of monitoring that convergence to establish whether I have taken enough draws.

Thanks

Jules

[TomDoan]

@MCMCPOSTPROC has a number of Markov Chain diagnostics, the use of which is covered in the Bayesian Econometrics e-course. However, unless your prior is seriously in conflict with the data, there are usually no convergence issues with that (switching between sigma|beta and beta|sigma). It's structural VAR's that can have major difficulties.

[Jules89]

Dear Tom,

the mcmcpostproc only gives the CD measure of Geweke and a Between-Within measure to check correlation of the chain.
Is there also a procedure or a code example which performs the convergence diagnostic of "Gelman and Rubin (1992) - Inference from Interative Simulaion Using Multiple Sequences"?

Thank you in advance

Best

Jules

[TomDoan]

No, though almost all convergence diagnostics are relatively simple statistics on a sequence of numbers. Multiple sequences generally only provide additional information if the chain is (very) sensitive to initial conditions and you have a way to come up with different initial conditions which can test that.

[Jules89]

Dear Tom,

I have rewritten the example code for the Bayesian VAR with Minnesota Prior for impulse response analysis. Furthermore I added some code to check for convergence of each VAR coefficient. The Code is the following:

Code: Select all

*** Data

open data haversample.rat
cal(q) 1959
data(format=rats) 1959:1 2006:4 ftb3 gdph ih cbhm

set loggdp = log(gdph)
set loginv = log(ih)
set logc   = log(cbhm)


*** VAR Setup and Storage Matrices

compute lags=6			                  ;*Number of lags
compute nstep=40                    ;*Number of response steps
compute ndraws=100000               ;*Number of keeper draws
compute nburn=95000   	             ;*Number of burn-in draws
compute keep = ndraws - nburn       ;*Number of kept draws
compute nvari = 4                   ;*Number of variables
compute ncoef = nvari*lags+1        ;*Number of coefficients per equation


compute tight=.1
compute other=.5


system(model=varmodel)
   variables loggdp loginv logc ftb3
   lags 1 to lags
   det constant
end(system)

estimate(ols,sigma)  ;*initial set of estimations
impulse(noprint,model=varmodel,steps=nstep,results=impulses)      ;*initial set of of impulse responses


compute nvar=%nvar                                                ;* Number of equations estimated
dec vect olssee(nvar)                                             ;* vector of variances
ewise olssee(i)=sqrt(%sigma(i,i))                                 ;* standard deviations for rescaling

declare vect[rect] %%responses(keep)                              ;* Storage Vector of Rectanculars for Impulse Response Draws
declare series[rect] bgibbs                                       ;* Storage Series of Rectangulars to save the VAR Coefficient Draws

gset bgibbs 1 keep = %zeros(ncoef,nvari)                          ;* To dimension the storage matrix inside the series



***********************************
********* Minnesota Prior *********
***********************************

source BVARBuildPriorMN.src

@BVARBuildPriorMN(model=varmodel,tight=tight,other=other) HBPRIORS HPRIORS
@BVARFinishPrior hbpriors hpriors bprior hprior

compute sigmad=%sigma
cmom(model=varmodel)



*********************************
********* Gibbs Sampler *********
*********************************



infobox(action=define,progress,lower=1,upper=ndraws) "Gibbs Sampler"
do draw= 1,ndraws
   infobox(current=draw)
   *
   * Draw b given sigma
   *
   compute bdraw =%ranmvkroncmom(%cmom,inv(sigmad),hprior,bprior)
   compute rssmat=%sigmacmom(%cmom,bdraw)
   *
   * Draw sigma given b
   *
   compute sigmad=%ranwisharti(%decomp(inv(rssmat)),%nobs)
   if draw<=nburn
      next

   compute %modelsetcoeffs(varmodel,bdraw)
   impulse(noprint,model=varmodel,results=impulses,steps=nstep,flatten=%%responses(draw-nburn))

   * Bookkeeping:

   compute bgibbs(draw-nburn) = bdraw                    ;* Saves the Draws of the VAR Coefficients

end do draw
infobox(action=remove)


display bgibbs(1)(1,4)

*******************************************
********* Convergence Diagnostics *********
*******************************************

declare series[vect] a
gset a 1 keep = %zeros(ncoef,1)

declare series b
gset b 1 ncoef = 0

declare series c
gset c 1 ncoef*nvari = 0

declare series d
gset d 1 ncoef*nvari = 0


do i=1,4
   do k = 1, keep
      do j = 1, ncoef

         compute a(k)(j) = bgibbs(k)(j,i)

      end do
   end do

      @MCMCPOSTPROC(ndraws=keep,cd=cd) a
      display cd

      set up 1 ncoef*nvari = 2
      set down 1 ncoef*nvari = -2

      do l = 1, ncoef

         compute b(l) = abs(cd(l))
         compute c(l+(i-1)*ncoef) = abs(cd(l))

      end do

end do


do i = 1,nvari*ncoef

   if c(i)>2
   compute d(i) = 1

end do
display "Number of non-converging coefficients" %sum(d)

graph(nodates,style=lines,overlay=spike,ovsamescale) 2
# up / 2
# c /11

My first question is: Is the code for the convergence diagnostic correct?
Second: If the code is correct, then you can see that a lot of coefficients are not converging (the number is 14 for 100000 draws with 95000 burn in replications). Even If I increase the number to 1000000 draws with 950000 burn-in samples the number of non-converging coefficients is still large.
Is it in general necessary that all coefficients converge? Because the number of coefficients in a VAR can be quite large and maybe single parameters don't have that much of an impact.
Do you know any scientific paper that deals with that problem?

Best Jules

[TomDoan]

That looks like a situation where if you look for something odd in 100 series you'll find it. If you comment out the recalculation of the sigmad, the draws will be serially independent, but you'll probably see exactly the same phenomenon. (Note that they're different elements of the VAR coefficients each time). I don't know if anyone has ever looked carefully at the behavior of the individual lag coefficients in a Gibbs sampler for a relatively large VAR. At any rate, the B/W measures don't point to any problems, and the draws are close to being uncorrelated, so all indications are that the sampler is working fine.

[Jules89]

In the description of the @mcmcpostproc procedure is written that the bw measure should follow an F distribution. What are the related degrees of freedom?


Thank you

Best Jules

[TomDoan]

BLOCKS-1,NDRAWS-BLOCKS

[Jules89]

Dear Tom,

is it possible to test the series of accepted Gibbs draws for convergence by testing every series for stationarity and serial indepence?
As far as I understand the MCMC methods correct, the aim is to take enough draws from the posterior such that the markov chain converges to a stationary distribution.
But when the chain converged, then the draws should be serially uncorrelated and the resulting series should be stationary.
So if I find the daws to be more or less fluctuating around a constant mean that should indicate that my Gibbs sampler converged, right?

Thanks

Best Jules

[TomDoan]

No. Gibbs draws will always have at least a mild serial correlation. Convergence of the chain means that the effects of the initial conditions have worn off, and you're getting draws from the unconditional distribution. Now with the lag coefficients and sigma's in a BVAR, the effects of the initial conditions are generally gone by about sweep #2. The only reason you seem to be having problems is that you're looking at 100 sequences and the CD measure probably has a tail which is quite a bit thicker than a Normal (probably because of the imprecision in the estimation of the subsample NSE).