机器学习系列:隐马尔可夫模型 (3)
接上一篇机器学习序列隐马尔可夫模型(2)本篇主要讲解隐马尔可夫模型的解码问题.五、解码问题解码问题或预测问题是在HMM学习问题已经解决的情况下即模型参数已确定对于某个观测序列O求概率最大的那个状态序列。解码问题一般可以用两种算法处理即近似算法和Viterbi算法。1. 近似算法近似算法是通过比较t时刻的概率值获取最有可能出现的状态从而构成一个状态序列把它作为预测的结果。时刻t最有可能的状态从而得到预测的状态序列。近似算法计算过程简单但不能保证它是整体上的最优状态序列而且其预测的序列可能有实际不发生的部分即有可能出现转移概率为0的相邻状态。下面给出的Viterbi算法可以获得整体上的最优状态序列。2.Viterbi算法Viterbi算法维特比算法是通过动态规划方法在各个状态节点中求出一条概率最大的路径其中的每一条路径即对应着一个状态序列。首先将各个时刻上的所有可能的状态结果用节点表示于是可以绘制出如下图所示的网状图其中每个时刻隐藏状态有N中选择也就是有N个对应的节点。解码问题的目标是在给定观测序列的情况下要在这个网格图中找到一条最佳路径使得概率最大。问题详述如下问题 已知模型参数给定观测序列, 求最优隐藏状态序列, 即因为所以问题等价于:为用动态规划方法求解问题这里先对每个节点定义阶段目标函数即定义在时刻t状态为的所有单个路径中概率最大值:如上图所示这里的表示t时刻的第i个节点的目标函数。因为所以就是最优路径中最后时刻T处选择的节点记以下的任务是依次求出, 为此先推导阶段目标函数的递推过程。由的定义有于是得到递推式,在前面给出的网状图中若在最优路径中已确定t时刻的状态是 欲回溯上一时刻t-1处最优路径中的状态等效地考虑如下目标函数于是定义一个时刻t处状态为时其上一时刻(t-1)处的最优状态序号函数最后可以得到Viterbi算法:----------------------------------------------Viterbi算法-------------------------------------------------------------------输入模型模型 观测序列输出最优状态序列1. 初始化2. 对于 t 2, 3, ..., T 依次递推 3. 最终路径回溯1 初始2 对t T-1, T-2, ..., 1 依次递推4. 求出最优状态序列--------------------------------------------------------------------------------------------------------------------------------以下给出Viterbi算法的MATLAB实现function predI myHMMviterbi(O, PIest,Aest, Best) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % HMM模型之解码问题 % % myHMMviterbi: 使用维特比算法解码找出最可能的隐藏状态序列 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %输入 % O: 1*T, 1个长度是T的观测序列矩阵 % PIest: N*1,初始状态概率的估计 % Aest N*N, A的估计 % Best N*M, B的估计 %输出 % predI: 最优状态序列 %2026.6.11 MiaoZhh %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% [s,T]size(O); [numStates,numObserv]size(Best); deltazeros(T,numStates); psizeros(T,numStates); %标记概率最大的那个状态下标 %% 各节点最大概率计算 %初始值计算 for i1:numStates delta(1,i)(Best(i,O(1,1)).*PIest(i,1)); psi(1,i)0; %表示无前驱 end %递推 for t1:T-1 AtAest.*delta(t,:); for i1:numStates [Max,idx]max(At(:,i)); delta(t1,i)Best(i,O(1,t1))*Max; psi(t1,i)idx; end end %% 最优路径回溯 predIzeros(1,T); %初始 tT [Max,predI(:,T)]max(delta(T,:)); %反推 for tT:(-1):2 predI(t-1)psi(t,predI(t)); end end六、验证以下在MATLAB中通过一个具体的例子来验证前面给出的学习问题和解码问题的算法。直接给出主函数如下function main() %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %验证HMM模型之学习问题的主函数 %其中调用的函数HMMgenerate用来生成在给定参数下的一组观测序列 %调用的函数myHMMtrain用来训练HMM模型以获得参数的估计值PIest Aest, Best %2026.8.14 MiaoZhh %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear all clc % 定义状态转移矩阵,发射矩阵和初始状态概率向量 A [0.9 0.1; 0.05 0.95]; % 2状态转移概率 B [1/6,1/6,1/6,1/6,1/6,1/6; % 状态1发射概率均匀分布 7/12,1/12,1/12,1/12,1/12,1/12]; % 状态2发射概率偏向1 PI[0.9;0.1]; % 初始状态概率向量偏向2 % 生成观测序列和真实状态序列 T1000; %序列长度 [O, I] myHMMgenerate(T, PI, A, B); % 初始猜测矩阵需满足概率分布 A_GUESS [0.85 0.15; 0.1 0.9]; B_GUESS [0.17 0.16 0.17 0.16 0.17 0.17; % 状态1 0.6 0.08 0.08 0.08 0.08 0.08]; % 状态2 PI_GUESS[0.5;0.5]; % 执行Baum-Welch算法 [PIest Aest, Best] myHMMtrain(O, PI_GUESS, A_GUESS, B_GUESS); % 对比真实与估计参数 disp(Baum-Welch估计的初始状态概率向量:); disp(PIest); disp(Baum-Welch估计的转移矩阵:); disp(Aest); disp(Baum-Welch估计的发射矩阵:); disp(Best); %% 预测 %用Viterbi算法预测 [Os, Is] hmmgenerate(100,Aest, Best); predI myHMMviterbi(Os, PIest,Aest, Best); %计算准确率如果有真实状态对比 trueStatesIs; likelyStatespredI; accuracy sum(likelyStates trueStates)/ numel(trueStates); fprintf(解码准确率: %.2f%%\n, accuracy * 100); end运行效果如下其中的myHMMgenerate函数的功能是用来产生一组给定模型参数情况下定长的观测序列和状态序列。这里我们模仿了MATLAB工具箱中的函数hmmgenerate. 以下给出其实现部分function [O, I] myHMMgenerate(T, PI, A, B) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %myHMMgenerate: 函数生成 一个HMM模型的状态序列和观测序列 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %输入 % T: 1*1, 序列长度 % PI: N*1参数PI的初步状态概率向量 % A: N*N, 状态转移概率矩阵 % B: N*M观测概率矩阵(发射概率矩阵) %输出 % O: 1*T, 长度是T的观测序列向量 % I: 1*T, 长度是T的观测序列向量 %2026.8.14 MiaoZhh %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% O zeros(1,T); I zeros(1,T); [numStates,numEmissions] size(B); % create three random sequences, one for initial state, one for state changes, one for emission initstatechange rand(1); statechange rand(1,T); randvals rand(1,T); % calculate cumulative probabilities PIc cumsum(PI,1); Ac cumsum(A,2); Bc cumsum(B,2); % normalize these just in case they dont sum to 1. PIcPIc./repmat(PIc(end,:),numStates,1); Ac Ac./repmat(Ac(:,end),1,numStates); Bc Bc./repmat(Bc(:,end),1,numEmissions); %PI initstate 1; for innerState numStates-1:-1:1 if initstatechange PIc(innerState) initstate innerState 1; break; end end % Assume that we start in state initstate. currentstate initstate; % main loop for count 1:T % calculate state transition stateVal statechange(count); state 1; for innerState numStates-1:-1:1 if stateVal Ac(currentstate,innerState) state innerState 1; break; end end % calculate emission val randvals(count); emit 1; for inner numEmissions-1:-1:1 if val Bc(state,inner) emit inner 1; break; end end % add values and states to output O(count) emit; I(count) state; currentstate state; end end《机器学习序列隐马尔可夫模型 》部分到此结束谢谢浏览。
