隐马尔可夫模型之python实战

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
Hidden_states = ("box 1", "box 2", "box 3") # 隐状态集合
Observations_states = ('red', 'white', 'red') # 观测状态集合

Start_probability = {'box 1': 0.2, 'box 2': 0.4, 'box 3': 0.4} # 初始状态

Hidden_transition_probability = { # 隐马尔可夫链中从当前盒子转移到其他盒子的概率
    'box 1' : {'box 1': 0.5, 'box 2': 0.2, 'box 3': 0.3},
    'box 2' : {'box 1': 0.3, 'box 2': 0.5, 'box 3': 0.2},
    'box 3' : {'box 1': 0.2, 'box 2': 0.3, 'box 3': 0.5},
   }

Hidden_observations_probability = {    # 原来叫emission_probability。这里表示每个盒子里面红白球的数量
    'box 1' : {'red': 0.5, 'white': 0.5},
    'box 2' : {'red': 0.4, 'white': 0.6},
    'box 3' : {'red': 0.7, 'white': 0.3},
   }

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
def viterbi(obs, states, start_p, trans_p, h2o_p): # Viterbi算法
    '''
    input:
        obs:观测状态集合O={红，白，红}
        states：隐状态集合Q={盒子1，盒子2，盒子3}
        start_p：初始状态Π=(0.2,0.4,0.4)T
        trans_p：转移状态矩阵A=【0.5 0.3 0.2
                              0.2 0.5 0.3
                              0.3 0.2 0.5】
        h2o_p：观测状态概率矩阵为：B=【0.5 0.4 
                                   0.7 0.5 
                                   0.6 0.3】
    output:
        prob:最优概率
        path[state]：最优路径
    '''
    V = [{}]
    path = {}

    # Initialize base cases (t == 0)
    for y in states:
        #t=1时刻，状态为i观测o1为红的概率，即δ1(i)
        V[0][y] = start_p[y] * h2o_p[y][obs[0]]
        # 初始状态，由start的概率，对应乘上发射概率，即由隐状态到观测状态的可能性
        path[y] = [y]

    # 开始遍历，t=2和t=3时刻
    for t in range(1,len(obs)):
        V.append({})
        newpath = {}

        for y in states:
            # 对于每个箱子，计算其由前一个的各个状态，到现在箱子的概率大小，取最大值。即求出最有可能到达现在箱子的路径
            # 前一个箱子转移到现在箱子的每个状态对应的路径大小都计算了，取最大值。作为V[t][y]的值，并更新路径
            (prob, state) = max([(V[t-1][y0] * trans_p[y0][y] * h2o_p[y][obs[t]], y0) for y0 in states])
            V[t][y] = prob
            newpath[y] = path[state] + [y]

        # Don't need to remember the old paths
        path = newpath

    print_dptable(V)
    (prob, state) = max([(V[len(obs) - 1][y], y) for y in states])
    return (prob, path[state])

1
2
3
4
5
6
viterbi(Observations_states,
                   Hidden_states,
                   Start_probability,
                   Hidden_transition_probability,
                   Hidden_observations_probability)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
states = ["box 1", "box 2", "box3"]
n_states = len(states)

observations = ["red", "white"]
n_observations = len(observations)

start_probability = np.array([0.2, 0.4, 0.4])

transition_probability = np.array([
  [0.5, 0.2, 0.3],
  [0.3, 0.5, 0.2],
  [0.2, 0.3, 0.5]
])

emission_probability = np.array([
  [0.5, 0.5],
  [0.4, 0.6],
  [0.7, 0.3]
])

1
2
3
4
5
6
7
8
model = hmm.MultinomialHMM(n_components=n_states)
model.startprob_=start_probability
model.transmat_=transition_probability
model.emissionprob_=emission_probability

seen = np.array([[0,1,0]]).T
logprob, box = model.decode(seen, algorithm="viterbi")

1
2
3
print("The ball picked:【", ", ".join(map(lambda x: observations[x[0]], seen)),"】")
print("The hidden box【", ", ".join(map(lambda x: states[x], box)),"】")