Examples



In [1]:

    
from stochmatrix import StochMatrix

Simple example



In [2]:

    
import numpy as np



In [3]:

    
# Taken from www.math.wustl.edu/~feres/Math450Lect04.pdf
P = np.zeros((10, 10))
P[[0, 0], [0, 2]] = 1./2
P[1, 1], P[1, 6] = 1./3, 2./3
P[2, 0] = 1
P[3, 4] = 1
P[[4, 4, 4], [3, 4, 8]] = 1./3
P[5, 5] = 1
P[6, 6], P[6, 8] = 1./4, 3./4
P[[7, 7, 7, 7], [2, 3, 7, 9]] = 1./4
P[8, 1] = 1
P[[9, 9, 9], [1, 4, 9]] = 1./3



In [4]:

    
np.get_printoptions()









    Out[4]:





{'edgeitems': 3,
 'formatter': None,
 'infstr': 'inf',
 'linewidth': 75,
 'nanstr': 'nan',
 'precision': 8,
 'suppress': False,
 'threshold': 1000}



In [5]:

    
np.set_printoptions(precision=3)



In [6]:

    
print P









    



[[ 0.5    0.     0.5    0.     0.     0.     0.     0.     0.     0.   ]
 [ 0.     0.333  0.     0.     0.     0.     0.667  0.     0.     0.   ]
 [ 1.     0.     0.     0.     0.     0.     0.     0.     0.     0.   ]
 [ 0.     0.     0.     0.     1.     0.     0.     0.     0.     0.   ]
 [ 0.     0.     0.     0.333  0.333  0.     0.     0.     0.333  0.   ]
 [ 0.     0.     0.     0.     0.     1.     0.     0.     0.     0.   ]
 [ 0.     0.     0.     0.     0.     0.     0.25   0.     0.75   0.   ]
 [ 0.     0.     0.25   0.25   0.     0.     0.     0.25   0.     0.25 ]
 [ 0.     1.     0.     0.     0.     0.     0.     0.     0.     0.   ]
 [ 0.     0.333  0.     0.     0.333  0.     0.     0.     0.     0.333]]

StochMatrix class



In [7]:

    
P = StochMatrix(P)

Whether P is irreducible



In [8]:

    
P.is_irreducible









    Out[8]:





False

Communication classes (strongly connected components) of P



In [9]:

    
P.comm_classes()









    Out[9]:





[[0, 2], [1, 6, 8], [3, 4], [5], [9], [7]]

Recurrent classes of P



In [10]:

    
P.rec_classes()









    Out[10]:





[[0, 2], [1, 6, 8], [5]]

Recurrent states of P



In [11]:

    
rec_states = [i for rec_class in P.rec_classes() for i in rec_class]



In [12]:

    
print rec_states









    



[0, 2, 1, 6, 8, 5]

Transient states of P



In [13]:

    
trans_states = [i for i in range(10) if i not in rec_states]



In [14]:

    
print trans_states









    



[3, 4, 7, 9]

Canonical form of P



In [15]:

    
permutation = rec_states + trans_states



In [16]:

    
permutation









    Out[16]:





[0, 2, 1, 6, 8, 5, 3, 4, 7, 9]



In [17]:

    
print P[permutation, :][:, permutation]









    



[[ 0.5    0.5    0.     0.     0.     0.     0.     0.     0.     0.   ]
 [ 1.     0.     0.     0.     0.     0.     0.     0.     0.     0.   ]
 [ 0.     0.     0.333  0.667  0.     0.     0.     0.     0.     0.   ]
 [ 0.     0.     0.     0.25   0.75   0.     0.     0.     0.     0.   ]
 [ 0.     0.     1.     0.     0.     0.     0.     0.     0.     0.   ]
 [ 0.     0.     0.     0.     0.     1.     0.     0.     0.     0.   ]
 [ 0.     0.     0.     0.     0.     0.     0.     1.     0.     0.   ]
 [ 0.     0.     0.     0.     0.333  0.     0.333  0.333  0.     0.   ]
 [ 0.     0.25   0.     0.     0.     0.     0.25   0.     0.25   0.25 ]
 [ 0.     0.     0.333  0.     0.     0.     0.     0.333  0.     0.333]]

Decomposing P to irreducible submatrices



In [18]:

    
for i, rec_class in enumerate(P.rec_classes()):
    print 'P_{0} ='.format(i)
    print P[rec_class, :][:, rec_class]









    



P_0 =
[[ 0.5  0.5]
 [ 1.   0. ]]
P_1 =
[[ 0.333  0.667  0.   ]
 [ 0.     0.25   0.75 ]
 [ 1.     0.     0.   ]]
P_2 =
[[ 1.]]



In [18]: