@MCMCPOSTPROC is a "post-processor" for Markov Chain Monte Carlo statistics. The information from the draws needs to be saved into a SERIES[VECTOR]. Each entry in the stats SERIES is a VECTOR of desired statistics generated from a draw. The job of @MCMCPOSTPROC is to extract sample information for each of those desired statistics.

 

Note that if the number of draws is large (more than 10000), it can take quite a while to do the diagnostics and the percentiles.

 

@MCMCPostProc(NDRAWS=# of draws,other options) stats

Parameters

Options

NDRAWS=# of draws (length of stats) (REQUIRED)

 

MEAN=(output) [VECTOR] averages of each component

STDERRS=(output) [VECTOR] standard errors of each component

CV=(output) [SYMMETRIC] estimated covariance matrix

NSE=(output) [VECTOR] numerical standard errors of each component

The MEAN and STDERRS are the sample mean and standard errors for each element of the VECTOR across the draws, so if you do MEAN=BMEANS, BMEANS(1) will be the sample mean across draws for the first statistic. CV is the full sample covariance matrix of the saved statistics. NSE are corresponding estimates of the standard errors of the MEAN component, which are computed using HAC estimators since the Markov Chain estimates in general are at least somewhat autocorrelated. The NSE should decrease to zero as the number of draws increases.

 

PERCENTILES=(input) ||desired percentiles|| [not used]

QUANTILES=(output) VECT[VECTOR] with percentiles

If you want percentiles of the distribution, use the PERCENTILES to make the request and QUANTILES (which is a VECTOR[VECTOR]) for the output. Each outer element of QUANTILES is for a statistic, and the inner elements are the computed percentiles in the same order as the request. So PERCENTILES=||.16,.50,.84||,QUANTILES=Q will mean that Q(1) will be a VECTOR with the .16,.50 and .84 quantiles of the first saved statistic, similarly Q(2) for the second statistic.

 

CD=(output) [VECTOR] Geweke CD measures

BW=(output) [VECTOR] Between-within measures

BLOCKS=number of blocks for BW [10]

CD and BW are diagnostics (see Geweke(1992) for a description). CD is a comparison of early draws against late draws and is mainly a measure of how well the chain has been "burned-in". Asymptotically, it's N(0,1) component by component with neither sign expected over the other. BW computes the ratio of between to within variation within the chain, using the number of blocks given by the BLOCKS option, that is, if you have 25000 draws and the default BLOCKS=10, it will look at 10 blocks of 2500. If the draws are independent, this will have an F distribution. A large value indicates that the draws tend to be highly correlated for long stretches.

Examples

@MCMCPOSTPROC produces no displayed output—it's up to you to decide what statistics you want and how to display them. This is for Gibbs sampling on a regression. BGIBBS has the regression coefficients for each draw. @MCMCPOSTPROC computes the mean, standard error, numerical standard error (which is an estimate of the precision with which the mean is estimated) and the CD measure. These are displayed in a REPORT, with each row labeled using the regressor label from the regression.

 

@mcmcpostproc(ndraws=ndraws,mean=bmean,$

stderrs=bstderrs,cd=bcd,nse=bnse) bgibbs

report(action=define)

report(atrow=1,atcol=1,align=center) "Variable" "Coeff" $

"Std Error" "NSE" "CD"

do i=1,%nreg

   report(row=new,atcol=1) %eqnreglabels(0)(i) bmean(i) $

     bstderrs(i) bnse(i) bcd(i)

end do i

report(action=format,atcol=2,tocol=3,picture="*.###")

report(action=format,atcol=4,picture="*.##")

report(action=show)

 

This is from Gibbs sampling on a two-regime Markov switching model.

@mcmcpostproc(ndraws=ndraws,mean=bmeans,stderrs=bstderrs) bgibbs

report(action=define)

report(atrow=1,atcol=2) "Mean" "S.D"

report(atrow=2,atcol=1) "mu(1)" bmeans(1) bstderrs(1)

report(atrow=3,atcol=1) "mu(2)" bmeans(2) bstderrs(2)

report(atrow=4,atcol=1) "phi(1)" bmeans(3) bstderrs(3)

report(atrow=5,atcol=1) "phi(2)" bmeans(4) bstderrs(4)

report(atrow=6,atcol=1) "sigmasq" bmeans(5) bstderrs(5)

report(atrow=7,atcol=1) "p" bmeans(6) bstderrs(6)

report(atrow=8,atcol=1) "q" bmeans(9) bstderrs(9)

report(action=format,picture="###.###")

report(action=show)

 

Copyright © 2026 Thomas A. Doan