Markov Chain 分析（一）

情境說明：

• 中央監控室有三個監控螢幕，分別顯示不同環境變化狀況。
• 當環境信號出現在各監控螢幕時，會持續顯示後續訊息在同一監控螢幕的機會有 50%。
• 當環境信號先出現在一號或二號監控螢幕時，緊接著顯示後續訊息在三號監控螢幕的機會有 50%。
• 當環境信號先出現在三號監控螢幕時，緊接著顯示後續訊息在一號或二號監控螢幕的機會有 25%。

探討：

1. 分別探討第一條訊息出現在各個監控螢幕時，緊接著出現 10 次訊息在各個監控螢幕的機會如何？
2. 是否具備穩態性質？假使有，各個監控螢幕出現訊息的穩態機率為何？需要幾次才能達到穩態？

圖解：

主程式碼：

1 # coding=utf8 2 3 from MarkovChain import MarkovChain 4 import numpy 5 6 7 if __name__ == '__main__': 8 worker = MarkovChain() 9 10 data = None 11 cases = [True, False] 12 13 if cases[0]: 14 # 範例 1: 15 data = numpy.array([ 16 [0.5, 0, 0.5], 17 [0, 0.5, 0.5], 18 [0.25, 0.25, 0.5]], 19 dtype=numpy.float64) 20 else: 21 # 範例 2: 22 data = numpy.array([ 23 [0.5, 0.25, 0.25], 24 [0, 1.0/3.0, 2.0/3.0], 25 [0.5, 0.5, 0]], 26 dtype=numpy.float64) 27 28 taskControl = [True, True] 29 30 if taskControl[0]: 31 # 直接計算穩態，有可能無解，表示無穩態。 32 33 result = worker.stateSolve(data) 34 worker.outputSateSolution(result) 35 36 if taskControl[1]: 37 # 在嘗試 n 次內，重複計算求穩態機率。 38 39 initData = numpy.array([[1, 0, 0], [0, 1, 0], [0, 0, 1]], dtype=numpy.float64) 40 41 for i in range(len(initData)): 42 print initData[i] 43 44 timesMax, result = worker.stateMatrix(initData[i], data, timesMax=30) 45 worker.outputSateMatrix(result, timesMax)

程式館碼：

1 # coding=utf8 2 3 import numpy 4 import math 5 6 7 class MarkovChain(): 8 def __init__(self): 9 pass 10 11 def dotMatrix(self, data1, data2): 12 result = numpy.dot(data1, data2) # 矩陣相乘 13 return result 14 15 def outputMatrix(self, data, timesMax): 16 print 'Time = %d' % (timesMax) 17 for row in data: 18 for col in row: 19 print '\t%.4f\t' % (col), 20 21 print 22 23 print 24 25 def stateMatrix(self, initData, data, timesMax=1): 26 result = data 27 28 for i in range(1, timesMax - 1): 29 prevResult = result 30 result = self.dotMatrix(data, result) 31 32 if numpy.allclose(result, prevResult, 0.0001, 0.0001): # 假使與前次計算接近 33 break 34 35 self.outputMatrix(result, i) 36 37 result = self.dotMatrix(initData, result) 38 return (i, result) 39 40 def stateSolve(self, data): 41 p = numpy.transpose(data) # 轉置矩陣 42 p2 = p - numpy.diag(numpy.ones(numpy.shape(p)[0])) # 矩陣相減，對角為 1 之矩陣 43 p2[numpy.shape(p)[0] - 1:] = [1] * numpy.shape(p)[0] # 替換矩陣末列為 0 之值為 1 44 p3 = numpy.zeros(numpy.shape(p)[0]) 45 p3[-1] = 1 46 result = numpy.linalg.solve(p2, p3) # 解方程式，各解即為穩態機率，倒數即為需要經過之次數 47 48 return result 49 50 def outputSateMatrix(self, data, timesMax): 51 print 'Max = %d' % (timesMax) 52 for row in data: 53 print '\t%.4f\t' % (row), 54 55 print '\n\n' 56 57 def outputSateSolution(self, data): 58 i = 1 59 print '%s\t%s\t%s' % ('Node', 'Probability', 'Round') 60 61 for e in data: 62 print '%4d\t%11.4f\t%5d' % (i, e, math.ceil(1/e)) 63 i += 1 64 65 print '\n\n'

程式輸出：

監控螢幕 穩態機率 輪轉次數 1 0.2500 4 2 0.2500 4 3 0.5000 2

分析結果：