Markov! Where are my umbrellas?

Today I learned about Markov Chains. The toy problem that we tackled was so challenging and figuring out the solution was so satisfying, I felt the urge to share immediately (which expedited my inaugural blog post).

Problem: I have four umbrellas, some at home, some in the office. I keep moving between home and office. I take an umbrella with me only if it rains. If it does not rain, I leave the umbrella behind. If there are no umbrellas where I am, and it's raining, I get wet. The probability that it rains on any given day is \(0.6\). What is the probability that I get wet?

What?! First of all, if it rains \(60\%\) of the time, I'm moving.

Ok so the (first) trick here is knowing that home/office doesn't matter. It only matters how many umbrellas I have with me in my current location. With that in mind, the random variable, say \(U\), can take the values \(\lbrace 0, 1, 2, 3, 4 \rbrace\).

For \(U = 0\), I will by default transition into a state with \(U = 4\) with \(100\%\) probability. This means I have zero umbrellas with me, and I will leave empty-handed. At my destination are four umbrellas.

For \(U = 1\), there is a \(60\%\) probability that I bring the umbrella with me, and transition to \(U = 4\). There is also a \(40\%\) probability that I leave the umbrella behind, and transition to \(U = 3\).

For \(U = 2\), there is a \(60\%\) probability that I transition to \(U = 3\), and \(40\%\) probability that I stay at \(U = 2\). Similar arguments can be made for \(U = 3\) and \(U = 4\).

With those transitions in mind, we can draw a Markov Chain as such:

The arrows and probabilities denote transitions and corresponding probabilities, respectively. How does this help? From here, we can populate a transition matrix, say \(M\), for Markov.

$$ M = \begin{bmatrix} 0.0 & 0.0 & 0.0 & 0.0 & 1.0 \\ 0.0 & 0.0 & 0.0 & 0.4 & 0.6 \\ 0.0 & 0.0 & 0.4 & 0.6 & 0.0 \\ 0.0 & 0.4 & 0.6 & 0.0 & 0.0 \\ 0.4 & 0.6 & 0.0 & 0.0 & 0.0 \end{bmatrix} $$

To understand this matrix, rows denote starting state and columns denote ending state. For example, row 2 column 3 (with zero indexing) means there is a \(60\%\) probability of transitioning from Node 2 to Node 3 (or \(U = 2\) to \(U = 3\)). Now, let us assign a state vector as our initial condition:

$$ s_0 = \begin{bmatrix} 0.0 \\ 0.0 \\ 0.0 \\ 1.0 \\ 0.0 \end{bmatrix} $$

The values in this vector can be arbitrary since we are seeking the long time steady state (provided they sum to \(100\%\)). This vector tells us the probability of being in each state. As we have it filled in, we gave \(U = 3\) a \(100\%\) probability (zero indexing still), i.e., we have three umbrellas with us currently. Now, to transition into the next state (i.e., march forward one time step), we just let the matrix do all the heavy lifting. By hitting the vector on the left with the transpose of \(M\), we get the new state condition:


$$ s_1 = M^T s_0 $$

We can repeat this process infinite times to reach the steady state condition. Note that in general, a matrix raised to the infinite power might not converge, particularly if its values are greater than \(1\). Fortunately for us, it converges. We can diagonalize \(M\) as such:


$$ M = PDP^{-1} $$

Here, \(P\) is a block matrix comprised of the right eigenvectors of \(M\). \(D\) is the diagonal matrix of \(M\), comprised of the eigenvalues along the diagonal. This allows us to multiply as many times as we like, e.g.:

$$ \begin{aligned} \lim_{n\rightarrow\infty}M^n & = \lim_{n\rightarrow\infty} \overbrace{PDP^{-1} \times PDP^{-1} \times ... \times PDP^{-1}}^{n~\text{times}} \\ & = \lim_{n\rightarrow\infty} \overbrace{PD(P^{-1}P)D(P^{-1}P)...DP^{-1}}^{n~\text{times}} \\ & = \lim_{n\rightarrow\infty} PD^nP^{-1} \end{aligned} $$

Note that the internal \(P\) and \(P^{-1}\) matrices cancel out. Ultimately, with some plugging and chugging, we get:


$$ \begin{aligned} s_{\infty} & = \lim_{n\rightarrow\infty} (M^T)^{n}s_o \\ & = \lim_{n\rightarrow\infty}(P^T)^{-1}D^nP^Ts_o \\ & = \begin{bmatrix} 0.09091 \\ 0.22727 \\ 0.22727 \\ 0.22727 \\ 0.22727 \end{bmatrix} \end{aligned} $$

The probability that I get wet will be the probability that it rains and \(U = 0\). In other words, the probability of \(U = 0\) times \(0.6\):


$$ \begin{aligned} p(\text{wet}) & = 0.09091 \times 0.6 \\ & = 0.05455 \end{aligned} $$

Update:
Let me be a little more explicit, as this is quite useful in interpreting the math. The \(P\) and \(D\) matrices for this case are:


$$ P = \begin{bmatrix} 0.597 & -0.583 & -0.776 & -0.447 & 0.594 \\ 0.480 & 0.148 & 0.511 & -0.447 & 0.183 \\ 0.128 & 0.414 & -0.320 & -0.447 & -0.545 \\ -0.286 & -0.619 & 0.175 & -0.447 & -0.331 \\ -0.561 & 0.290 & -0.056 & -0.447 & 0.455 \end{bmatrix} $$

and

$$ D = \begin{bmatrix} -0.939 & 0.0 & 0.0 & 0.0 & 0.0 \\ 0.0 & -0.497 & 0.0 & 0.0 & 0.0 \\ 0.0 & 0.0 & 0.072 & 0.0 & 0.0 \\ 0.0 & 0.0 & 0.0 & 1.0 & 0.0 \\ 0.0 & 0.0 & 0.0 & 0.0 & 0.765 \end{bmatrix} $$

Raising \(D\) to the \(n^{\text{th}}\) power is the same as raising each eigenvalue to the \(n^{\text{th}}\) power. Since four of the five values are fractions, they quickly become zero. The eigenvalue with value unity remains the same. Then, as \(M\) is taken to the \(n^{\text{th}}\) power, we obtain


$$ \lim_{n\rightarrow\infty}(M^T)^n = \begin{bmatrix} 0.091 & 0.091 & 0.091 & 0.091 & 0.091 \\ 0.227 & 0.227 & 0.227 & 0.227 & 0.227 \\ 0.227 & 0.227 & 0.227 & 0.227 & 0.227 \\ 0.227 & 0.227 & 0.227 & 0.227 & 0.227 \\ 0.227 & 0.227 & 0.227 & 0.227 & 0.227 \end{bmatrix} $$

Notice that this is precisely the steady state vector repeated as columns. This ensures that we will recover the same steady state vector every time, as long as the initial vector sums to 100%.

However, it has also been pointed out to me that the end can be solved more analytically. So let's take a couple steps back. Rather than multiplying an initial state vector \(n\)-times, all we have to do is multiply the final state vector once. Since it is the final state, it will remain unchanged under the transformation. Specifically:


$$ \begin{aligned} s_{\infty+1} = s_{\infty} & = M^T s_{\infty} \\ 0 & = M^T s_{\infty} - s_{\infty} \\ 0 & = (M^T - I) s_{\infty} \end{aligned} $$

It turns out \((M^T-I)\) only has rank \(4\). Another equation can be added to one of the rows, conditioning that the elements of the state vector sum to \(1\), i.e.,


$$ 1 = \pi_0 + \pi_1 + \pi_2 + \pi_3 + \pi_4 $$

Here, the \(\pi\)'s denote the steady state probabilities of each node. Explicitly, we have:


$$ \begin{bmatrix} 0 \\ 0 \\ 0 \\ 0 \\ 1 \end{bmatrix} = \begin{bmatrix} -1.0 & 0.0 & 0.0 & 0.0 & 0.4 \\ 0.0 & -1.0 & 0.0 & 0.4 & 0.6 \\ 0.0 & 0.0 & -0.6 & 0.6 & 0.0 \\ 0.0 & 0.4 & 0.6 & -1.0 & 0.0 \\ 2.0 & 1.6 & 1.0 & 1.0 & 0.0 \end{bmatrix} \begin{bmatrix} \pi_0 \\ \pi_1 \\ \pi_2 \\ \pi_3 \\ \pi_4 \end{bmatrix} $$

Finally, the modified matrix equation can be solved with matrix multiplication to get:


$$ s_{\infty} = \begin{bmatrix} \pi_0 \\ \pi_1 \\ \pi_2 \\ \pi_3 \\ \pi_4 \end{bmatrix} = \begin{bmatrix} 0.09091 \\ 0.22727 \\ 0.22727 \\ 0.22727 \\ 0.22727 \end{bmatrix} $$.