annotate toolboxes/FullBNT-1.0.7/bnt/examples/static/sprinkler1.m @ 0:cc4b1211e677 tip

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