I have some questions regarding convergence diagnostics, especially how to read the univariate and multivariate plots.
In the attached file, you may find six plots. “convergence_udiag1_short” and “convergence_udiag2_short” are for univariate measures and “convergence_mdiag_short” is for multivariate measures. These are based 400000 draws. According to the relevant documents, for example “Pfeifer (2014): An Introduction to Graphs in Dynare at sites.google.com/site/pfeiferecon/dynare”, the parameters in “convergence_udiag1_short” seem to converge because the two lines seem to stabilize horizontally and be close to each other. While in “convergence_udiag2_short”, the parameters seem to be problematic as the two lines have not stabilized horizontally after 400000 draws, although the two lines are quite close to each other. In “convergence_mdiag_short”, the multivariate measures seem to behave well. In order to see if the problem arises because the chain is not long enough, I run MCMC with “load_mh_file” and the chain now contains 1 million draws. The plots are “convergence_udiag1_long”, “convergence_udiag2_long” and “convergence_mdiag_long”. As for the short chain, “convergence_udiag1_long” and “convergence_mdiag_long” seem OK. But for “convergence_udiag2_long”, I am not quite sure. It seems that the parameters in this plot mix very slowly.
My first question is how to understand these plots and more generally, how to judge the convergence of MCMC. In this case, the multivariate plot seems ok but the univariate plots do not. Do I need even more draws to ensure all the chains to converge?
To provide further evidence of the model, I also included two mode_check plots and two prior_posterior plots. I found that two persistence parameters of shock processes are very high, e.g. 0.998 for rho_d and rho_g. So the model itself seems not to absorb the persistence from the data and the shock processes do that.
My second question is that if the slow or non convergence is related to the potentially problematic parameter estimates?
See the original paper Brooks/Gelman for the interpretation. You want convergence for every single parameter. Convergence of the multivariate statistic is necessary but obviously not sufficient (it is basically just a summary that takes together all univariate diagnostics, but one where jumps in two parameters at the same time that leave the likelihood almost unaffacted are not detected). If you find this too adhoc, use just one chain and the Geweke diagnostics. Here, you get a p-value for a test statistic.
High autocorrelation often is an indication of problems in the data treatment (most often forgotten constants in the observation equation). However, there is usually not a direct connection between autocorrelation in the model and the one of the MCMC. A problem rather occurs if you are close to the upper or lower bound of the parameters. In this case, depending on the Hessian you might get many draws that are rejected. Mixing will be poor due to a too low acceptance rate.
Regarding your second point, I think my observation equation and data treatment should be OK. However, what is problematic in my model, as you said, is that autocorrelation coefficient of shock process is too high, which is close to the upper bound of the parameter. Specifically, I have three shocks, a preference shock (d shock), a government spending shock (g shock) and a TFP growth shock (z shock). The AR coefficients of d shock and g shock are 0.9987 and 0.9979 at the mode. The standard deviation of these two parameters at the mode are almost 0. If looking at the mode_check plot in the attachment, rho_d is almost flat around the mode. If looking at the prior_posterior plot, these two parameters are spikes around the mode. I also attach the trace plots for these two parameters, which are drifting all the time among the 1 million draws and there is no signs of mixing.
I am wondering if all these evidence are related to each other and how to deal with it. I am afraid this is not merely a numerical problem. It may be related to the model specification. The current model I am estimating is close to a small-scale DSGE model as in An and Schorfheide (2007). The model only has price rigidity. The other popular frictions are not included, such as habit formation, wage rigidity and etc. Without enough frictions, it seems that the shock processes are required to be very persistent (close to unit root) to match the data. I stick to this model because it is possible and easy to derive analytical form of boundary between determinacy and indeterminacy regions. Within this model, do you think there is any possible ways to deal with the current problem?
You need to find out why the model assigns close to a unit root to the data. This is very uncommon, and as I said, usually hints at a problem with the data specification or the model. If your data for example is in first differences, it would be strange for the model to imply a unit root as this would imply the data to be I(2).