diff options
Diffstat (limited to 'sourcecodes/bnt-master/BNT/examples/static/Belprop')
13 files changed, 431 insertions, 0 deletions
diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/CVS/Entries b/sourcecodes/bnt-master/BNT/examples/static/Belprop/CVS/Entries new file mode 100644 index 00000000..6550caa0 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/CVS/Entries @@ -0,0 +1,11 @@ +/belprop_loop1_discrete.m/1.1.1.1/Wed May 29 15:59:54 2002// +/belprop_loop1_gauss.m/1.1.1.1/Wed May 29 15:59:54 2002// +/belprop_loopy_cg.m/1.1.1.1/Wed May 29 15:59:54 2002// +/belprop_loopy_discrete.m/1.1.1.1/Wed May 29 15:59:54 2002// +/belprop_loopy_gauss.m/1.1.1.1/Wed May 29 15:59:54 2002// +/belprop_polytree_cg.m/1.1.1.1/Wed May 29 15:59:54 2002// +/belprop_polytree_discrete.m/1.1.1.1/Tue Oct 1 18:21:26 2002// +/belprop_polytree_gauss.m/1.1.1.1/Wed May 29 15:59:54 2002// +/bp1.m/1.1.1.1/Wed May 29 15:59:54 2002// +/gmux1.m/1.1.1.1/Wed May 29 15:59:54 2002// +D diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/CVS/Repository b/sourcecodes/bnt-master/BNT/examples/static/Belprop/CVS/Repository new file mode 100644 index 00000000..f3d573bd --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/examples/static/Belprop diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/CVS/Root b/sourcecodes/bnt-master/BNT/examples/static/Belprop/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loop1_discrete.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loop1_discrete.m new file mode 100644 index 00000000..20faed5e --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loop1_discrete.m @@ -0,0 +1,26 @@ +% Compare different loopy belief propagation algorithms on a graph with a single loop. +% LBP should give exact results if it converges. + +N = 4; +dag = zeros(N,N); +C = 1; S = 2; R = 3; W = 4; +dag(C,[R S]) = 1; +dag(R,W) = 1; +dag(S,W)=1; +ns = 2*ones(1,N); +bnet = mk_bnet(dag, ns); +for i=1:N + bnet.CPD{i} = tabular_CPD(bnet, i); +end + +engines = {}; +engines{end+1} = jtree_inf_engine(bnet); +engines{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel'); +engines{end+1} = belprop_fg_inf_engine(bnet_to_fgraph(bnet)); +engines{end+1} = belprop_inf_engine(bnet, 'protocol', 'parallel'); + +% belprop_fg does not support marginal_family +% belprop_fg and belprop do not support loglik even on discrete +[time, engines] = cmp_inference_static(bnet, engines, 'maximize', 0, 'exact', 1, 'observed', 2, ... + 'check_ll', 0, 'singletons_only', 1, 'check_converged', 2:4); + diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loop1_gauss.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loop1_gauss.m new file mode 100644 index 00000000..e547e44a --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loop1_gauss.m @@ -0,0 +1,25 @@ +% Compare different loopy belief propagation algorithms on a graph with a single loop. +% LBP should give exact results if it converges. + +N = 4; +dag = zeros(N,N); +C = 1; S = 2; R = 3; W = 4; +dag(C,[R S]) = 1; +dag(R,W) = 1; +dag(S,W)=1; +ns = 2*ones(1,N); +bnet = mk_bnet(dag, ns, 'discrete', []); +for i=1:N + bnet.CPD{i} = gaussian_CPD(bnet, i); +end + +engines = {}; +engines{end+1} = jtree_inf_engine(bnet); +engines{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel', 'max_iter', 20); +%engines{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel', 'max_iter', 20, 'filename', ... +% '/home/eecs/murphyk/matlab/gausspearl.txt', 'tol', 1e-5); + +% pearl gaussian does not compute loglik +[time, engines] = cmp_inference_static(bnet, engines, 'maximize', 0, 'exact', 1, 'observed', [2], ... + 'check_ll', 0, 'singletons_only', 0, 'check_converged', [2]); + diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loopy_cg.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loopy_cg.m new file mode 100644 index 00000000..01f36b03 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loopy_cg.m @@ -0,0 +1,22 @@ +% Same as cg1, except we assume all discretes are observed, +% and use loopy for approximate inference. + +ns = 2*ones(1,9); +F = 1; W = 2; E = 3; B = 4; C = 5; D = 6; Min = 7; Mout = 8; L = 9; +n = 9; +dnodes = [B F W]; +cnodes = mysetdiff(1:n, dnodes); + +%bnet = mk_incinerator_bnet(ns); +bnet = mk_incinerator_bnet; + +bnet.observed = [dnodes E]; + +engines = {}; +engines{end+1} = jtree_inf_engine(bnet); +engines{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel'); +nengines = length(engines); + + +[time, engines] = cmp_inference_static(bnet, engines, 'maximize', 0, 'check_ll', 0, ... + 'singletons_only', 0, 'exact', 1, 'check_converged', 2); diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loopy_discrete.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loopy_discrete.m new file mode 100644 index 00000000..c46e6d02 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loopy_discrete.m @@ -0,0 +1,13 @@ +% Compare different loopy belief propagation algorithms on a graph with many loops + +bnet = mk_asia_bnet('orig'); + +engines = {}; +engines{end+1} = jtree_inf_engine(bnet); +engines{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel'); +engines{end+1} = belprop_fg_inf_engine(bnet_to_fgraph(bnet)); +engines{end+1} = belprop_inf_engine(bnet, 'protocol', 'parallel'); + +[time, engines] = cmp_inference_static(bnet, engines, 'maximize', 0, 'exact', 1, 'observed', [1 3 5], ... + 'check_ll', 0, 'singletons_only', 1, 'check_converged', 2:4); + diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loopy_gauss.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loopy_gauss.m new file mode 100644 index 00000000..a9924aed --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_loopy_gauss.m @@ -0,0 +1,12 @@ +% Compare different loopy belief propagation algorithms on a graph with many loops +% If LBP converges, the means should be exact + +bnet = mk_asia_bnet('gauss'); + +engines = {}; +engines{end+1} = jtree_inf_engine(bnet); +engines{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel'); + +[time, engines] = cmp_inference_static(bnet, engines, 'maximize', 0, 'exact', 1, 'observed', [1 3 5], ... + 'check_ll', 0, 'singletons_only', 0, 'check_converged', 2); + diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_polytree_cg.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_polytree_cg.m new file mode 100644 index 00000000..70aa03da --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_polytree_cg.m @@ -0,0 +1,33 @@ +% Inference on a conditional Gaussian model + +% Make the following polytree, where all arcs point down + +% 1 2 +% \ / +% 3 +% / \ +% 4 5 + +N = 5; +dag = zeros(N,N); +dag(1,3) = 1; +dag(2,3) = 1; +dag(3, [4 5]) = 1; + +ns = [2 1 2 1 2]; + +dnodes = 1; +%onodes = [1 5]; +bnet = mk_bnet(dag, ns, 'discrete', dnodes, 'observed', dnodes); + +bnet.CPD{1} = tabular_CPD(bnet, 1); +for i=2:N + bnet.CPD{i} = gaussian_CPD(bnet, i); +end + +engine = {}; +engine{end+1} = jtree_inf_engine(bnet); +engine{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel'); + +[time, engine] = cmp_inference_static(bnet, engine, 'maximize', 0, 'check_ll', 0, ... + 'singletons_only', 0, 'observed', [1 3]); diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_polytree_discrete.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_polytree_discrete.m new file mode 100644 index 00000000..d8a3a229 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_polytree_discrete.m @@ -0,0 +1,38 @@ +% Make the following polytree, where all arcs point down + +% 1 2 +% \ / +% 3 +% / \ +% 4 5 + +N = 5; +dag = zeros(N,N); +dag(1,3) = 1; +dag(2,3) = 1; +dag(3, [4 5]) = 1; + +ns = 2*ones(1,N); % binary nodes + +onodes = [1 5]; + +bnet = mk_bnet(dag, ns, 'observed', onodes); + +if 0 +seed = 0; +rand('state', seed); +randn('state', seed); +end + +for i=1:N + %bnet.CPD{i} = tabular_CPD(bnet, i); + bnet.CPD{i} = noisyor_CPD(bnet, i); +end + +engine = {}; +engine{end+1} = jtree_inf_engine(bnet); +engine{end+1} = pearl_inf_engine(bnet, 'protocol', 'tree'); +engine{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel'); + +[err, time] = cmp_inference_static(bnet, engine, 'maximize', 0, 'check_ll', 1); + diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_polytree_gauss.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_polytree_gauss.m new file mode 100644 index 00000000..1823c107 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/belprop_polytree_gauss.m @@ -0,0 +1,135 @@ +% Do the example from Satnam Alag's PhD thesis, UCB ME dept 1996 p46 + +% Make the following polytree, where all arcs point down + +% 1 2 +% \ / +% 3 +% / \ +% 4 5 + +N = 5; +dag = zeros(N,N); +dag(1,3) = 1; +dag(2,3) = 1; +dag(3, [4 5]) = 1; + +ns = [2 1 2 1 2]; + +bnet = mk_bnet(dag, ns, 'discrete', []); + +bnet.CPD{1} = gaussian_CPD(bnet, 1, 'mean', [1 0]', 'cov', [4 1; 1 4]); +bnet.CPD{2} = gaussian_CPD(bnet, 2, 'mean', 1, 'cov', 1); +B1 = [1 2; 1 0]; B2 = [2 1]'; +bnet.CPD{3} = gaussian_CPD(bnet, 3, 'mean', [0 0]', 'cov', [2 1; 1 1], ... + 'weights', [B1 B2]); +H1 = [1 1]; +bnet.CPD{4} = gaussian_CPD(bnet, 4, 'mean', 0, 'cov', 1, 'weights', H1); +H2 = [1 0; 1 1]; +bnet.CPD{5} = gaussian_CPD(bnet, 5, 'mean', [0 0]', 'cov', eye(2), 'weights', H2); + +engine = {}; +engine{end+1} = jtree_inf_engine(bnet); +engine{end+1} = pearl_inf_engine(bnet, 'protocol', 'tree'); +engine{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel'); +E = length(engine); + +if 1 +% no evidence +evidence = cell(1,N); +ll = zeros(1,E); +for e=1:E + [engine{e}, ll(e)] = enter_evidence(engine{e}, evidence); + add_ev = 1; + m = marginal_nodes(engine{e}, 3, add_ev); + assert(approxeq(m.mu, [3 2]')) + assert(approxeq(m.Sigma, [30 9; 9 6])) + + m = marginal_nodes(engine{e}, 4, add_ev); + assert(approxeq(m.mu, 5)) + assert(approxeq(m.Sigma, 55)) + + m = marginal_nodes(engine{e}, 5, add_ev); + assert(approxeq(m.mu, [3 5]')) + assert(approxeq(m.Sigma, [31 39; 39 55])) +end +end + +if 1 +% evidence on leaf 5 +evidence = cell(1,N); +evidence{5} = [5 5]'; +for e=1:E + [engine{e}, ll(e)] = enter_evidence(engine{e}, evidence); + add_ev = 1; + m = marginal_nodes(engine{e}, 3, add_ev); + assert(approxeq(m.mu, [4.4022 1.0217]')) + assert(approxeq(m.Sigma, [0.7011 -0.4891; -0.4891 1.1087])) + + m = marginal_nodes(engine{e}, 4, add_ev); + assert(approxeq(m.mu, 5.4239)) + assert(approxeq(m.Sigma, 1.8315)) + + m = marginal_nodes(engine{e}, 1, add_ev); + assert(approxeq(m.mu, [0.3478 1.1413]')) + assert(approxeq(m.Sigma, [1.8261 -0.1957; -0.1957 1.0924])) + + m = marginal_nodes(engine{e}, 2, add_ev); + assert(approxeq(m.mu, 0.9239)) + assert(approxeq(m.Sigma, 0.8315)) + + m = marginal_nodes(engine{e}, 5, add_ev); + assert(approxeq(m.mu, evidence{5})) + assert(approxeq(m.Sigma, zeros(2))) +end +end + +if 1 +% evidence on leaf 4 (non-info-state version is uninvertible) +evidence = cell(1,N); +evidence{4} = 10; +for e=1:E + [engine{e}, ll(e)] = enter_evidence(engine{e}, evidence); + add_ev = 1; + m = marginal_nodes(engine{e}, 3, add_ev); + assert(approxeq(m.mu, [6.5455 3.3636]')) + assert(approxeq(m.Sigma, [2.3455 -1.6364; -1.6364 1.9091])) + + m = marginal_nodes(engine{e}, 5, add_ev); + assert(approxeq(m.mu, [6.5455 9.9091]')) + assert(approxeq(m.Sigma, [3.3455 0.7091; 0.7091 1.9818])) + + m = marginal_nodes(engine{e}, 1, add_ev); + assert(approxeq(m.mu, [1.9091 0.9091]')) + assert(approxeq(m.Sigma, [2.1818 -0.8182; -0.8182 2.1818])) + + m = marginal_nodes(engine{e}, 2, add_ev); + assert(approxeq(m.mu, 1.2727)) + assert(approxeq(m.Sigma, 0.8364)) +end +end + + +if 1 +% evidence on leaves 4,5 and root 2 +evidence = cell(1,N); +evidence{2} = 0; +evidence{4} = 10; +evidence{5} = [5 5]'; +for e=1:E + [engine{e}, ll(e)] = enter_evidence(engine{e}, evidence); + add_ev = 1; + m = marginal_nodes(engine{e}, 3, add_ev); + assert(approxeq(m.mu, [4.9964 2.4444]')); + assert(approxeq(m.Sigma, [0.6738 -0.5556; -0.5556 0.8889])); + + m = marginal_nodes(engine{e}, 1, add_ev); + assert(approxeq(m.mu, [2.2043 1.2151]')); + assert(approxeq(m.Sigma, [1.2903 -0.4839; -0.4839 0.8065])); +end +end + +if 1 + [time, engine] = cmp_inference_static(bnet, engine, 'maximize', 0, 'check_ll', 0, ... + 'singletons_only', 0, 'observed', [1 3 5]); +end diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/bp1.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/bp1.m new file mode 100644 index 00000000..93cba443 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/bp1.m @@ -0,0 +1,24 @@ +% Compare different loopy belief propagation algorithms on a graph with a single loop. +% LBP should give exact results if it converges. + +seed = 0; +rand('state', seed); +randn('state', seed); + +N = 2; +dag = zeros(N,N); +dag(1,2)=1; +ns = ones(1,N); +bnet = mk_bnet(dag, ns, 'discrete', []); +for i=1:N + %bnet.CPD{i} = gaussian_CPD(bnet, i, 'mean', 0); + bnet.CPD{i} = gaussian_CPD(bnet, i); +end + +engines = {}; +engines{end+1} = jtree_inf_engine(bnet); +engines{end+1} = pearl_inf_engine(bnet, 'protocol', 'tree'); + +[time, engines] = cmp_inference_static(bnet, engines, 'maximize', 0, 'exact', 1:2, 'observed', [2], ... + 'check_ll', 0, 'singletons_only', 1, 'check_converged', []); + diff --git a/sourcecodes/bnt-master/BNT/examples/static/Belprop/gmux1.m b/sourcecodes/bnt-master/BNT/examples/static/Belprop/gmux1.m new file mode 100644 index 00000000..21846b68 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/examples/static/Belprop/gmux1.m @@ -0,0 +1,90 @@ +% Test gmux. +% The following model, where Y is a gmux node, +% and M is set to 1, should be equivalent to X1 -> Y +% +% X1 Xn M +% \ | / +% Y + +n = 3; +N = n+2; +Xs = 1:n; +M = n+1; +Y = n+2; +dag = zeros(N,N); +dag([Xs M], Y)=1; + +dnodes = M; +ns = zeros(1, N); +sz = 2; +ns(Xs) = sz; +ns(M) = n; +ns(Y) = sz; + +bnet = mk_bnet(dag, ns, 'discrete', M, 'observed', [M Y]); + +psz = ns(Xs(1)); +selfsz = ns(Y); + +W = randn(selfsz, psz); +mu = randn(selfsz, 1); +Sigma = eye(selfsz, selfsz); + +bnet.CPD{M} = root_CPD(bnet, M); +for i=Xs(:)' + bnet.CPD{i} = gaussian_CPD(bnet, i, 'mean', zeros(psz, 1), 'cov', eye(psz, psz)); +end +bnet.CPD{Y} = gmux_CPD(bnet, Y, 'mean', mu, 'weights', W, 'cov', Sigma); + +evidence = cell(1,N); +yval = randn(selfsz, 1); +evidence{Y} = yval; +m = 2; +%notm = not(m-1)+1; % only valid for n=2 +notm = mysetdiff(1:n, m); +evidence{M} = m; + +engines = {}; +engines{end+1} = jtree_inf_engine(bnet); +engines{end+1} = pearl_inf_engine(bnet, 'protocol', 'parallel'); + +for e=1:length(engines) + engines{e} = enter_evidence(engines{e}, evidence); + mXm{e} = marginal_nodes(engines{e}, Xs(m)); + + % Since M=m, only Xm was updated. + % Hence the posterior on Xnotm should equal the prior. + for i=notm(:)' + mXnotm = marginal_nodes(engines{e}, Xs(i)); + assert(approxeq(mXnotm.mu, zeros(psz,1))) + assert(approxeq(mXnotm.Sigma, eye(psz, psz))) + end +end + +% Check that all engines give the same posterior +for e=2:length(engines) + assert(approxeq(mXm{e}.mu, mXm{1}.mu)) + assert(approxeq(mXm{e}.Sigma, mXm{1}.Sigma)) +end + + +% Compute the correct posterior by building Xm -> Y + +N = 2; +dag = zeros(N,N); +dag(1, 2)=1; +ns = [psz selfsz]; +bnet = mk_bnet(dag, ns, 'discrete', [], 'observed', 2); + +bnet.CPD{1} = gaussian_CPD(bnet, 1, 'mean', zeros(psz, 1), 'cov', eye(psz, psz)); +bnet.CPD{2} = gaussian_CPD(bnet, 2, 'mean', mu, 'cov', Sigma, 'weights', W); + +jengine = jtree_inf_engine(bnet); +evidence = {[], yval}; +jengine = enter_evidence(jengine, evidence); % apply Bayes rule to invert the arc +mX = marginal_nodes(jengine, 1); + +for e=1:length(engines) + assert(approxeq(mX.mu, mXm{e}.mu)) + assert(approxeq(mX.Sigma, mXm{e}.Sigma)) +end |
