How to Find Stationary Distributions of a Squared Transition Matrix?

Question

Given a Markov chain with a transition matrix $p$ \begin{pmatrix} 0 & 0 & .1 & .9 \\ 0 & 0 & .6 & .4 \\ \ .8 & .2 & 0 & 0 \\ .4 & .6 & 0 & 0 \end{pmatrix} I can easily determine the stationary distribution. However, if I want to find the stationary distributions of $p^2$ \begin{pmatrix} .44 & .56 & 0 & 0 \\ .64 & .36 & 0 & 0 \\ \ 0 & 0 & .2 & .8 \\ 0 & 0 & .4 & .6 \end{pmatrix} I have much more difficulty.

My first step is to set up a system of equations:

$.44π_1+.64π_2+0π_3+0π_4 = π_1$

$.56π_1+.36π_2+0π_3+0π_4 = π_2$

$0π_1+0π_2+.2π_3+.4π_4 = π_3$

$0π_1+0π_2+.8π_3+.6π_4 = π_4$

Since $π_1+π_2+π_3+π_4=1$, I replace any equation in the system with this equation and subtract $π_1, π_2, π_3$ from the other equations:

$-.56π_1+.64π_2+0π_3+0π_4 = 0$

$.56π_1-.64π_2+0π_3+0π_4 = 0$

$0π_1+0π_2-.8π_3+.4π_4 = 0$

$π_1+π_2+π_3+π_4 = 1$

I now have: \begin{equation} \begin{pmatrix} π_1 & π_2 & π_3 & π_4 \\ \end{pmatrix} * \begin{pmatrix} -.56 & .56 & 0 & 1 \\ .64 & -.64 & 0 & 1 \\ \ 0 & 0 & -.8 & 1 \\ 0 & 0 & .4 & 1 \end{pmatrix} = \begin{pmatrix} 0&0&0&1\\ \end{pmatrix} \end{equation}

But since this matrix is not invertible, I cannot solve for \begin{pmatrix} π_1 & π_2 & π_3 & π_4 \\ \end{pmatrix} by multiplying \begin{equation} \begin{pmatrix} 0&0&0&1\\ \end{pmatrix} * \begin{pmatrix} -.56 & .56 & 0 & 1 \\ .64 & -.64 & 0 & 1 \\ \ 0 & 0 & -.8 & 1 \\ 0 & 0 & .4 & 1 \end{pmatrix}^{-1} \end{equation}

Can anyone help me find the stationary distributions for $p^2$?

Answer 1

Your method is rather unusual. There’s no guarantee that the matrix that you form will be invertible. You’ve been lucky so far in only encountering Markov chains that don’t separate into noncommunicating classes, like this one. Each of the classes has its own stationary distribution, and the general solution is a convex combination of them.

A more conventional way to do this is to directly solve the equation $\mathbf\pi P=\mathbf\pi$, or $\mathbf\pi(P-I)=0$. This can be done by computing the null space of $(P-I)^T$, which in turn can be done via Gaussian elimination. However, since the chain decomposes into two noncommunicating classes, it might be simpler to work out the two individual stationary distributions, call them $\mathbf\pi_A$ and $\mathbf\pi_B$, which you should be able to do using your method on the $2\times2$ submatrices. Any convex combination $(1-\lambda)\mathbf\pi_A+\lambda\mathbf\pi_B$, with $0\le\lambda\le1$, of these two stationary distributions will also be a stationary distribution.

To continue with your method, proceed as you might for any system of linear equations. First, move the variables to the other of the product by transposing everything, then augment the coefficient matrix with the vector on the right-hand side of the equation, producing $$\pmatrix{-.56&.64&0&0&0\\.56&-.64&0&0&0\\0&0&-.8&.4&0\\1&1&1&1&1}.$$ Perform Gaussian elimination (a.k.a. row.reduction) on this matrix. (Moving the bottom row to the top before doing anything else will make this easier.) There will be an infinite number of solutions, corresponding to the fact mentioned above that there is no unique stationary distribution for this particular chain.

Answer 2

You can use Gaussian elimination on this linear system. Since your coefficient matrix has rank $3$, there will be a free variable, and the solution will be nonunique.

Answer 3

Hint Suppose we have a Markov chain with a block diagonal transition matrix $$Q = \pmatrix{Q_1&\\&Q_2}$$ like $p^2$. For simplicity we'll treat the $2$-block case, but the general case is similar.

If $\Pi$ is a vector representing a stationary distribution of $Q$, we can decompose $\Pi = \pmatrix{\Pi_1&\Pi_2}$ into vectors $\Pi_1, \Pi_2$ of appropriate sizes, and then the equation $\Pi Q = \Pi$, that is,

$$\pmatrix{\Pi_1&\Pi_2} \pmatrix{Q_1&\\&Q_2} = \pmatrix{\Pi_1&\Pi_2} ,$$ decomposes into the independent equations $$\Pi_1 Q_1 = \Pi_1, \qquad \Pi_2 Q_2 = \Pi_2 .$$

In particular, if $\Sigma_1$ is a stationary distribution for $Q_1$ and $\Sigma_2$ is a stationary distribution for $Q_2$, we can check that $$\pmatrix{\lambda_1 \Sigma_1 & \lambda_2 \Sigma_2}$$ is a stationary distribution for $Q$ for any $\lambda_1, \lambda_2 > 0$ such that $\lambda_1 + \lambda_2 = 1$. Put another way, if we identify $\Sigma_1$ and $\Sigma_2$ respectively with the stationary distributions $\pmatrix{\Sigma_1 & 0}$ and $\pmatrix{0 & \Sigma_2}$ of $Q$, any convex combination of those distributions is again a stationary distribution.

In any case, the free parameter here shows that there is not a unique stationary parameter. To find a stationary distribution in general, we can solve $\Pi Q = \Pi$, or just as well $\Pi (Q - I) = 0$ subject to $\Pi \pmatrix{1 \\ \cdots \\ 1} = 1$.

It might be instructive to consider the extreme case $Q = I$.