问题

现要求通过吉布斯采样方法，利用该网络进行概率推理（计算 P(R=T|S=F, W=T)、P2(C=F|W=T)的概率值）。

原理

吉布斯采样的核心思想为一维一维地进行采样，采某一个维度的时候固定其他的维度，在本次实验中，假定上一个采样的样本为<C(True)、S(False)、R(True)、W(False)>，此时对C维度进行采样，吉布斯采样将会利用 P(C|S=False,R=True,W=False)的分布得到一个新的 C的值（False），并将该值代替原先的值产生一个新的样本<C(False)、S(False)、R(True)、W(False)>。

给定分布π(C,S,R,W)
step 1. t=0 时刻产生一个初始状态<C0,S0,R0,W0>,count = 0,采样序列 x=[<C0,S0,R0,W0>]
step 2. 从条件概率分布 P(C|S0,R0)采样得到<C1,S0,R0,W0>，加入到 x[C 已知跳转下一步]
step 3. 从条件概率分布 P(S|C1,R0,W0)采样得到<C1,S1,R0,W0>，加入到 x[S 已知跳转下一步]
step 4. 从条件概率分布 P(R|C1,S1,W0)采样得到<C1,S1,R1,W0>，加入到 x[R 已知跳转下一步]
step 5. 从条件概率分布 P(W|S1,R1)采样得到<C1,S1,R1,W1>，加入到 x[W 已知跳转下一步]
step 6. count = count + 1,如果 count < 指定采样次数跳转到 step2
step 7. 在采样序列中统计满足条件的样本数量，除以总采样数即为所求。
数据结构使用一个 4*2*2*2*2 的矩阵 M 用于存储条件概率分布。如 M[0,:,0,0,0]即表示
P(C|S=False,R=False,W=False)的分布，M[0,:,1,0,0]表示P(C|S=True,R=False,W=False)的分布。该矩阵可通过给定的贝叶斯网络进行构建。

解答

# -*- coding:utf-8 -*-# Gibbs samplingimport numpy as np
import copyclass Gibbs:def __init__(self,query_id=1):self.x = []self.query_id = query_idassert query_id == 1 or query_id ==2self.tran_matrix = np.zeros((4,2,2,2,2))# 0 : C Cloudy# 1 : S Sprinkler# 2 : R Rain# 3 : W Wet grass# 计算条件概率分布# P(C) = 0.5self.tran_matrix[0] = 0.5 # P(C) = 0.5# P(S|C=T) = 0.1,P(S|C=F) = 0.5self.tran_matrix[1,1,0,:,:] = self.tran_matrix[0,1,0,:,:] * (1-0.1)self.tran_matrix[1,1,1,:,:] = self.tran_matrix[0,1,1,:,:] * 0.1self.tran_matrix[1,0,0,:,:] = self.tran_matrix[0,0,0,:,:] * (1-0.5)self.tran_matrix[1,0,1,:,:] = self.tran_matrix[0,0,1,:,:] * 0.5# P(R|C=T) = 0.8,P(R|C=F) = 0.2self.tran_matrix[2,1,:,0,:] = self.tran_matrix[1,1,:,0,:] * (1-0.8)self.tran_matrix[2,1,:,1,:] = self.tran_matrix[1,1,:,1,:] * 0.8self.tran_matrix[2,0,:,0,:] = self.tran_matrix[1,0,:,0,:] * (1-0.2)self.tran_matrix[2,0,:,1,:] = self.tran_matrix[1,0,:,1,:] * 0.2# P(W|S=T,R=T) = 0.99, P(W|S=T,R=F) = 0.9# P(W|S=F,R=T) = 0.9, P(W|S=F,R=F) = 0self.tran_matrix[3,:,1,1,0] = self.tran_matrix[2,:,1,1,0] * (1-0.99)self.tran_matrix[3,:,1,1,1] = self.tran_matrix[2,:,1,1,1] * 0.99self.tran_matrix[3,:,1,0,0] = self.tran_matrix[2,:,1,0,0] * (1-0.9)self.tran_matrix[3,:,1,0,1] = self.tran_matrix[2,:,1,0,1] * 0.9self.tran_matrix[3,:,0,1,0] = self.tran_matrix[2,:,0,1,0] * (1-0.9)self.tran_matrix[3,:,0,1,1] = self.tran_matrix[2,:,0,1,1] * 0.9self.tran_matrix[3,:,0,0,0] = self.tran_matrix[2,:,0,0,0] * (1-0)self.tran_matrix[3,:,0,0,1] = self.tran_matrix[2,:,0,0,1] * 0self.tran_matrix[0] = self.tran_matrix[3] / (self.tran_matrix[3].sum(axis=0,keepdims=True) + 1e-9)self.tran_matrix[1] = self.tran_matrix[3] / (self.tran_matrix[3].sum(axis=1,keepdims=True) + 1e-9)self.tran_matrix[2] = self.tran_matrix[3] / (self.tran_matrix[3].sum(axis=2,keepdims=True) + 1e-9)self.tran_matrix[3] = self.tran_matrix[3] / (self.tran_matrix[3].sum(axis=3,keepdims=True) + 1e-9)# 初始化样本if self.query_id == 1:# P(R=T|S=F,W=T)self.ignore_var_idx = [1,3] # S=F,W=Tself.x.append([True,False,True,True])else:# P(C=F|W=T)self.ignore_var_idx = [3] # W=Tself.x.append([True,False,True,True])self._sample_axis = 0self._var_num = 4def sample(self,sample_num:int):for _ in range(sample_num * (self._var_num - len(self.ignore_var_idx))):last_x = copy.copy(self.x[-1])sample_axis = self._next_sample_axis()last_x[sample_axis]=Truesample_prob = self.tran_matrix[sample_axis,int(last_x[0]),int(last_x[1]),\int(last_x[2]),int(last_x[3])]if np.random.rand() < sample_prob:last_x[sample_axis] = Trueelse:last_x[sample_axis] = Falseself.x.append(last_x)self.x = self.x[::self._var_num - len(self.ignore_var_idx)]def _next_sample_axis(self):self._sample_axis += 1self._sample_axis %= self._var_numwhile self._sample_axis in self.ignore_var_idx:self._sample_axis += 1self._sample_axis %= self._var_numreturn self._sample_axisdef calculate_ans(self):count = 0for x in self.x:if x[2] and self.query_id==1: count += 1if not x[0] and self.query_id==2: count += 1return count / len(self.x)def reset(self):self.x = self.x[:1]gibbs = Gibbs(1)
gibbs.sample(100)
ans = gibbs.calculate_ans()
gibbs.reset()
print('P(R=T|S=F,W=T)采样100次，结果为:',ans)
gibbs.sample(500)
ans = gibbs.calculate_ans()
gibbs.reset()
print('P(R=T|S=F,W=T)采样500次，结果为:',ans)
gibbs.sample(1000)
ans = gibbs.calculate_ans()
gibbs.reset()
print('P(R=T|S=F,W=T)采样1000次，结果为:',ans)gibbs = Gibbs(2)
gibbs.sample(100)
ans = gibbs.calculate_ans()
gibbs.reset()
print('P(C=F|W=T)采样100次，结果为:',ans)
gibbs.sample(500)
ans = gibbs.calculate_ans()
gibbs.reset()
print('P(C=F|W=T)采样500次，结果为:',ans)
gibbs.sample(1000)
ans = gibbs.calculate_ans()
gibbs.reset()
print('P(C=F|W=T)采样1000次，结果为:',ans)

运行结果

P(R=T|S=F, W=T)≈1
P(C=F|W=T)≈0.41

【人工智能】对贝叶斯网络进行吉布斯采样相关推荐

【人工智能】— 贝叶斯网络、概率图模型、全局语义、因果链、朴素贝叶斯模型、枚举推理、变量消元
[人工智能]- 贝叶斯网络频率学派 vs. 贝叶斯学派贝叶斯学派 Probability(概率): 独立性/条件独立性: Probability Theory(概率论): Graphical mo ...
电子科技大学人工智能期末复习笔记（四）：概率与贝叶斯网络
目录前言概率概率公式贝叶斯公式链式条件概率例题 1. 求联合概率分布/边缘概率分布/条件概率分布 2. 灵活运用贝叶斯公式概率总结贝叶斯网络判断独立性两个事件独立的判断条件独立性 ...
贝叶斯网络之父Judea Pearl：新因果科学与数据科学、人工智能的思考
来源:AI科技评论本文约6000字,建议阅读10分钟本文介绍贝叶斯网络之父 Judea Pearl <新因果科学与数据科学.人工智能的思考>的报告. 标签:人工智能 6月21日,图灵奖 ...
贝叶斯网络之父Judea Pearl：要建立真正的人工智能，少不了因果推理
来源:专知参与 | Yingying, Xiaowen, Sanglei 2011年图灵奖得主,贝叶斯网络之父Judea Pearl认为现在人工智能的发展进入的新的瓶颈期,各种新的成果不过本质上不过 ...
人工智能学习（十）：什么是贝叶斯网络——伯克利版
目录 10.1 概率建模 10.1.1 独立性 10.1.2 条件独立 10.1.2.1 条件独立和链式法则 10.2 贝叶斯网络 10.2.1 图形化的模型符号 10.2.2 贝叶斯网络的构建 10 ...
[人工智能AI]之贝叶斯网络
[人工智能AI]之贝叶斯网络(Bayesian network) 部分图片和来源自: NJU-人工智能-高阳教授的课件通俗地讲,贝叶斯网络就是用一组有向无环图,表示多个事件的因果依赖关系,并借此完 ...
在Python中使用贝叶斯网络的实例
我们在之前的文章中(请见文末给出的参考资料[1])已经介绍了贝叶斯网络的基本原理,以及基于贝叶斯网络进行概率推断(Exact Inference)的消去法.本文将结合一个具体的例子来演示在Python ...
基于在软件工程中对贝叶斯网络的循证决策
基于在软件工程中对贝叶斯网络的循证决策摘要:在软件工程中的推荐系统应该设计成集成依据并成为从业人员的经验.贝叶斯网络为以证据为基础的决策提供了自然统计框架,通过结合现有证据的综合概要与相关的不确定性 ...
贝叶斯网络之父Judea Pearl力荐、LeCun点赞，这篇长论文全面解读机器学习中的因果关系...
来源:机器之心作者:Bernhard Schölkopf 图灵奖得主.贝叶斯网络之父 Judea Pearl 曾自嘲自己是「AI 社区的反叛者」,因为他对人工智能发展方向的观点与主流趋势相反.Pea ...

【人工智能】对贝叶斯网络进行吉布斯采样

问题

原理

解答

运行结果

【人工智能】对贝叶斯网络进行吉布斯采样相关推荐

最新文章

热门文章