computing-spd.py

# -*- coding: utf-8 -*- 
# <nbformat>3.0</nbformat> 
 
# <codecell> 
 
""" 
__author__ = Aziz Altowayan 
__date__   = April, 09 2014 
Title:  
    Self-Organizing Data Structures with MTF and Transposition heuristics. 
Description: 
    This listing uses the transition matrices (with case n=3) of Move-To-Front rule and Transposition rule to: 
        1. Compute the Stationary Probability Distributions SPD Prob(π). 
        2. Compute the expected cost of both heuristics ec(A). 
 
# Here n = 3 is demonstrated as a file of 3 records: 
    R1 = A, R2 = B, and R3 = C with respective probabilities of a, b, and c. 
    And we have n! permutations thus 6x6 matrix. 
 
# Computing the expected cost of a heuristic(A), as: 
ec(A) = sum(Prob(π)) * sum(cost(π)), where: 
    # Prob(π): πj = sum(πi * pij)   "SPD, elements of the eigenvector of the matrix". 
    # cost(π) = sum(i * Pπ(i))      "Avg. number of comparisons required to search a record". 
""" 
 
import numpy as np 
from numpy.linalg import eig 
 
# Transition probabilities values "searching a letter". 
a, b, c = 0.6, 0.3, 0.1 
 
# Move-To-Front transition matrix 
mtf = np.matrix([ 
        [a, b, c, 0, 0, 0], 
        [a, b, 0, c, 0, 0], 
        [0, 0, c, 0, a, b], 
        [0, 0, 0, c, a, b], 
        [0, b, c, 0, a, 0], 
        [a, 0, 0, c, 0, b]]) 
 
# ec(Move-To-Front) 
S, U = eig(mtf.T)   # S: eigenvalues, U: the normalized (unit "length") eigenvectors (i.e. col v[:,i] is eigenvector of w[i]) 
stationary = np.array(U[:,np.where(np.abs(S-1.) < 1e-8)[0][0]].flat) 
stationary = stationary / np.sum(stationary) 
 
# cost(π) = sum(Pπ(i) * i): 
print 'Move-To-Front matrix:\n', mtf, '\n\ncost(π):' 
for row in range(len(mtf)): 
    print ' ',sum( [mtf[row,[i]] * (i + 1) for i in range(len(mtf))] ), 
 
# Prob(π) SPD 
print '\nProb(π):' 
for prob in stationary: 
    print ' ', prob, 
 
# Transposition Rule transition matrix 
transpos = np.matrix([ 
        [a, b, c, 0, 0, 0], 
        [a, b, 0, c, 0, 0], 
        [b, 0, a, 0, c, 0], 
        [0, a, 0, b, 0, c], 
        [0, 0, a, 0, c, b], 
        [0, 0, 0, b, a, c]]) 
 
# ec(Transposition) 
S, U = eig(transpos.T) 
stationary = np.array(U[:,np.where(np.abs(S-1.) < 1e-8)[0][0]].flat) 
stationary = stationary / np.sum(stationary) 
 
# cost(π) = sum(Pπ(i) * i): 
print '\n\nTransposition matrix:\n', transpos, '\n\ncost(π):' 
for row in range(len(transpos)): 
    print ' ', sum( [transpos[row,[i]] * (i + 1) for i in range(len(transpos))] ), 
 
# Prob(π) SPD 
print '\nProb(π):' 
for prob in stationary: 
    print ' ', prob, 
 
# <codecell>