每日算法 — 使用java实现隐马尔可夫模型：前向算法与维特比解码

引言

隐马尔可夫模型（Hidden Markov Model, HMM）是序列数据建模的经典概率图模型，广泛应用于语音识别、自然语言处理、生物信息学等领域。HMM假设系统存在一组不可直接观测的隐状态，这些隐状态按马尔可夫链演化，并生成可观测的输出。本文将用Java实现HMM的核心算法——前向算法（概率计算）、维特比算法（最优路径解码）与Baum-Welch参数估计，帮助读者掌握序列标注与模式识别的底层原理。

核心原理

HMM的形式化定义

一个HMM由五元组 λ = (Q, V, A, B, π) 定义：

Q = {q₁, q₂, …, q_N}: 隐状态集合，共N个状态
V = {v₁, v₂, …, v_M}: 观测符号集合，共M种观测
A = [a_{ij}]: 状态转移概率矩阵，a_{ij} = P(q_j | q_i)
B = [b_j(k)]: 观测概率矩阵，b_j(k) = P(v_k | q_j)
π = [π_i]: 初始状态概率分布

三大基本问题

评估问题：给定模型λ和观测序列O，计算P(O|λ)——用前向算法解决
解码问题：给定模型λ和观测序列O，求最可能的隐状态序列——用维特比算法解决
学习问题：给定观测序列O，估计模型参数λ——用Baum-Welch算法解决

前向算法

定义前向概率 α_t(i) = P(O₁, O₂, …, O_t, q_t = i | λ)，表示到时刻t观测到部分序列且状态为i的概率。

α_{t+1}(j) = [Σ_i α_t(i) · a_{ij}] · b_j(O_{t+1})

最终 P(O|λ) = Σ_i α_T(i)

维特比算法

定义 δ_t(i) 为到时刻t以状态i结尾的所有路径的最大概率：

δ_{t+1}(j) = max_i [δ_t(i) · a_{ij}] · b_j(O_{t+1})

通过记录每个时刻的最优前驱状态，最终回溯得到全局最优路径。

Java实现

项目结构

src/
├── HMM.java            // HMM模型与核心算法
├── HMMTrainer.java     // Baum-Welch训练器
└── HMMDemo.java        // 演示主类

HMM核心类

import java.util.Arrays;

/**
 * 隐马尔可夫模型实现
 */
public class HMM {
    public final int numStates;      // 隐状态数 N
    public final int numObservations; // 观测符号数 M

    public double[][] A;  // 状态转移矩阵 [N][N]
    public double[][] B;  // 观测概率矩阵 [N][M]
    public double[] pi;   // 初始概率 [N]

    public HMM(int numStates, int numObservations) {
        this.numStates = numStates;
        this.numObservations = numObservations;
        this.A = new double[numStates][numStates];
        this.B = new double[numStates][numObservations];
        this.pi = new double[numStates];
    }

    /**
     * 前向算法: 计算观测序列的概率 P(O|lambda)
     */
    public double forward(int[] observations) {
        int T = observations.length;
        double[][] alpha = new double[T][numStates];

        // 初始化: t=0
        for (int i = 0; i < numStates; i++) {
            alpha[0][i] = pi[i] * B[i][observations[0]];
        }

        // 递推: t=1 to T-1
        for (int t = 1; t < T; t++) {
            for (int j = 0; j < numStates; j++) {
                double sum = 0;
                for (int i = 0; i < numStates; i++) {
                    sum += alpha[t-1][i] * A[i][j];
                }
                alpha[t][j] = sum * B[j][observations[t]];
            }
        }

        // 终止: 求和
        double prob = 0;
        for (int i = 0; i < numStates; i++) {
            prob += alpha[T-1][i];
        }
        return prob;
    }

    /**
     * 维特比算法: 求最优隐状态序列
     * @return [最优路径概率, 最优状态序列]
     */
    public Object[] viterbi(int[] observations) {
        int T = observations.length;
        double[][] delta = new double[T][numStates];
        int[][] psi = new int[T][numStates]; // 记录前驱状态

        // 初始化
        for (int i = 0; i < numStates; i++) {
            delta[0][i] = pi[i] * B[i][observations[0]];
            psi[0][i] = 0;
        }

        // 递推
        for (int t = 1; t < T; t++) {
            for (int j = 0; j < numStates; j++) {
                double maxProb = -1;
                int maxState = 0;
                for (int i = 0; i < numStates; i++) {
                    double prob = delta[t-1][i] * A[i][j];
                    if (prob > maxProb) {
                        maxProb = prob;
                        maxState = i;
                    }
                }
                delta[t][j] = maxProb * B[j][observations[t]];
                psi[t][j] = maxState;
            }
        }

        // 终止: 找最终最优状态
        double maxFinalProb = -1;
        int finalState = 0;
        for (int i = 0; i < numStates; i++) {
            if (delta[T-1][i] > maxFinalProb) {
                maxFinalProb = delta[T-1][i];
                finalState = i;
            }
        }

        // 回溯最优路径
        int[] path = new int[T];
        path[T-1] = finalState;
        for (int t = T - 2; t >= 0; t--) {
            path[t] = psi[t+1][path[t+1]];
        }

        return new Object[]{maxFinalProb, path};
    }

    /**
     * 后向算法: 计算后向概率 beta_t(i)
     */
    public double[][] backward(int[] observations) {
        int T = observations.length;
        double[][] beta = new double[T][numStates];

        // 初始化: t=T-1
        Arrays.fill(beta[T-1], 1.0);

        // 递推: t=T-2 downto 0
        for (int t = T - 2; t >= 0; t--) {
            for (int i = 0; i < numStates; i++) {
                double sum = 0;
                for (int j = 0; j < numStates; j++) {
                    sum += A[i][j] * B[j][observations[t+1]] * beta[t+1][j];
                }
                beta[t][i] = sum;
            }
        }
        return beta;
    }

    /**
     * Baum-Welch单步更新 (EM算法的一次迭代)
     */
    public void baumWelchStep(int[] observations) {
        int T = observations.length;
        double[][] alpha = computeAlpha(observations);
        double[][] beta = backward(observations);

        // 计算gamma_t(i): 时刻t处于状态i的概率
        double[][] gamma = new double[T][numStates];
        for (int t = 0; t < T; t++) {
            double denom = 0;
            for (int i = 0; i < numStates; i++) {
                gamma[t][i] = alpha[t][i] * beta[t][i];
                denom += gamma[t][i];
            }
            for (int i = 0; i < numStates; i++) {
                gamma[t][i] /= denom;
            }
        }

        // 计算xi_t(i,j): 时刻t从状态i转移到j的概率
        double[][][] xi = new double[T-1][numStates][numStates];
        for (int t = 0; t < T - 1; t++) {
            double denom = 0;
            for (int i = 0; i < numStates; i++) {
                for (int j = 0; j < numStates; j++) {
                    xi[t][i][j] = alpha[t][i] * A[i][j] * B[j][observations[t+1]] * beta[t+1][j];
                    denom += xi[t][i][j];
                }
            }
            for (int i = 0; i < numStates; i++) {
                for (int j = 0; j < numStates; j++) {
                    xi[t][i][j] /= denom;
                }
            }
        }

        // 重估计参数
        // 更新pi
        for (int i = 0; i < numStates; i++) {
            pi[i] = gamma[0][i];
        }

        // 更新A
        for (int i = 0; i < numStates; i++) {
            double gammaSum = 0;
            for (int t = 0; t < T - 1; t++) gammaSum += gamma[t][i];
            for (int j = 0; j < numStates; j++) {
                double xiSum = 0;
                for (int t = 0; t < T - 1; t++) xiSum += xi[t][i][j];
                A[i][j] = gammaSum > 0 ? xiSum / gammaSum : 0;
            }
        }

        // 更新B
        for (int j = 0; j < numStates; j++) {
            double gammaSum = 0;
            for (int t = 0; t < T; t++) gammaSum += gamma[t][j];
            for (int k = 0; k < numObservations; k++) {
                double obsSum = 0;
                for (int t = 0; t < T; t++) {
                    if (observations[t] == k) obsSum += gamma[t][j];
                }
                B[j][k] = gammaSum > 0 ? obsSum / gammaSum : 0;
            }
        }
    }

    private double[][] computeAlpha(int[] observations) {
        int T = observations.length;
        double[][] alpha = new double[T][numStates];
        for (int i = 0; i < numStates; i++) {
            alpha[0][i] = pi[i] * B[i][observations[0]];
        }
        for (int t = 1; t < T; t++) {
            for (int j = 0; j < numStates; j++) {
                double sum = 0;
                for (int i = 0; i < numStates; i++) {
                    sum += alpha[t-1][i] * A[i][j];
                }
                alpha[t][j] = sum * B[j][observations[t]];
            }
        }
        return alpha;
    }
}

演示主类

/**
 * HMM演示: 天气预测模型
 * 隐状态: 晴天(0), 雨天(1)
 * 观测: 散步(0), 购物(1), 打扫(2)
 */
public class HMMDemo {
    public static void main(String[] args) {
        HMM hmm = new HMM(2, 3);

        // 初始化模型参数
        hmm.pi = new double[]{0.6, 0.4};
        hmm.A = new double[][]{
            {0.7, 0.3},  // 晴->晴 0.7, 晴->雨 0.3
            {0.4, 0.6}   // 雨->晴 0.4, 雨->雨 0.6
        };
        hmm.B = new double[][]{
            {0.6, 0.3, 0.1},  // 晴天: 散步0.6, 购物0.3, 打扫0.1
            {0.1, 0.4, 0.5}   // 雨天: 散步0.1, 购物0.4, 打扫0.5
        };

        // 观测序列: 散步->购物->打扫
        int[] observations = {0, 1, 2};

        // 1. 前向算法: 计算观测序列概率
        double prob = hmm.forward(observations);
        System.out.println("观测序列概率 P(O|lambda): " + String.format("%.6f", prob));

        // 2. 维特比算法: 解码最优隐状态序列
        Object[] result = hmm.viterbi(observations);
        double maxProb = (double) result[0];
        int[] path = (int[]) result[1];
        String[] states = {"晴天", "雨天"};

        System.out.println("最优路径概率: " + String.format("%.6f", maxProb));
        System.out.print("最优隐状态序列: ");
        for (int s : path) {
            System.out.print(states[s] + " ");
        }
        System.out.println();

        // 3. Baum-Welch训练
        System.out.println("\n=== Baum-Welch参数重估计 ===");
        hmm.baumWelchStep(observations);
        System.out.println("更新后的转移矩阵A:");
        for (double[] row : hmm.A) {
            System.out.println(java.util.Arrays.toString(row));
        }
    }
}

复杂度分析

算法	时间复杂度	空间复杂度
前向算法	O(T · N²)	O(T · N)
维特比算法	O(T · N²)	O(T · N)
Baum-Welch(单次迭代)	O(T · N² · M)	O(T · N²)

其中 T 为序列长度，N 为隐状态数，M 为观测符号数。HMM的计算效率使其在实时语音识别等场景中仍具实用价值。

总结

本文用Java完整实现了隐马尔可夫模型的三大核心算法：前向算法高效计算序列概率、维特比算法通过动态规划解码最优路径、Baum-Welch算法基于EM框架迭代优化模型参数。HMM是理解更复杂序列模型（如CRF、RNN、Transformer）的重要基础，读者可将其应用于中文分词、词性标注、异常检测等实际任务。

Tagged with: HMM 前向算法动态规划维特比算法隐马尔可夫模型