1、维特比算法(Viterbi Algorithm)
维特比算法是一种动态规划算法,用于在给定观察序列的情况下,找到最可能的状态序列(即解码问题)。简而言之,它解决的是给定观测序列和模型参数时,确定最有可能产生这些观测的隐藏状态序列的问题。维特比算法广泛应用于语音识别、生物信息学(如DNA序列分析)和其他领域,其中需要从一系列观测数据中推断出最可能的隐含状态序列。
算法通过构建一个路径概率图(通常是一个表格),利用动态规划逐步计算和记录每个状态的最大概率路径。对于每个时间点,它记录下达到该状态的最大概率以及达到该概率的路径。最终,算法回溯这些记录,找到概率最高的路径作为最可能的状态序列。
import numpy as np # 定义模型参数 A = np.array([[0.7, 0.3], [0.4, 0.6]]) # 状态转移概率矩阵 B = np.array([[0.1, 0.9], [0.8, 0.2]]) # 观察概率矩阵 pi = np.array([0.5, 0.5]) # 初始状态概率分布 observations = [0, 1, 0] # 观察序列,0: 携带伞, 1: 不携带伞 # 实现维特比算法 def viterbi(A, B, pi, observations): N = A.shape[0] # 状态数量 T = len(observations) # 观察序列长度 V = np.zeros((N, T)) # 存储每一步的最大概率 path = np.zeros((N, T), dtype=int) # 存储路径 # 初始化 V[:, 0] = pi * B[:, observations[0]] # 动态规划 for t in range(1, T): for n in range(N): seq_probs = V[:, t-1] * A[:, n] * B[n, observations[t]] V[n, t] = np.max(seq_probs) path[n, t] = np.argmax(seq_probs) # 回溯找到最优路径 best_path = np.zeros(T, dtype=int) best_path[-1] = np.argmax(V[:, -1]) for t in range(T-2, -1, -1): best_path[t] = path[best_path[t+1], t+1] return best_path, np.max(V[:, -1]) # 应用维特比算法 best_path, best_prob = viterbi(A, B, pi, observations) # 转换状态序列为文字 state_names = ['晴天', '雨天'] best_path_names = [state_names[state] for state in best_path] print(best_path_names, best_prob)
2、鲍姆-韦尔奇算法(Baum-Welch Algorithm)
鲍姆-韦尔奇算法是一种特殊的期望最大化(Expectation-Maximization, EM)算法,用于未知参数的统计模型中参数的估计,特别是在HMM中。它主要用于学习给定观察序列时HMM的参数(即训练HMM)。鲍姆-韦尔奇算法被用于从观测数据中学习HMM的模型参数,例如状态转移概率、观测概率和初始状态概率。它在语言模型、遗传序列分析等领域有广泛应用。
算法交替执行两个步骤,期望步骤(E步)和最大化步骤(M步)。在E步,算法使用当前的模型参数估计隐藏状态的概率分布;在M步,算法更新模型参数以最大化在E步中计算得到的期望对数似然。这个过程反复迭代,直到模型参数收敛。
import numpy as np # 随机初始化模型参数 np.random.seed(42) # 确保示例的可重复性 A_init = np.random.rand(2, 2) B_init = np.random.rand(2, 2) pi_init = np.random.rand(2) # 归一化 A_init /= A_init.sum(axis=1)[:, np.newaxis] B_init /= B_init.sum(axis=1)[:, np.newaxis] pi_init /= pi_init.sum() # 假设的观察序列(0: X, 1: Y) observations = np.random.choice(2, 10) # 鲍姆-韦尔奇算法(Baum-Welch Algorithm) def baum_welch(observations, A_init, B_init, pi_init, iterations=100): N = A_init.shape[0] # 状态数量 M = B_init.shape[1] # 观察符号数量 T = len(observations) # 观察序列长度 for _ in range(iterations): # E步骤:使用当前参数计算前向概率和后向概率 alpha = np.zeros((N, T)) beta = np.zeros((N, T)) alpha[:, 0] = pi_init * B_init[:, observations[0]] beta[:, -1] = 1 # 计算alpha for t in range(1, T): for j in range(N): alpha[j, t] = np.sum(alpha[:, t-1] * A_init[:, j]) * B_init[j, observations[t]] # 计算beta for t in range(T-2, -1, -1): for i in range(N): beta[i, t] = np.sum(A_init[i, :] * B_init[:, observations[t+1]] * beta[:, t+1]) # 计算gamma和xi(辅助变量) gamma = np.zeros((N, T)) xi = np.zeros((N, N, T-1)) for t in range(T-1): denom = np.sum([alpha[i, t] * beta[i, t] for i in range(N)]) for i in range(N): gamma[i, t] = (alpha[i, t] * beta[i, t]) / denom for j in range(N): xi[i, j, t] = (alpha[i, t] * A_init[i, j] * B_init[j, observations[t+1]] * beta[j, t+1]) / denom gamma[:, -1] = alpha[:, -1] * beta[:, -1] / np.sum(alpha[:, -1] * beta[:, -1]) # M步骤:更新参数 pi_init = gamma[:, 0] A_init = np.sum(xi, 2) / np.sum(gamma[:, :-1], axis=1)[:, np.newaxis] for k in range(M): mask = observations == k B_init[:, k] = np.sum(gamma[:, mask], axis=1) / np.sum(gamma, axis=1) return A_init, B_init, pi_init # 应用鲍姆-韦尔奇算法 A_final, B_final, pi_final = baum_welch(observations, A_init, B_init, pi_init) print(A_final, B_final, pi_final)
3、应用示例
Python中应用维特比算法和鲍姆-韦尔奇算法通常涉及使用隐马尔可夫模型(HMM)。这两种算法在机器学习领域中有广泛应用,例如在语音识别、自然语言处理、生物信息学等领域。
1)维特比算法应用
如我们有一个简单的天气模型,状态包括晴朗和雨天,观测包括人们穿的是T恤、夹克或雨衣。我们想要根据观测序列来预测最可能的天气序列。
import numpy as np # 状态集合 states = ['晴朗', '雨天'] # 观测序列 observations = ['夹克', 'T恤', '雨衣'] # 观测集合 obs_dict = {'T恤': 0, '夹克': 1, '雨衣': 2} # 状态转移概率 trans_prob = {'晴朗': {'晴朗': 0.7, '雨天': 0.3}, '雨天': {'晴朗': 0.4, '雨天': 0.6}} # 观测概率 emit_prob = {'晴朗': {'T恤': 0.6, '夹克': 0.3, '雨衣': 0.1}, '雨天': {'T恤': 0.1, '夹克': 0.4, '雨衣': 0.5}} # 初始状态概率 start_prob = {'晴朗': 0.5, '雨天': 0.5} def viterbi(obs, states, start_p, trans_p, emit_p): V = [{}] for st in states: V[0][st] = {"prob": start_p[st] * emit_p[st][obs[0]], "prev": None} for t in range(1, len(obs)): V.append({}) for st in states: max_tr_prob = max(V[t-1][prev_st]["prob"]*trans_p[prev_st][st] for prev_st in states) for prev_st in states: if V[t-1][prev_st]["prob"] * trans_p[prev_st][st] == max_tr_prob: max_prob = max_tr_prob * emit_p[st][obs[t]] V[t][st] = {"prob": max_prob, "prev": prev_st} break opt = [] max_prob = max(value["prob"] for value in V[-1].values()) previous = None for st, data in V[-1].items(): if data["prob"] == max_prob: opt.append(st) previous = st break for t in range(len(V) - 2, -1, -1): opt.insert(0, V[t + 1][previous]["prev"]) previous = V[t + 1][previous]["prev"] print('最可能的天气序列是: ' + ' '.join(opt)) print('序列概率: {:.5f}'.format(max_prob)) # 应用Viterbi算法 viterbi(observations, states, start_prob, trans_prob, emit_prob)
2)鲍姆-韦尔奇算法应用
如有一个系统,它在两种状态下运行(例如“高效”和“低效”),并且能够观察到的输出是某种性能指标(例如“高”、“中”、“低”)。我们的目标是基于观察到的性能指标序列来学习状态转移概率、观测概率和初始状态概率。
import numpy as np # 定义状态空间 states = ["高效", "低效"] # 定义观测空间 observations = ["高", "中", "低"] # 定义初始状态概率 startprob = np.array([0.6, 0.4]) # 定义状态转移概率矩阵 transmat = np.array([[0.7, 0.3], [0.2, 0.8]]) # 定义观测概率矩阵 emissionprob = np.array([[0.4, 0.3, 0.3], [0.2, 0.5, 0.3]]) # 观测序列 observed_sequence = ["高", "中", "低", "高", "中"] # 使用鲍姆-韦尔奇算法训练模型 def baum_welch(observations, startprob, transmat, emissionprob): n_states = len(states) n_observations = len(observations) # 初始化 alpha = np.zeros((n_states, n_observations)) beta = np.zeros((n_states, n_observations)) gamma = np.zeros((n_states, n_observations)) xi = np.zeros((n_states, n_states, n_observations - 1)) # 前向算法 alpha[:, 0] = startprob * emissionprob[:, observations[0]] for t in range(1, n_observations): for i in range(n_states): alpha[i, t] = (alpha[:, t-1] * transmat[i, :]) * emissionprob[i, observations[t]] # 后向算法 beta[:, n_observations-1] = np.ones(n_states) for t in range(n_observations-2, -1, -1): for i in range(n_states): beta[i, t] = (transmat[:, i] * emissionprob[:, observations[t+1]] * beta[:, t+1]).sum() # 计算gamma和xi for t in range(n_observations): for i in range(n_states): gamma[i, t] = alpha[i, t] * beta[i, t] / alpha[:, t].sum() for j in range(n_states): xi[i, j, t] = alpha[i, t] * transmat[i, j] * emissionprob[j, observations[t+1]] * beta[j, t+1] / alpha[:, t].sum() # 更新参数 startprob = gamma[:, 0] transmat = xi.sum(axis=2) / gamma[:, :-1].sum(axis=1)[:, None] emissionprob = (gamma[:, None, :] * observations).sum(axis=2) / gamma.sum(axis=1)[:, None] return startprob, transmat, emissionprob # 训练模型 startprob, transmat, emissionprob = baum_welch(observed_sequence, startprob, transmat, emissionprob) # 打印结果 print("初始状态概率:", startprob) print("状态转移概率矩阵:", transmat) print("观测概率矩阵:", emissionprob)