hmm(二）

最新推荐文章于 2021-07-15 23:14:51 发布
原创 最新推荐文章于 2021-07-15 23:14:51 发布 · 617 阅读 · 0 ·
CC 4.0 BY-SA版权
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。

参数估计

根据观察到的序列集来找到一个最有可能的 HMM


代码

#coding:utf-8

class Hmm:
	def __init__(self, description=None):
		self.description = description

	@property
	def hidden_states(self):
	    return self._hidden_states
	@hidden_states.setter
	def hidden_states(self, vals):
		self._hidden_states = vals

	@property
	def observed_states(self):
	    return self._observed_states
	@observed_states.setter
	def observed_states(self, vals):
		self._observed_states = vals

	@property
	def transition_matrix(self):
	    return self._transition_matrix
	@transition_matrix.setter
	def transition_matrix(self, vals):
		self._transition_matrix = vals

	@property
	def emission_matrix(self):
	    return self._emission_matrix
	@emission_matrix.setter
	def emission_matrix(self, vals):
		self._emission_matrix = vals

	@property
	def start_prob(self):
		return self._start_prob
	@start_prob.setter
	def start_prob(self, val):
		self._start_prob = val

	@property
	def final_prob(self):
		return self._final_prob
	@final_prob.setter
	def final_prob(self, val):
		self._final_prob = val

	@property
	def num_states(self):
	    return self._num_states
	@num_states.setter
	def num_states(self, val):
		self._num_states = val

	@property
	def num_ob_states(self):
	    return self._num_ob_states
	@num_states.setter
	def num_ob_states(self, val):
		self._num_ob_states = val
	
	


def hmmfunc_wrapper(func):
	def inner_wrapper(*args, **kwargs):
		print '*** before called method=[%s]***' % func.__name__
		func(args, kwargs)
		print '*** after called method=[%s]***' % func.__name__
	return inner_wrapper

import copy
import numpy as np

class HmmFunc:
	def __init__(self, description):
		self.description = description

	def _find_start_state(self,state, src_states):
		for i , val in enumerate(src_states):
			if val == state:
				return i
		return 0

	def cal_likehood(self, observerd_datas, hmm):

		for j,state in enumerate(observerd_datas):
			if j == 0:
				#alpha(t,j) = a0(j) * emission(j,0)
				start_idx = self._find_start_state(observerd_datas[0], hmm.observed_states)
				alpha = [1.0 * val * hmm.emission_matrix[i][start_idx] for i , val in enumerate(hmm.start_prob)]
				continue

			tmp = [0.0] * len(alpha)
			# Transition : N * N
			# Emission   : N * Ob
			# alpha(t,j) = sum_k( alpha(t-1,k)* transition(k,m) )* emission(m,j)
			idx = self._find_start_state(state,hmm.observed_states)
			num_states = hmm.num_states
			for k in range(num_states):
				sum = 0
				for m in range(num_states):
					sum += alpha[m] * hmm.transition_matrix[m][k]
				#endfor m
				tmp[k] = sum * hmm.emission_matrix[k][idx]
			#endfor k
			alpha = [ v for v in tmp]
		return np.sum(alpha)

	def cal_most_prob_seq(self, observerd_datas, hmm):
		#max_probs_end in step j
		seq_len = len(observerd_datas)
		num_states = hmm.num_states
		backwards = [(num_states * [1.0])] * seq_len
		delta = [(num_states * [1.0])] * seq_len
		for j , state in enumerate(observerd_datas):
			if j == 0:
				start_idx = self._find_start_state(observerd_datas[0], hmm.observed_states)
				delta[0] = [1.0 * val * hmm.emission_matrix[i][start_idx] for i , val in enumerate(hmm.start_prob)]
				continue
			idx = self._find_start_state(state, hmm.observed_states)
			for k in range(num_states):
				max_val = -1000
				maxid = 0
				for m in range(num_states):
					 tmp = delta[j-1][m] * hmm.transition_matrix[m][k]
					 if tmp > max_val:
					 	max_val = tmp
					 	maxid = m
				delta[j][k] = max_val * hmm.emission_matrix[k][idx]
				backwards[j][k] = maxid
		#self.backwards = backwards
		#self.delta = delta
		idx = self._tranverse(delta, backwards, seq_len, num_states)
		return self._get_most_prob(idx, hmm.hidden_states)

	def _tranverse(self, delta, backwards, seq_len, num_states):
		max_val = delta[seq_len-1][0]
		max_id  = 0
		for key ,val in enumerate(delta[seq_len -1]):
			if val > max_val:
				val = max_val
				max_id = key - 1

		ret = []
		ret.append(max_id)
		for trs in range(seq_len -1 , 0, -1):
			ret.append(backwards[trs][max_id])
			max_id = backwards[trs][max_id]
		return ret

	def _get_most_prob(self, idx, hidden_states):
		a = [ hidden_states[i] for i in idx]
		#return a.reverse()
		a.reverse()
		return a

	def do_inference(self):
		pass

def normnize_arr(arr):
	sumv = sum(arr)
	return [ v / sumv for v in arr]



class HmmFucEm(HmmFunc):
	def __init__(self, hmm=None):
		self.description = None
		self.hmm = hmm
		HmmFunc.__init__(self,None)

	def cal_step_forward(self, observerd_datas):
		hmm = self.hmm
		ret = []
		for j,state in enumerate(observerd_datas):
			if j == 0:
				#alpha(t,j) = a0(j) * emission(j,0)
				start_idx = self._find_start_state(observerd_datas[0], hmm.observed_states)
				alpha = [1.0 * val * hmm.emission_matrix[i][start_idx] for i , val in enumerate(hmm.start_prob)]
				ret.append(copy.copy(alpha))
				continue

			tmp = [0.0] * len(alpha)
			# Transition : N * N
			# Emission   : N * Ob
			# alpha(t,j) = sum_k( alpha(t-1,k)* transition(k,m) )* emission(m,j)
			idx = self._find_start_state(state,hmm.observed_states)
			num_states = hmm.num_states
			for k in range(num_states):
				sum = 0
				for m in range(num_states):
					sum += alpha[m] * hmm.transition_matrix[m][k]
				#endfor m
				tmp[k] = sum * hmm.emission_matrix[k][idx]
			#endfor k
			alpha = [ v for v in tmp]
			ret.append(copy.copy(alpha))

			self.forward = copy.copy(ret)
		return np.sum(self.forward[num_states-1])

	def cal_step_backward(self, observerd_datas):
		hmm = self.hmm
		ret = []
		num_states = hmm.num_states
		ob_len = len(observerd_datas)
		# init set all one's
		beta = [1.0] * num_states
		ret.append(copy.copy(beta))
		for j in range(ob_len-1, 0,-1):
			state = observerd_datas[j] 
			tmp = [0.0] * len(beta)
			# Transition : N * N
			# Emission   : N * Ob
			# beta(t, k) =  sum_m( beta(t+1,m) * transition( k, m) * emision(m, obj(t+1)) )
			idx = self._find_start_state(state,hmm.observed_states)
			# k(n+1) -- j
			for k in range(num_states):
				sumv = 0
				for m in range(num_states):
					sumv +=  beta[m] * hmm.transition_matrix[k][m] * hmm.emission_matrix[m][idx]
				#endfor m
				tmp[k] = sumv
			#endfor k
			beta = [ v for v in tmp]
			ret.append(copy.copy(beta))
			tp = []
			for v in range(len(ret)-1,-1,-1):
				tp.append(ret[v])
		#return prob j ==0

		self.backwards = copy.copy(tp)
		state = observerd_datas[0]
		idx = self._find_start_state(state,hmm.observed_states)
		prob = 0.0
		for v in range(num_states):
			prob += beta[v] * hmm.start_prob[v] * hmm.emission_matrix[v][idx]
		return prob

	def _cal_gamma_t_expectation_maximization_e(self, num_states, ob_len):
		ret=[]
		for t in range(ob_len):
			v = [0.0] * num_states
			sumall = 0.0
			for j in range(num_states):
				v[j] = self.forward[t][j] * self.backwards[t][j]
				sumall += v[j]
			# for j in range(num_states):
			# 	v[j] /= sumall
			v = np.array(v) / sumall
			ret.append(v)
		return ret

	def _cal_gamma_ij_t_expectation_maximization_e(self, num_states, ob_len, obs):
		tmp = []
		for t in range(ob_len - 1):
			idx = self._find_start_state( obs[t+1], self.hmm.observed_states)
			gamma_ij_t = [[0.0] * num_states] * num_states
			#gamma_ij_t = np.array(gamma_ij_t)
			sum_ij = 0.0
			for i in range(num_states):
				for j in range(num_states):
					gamma_ij_t[i][j] = self.forward[t][i] 
					gamma_ij_t[i][j] *= self.hmm.transition_matrix[i][j] 
					gamma_ij_t[i][j] *= self.hmm.emission_matrix[j][idx] 
					gamma_ij_t[i][j] *= self.backwards[t+1][j]
					sum_ij += gamma_ij_t[i][j]
			#endfor_i
			# for i in range(num_states):
			# 	for j in range(num_states):
			# 		gamma_ij_t[i][j] /=1
			gamma_ij_t = np.array(gamma_ij_t) / sum_ij
			tmp.append(gamma_ij_t)
		#endfor_t
		return tmp
	def _update_params_expectation_maximization_m(self, obs):
		num_states = self.hmm.num_states
		ob_len = len(obs)
		gamma_t = self._cal_gamma_t_expectation_maximization_e(num_states, ob_len)
		gamma_ij_t = self._cal_gamma_ij_t_expectation_maximization_e(num_states, ob_len, obs)
		#update params
		self.hmm.start_prob = [v for v in gamma_t[0]]

		#update transition_matrix
		for i in range(num_states):
			for j in range(num_states):
				denum = 0.0
				num = 0.0
				for t in range(ob_len-1):
					denum += gamma_ij_t[t][i][j]
					num += gamma_t[t][i]
				self.hmm.transition_matrix[i][j] = 1.0 * denum/ num
			self.hmm.transition_matrix[i] = normnize_arr(self.hmm.transition_matrix[i])
		#endfor transition_matrix

		#update emission_matrix
		for k in range(self.hmm.num_ob_states):
			for j in range(num_states):
				self.hmm.emission_matrix[j][k] = 0.0

		for k in range(self.hmm.num_ob_states):	
			sum_gamma_t_j = 0.0
			for j in range(num_states):
				for t in range(ob_len):
					sum_gamma_t_j += gamma_t[t][j]
					idx = self._find_start_state(obs[t], self.hmm.observed_states)
					if idx == k:
						self.hmm.emission_matrix[j][k] += gamma_t[t][j]
			for j in range(num_states):
				self.hmm.emission_matrix[j][k] /= sum_gamma_t_j
		for i in range(self.hmm.num_states):
			self.hmm.emission_matrix[i] = normnize_arr(self.hmm.emission_matrix[i])

		#end_all

	def _do_safe_mathlog(self):

		num_states = self.hmm.num_states
		num_ob_states = self.hmm.num_ob_states

		for k in range(num_states):
			self.hmm.start_prob[k] = 0.00001 + 0.99999 * self.hmm.start_prob[k]

		for k in range(num_states):
			for j in range(num_states):
				self.hmm.transition_matrix[k][j] = 0.00001 + 0.99999 * self.hmm.transition_matrix[k][j]


		for k in range(num_states):
			for s in range(num_ob_states):
				self.hmm.emission_matrix[k][s] = 0.00001 + 0.99999 * self.hmm.emission_matrix[k][s]



	def _cal_prob_em(self, obs):

		from math import log

		num_states = self.hmm.num_states
		num_ob_states = self.hmm.num_ob_states
		ob_len = len(obs)
		gamma_t = self._cal_gamma_t_expectation_maximization_e(num_states, ob_len)
		gamma_ij_t = self._cal_gamma_ij_t_expectation_maximization_e(num_states, ob_len, obs)

		Qsum = 0.0

		self._do_safe_mathlog()

		for k in range(num_states):
			Qsum += gamma_t[0][k] * log(self.hmm.start_prob[k])

		gamma_ij_t_sum = [[1.0] * num_states] * num_states
		for t in range(ob_len-1):
			for i in range(num_states):
				for j in range(num_states):
					gamma_ij_t_sum[i][j] += gamma_ij_t[t][i][j]

		for i in range(num_states):
			for j in range(num_states):
				Qsum += gamma_ij_t_sum[i][j] * log(self.hmm.transition_matrix[i][j])


		for i in range(num_states):
			for j in range(num_ob_states):
				for t in range(ob_len):
					Qsum += gamma_t[t][i] * log(self.hmm.emission_matrix[i][j])

		return Qsum




	def do_inference(self, obs):
		QsumOld = self._cal_prob_em(obs)
		self._update_params_expectation_maximization_m(obs)
		QsumNew = self._cal_prob_em(obs)

		from math import fabs

		while fabs(QsumNew - QsumOld) > 0.0000001:
			QsumOld = QsumNew
			self._update_params_expectation_maximization_m(obs)
			QsumNew = self._cal_prob_em(obs)

			print 'QsumNew=%s,QsumOld=%s' % (QsumNew,QsumOld)

		return QsumNew


if __name__ == '__main__':
	hmm = Hmm()
	hmm.hidden_states = ['Sunny','Cloudy','Rainy']
	hmm.observed_states = ['Damp','Soggy','Dry','Dryish']
	hmm.num_states = 3
	#hmm.num_ob_states = len(hmm.observed_states) 
	hmm.num_ob_states = 4
	hmm.start_prob =[0.63,0.17,0.20]
	hmm.transition_matrix=[[0.5,0.25,0.25], [0.375,0.125,0.375],[0.125,0.675,0.375]]
	hmm.emission_matrix=[[0.15,0.05,0.6,0.20],[0.25,0.25,0.25,0.25],[0.35,0.5,0.05,0.10]]
	
      
	func = HmmFunc('cal')
	#given transition matrix, emission_matrix, init_prob , observerd_sequences to calculate probility
	prob = func.cal_likehood(['Damp','Soggy','Dry'], hmm)
	print 'cal_likehood array=%s with probility=%s' % ( ['Damp','Soggy','Dry'], prob)
	
	#given transition matrix, emission_matrix, init_prob, observerd_sequences to calculate most prob hidden states
	most_prob = func.cal_most_prob_seq(['Damp','Soggy','Dry'], hmm)
	print 'cal_most_prob_seq=%s ' % (most_prob)



	emFunc = HmmFucEm(hmm=hmm)
	forward = emFunc.cal_step_forward(['Damp','Soggy','Dry'])
	#print forward
	print forward, emFunc.forward
	backwards = emFunc.cal_step_backward(['Damp','Soggy','Dry'])
	print backwards, emFunc.backwards


	pp = 0.0
	for i , j in zip(emFunc.forward[2], emFunc.backwards[2]):
		pp += i * j
	print 'pp=%s' % pp

	#do inference
	print '++++++++++++++==do_inference_++++++++++++='
	inferenceFunc = HmmFucEm(hmm=hmm)
	obs = ['Damp','Soggy','Dry']
	inferenceFunc.cal_step_forward(obs)
	inferenceFunc.cal_step_backward(obs)
	inferenceFunc.do_inference(obs)
	# print inferenceFunc.hmm.start_prob

	# import CheckFunc as checkFunc

	# print 'init_prob=%s' % checkFunc.check(inferenceFunc.hmm.start_prob)


	# print inferenceFunc.hmm.transition_matrix
	# print inferenceFunc.hmm.emission_matrix

标注


1. 只用了一个sample

2.参考《统计学习方法-李航》


一文搞懂HMM（隐马尔可夫模型）
07-16 5094
什么是熵(Entropy) 简单来说，熵是表示物质系统状态的一种度量，用它老表征系统的无序程度。熵越大，系统越无序，意味着系统结构和运动的不确定和无规则；反之，，熵越小，系统越有序，意味着具有确定和有规则的运动状态。熵的中文意思是热量被温度除的商。负熵是物质系统有序化，组织化，复杂化状态的一种度量。 熵最早来原于物理学. 德国物理学家鲁道夫·克劳修斯首次提出熵的概念，用来表示任何一种...
hmm
whuawell的博客
04-16 462
介绍 代码 #coding:utf-8 class Hmm: def __init__(self, description=None): self.description = description @property def hidden_states(self): return self._hidden_states @hidden_states.sett
HMM（1）
u014168855的博客
03-26 467
1.HMM介绍    隐马尔可夫模型（HMM，Hidden Markov Model）可用于标注问题，在语音识别,NLP,生物信息，模式识别领域被证明是有效的算法。     HMM是关于时序的概率模型，描述由一个隐藏的马尔可夫链生成不可观测的状态随机序列，再由各个状态生成观测随机序列的过程。     隐马尔可夫模型随机生成的状态随机序列，称为状态序列；每个状态生成一个观测，由此生成的观测随机序列，...
HMM详解
Flying Dreams
10-12 2605
niufen.zip_HMM 图像_HMM 聚类_HMM聚类_二维 HMM_聚类 系统
09-19
标题中的“niufen.zip_HMM图像_HMM聚类_HMM聚类_二维HMM_聚类系统”揭示了这个压缩包文件的主要内容，涉及到的是基于隐马尔科夫模型（Hidden Markov Model, HMM）的图像处理和聚类算法。HMM是一种统计建模方法，常...
HMM.rar_HMM_HMM ppt
09-21
二、HMM的三个基本问题 1. 参数估计（Learning）：给定观测序列，如何估计出模型的参数，包括初始状态概率分布π和状态转移概率矩阵A。 2. 前向后向算法（Forward-Backward Algorithm）：用于计算在给定模型下，观测...
HMM天气预测.rar_GUI界面_HMM predict_gui hmm_天气gui_天气预测 hmm
07-15
用隐马尔科夫模型进行未来天气的预测，并设计了一个GUI的交互界面
HMM.rar_HMM
09-24
二、HMM的模型结构 HMM模型通常用一个有向图来表示，图中的节点代表状态，边代表状态间的转移，节点上的箭头表示状态转移的方向。每个状态有向相邻的状态数量不限，但每个状态至多有一个前驱和一个后继状态。同时，...
HMM.zip_hmm c语言_site:www.pudn.com
09-24
在C语言中，这通常涉及二维数组的初始化、更新和遍历。 3. 参数估计：实现Baum-Welch算法或最大似然估计来更新模型参数。这可能涉及到迭代过程，每次迭代都要计算新的概率并更新模型。 4. 序列匹配和识别：利用...
python实现hmm
10-16
用python实现的Hmm，有详细的代码注解，对Hmm感兴趣的同学值得下载学习
HMM初探
shiqi,bao的博客
07-07 253
马尔科夫性质： 无记忆。下一状态的概率分布只能由当前状态决定，在时间序列中它前面的时间均与之无关。 例如天气，如果我们假定天气是马尔科夫的，那么今天的天气仅仅与昨天的天气存在关联，而与前天之前的天气没有关系。
HMM
zzx3163967592的博客
12-23 320
https://zhuanlan.zhihu.com/p/62529131 一、知识回顾： HMM(隐马尔可夫模型)的完整参数集合可以用一个五元组来表示。式中各个参数含义如下： N为HMM的隐状态个数。令符号为t时刻的状态，则隐含状态序列可表示为，。其中，为隐状态集。 M为HMM的可观测状态个数。令符号为t时刻的观测，则观测序列可表示为，。其中，为观测集。 ...
ＨＭＭ
代号1998的专栏
02-18 241
HMM（隐马尔可夫模型）
hanwanbing的博客
07-15 4712
原文：https://www.cnblogs.com/skyme/p/4651331.html 笔记： 一、熵(Entropy) 熵是表示物质系统状态的一种度量，用它表征系统的无序程度。熵越大，系统越无序，意味着系统结构和运动的不确定和无规则；反之，熵越小，系统越有序，意味着具有确定和有规则的运动状态。熵的中文意思是热量被温度除的商。负熵是物质系统有序化，组织化，复杂化状态的一种度量。 二、最大熵模型 这个“熵“并不是指热力学上熵的概念，而是由信息论男神克劳德·艾尔伍德...
HMM在二维形状识别中的应用探索
"本文介绍了HMM（隐马尔科夫模型）在二维形状识别中的应用，适合初学者了解和学习。" 在模式识别领域，尤其是在计算机视觉中，二维形状识别是一项基础且关键的任务。人类视觉系统识别物体时，形状通常是最先被注意...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值