专业的网站建设费用,国外ip地址怎么弄,零基础室内设计难学吗,哪里建设网站不需要备案问题
现要求通过吉布斯采样方法#xff0c;利用该网络进行概率推理#xff08;计算 P(RT|SF, WT)、P2(CF|WT)的概率值#xff09;。
原理
吉布斯采样的核心思想为一维一维地进行采样#xff0c;采某一个维度的时候固定其他的维度#xff0c;在本次实验中#xff0c;假…
问题
现要求通过吉布斯采样方法利用该网络进行概率推理计算 P(RT|SF, WT)、P2(CF|WT)的概率值。
原理
吉布斯采样的核心思想为一维一维地进行采样采某一个维度的时候固定其他的维度在本次实验中假定上一个采样的样本为C(True)、S(False)、R(True)、W(False)此时对C维度进行采样吉布斯采样将会利用 P(C|SFalse,RTrue,WFalse)的分布得到一个新的 C的值False并将该值代替原先的值产生一个新的样本C(False)、S(False)、R(True)、W(False)。
给定分布π(C,S,R,W) step 1. t0 时刻产生一个初始状态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|SFalse,RFalse,WFalse)的分布M[0,:,1,0,0]表示P(C|STrue,RFalse,WFalse)的分布。该矩阵可通过给定的贝叶斯网络进行构建。
解答
# -*- coding:utf-8 -*-# Gibbs samplingimport numpy as np
import copyclass Gibbs:def __init__(self,query_id1):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|CT) 0.1,P(S|CF) 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|CT) 0.8,P(R|CF) 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|ST,RT) 0.99, P(W|ST,RF) 0.9# P(W|SF,RT) 0.9, P(W|SF,RF) 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(axis0,keepdimsTrue) 1e-9)self.tran_matrix[1] self.tran_matrix[3] / (self.tran_matrix[3].sum(axis1,keepdimsTrue) 1e-9)self.tran_matrix[2] self.tran_matrix[3] / (self.tran_matrix[3].sum(axis2,keepdimsTrue) 1e-9)self.tran_matrix[3] self.tran_matrix[3] / (self.tran_matrix[3].sum(axis3,keepdimsTrue) 1e-9)# 初始化样本if self.query_id 1:# P(RT|SF,WT)self.ignore_var_idx [1,3] # SF,WTself.x.append([True,False,True,True])else:# P(CF|WT)self.ignore_var_idx [3] # WTself.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_id1: count 1if not x[0] and self.query_id2: 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(RT|SF,WT)采样100次结果为:,ans)
gibbs.sample(500)
ans gibbs.calculate_ans()
gibbs.reset()
print(P(RT|SF,WT)采样500次结果为:,ans)
gibbs.sample(1000)
ans gibbs.calculate_ans()
gibbs.reset()
print(P(RT|SF,WT)采样1000次结果为:,ans)gibbs Gibbs(2)
gibbs.sample(100)
ans gibbs.calculate_ans()
gibbs.reset()
print(P(CF|WT)采样100次结果为:,ans)
gibbs.sample(500)
ans gibbs.calculate_ans()
gibbs.reset()
print(P(CF|WT)采样500次结果为:,ans)
gibbs.sample(1000)
ans gibbs.calculate_ans()
gibbs.reset()
print(P(CF|WT)采样1000次结果为:,ans)
运行结果 P(RT|SF, WT)≈1 P(CF|WT)≈0.41