Baum-Welch算法用来对隐马尔可夫模型的参数进行学习，Baum-Welch算法是EM算法的一种特例，属于非监督学习算法，下边第一部分程序的迭代次数为1的算法代码，其中调用了计算中所需要的两个变量Gamma和Xi的计算代码，一并在此给出，第二部分为迭代n次的代码，代码中已经给出了示例。
1，迭代次数为1的Baum-Welch算法代码
function [A_1,B_1,Pi_1] = BaumWelchAlgo(A,B,Pi,O)% 函数功能：实现迭代次数为1的Baum-Welch算法%% 参考文献：1，李航《统计学习方法》%2，Bilmes L,et al. A Gentle Tutorial of EM Algorithm and its Application to%Parameter Estimation for Gaussian Mixture and Hidden Markov Models%3，Rabiner L. A Tutorial on Hidden Markov Models and Selected%Applications in Speech Recognition%% 思路：% pi_i = Gamma_1(i);% a_ij = sum_t-from-1-to-(T-1)(Xi_t(i,j)) / sum_t-from-1-to-T(Gamma_t(i));% b_j(k) = sum_t-from-1-to-T_ot-=-vk(Gamma_t(i)) / sum_t-from-1-to-T(Gamma_t(j));% 其中：% sum_i-from-1-to-N(pi_i) = 1% sum_j-from-1-to-N(a_ij) = 1% sum_k-from-1-to-M(b_j(k)) = 1% % 输入：%第一次选取的模型参数A,B,Pi和观测序列O。% 输出：%经过一次迭代Baum-Welch算法之后得到的模型参数A_1,B_1,Pi_1%% 示例：% 输入% A =%0.5000 0.2000 0.3000%0.3000 0.5000 0.2000%0.2000 0.3000 0.5000% B =%0.5000 0.5000%0.4000 0.6000%0.7000 0.3000% Pi =%0.2000%0.4000%0.4000% O =%1%2%1%2% 输出：% >>[a,b,c] = BaumWelchAlgo(A,B,Pi,O)% % a =% 0.5151 0.2212 0.2637% 0.2987 0.5231 0.1783% 0.2210 0.3674 0.4116% b =%0.4325 0.5675%0.4257 0.5743%0.6389 0.3611% c =%0.1872%0.3241%0.4887A_size = size(A);B_size = size(B);O_size = size(O);N = A_size(1,1);%状态集个数M = A_size(1,2);Y = B_size(1,2);%B的列数，即观测集个数K = O_size(1,1);% ----调用GammaXiAlgo函数，计算Gamma矩阵和Xi矩阵----[Gamma,Xi] = GammaXiAlgo(A,B,Pi,O);% ---------------计算转移概率矩阵A_1---------------% 转移概率a_ij = 在观测O下由状态i转移到状态j的期望值 / 在观测O下由状态i转移的期望值% = expected number of transitions from state Si to state Sj /% excepted number of transitions from state SiA_1 = zeros();f = sum(Gamma,2);%计算Gamma的每一行之和，表示在观测时间序列的所有时刻都处于某一个状态的概率和for i = 1:Mfor j = 1:NA_1(i,j) = sum(sum(Xi(i,j,:))) / (f(i,1) - Gamma(i,K));endend% ---------------计算观测概率矩阵B_1---------------% 观测概率b_j(k) = 满足条件(ot = vk)的观测O下由状态i转移的期望值 / 观测O下由状态i转移的期望值%= excepted number of times in state j and observing symbol%vk / excepted number of times in state jB_1 = zeros();G = zeros();c = 0;for y = 1:Yfor j = 1:Nfor t = 1:Kif O(t,1) == y % 用于比较输出的观测值是否由对应的状态值输出，即条件ot = vk;c = c + 1;G(c,1) = Gamma(j,t);endendB_1(j,y) = sum(G) / f(j,1);c = 0;G = zeros();endend% ---------------计算初始概率pi_1---------------Pi_1 = zeros();for i = 1:MPi_1(i,1) = Gamma(i,1);end
2，变量Gamma和Xi的计算代码
function [Gamma,Xi] = GammaXiAlgo(A,B,Pi,O)% 函数功能：计算Baum-Welch算法公式中的gamma概率和xi概率%% 参考文献：李航《统计学习方法》%% 思路：% gamma_t(i) = P(it = qi|O,lambda) = P(it = qi,O|lambda) / P(O|lambda);% xi_t(i) = P(i_t = qi,i_t+1 = qj|O,lambda) = P(i_t = qi,i_t+1 = qj,O|lambda) / P(O|lambda)% gamma_t(i)表示在时刻t处于状态qi的概率% xi_t(i,j)表示在时刻t处于状态qi，且在时刻t+1处于状态qj的概率%其中t的范围是观测时间段长度，即序列O的长度K%qi的范围是状态集合中状态的个数，即转移矩阵A的维数N% % 输入：隐马尔可夫模型lambda(A,B,Pi)，观测序列O；%其中A为转移概率，B为观测概率（行数表示状态集合中的元的个数，列数表示观测集中的元的个数）；%Pi为初始概率，初始概率Pi和观测序列O须为列矩阵。% 输出：矩阵Gamma：存放不同时刻t处于不同状态qi的gamma概率值，为N*K阶% 矩阵Xi：存放不同时刻t处于不同状态qi且t+1处于不同状态qj的概率值，为N*N*K阶% 其中N*N维矩阵中，列对应t时刻，行对应t+1时刻，第K维表示的是时间序列% 示例：% y(:,:,2) =% 0.1474 0.0523 0.1157% 0.1268 0.1874 0.1106% 0.0464 0.0617 0.1518% 那么，y(:,:,2)中的2表示此时处在时刻t=2；% 0.1474表示在时刻2处于状态1且在时刻3处于状态1的概率；% 0.1268表示在时刻2处于状态2且在时刻3处于状态1的概率；% 0.0523表示在时刻2处于状态1且在时刻3处于状态2的概率；% 以此类推。%% 示例：% 输入：% A =%0.5000 0.2000 0.3000%0.3000 0.5000 0.2000%0.2000 0.3000 0.5000% B =%0.5000 0.5000%0.4000 0.6000%0.7000 0.3000% Pi =%0.2000%0.4000%0.4000% O2 =%1%2%1%2% 输出：% >>[x,y] = GammaXiAlgo(A,B,Pi,O2)% x =%0.1872 0.3154 0.3205 0.3508%0.3241 0.4248 0.3014 0.4192%0.4887 0.2599 0.3781 0.2301% y(:,:,1) =%0.1024 0.0462 0.0386%0.0983 0.1847 0.0412%0.1147 0.1939 0.1801% y(:,:,2) =%0.1474 0.0523 0.1157%0.1268 0.1874 0.1106%0.0464 0.0617 0.1518% y(:,:,3) =%0.1742 0.0836 0.0627%0.0886 0.1773 0.0355%0.0879 0.1583 0.1319P1 = zeros();A_size = size(A);O_size = size(O);N = A_size(1,1);M = A_size(1,2);K = O_size(1,1);% 调用前向后向算法函数ForwardBackwardAlgo_1，计算得到alpha矩阵P1，beta矩阵P2和观测序列O的概率p1[P1,P2,p1] = ForwardBackwardAlgo(A,B,Pi,O);% t表示观测序列时间长度数，i表示状态数；% -------------计算单个状态概率gamma = Gamma--------------Gamma = zeros();for t = 1:Kfor i = 1:MGamma(i,t) = P1(i,t) * P2(i,t) / p1;endend% --------------计算两个状态概率xi = Xi-----------------Xi = zeros();for t = 1:K-1for i = 1:Mfor j = 1:NXi(i,j,t) = P1(i,t) * A(i,j) * B(j,O(t+1,1)) * P2(j,t+1) / p1;endendend
3，迭代次数为n的Baum-Welch算法
function [A_n,B_n,Pi_n] = BaumWelchAlgo_n(A,B,Pi,O,n)% 函数功能：n次迭代的BaumWelchAlgo_n算法for r = 1:n[A,B,Pi] = BaumWelchAlgo(A,B,Pi,O);endA_n = A;B_n = B;Pi_n = Pi;
下一篇介绍HMM解码程序-维特比算法。
