Markov chain simulation

Question

I'm wondering whether there is an algorithm to simulate a discrete Markov chain with a specific number of occurrences of state knowing the transition matrix way.

For example, how to simulate in R a Markov chain of length $n$ with $p$ occurrences ($p < n$)

TransitionMatrix<- matrix(c(0.7, 0.3, 0.4, 0.6),byrow=TRUE, nrow=2)
colnames(TransitionMatrix) <- c('0','1') row.names(TransitionMatrix) <- c('0','1')

Answer 1

Suppose you have transition probabilities $p_1 \leq p_2 \dots \leq p_k$, $\sum p_k = 1$. Define the distribution function in pseudo code:

function cdf(state) { a = math.random() % randomly generated number in U(0,1)
if 0 <= a < p_1, return 1
if p_1 <= a < p_1 + p_2, return 2
...
if p_1 + ... + p_{k-1} <= a < 1, return k

In the above, I have assumed that the transition probabilities are all the same, i.e. independent of the state. Otherwise, the $p_i$ will have to depend on state.

The above code is just saying that if I have numbers $p_1$, ..., $p_n$ such that $\sum p_i = 1$, then I can partition the unit interval $[0,1]$ with these numbers. The area of the rectangle [$p_1 + ... + p_{k-1}, p_1 + ... + p_k$] with unit height is $p_k$, and that's the probability we want for the $k^{th}$ transition. This is the standard way of generating probability distributions, since if $X$ has cdf $F$, then $F(X) \sim U(0,1)$.

Now pick an initial state (say $s = 1$) and a state you want to measure, say $r = 4$.

count = 0
for i = 1 : n
   s = cdf(s) // transition with the appropriate probabilities
   if s == 4, count++
   end
end

I should mention this algorithm is the worst.

Answer 2

First a reformulation of the question, then a pseudo-algorithm to solve it.

Reformulation: Let $n\geqslant1$ and assume that $X=(X_t)_{1\leqslant t\leqslant n}$ is a Markov chain of length $n$ on the state space $\{0,1\}$ with transition probability matrix $\begin{pmatrix}1-a & a\\ b& 1-b\end{pmatrix}$ for some $a$ and $b$ in $(0,1)$.

Fix $k\leqslant n$ and consider the set $\Gamma_k^n$ of paths $\gamma=(x_t)_{1\leqslant t\leqslant n}$ of length $n$ visiting $k$ times the state $1$. For every $x$ and $y$ in $\{0,1\}$ and every path $\gamma$ in $\Gamma_k^n$, let $N_\gamma(xy)$ denote the number of times $t\leqslant n$ such that $(x_t,x_{t+1})=(x,y)$. Finally, let $M(\gamma)=N_\gamma(01)$.

Then, forgetting a possible discrepancy of $\pm1$ here and there, one sees that $N_\gamma(10)=M(\gamma)$, $N_\gamma(00)=n-k-M(\gamma)$ and $N_\gamma(11)=k-M(\gamma)$, hence, for every $\gamma$ in $\Gamma_k^n$, $$ P(X=\gamma)= (1-a)^{N_\gamma(00)}a^{N_\gamma(01)}b^{N_\gamma(10)}(1-b)^{N_\gamma(11)}, $$ is also, since $n$ and $k$ are fixed, $$ P(X=\gamma)= (1-a)^{n-k-M(\gamma)}a^{M(\gamma)}b^{M(\gamma)}(1-b)^{k-M(\gamma)}\propto c^{M(\gamma)}, $$ where $$ c=\frac{ab}{(1-a)(1-b)}. $$ Thus, there exists some $A_k^n$ independent of $\gamma$ such that, for every $\gamma$ in $\Gamma_k^n$, $$ P(X=\gamma\mid X\ \text{visits}\ 1\ \text{exactly}\ k\ \text{times})=A_k^n\cdot c^{M(\gamma)}, $$ and the question is to simulate a random path in $\Gamma_k^n$ following this distribution.

Pseudo-algorithm: This assumes one is able to generate uniformly subsets $T$ of $\{1,2,\ldots,n\}$ of size $k$.

• Choose a path $\gamma$ uniformly in $\Gamma_k^n$, that is, choose uniformly a subset $T\subseteq\{1,2,\ldots,n\}$ of size $k$ and define $\gamma=(x_t)$ by $x_t=1$ if $t$ is in $T$ and $x_t=0$ otherwise.
• Compute $M(\gamma)$.
• Accept $\gamma$ with probability proportional to $c^{M(\gamma)}$.
• Repeat until a path $\gamma$ is accepted. Return $\gamma$.

To generate $T$, the following procedure might prove useful:

• Let $T=\varnothing$.
• Choose $x$ uniformly in $\{1,2,\ldots,n\}$. If $x$ is not already in $T$, add $x$ to $T$.
• Repeat until $T$ has size $k$. Return $T$.