# Finding steady state probabilities by solving equation system

(I know that there are numerous questions on this, but my problem is in actually solving the equations, which isn't the problem in other questions.)

I'm trying to figure out the steady state probabilities for a Markov Chain, but I'm having problems with actually solving the equations that arise. So,

$$V = V\times \left[ { \begin{array}{ccc} 0 & 0.4 & 0 & 0.6\\ 0 & 0.5 & 0.5 & 0 \\ 0.4 & 0 & 0.2 & 0.4\\ 0 & 0 & 1 & 0 \end{array} }\right]$$ Where $V = \left[ \begin{array}{cccc}a & b & c & d \end{array} \right]$

So typing out the equations that I get from multiplying V with the transition matrix.

\begin{align} &0a &+ &0b &+ &0.4c &+ &0d&= a\\ &0.4a &+ &0.5b &+ &0c &+ &0d &= b\\ &0a &+ &0.5b &+ &0.2c &+ &d &=c\\ &0.6a &+ &0b &+ &0.4c &+ &0d &=d \end{align}

How do I solve these equations?

I know the answer for a, b, c, d but I get lost in trying to get to the answer from the equations.

## 3 Answers

Start with the good equations:

\begin{align} &0a &+ &0b &+ &0.4c &+ &0d&= a\\ &0.4a &+ &0.5b &+ &0c &+ &0d &= b\\ &0a &+ &0.5b &+ &0.2c &+ &d &=c\\ &0.6a &+ &0b &+ &0.4c &+ &0d &=d \end{align}

now you can just substitute:

\begin{align} &-a & & & &0.4c & &&= 0\\ &0.4a &- &0.5b & & & &&= 0\\ & & &0.5b &- &0.8c &+ &d &=0\\ &0.6a & & &+ &0.4c &- &d &=0 \end{align}

substitute $a = 0.4c$: \begin{align} -&0.5b &+ &0.16c & &&= 0\\ &0.5b &- &0.8c &+ &d &=0\\ && &0.64c &- &d &=0 \end{align}

substitute $d = 0.64c$: \begin{align} -&0.5b&+&0.16c =0\\ &0.5b&-&0.16c =0 \end{align} you get a consistent system.

This confirms that the system has solutions. The general solution is

$$(0.4c, 0.32c, c, 0.64c)$$

Now choose $c$ such as $a+b+c+d = 1$ and check that $a>0, b>0, c>.0, d>0$.

• @D.Singh look at my edit. – mookid Nov 13 '14 at 21:34
• i change every $a$ to $0.4c$. This is a substitution. – mookid Nov 13 '14 at 22:00
• Of course. And now it adds up perfectly! Thanks a lot. – DSF Nov 13 '14 at 22:02
• well, $a+b+c+d = 2.36 c$ so $c = 1/2.36$ – mookid Nov 13 '14 at 23:17

Another thing you can do, which makes the procedure completely general - in addition to the four equations you get from

$$v(P-I) = 0$$

Where $P$ is your transition probabilities matrix above, you also need -

$$v_1 + v_2 +v_3 + \dots = 1$$ So, just add a column to your $P-I$ matrix with all ones. This new matrix has n rows and n+1 columns. Let's call it $Q$. And your right hand side vector becomes $b = [0, 0, \dots, 1]$, meaning repeat zeros until n and then add a 1. Now you get -

$$vQ = b$$ Post multiply both sides by $Q^T$

$$vQQ^T = bQ^T$$

Which leads to -

$$v = bQ^T(QQ^T)^{-1}$$

Here, $v$ is a column vector. Most linear algebra solvers expect the inverse matrix first multiplied by the vector, which is not the form of the result above. To convert it into this format (so it can be easily passed to these solvers), simply take the transpose.

$$v^T = (Q^T Q)^{-1}(Qb^T)$$

However, note that $Qb^T$ is simply a vector of ones (since all but the last column of $b$ is zero and the last row of $Q^T$ is all ones). So, no need to explicitly calculate $Qb^T$ it is simply a vector of ones. Here is a python method that does this for you:

import numpy as np

def steady_state_prop(
p=np.matrix([
[0,.2,.4,.4],
[.2,0,.4,.4],
[.2,.3,0,.5],
[.3,.4,.3,0]
])):
dim = p.shape[0]
q = (p-np.eye(dim))
ones = np.ones(dim)
q = np.c_[q,ones]
QTQ = np.dot(q, q.T)
bQT = np.ones(dim)
return np.linalg.solve(QTQ,bQT)


The equation you want to solve is $[a\ b\ c\ d](V-I)=0$, where $I$ is the identity matrix which has $1$'s on the diagonal and $0$'s elsewhere. To solve the equation, use Gaussian elimination.

$$V-I=\left[ { \begin{array}{ccc} -1 & 0.4 & 0 & 0.6\\ 0 & -0.5 & 0.5 & 0 \\ 0.4 & 0 & -0.8 & 0.4\\ 0 & 0 & 1 & -1 \end{array} }\right]$$

• So I have to apply Gaussian Elimination to the matrix V-I, that you have written out? – DSF Nov 13 '14 at 21:23
• Yes. Once you do that you will get a matrix in row echelon form, which you then multiply by $[a\ b\ c\ d]$ to see what the relationship between $a$, $b$, $c$, and $d$ has to be. Then you have to make it so that $a+b+c+d=1$ and you'll get your probability distribution. – Matt Samuel Nov 13 '14 at 21:26