一、隐马尔科夫HMM假设:
有且只有3种天气:0晴天,1阴天,2雨天
各种天气间的隔天转化概率mp:
mp[3][3] | 晴天 | 阴天 | 雨天 |
晴天 | 0.33333 | 0.33333 | 0.33333 |
阴天 | 0.33333 | 0.33333 | 0.33333 |
雨天 | 0.33333 | 0.33333 | 0.33333 |
有2种活动: 0去公园,1不去公园
各种天气下进行各种活动的概率:
w2a[3][2] | 去公园 | 不去公园 |
晴天 | 0.75 | 0.25 |
阴天 | 0.4 | 0.6 |
雨天 | 0.25 | 0.75 |
观察5天的活动序列:0 0 1 0 1 ;(0去公园,1不去公园)
5天的动作观察序列 O[5] | ||||
1天 | 2天 | 3天 | 4天 | 5天 |
去公园 | 去公园 | 不去 | 去公园 | 不去 |
0 | 0 | 1 | 0 | 1 |
第0天的天气概率pi:
pi[3] | 第0天天气概率 |
晴天 | 0.5 |
阴天 | 0.3 |
雨天 | 0.2 |
定义几个宏及变量表示HMM:
#define T 5 //观察N天
#define M 3 //每天可能有M种天气
#define N 2 //动作种类,去公园+不去公园
float mp[M][M]; //相邻两天的天气转换概率
float w2a[M][N]; //weather to action,各天气下采取各动作的概率
int O[T]; //N天观察到的天气序列
int bestPath[T][M];//bestPath[t][i]表示 (第t天处于第i种天气状态且出现“O0-Ot”观察序列的概率最大的)路径在时间t-1时刻的天气状态
float pi[M]; //一个M维向量,表示第一天各种天气出现的概率
二、前向算法--出现观察动作序列的概率
即求出现:“1去公园,2去公园,3不去公园,4去公园,5不去公园” 动作序列的概率。
这是一个求和问题,用DP思想的前向算法解决。
1、构造矩阵float sumP[T][M];
sumP[t][i]表示:在时间t,处于天气i,且出现观察动作序列O0-Ot的概率
2、填表求解过程:
//初始表
for(i=0; i<M; i++)
{
sumP[0][i] = pi[i]*w2a[ i ][ O[0] ];
}
//填表
for(t=1; t<T; t++)
{
for (i=0; i<M; i++)
{
sumP[t][i]=0;
for (j=0; j<M; j++)
{
sumP[t][i] += sumP[t-1][j]*mp[j][i]*w2a[ i ][ O[t] ];
}
}
}
求观察序列出现概率:
float sumPP=0;
for(i=1;i<M;i++)
{
sumPP += sumP[T-1][i];
}
cout<<"使用前向算法算出,观察序列的出现概率为"<<sumPP<<endl;
三、维特比解码-哪种天气下出现观察动作序列的概率最大
在哪种天气序列下,出现观察序列的概率最大,并求 (出现该天气序列和观察动作序列事件)的最大概率。
这是最优化问题,DP思想求最大,采用维特比解码
1、构造矩阵maxP[T][M]
maxP[t][i], 算出第t天处于第i状态且出现观察序列O0~Ot的路径中,(出现天气路径和观察路径)概率最大的路径的概率。
bestPath[t][i],表示是第t天,天气状态为i状态出现观察序列O0~Ot的最大概率路径下(即取得maxP[t][i]),第t-1天的天气状态。
2、填表过程:
//初始表
float maxpp=0;
int maxPre=0;
for(i=0;i<M;i++)
{
maxP[0][i] = pi[i]*w2a[ i ][ O[0] ];
if(maxP[0][i]>maxpp)
{
bestPath[0][i]=-1;
maxpp=maxP[0][i];
}
}
//后填表
for(t=1;t<T;t++) //每一天
{
for(i=0;i<M;i++)//maxP[t][j],要算出第t天处于第i状态且出现观察序列O0~Ot的路径中,(出现天气路径和观察路径)概率最大的路径的概率
{
maxpp=0;
maxPre=0;
for(j=0;j<M;j++)//第t-1天,天气为第j状态
{
float temp = maxP[t-1][j]*mp[j][i]*w2a[i][O[t]];
if(temp > maxpp)
{
maxpp = temp;
maxPre = j;
}
}
maxP[t][i] = maxpp;
bestPath[t][i] = maxPre;
}
}
输出最大概率及其天气路径:
float maxEndP=0;
int lastChoice;
for(i=0; i<M; i++)
{
if(maxP[T-1][i] > maxEndP)
{
maxEndP = maxP[T-1][i];
lastChoice = i;
}
}
cout<<"最大的概率为:"<<maxEndP<<endl;
cout<<"所有路径中,出现观察序列概率的最大的天气路径为:"<<endl;
cout<<lastChoice<<" ";
for(t=T-1; t>0; t--)
{
cout<<bestPath[t][lastChoice]<<" ";
lastChoice = bestPath[t][lastChoice];
}
cout<<endl;
四、代码:
#include<iostream>
using namespace std;
#define T 5 //观察N天
#define M 3 //每天可能有M种天气
#define N 2 //动作种类,去公园+不去公园
float mp[M][M]; //相邻两天的天气转换概率
float w2a[M][N]; //weather to action,各天气下采取各动作的概率
int O[T]; //N天观察到的天气序列
int bestPath[T][M];//bestPath[t][i]表示 (第t天处于第i种天气状态且出现“O0-Ot”观察序列的概率最大的)路径在时间t-1时刻的天气状态
float pi[M]; //一个M维向量,表示第一天各种天气出现的概率,【通常包含一个1,其余为0】
void initHMM()
{
int i,j,k;
//输入天气转移概率
for(i=0; i<M; i++)
{
for(j=0;j<M;j++)
{
cin>>mp[i][j];
}
}
//输入各天气下,采取不同动作的概率
for(i=0; i<M; i++)
{
for(j=0; j<N; j++)
{
cin>>w2a[i][j];
}
}
//输入pi
for(i=0; i<M; i++)
{
cin>>pi[i];
}
//输入观察序列
for(i=0; i<T; i++)
{
cin>>O[i];
}
}
/*维特比算法
已知HMM模型,转移概率(天气转换概率),初始状态(第0天的天气或天气概率),各天气下采取各动作的概率,观察序列(动作序列)
1.求出(使观察序列出现概率最大的)天气路径,
即在那种天气序列下,出现观察序列的概率最大
2.求出 (出现天气路径且出现观察序列的事件)的最大概率
这是一种最优化问题,求最大,DP思想
*/
float viterbi()
{
float maxP[T][M];
int i,j,t;
//初始表
float maxpp=0;
int maxPre=0;
for(i=0;i<M;i++)
{
maxP[0][i] = pi[i]*w2a[ i ][ O[0] ];
if(maxP[0][i]>maxpp)
{
bestPath[0][i]=-1;
maxpp=maxP[0][i];
}
}
//后填表
for(t=1;t<T;t++) //每一天
{
for(i=0;i<M;i++)//maxP[t][j],要算出第t天处于第i状态且出现观察序列O0~Ot的路径中,(出现天气路径和观察路径)概率最大的路径的概率
{
maxpp=0;
maxPre=0;
for(j=0;j<M;j++)//第t-1天,天气为第j状态
{
float temp = maxP[t-1][j]*mp[j][i]*w2a[i][O[t]];
if(temp > maxpp)
{
maxpp = temp;
maxPre = j;
}
}
maxP[t][i] = maxpp;
bestPath[t][i] = maxPre;
}
}
float maxEndP=0;
int lastChoice;
for(i=0; i<M; i++)
{
if(maxP[T-1][i] > maxEndP)
{
maxEndP = maxP[T-1][i];
lastChoice = i;
}
}
cout<<"最大的概率为:"<<maxEndP<<endl;
cout<<"所有路径中,出现观察序列概率的最大的天气路径(逆序)为:"<<endl;
cout<<lastChoice<<" ";
for(t=T-1; t>0; t--)
{
cout<<bestPath[t][lastChoice]<<" ";
lastChoice = bestPath[t][lastChoice];
}
cout<<endl;
return maxEndP;
}
/*
前向算法
已知HMM模型,转移概率(天气转换概率),初始状态(第0天的天气或天气概率),各天气下采取各动作的概率,观察序列(动作序列)
1.求出观察序列出现的概率为多大
即,在给定的初始天气,转移天气概率,及天气动作概率的条件下,出现观察到的动作概率是多少
这是一个加法求和问题,DP思想
*/
void forward()
{
float sumP[T][M];
//sumP[t][i]表示:在时间t,处于天气i,且出现观察动作序列O0-Ot的概率
int t,i,j,k;
//初始表
for(i=0; i<M; i++)
{
sumP[0][i] = pi[i]*w2a[ i ][ O[0] ];
}
//填表
for(t=1; t<T; t++)
{
for (i=0; i<M; i++)
{
sumP[t][i]=0;
for (j=0; j<M; j++)
{
sumP[t][i] += sumP[t-1][j]*mp[j][i]*w2a[ i ][ O[t] ];
}
}
}
float sumPP=0;
for(i=1;i<M;i++)
{
sumPP += sumP[T-1][i];
}
cout<<"使用前向算法算出,观察序列的出现概率为"<<sumPP<<endl;
}
int main()
{
initHMM();
viterbi();
forward();
system("pause");
}
/*
0.33333 0.33333 0.33333
0.33333 0.33333 0.33333
0.33333 0.33333 0.33333
0.75 0.25
0.4 0.6
0.25 0.75
0.5 0.3 0.2
0 0 1 0 1
*/
五、运行结果:
0.33333 0.33333 0.33333
0.33333 0.33333 0.33333
0.33333 0.33333 0.33333
0.75 0.25
0.4 0.6
0.25 0.75
0.5 0.3 0.2
0 0 1 0 1
最大的概率为:0.00146479
所有路径中,出现观察序列概率的最大的天气路径为:
2 0 2 0 0
使用前向算法算出,观察序列的出现概率为0.0284842
六、参数估计
在本例子中的参数为(lamda):
第一天的状态概率pi、隐状态(weather)隔天互转的概率表mp、各种天气下发生各种动作的概率表 w2a。
训练数据情况:
如果有一批样本,
1 同时给出每个样本的天气序列和动作序列,那么以上参数只要简单统计一下就能算出。属于有监督的学习。
2 如果训练样本不给出天气序列(隐状态序列)那么需要使用Baul-Welch方法计算,Baul-Welch本质是一个EM(最大期望)算法。属于无监督的学习。
自己的一个想法:初始化lamda为合理的随机值,有一批无隐状态信息的样本,通过维特比解码可以算出最可能出现的“隐状态序列”及其概率,这样就可以使用1方法重新计算lamda。迭代反复计算样lamda,直到lamda收敛。
实际的Baul-Welch算法,前向算法结合一个类似的后向过程算lamda。看的不太懂:http://blog.csdn.net/xueyingxue001/article/details/52396494