（1）给定一个观察序列O=O1O2...OT,和模型μ=（A,B,π），如何快速有效地选择在一定意义下“最优”的状态序列Q=q1q2...qT，使该状态最好地解释观察序列。

二：实例篇

（1）假设连续观察3天的海藻湿度为(Dry,Damp,Soggy),求这三天最可能的天气情况。天气只有三类(Sunny,Cloudy,Rainy)，而且海藻湿度和天气有一定的关系。

4. 状态转移矩阵

   Sunny        Cloudy   Rainy
Sunny 0.5         0.375  0.125
Cloudy 0.25 0.125 0.625
Rainy 0.25 0.375 0.375

Cloudy(昨天)->Sunny(今天)的概率是0.25;

Sunny(昨天)->Rainy(今天)的概率是0.125.

5. 混淆矩阵（海藻湿度与天气的相关性）：

（2）计算方法：

由HMM可知，Day2的天气仅取决于Day1；Day3的天气又只取决于Day2的天气。

step1：Day1由于是初始状态，我们分别求

P(Day1-Sunny)=0.63*0.6;
P(Day1-Cloudy)=0.17*0.25;
P(Day1-Rain)=0.20*0.05;
Choose max{ P(Day1-Sunny) , P(Day1-Cloudy),P(Day1-Rainy)}, 得到P(Day1-Sunny)最大，得出第1天Sunny的概率最大。

step2: Day2的天气又取决于Day1的天气状况，同时也受Day2观察的海藻情况影响。

P(Day2-Sunny)= max{ P(Day1-Sunny)*0.5, P(Day1-Cloudy)*0.25,  P(Day1-Rainy)*0.25} *0.15;
P(Day2-Cloudy)= max{ P(Day1-Sunny)*0.375,  P(Day1-Cloudy)*0.125, P(Day1-Rainy)*0.625} *0.25;
P(Day2-Rainy)= max{ P(Day1-Sunny)*0.125,  P(Day1-Cloudy)*0.625 , P(Day1-Rainy)*0.375} *0.35;
Choosemax{ P(Day2-Sunny) , P(Day2-Cloudy), P(Day2-Rainy)},得到P(Day2-Rainy)最大，得出第2天Rainy的概率最大。
故{Sunny,Rainy}是前两天最大可能的天气序列。

step3: Day3的天气又取决于Day2的天气状况，同时也受Day3观察的海藻情况影响。

P(Day3-Sunny)= max{ P(Day2-Sunny)*0.5, P(Day2-Cloudy)*0.25,  P(Day2-Rainy)*0.25} *0.05;
P(Day3-Cloudy)= max{ P(Day2-Sunny)*0.375,  P(Day2-Cloudy)*0.125, P(Day2-Rainy)*0.625} *0.25;
P(Day3-Rainy)= max{ P(Day2-Sunny)*0.125,  P(Day2-Cloudy)*0.625, P(Day2-Rainy)*0.375} *0. 05;
Choosemax{ P(Day3-Sunny) , P(Day3-Cloudy), P(Day3-Rainy)},得到 P(Day3-Rainy)最大，得出第3天Rainy的概率最大。

(3)代码实现

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
/* 题目描述：给定一个观察序列O=O1O2...OT,和模型μ=（A,B,π），
如何快速有效地选择在一定意义下“最优”的状态序列Q=q1q2...qT，
使该状态最好地解释观察序列。
下面a是状态转移矩阵； b是发射矩阵，list 是观察序列；
*/
const int HID_STATUSES = 3;  // 隐藏的状态个数
const int VIS_ST = 4;        // 观察到的实际所有可能的情况（状态）个数
const int N = 3;             // 观察的天数（序列长度）
int main()
{
    float a[HID_STATUSES][HID_STATUSES] = {{0.5,0.375,0.125},{0.25,0.125,0.625},{0.25,0.375,0.375}};//状态转移矩阵
    float b[HID_STATUSES][VIS_ST] = {{0.6,0.2,0.15,0.05},{0.25,0.25,0.25,0.25,},{0.05,0.10,0.35,0.50,}};//b是发射矩阵
    float ini[HID_STATUSES] = {0.63,0.17,0.20};// 初始的
    float result[N][HID_STATUSES];
    int list[N] = {0,2,3};// 测试用的观察到的状态序列,4天的观测值；由于就有2中观察状态，所以用
    // 0 和 1 区分，实际中应该是input的，以及观察的天数也是input的或者read file
    int max[N][HID_STATUSES];
    float tmp;
    int i,j,k;
    //step1:Initialization，其中{0.2,0.4,0.4} 是π
    for(i=0;i<HID_STATUSES;++i)
    {
        result[0][i] = ini[i]*b[i][list[0]];
    }

    int count = 1, max_node;// count 是list的下标，从1开始，第一天是step1已经初始化
    float max_v;
    //step2:归纳运算 ，i 从1开始，应为step1已经初始化了
    for (i=1; i<N; i++)
    {
        for(j=0; j<HID_STATUSES; j++)
        {
            tmp = result[i-1][0] * a[0][j] * b[j][list[count]];
            max[i][j] = 0;
            for(k=1; k<HID_STATUSES; k++)
            {
                if(result[i-1][k] * a[k][j] * b[j][list[count]] > tmp)
                {
                    tmp = result[i-1][k] * a[k][j]* b[j][list[count]];
                    max[i][j] = k;
                }
               result[i][j] = tmp;
            }
            max_v = result[N-1][0];
            max_node = 0;
            for (k=1; k<HID_STATUSES; k++)
            {
                if(result[N-1][k] > max_v)
                {
                    max_v = result[N-1][k];
                    max_node = k;
                }
            }
        }
        count += 1;
    }
    //step3:终结
   for (i=0; i<N; i++)
    {
        for(j=0; j<HID_STATUSES; j++)
        {
            printf("%d %d     %f\n",i+1,j+1,result[i][j]);

        }
    }
    printf("Pmax= %f\n", max_v);
    printf("step%d: %d   \n",N,max_node+1);
    //step4:回溯
    for(k=N-1; k>=1; k--)
    {
        printf("step%d: %d  \n",k, max[k][max_node]+1);
        max_node = max[k][max_node];
    }
    return 0;
}

大数据之道 HMM系列