|
wolffd@0
|
1 % Lawn sprinker example from Russell and Norvig p454
|
|
wolffd@0
|
2 % For a picture, see http://www.cs.berkeley.edu/~murphyk/Bayes/usage.html#basics
|
|
wolffd@0
|
3
|
|
wolffd@0
|
4 N = 4;
|
|
wolffd@0
|
5 dag = zeros(N,N);
|
|
wolffd@0
|
6 C = 1; S = 2; R = 3; W = 4;
|
|
wolffd@0
|
7 dag(C,[R S]) = 1;
|
|
wolffd@0
|
8 dag(R,W) = 1;
|
|
wolffd@0
|
9 dag(S,W)=1;
|
|
wolffd@0
|
10
|
|
wolffd@0
|
11 false = 1; true = 2;
|
|
wolffd@0
|
12 ns = 2*ones(1,N); % binary nodes
|
|
wolffd@0
|
13
|
|
wolffd@0
|
14 %bnet = mk_bnet(dag, ns);
|
|
wolffd@0
|
15 bnet = mk_bnet(dag, ns, 'names', {'cloudy','S','R','W'}, 'discrete', 1:4);
|
|
wolffd@0
|
16 names = bnet.names;
|
|
wolffd@0
|
17 %C = names{'cloudy'};
|
|
wolffd@0
|
18 bnet.CPD{C} = tabular_CPD(bnet, C, [0.5 0.5]);
|
|
wolffd@0
|
19 bnet.CPD{R} = tabular_CPD(bnet, R, [0.8 0.2 0.2 0.8]);
|
|
wolffd@0
|
20 bnet.CPD{S} = tabular_CPD(bnet, S, [0.5 0.9 0.5 0.1]);
|
|
wolffd@0
|
21 bnet.CPD{W} = tabular_CPD(bnet, W, [1 0.1 0.1 0.01 0 0.9 0.9 0.99]);
|
|
wolffd@0
|
22
|
|
wolffd@0
|
23
|
|
wolffd@0
|
24 CPD{C} = reshape([0.5 0.5], 2, 1);
|
|
wolffd@0
|
25 CPD{R} = reshape([0.8 0.2 0.2 0.8], 2, 2);
|
|
wolffd@0
|
26 CPD{S} = reshape([0.5 0.9 0.5 0.1], 2, 2);
|
|
wolffd@0
|
27 CPD{W} = reshape([1 0.1 0.1 0.01 0 0.9 0.9 0.99], 2, 2, 2);
|
|
wolffd@0
|
28 joint = zeros(2,2,2,2);
|
|
wolffd@0
|
29 for c=1:2
|
|
wolffd@0
|
30 for r=1:2
|
|
wolffd@0
|
31 for s=1:2
|
|
wolffd@0
|
32 for w=1:2
|
|
wolffd@0
|
33 joint(c,s,r,w) = CPD{C}(c) * CPD{S}(c,s) * CPD{R}(c,r) * ...
|
|
wolffd@0
|
34 CPD{W}(s,r,w);
|
|
wolffd@0
|
35 end
|
|
wolffd@0
|
36 end
|
|
wolffd@0
|
37 end
|
|
wolffd@0
|
38 end
|
|
wolffd@0
|
39
|
|
wolffd@0
|
40 joint2 = repmat(reshape(CPD{C}, [2 1 1 1]), [1 2 2 2]) .* ...
|
|
wolffd@0
|
41 repmat(reshape(CPD{S}, [2 2 1 1]), [1 1 2 2]) .* ...
|
|
wolffd@0
|
42 repmat(reshape(CPD{R}, [2 1 2 1]), [1 2 1 2]) .* ...
|
|
wolffd@0
|
43 repmat(reshape(CPD{W}, [1 2 2 2]), [2 1 1 1]);
|
|
wolffd@0
|
44
|
|
wolffd@0
|
45 assert(approxeq(joint, joint2));
|
|
wolffd@0
|
46
|
|
wolffd@0
|
47
|
|
wolffd@0
|
48 engine = jtree_inf_engine(bnet);
|
|
wolffd@0
|
49
|
|
wolffd@0
|
50 evidence = cell(1,N);
|
|
wolffd@0
|
51 evidence{W} = true;
|
|
wolffd@0
|
52
|
|
wolffd@0
|
53 [engine, ll] = enter_evidence(engine, evidence);
|
|
wolffd@0
|
54
|
|
wolffd@0
|
55 m = marginal_nodes(engine, S);
|
|
wolffd@0
|
56 p1 = m.T(true) % P(S=true|W=true) = 0.4298
|
|
wolffd@0
|
57 lik1 = exp(ll); % P(W=true) = 0.6471
|
|
wolffd@0
|
58 assert(approxeq(p1, 0.4298));
|
|
wolffd@0
|
59 assert(approxeq(lik1, 0.6471));
|
|
wolffd@0
|
60
|
|
wolffd@0
|
61 pSandW = sumv(joint(:,true,:,true), [C R]); % P(S,W) = sum_cr P(CSRW)
|
|
wolffd@0
|
62 pW = sumv(joint(:,:,:,true), [C S R]);
|
|
wolffd@0
|
63 pSgivenW = pSandW / pW; % P(S=t|W=t) = P(S=t,W=t)/P(W=t)
|
|
wolffd@0
|
64 assert(approxeq(pW, lik1))
|
|
wolffd@0
|
65 assert(approxeq(pSgivenW, p1))
|
|
wolffd@0
|
66
|
|
wolffd@0
|
67
|
|
wolffd@0
|
68 m = marginal_nodes(engine, R);
|
|
wolffd@0
|
69 p2 = m.T(true) % P(R=true|W=true) = 0.7079
|
|
wolffd@0
|
70
|
|
wolffd@0
|
71 pRandW = sumv(joint(:,:,true,true), [C S]); % P(R,W) = sum_cr P(CSRW)
|
|
wolffd@0
|
72 pRgivenW = pRandW / pW; % P(R=t|W=t) = P(R=t,W=t)/P(W=t)
|
|
wolffd@0
|
73 assert(approxeq(pRgivenW, p2))
|
|
wolffd@0
|
74
|
|
wolffd@0
|
75
|
|
wolffd@0
|
76 % Add extra evidence that R=true
|
|
wolffd@0
|
77 evidence{R} = true;
|
|
wolffd@0
|
78
|
|
wolffd@0
|
79 [engine, ll] = enter_evidence(engine, evidence);
|
|
wolffd@0
|
80
|
|
wolffd@0
|
81 m = marginal_nodes(engine, S);
|
|
wolffd@0
|
82 p3 = m.T(true) % P(S=true|W=true,R=true) = 0.1945
|
|
wolffd@0
|
83 assert(approxeq(p3, 0.1945))
|
|
wolffd@0
|
84
|
|
wolffd@0
|
85
|
|
wolffd@0
|
86 pSandRandW = sumv(joint(:,true,true,true), [C]); % P(S,R,W) = sum_c P(cSRW)
|
|
wolffd@0
|
87 pRandW = sumv(joint(:,:,true,true), [C S]); % P(R,W) = sum_cs P(cSRW)
|
|
wolffd@0
|
88 pSgivenWR = pSandRandW / pRandW; % P(S=t|W=t,R=t) = P(S=t,R=t,W=t)/P(W=t,R=t)
|
|
wolffd@0
|
89 assert(approxeq(pSgivenWR, p3))
|
|
wolffd@0
|
90
|
|
wolffd@0
|
91 % So the sprinkler is less likely to be on if we know that
|
|
wolffd@0
|
92 % it is raining, since the rain can "explain away" the fact
|
|
wolffd@0
|
93 % that the grass is wet.
|
|
wolffd@0
|
94
|
|
wolffd@0
|
95 lik3 = exp(ll); % P(W=true, R=true) = 0.4581
|
|
wolffd@0
|
96 % So the combined evidence is less likely (of course)
|
|
wolffd@0
|
97
|
|
wolffd@0
|
98
|
|
wolffd@0
|
99
|
|
wolffd@0
|
100
|
|
wolffd@0
|
101 % Joint distributions
|
|
wolffd@0
|
102
|
|
wolffd@0
|
103 evidence = cell(1,N);
|
|
wolffd@0
|
104 [engine, ll] = enter_evidence(engine, evidence);
|
|
wolffd@0
|
105 m = marginal_nodes(engine, [S R W]);
|
|
wolffd@0
|
106
|
|
wolffd@0
|
107 evidence{R} = 2;
|
|
wolffd@0
|
108 [engine, ll] = enter_evidence(engine, evidence);
|
|
wolffd@0
|
109 m = marginal_nodes(engine, [S R W]);
|
|
wolffd@0
|
110
|
|
wolffd@0
|
111
|
|
wolffd@0
|
112
|
