I'm currently studying the Perfect Sampling approach to Markov Chain Monte Carlo proposed by Propp and Wilson in 1996.Though having understood most important aspects, I'm still struggling to figure out their example used for explaining that one must reuse the random numbers generated in earlier iterations. The example reads as follows and is allegedly "simple to check".
It is sometimes desirable to view the process as an iterative one, in which one successively starts up $n$ copies of the chain at times $-1$, $-2$, etc., until one has gone sufficiently far back in the past to allow the different histories to coalesce by time $0$. However, when one adopts this point of view (and we will want to do this in the next subsection), it is important to bear in mind that the random bits that one uses in going from time $t$ to time $t+1$ must be the same for the many sweeps one might make through this time-step. If one ignores this requirement, then there will in general be bias in the samples that one generates. The curious reader may verify this by considering the Markov chain whose states are $0$, $1$, and $2$, and in which transitions are implemented using a fair coin, by the rule that one moves from state $i$ to state $\min(i+1,2)$ if the coin comes up heads and to state $\max(i-1,O)$ otherwise. It is simple to check that if one runs an incorrect version of our scheme in which entirely new random bits are used every time the chain gets restarted further into the past, the samples one gets will be biased in favor of the extreme states $0$ and $2$.
So obviously the transition matrix of the described Markov Chain is $$ P = \begin{pmatrix} 1/2 & 1/2 & 0 \\ 1/2 & 0 & 1/2 \\ 0 & 1/2 & 1/2 \end{pmatrix} $$ with state space $S = \lbrace 0,1,2 \rbrace$. So the stationary distribution is $ \pi = \begin{pmatrix} 1/3 & 1/3 & 1/3 \end{pmatrix}$ and the claim is that when applying Perfect Sampling wrongly( that is with regeneration of the random transitions) we would get states $0$ and $2$ with a higher frequency than $1/3$.
I tried modeling the dynamics of the coupled Markov-Chains started at each of the 3 different states. To do so I defined the following $2$-dimensional Markov-Chain $$ X_t = (X_t^1, X_t^2), $$ where $X_t^1$ denotes the number of the $3$-chains currently in state $0$ and $X_t^2$ denotes the number of the $3$-chains currently in state $1$. The number of chains in state $2$ would then follow from $3 = X_t^1 + X_t^2$. Then the state space of $X$ is $$ S = \lbrace (0,0), (0,1), (0,2), (0,3), (1,0), (1,1), (1,2), (2,0), (2,1), (3,0) \rbrace $$ and the respective transition matrix I came up with is (in R code)
m <- matrix(c(.5, 0, 0, .5, 0, 0, 0, 0, 0, 0,
.25, 0, .25, 0, .25, 0, .25, 0, 0, 0,
.25, .25, 0, 0, 0, 0, 0, .25, .25, 0,
.5, 0, 0, 0, 0, 0, 0, 0, 0, .5,
0, .25, 0, .25, .25, 0, .25, 0, 0, 0,
0, .125, .125, 0, .125, .25, .125, .125, .125, 0,
0, .25, 0, 0, .25, 0, 0, 0, .25, .25,
0, 0, .25, .25, 0, 0, 0, .25, .25, 0,
0, 0, .25, 0, 0, 0, .25, .25, 0, .25,
0, 0, 0, .5, 0, 0, 0, 0, 0, .5),
byrow=TRUE, nrow=10)
But using this transition matrix and computing the distributions for a start in the initial state $(1,1)$ suggests that the $3$ "coalescence states" $(3,0), (0,3) $ and $(0,0)$ have equal probability.