Gibbs sampler gets stuck in local mode I am very new to statistics and trying to implement a Gibbs sampler. However, according to wikipedia https://en.wikipedia.org/wiki/Gibbs_sampling and this discussion thread http://metaoptimize.com/qa/questions/7831/disadvatage-of-gibbs-sampler, it is possible for Gibbs sampling to get stuck in a local mode. When this is the case, how do people tend to resolve this issue?
I am attempting to sample a distribution in the context of bioinformatics, where the joint distribution is difficult to compute analytically. The reference paper is here http://bioinformatics.oxfordjournals.org/content/31/12/i97.long, though it might be too wordy to go through. 
Any general directions will be appreciated.
 A: This is certainly possible. It often happens when variables are strongly correlated. 
For simplicity, consider a two-parameter model. Because Gibbs sampling alters only one variable at a time, it can only move vertically or horizontally on the Cartesian plane. It will be unable to reach regions of high posterior probability that lie diagonally from one another. 
You may be able to resolve this problem by reparameterizing your model so the posterior is roughly spherical or by jointly updated blocks of correlated parameters. For an example of how to reparameterize, you could take a regression model and standardize the covariates. This makes is so that the slope won't have much effect on the fit of the intercept and vice versa. Without standardizing, data may sit far from the origin, and a small change in slope could demand a huge change in intercept for the regression line to fall across the data.
MCMC implementation is difficult and error prone, so you may want to consider tools such as STAN. STAN provides already-implemented state-of-the-art MCMC methods via a language similar in spirit to R's formula syntax. It has interfaces with many common scientific computing languages. (I am not affiliated with STAN -- merely a believer.)
A: Here's a simple case where Gibbs sampling gets stuck:
Imagine $(\theta_1,\theta_2)$ is a 50-50 mixture of two bivariate normals which each have independent components with variance 1.
The first mixture component is centered at $(10,-10)$. The second is centered at $(-10,10)$.

If you're in the top-left (green) mode, you'll stay there. If you're in the red mode, you'll stay there instead.
Each component is a local mode in the bivariate distribution. When you start, you'll right away get stuck in one mode or the other -- you can easily move around within either mode, but not between them.
In a situation like this one, there are a variety of strategies that might be used. For example, sometimes it is feasible to sample affected variables as a block. Sometimes you can integrate some variables out, which can reduce or even fix the problem. Sometimes you might be able to reparameterize (rewriting as $(\psi_1,\psi_2)=(\theta_1+\theta_2,\theta_1-\theta_2)$ would work for the abobe example). You may be able to add a variable (or more than one) that allow movement between modes that could not previously communicate. There are numerous other strategies; some more suited to particular kinds of multiple modes.
You might also move to a completely different approach, such as Hamiltonian Monte Carlo.
