Markov Chains in Python: a Simple Weather Model

Weather those days is so unstable and behaves strangely for the season. So, for a rainy Friday evening here is an interesting case for Python use.

Markov property is the memoryless property of a stochastic process, namely the property of future states to depend only upon the present state, not the past states. Markov chain is a process that exhibits Markov property.

Using Wikipedia example of a very simple weather model (https://en.wikipedia.org/wiki/Examples_of_Markov_chains#A_very_simple_weather_model), there is 90% probability for a sunny day tomorrow given a sunny day today and 10% probability for a rainy day tomorrow given sunny day today. Also, there is 50% probability for a sunny day tomorrow given a rainy day today and 50% probability for a rainy day tomorrow given a rainy day today. Briefly:

To solve for the weather tomorrow we need two numbers: (1) transition probabilities (as described above) and (2) initial state - what is the weather today (for instance, if today is sunny - we use (1,0), if it is rainy (0,1). Finally, we use matrix multiplication (bear in mind the difference between matrix multiplication and element-wise multiplication) to solve the problem.  The resulted state depends on previous state (changing initial state) and transition probabilities (the transition probabilities matrix remains unchanged).

Next day probability for a sunny day is 90% and probability for a rainy day is 10%. After certain number of predictions, eventually the steady state is reached and the weather prediction for the next day is independent of the initial weather. In our example at some point we reach a probability for a sunny day of 83%. and it barely changes (you can see the result below).

Here is a simple Python implementation of the weather model for a 20-day period: 

import numpy as np

transition=np.array([[0.9, 0.1],[0.5, 0.5]])

initial=np.array([1,0])

for i in range(20):

    n = 0

    counter = 0

    while counter <= n:

        counter+=1

        initial=np.dot(initial,transition)

    print (initial)

Result:

[ 0.9  0.1]

[ 0.86  0.14]

[ 0.844  0.156]

[ 0.8376  0.1624]

[ 0.83504  0.16496]

[ 0.834016  0.165984]

[ 0.8336064  0.1663936]

[ 0.83344256  0.16655744]

[ 0.83337702  0.16662298]

[ 0.83335081  0.16664919]

[ 0.83334032  0.16665968]

[ 0.83333613  0.16666387]

[ 0.83333445  0.16666555]

[ 0.83333378  0.16666622]

[ 0.83333351  0.16666649]

[ 0.8333334  0.1666666]

[ 0.83333336  0.16666664]

[ 0.83333334  0.16666666]

[ 0.83333334  0.16666666]

[ 0.83333334  0.16666666]