diff options
Diffstat (limited to 'sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine')
18 files changed, 840 insertions, 0 deletions
diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Entries b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Entries new file mode 100644 index 00000000..055aa4df --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Entries @@ -0,0 +1,11 @@ +/README/1.1.1.1/Sun May 11 15:39:50 2003// +/clq_containing_nodes.m/1.1.1.1/Wed May 29 11:59:46 2002// +/enter_evidence.m/1.1.1.1/Wed Mar 12 10:38:00 2003// +/marginal_difclq_nodes.m/1.1.1.1/Fri Feb 21 11:20:32 2003// +/marginal_nodes.m/1.1.1.1/Fri Feb 21 11:13:10 2003// +/marginal_singleclq_nodes.m/1.1.1.1/Wed Jan 29 11:23:58 2003// +/problems.txt/1.1.1.1/Wed May 29 11:59:46 2002// +/push.m/1.1.1.1/Mon Feb 10 15:38:04 2003// +/push_pot_toclique.m/1.1.1.1/Wed May 29 11:59:46 2002// +/stab_cond_gauss_inf_engine.m/1.1.1.1/Fri Mar 28 17:12:42 2003// +D diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Entries.Log b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Entries.Log new file mode 100644 index 00000000..24f16336 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Entries.Log @@ -0,0 +1 @@ +A D/Old//// diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Repository b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Repository new file mode 100644 index 00000000..849daef7 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/inference/static/@stab_cond_gauss_inf_engine diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Root b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/CVS/Entries b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/CVS/Entries new file mode 100644 index 00000000..ce0c4813 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/CVS/Entries @@ -0,0 +1,2 @@ +/initialize_engine.m/1.1.1.1/Wed May 29 11:59:46 2002// +D diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/CVS/Repository b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/CVS/Repository new file mode 100644 index 00000000..eb292815 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/inference/static/@stab_cond_gauss_inf_engine/Old diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/CVS/Root b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/initialize_engine.m b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/initialize_engine.m new file mode 100644 index 00000000..6fb51c2e --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/Old/initialize_engine.m @@ -0,0 +1,65 @@ +function [engine, loglik] = initialize_engine(engine) +%initialize +bnet = bnet_from_engine(engine); +ns = bnet.node_sizes(:); +N = length(bnet.dag); + +pot_type = 'scg' +check_for_cd_arcs([], bnet.cnodes, bnet.dag); + +% Evaluate CPDs with evidence, and convert to potentials +pot = cell(1, N); +C = length(engine.cliques); +inited = zeros(1, C); +clpot = cell(1, C); +evidence = cell(1, N); +for n=1:N + fam = family(bnet.dag, n); + e = bnet.equiv_class(n); + pot{n} = CPD_to_scgpot(bnet.CPD{e}, fam, ns, bnet.cnodes, evidence); + cindex = engine.clq_ass_to_node(n); + if inited(cindex) + %clpot{cindex} = direct_combine_pots(clpot{cindex}, pot{n}); + clpot{cindex} = direct_combine_pots(pot{n}, clpot{cindex}); + else + clpot{cindex} = pot{n}; + inited(cindex) = 1; + end +end + +for i=1:C + if inited(i) == 0 + clpot{i} = scgpot([], [], [], []); + end +end + +seppot = cell(C, C); +% separators are is not need to initialize + +% collect to root (node to parents) +for n=engine.postorder(1:end-1) + for p=parents(engine.jtree, n) + [margpot, comppot] = complement_pot(clpot{n}, engine.separator{p,n}); + margpot = marginalize_pot(clpot{n}, engine.separator{p,n}); + clpot{n} = comppot; + %seppot{p, n} = margpot; + clpot{p} = combine_pots(clpot{p}, margpot); + %clpot{p} = combine_pots(margpot, clpot{p}); + end +end + +temppot = clpot; +%temppot = clpot{engine.root}; +for n=engine.preorder + for c=children(engine.jtree, n) + seppot{n,c} = marginalize_pot(temppot{n}, engine.separator{n,c}); + %seppot{n,c} = marginalize_pot(clpot{n}, engine.separator{n,c}); + %clpot{c} = direct_combine_pots(clpot{c}, seppot{n,c}); + temppot{c} = direct_combine_pots(temppot{c}, seppot{n,c}); + end +end + +engine.clpot = clpot; +engine.seppot = seppot; + + diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/README b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/README new file mode 100644 index 00000000..e905e28c --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/README @@ -0,0 +1,12 @@ +% Stable conditional Gaussian inference +% Originally written by Huang, Shan <shan.huang@intel.com> 2001 +% Fixed by Rainer Deventer 2003 + + +@techreport{Lauritzen99, + author = "S. Lauritzen and F. Jensen", + title = "Stable Local Computation with Conditional {G}aussian Distributions", + year = 1999, + number = "R-99-2014", + institution = "Dept. Math. Sciences, Aalborg Univ." +} diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/clq_containing_nodes.m b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/clq_containing_nodes.m new file mode 100644 index 00000000..c64d2bff --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/clq_containing_nodes.m @@ -0,0 +1,24 @@ +function c = clq_containing_nodes(engine, nodes, fam) +% CLQ_CONTAINING_NODES Find the lightest clique (if any) that contains the set of nodes +% c = clq_containing_nodes(engine, nodes, family) +% +% If the optional 'family' argument is specified, it means nodes = family(nodes(end)). +% (This is useful since clq_ass_to_node is not accessible to outsiders.) +% Returns c=-1 if there is no such clique. + +if nargin < 3, fam = 0; else fam = 1; end + +if length(nodes)==1 + c = engine.clq_ass_to_node(nodes(1)); +elseif fam + c = engine.clq_ass_to_node(nodes(end)); +else + B = engine.cliques_bitv; + w = engine.clique_weight; + clqs = find(all(B(:,nodes), 2)); % all selected columns must be 1 + if isempty(clqs) + c = -1; + else + c = clqs(argmin(w(clqs))); + end +end diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/enter_evidence.m b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/enter_evidence.m new file mode 100644 index 00000000..4b02fc8f --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/enter_evidence.m @@ -0,0 +1,260 @@ +function [engine, loglik] = enter_evidence(engine, evidence, varargin) +% ENTER_EVIDENCE enter evidence to engine including discrete and continuous evidence +% [engine, ll] = enter_evidence(engine, evidence) +% +% ll is always 0, which is wrong. + +if ~isempty(engine.evidence) + bnet = bnet_from_engine(engine); + engine = stab_cond_gauss_inf_engine(bnet); + engine.evidence = evidence; +else + engine.evidence = evidence; + bnet = bnet_from_engine(engine); +end + +engine.evidence = evidence; +bnet = bnet_from_engine(engine); +ns = bnet.node_sizes(:); +observed = ~isemptycell(evidence); +onodes = find(observed); +hnodes = find(isemptycell(evidence)); +cobs = myintersect(bnet.cnodes, onodes); +dobs = myintersect(bnet.dnodes, onodes); + +engine = incorporate_dis_evidence(engine, dobs, evidence); +l = length(cobs); +for i = 1:l + node = cobs(i); + engine = incorporate_singleconts_evidence(engine, node, evidence); +end +clpot = engine.clpot; + +clq_num = length(engine.cliques); +for n=engine.postorder(1:end-1) + for p=parents(engine.jtree, n) + [margpot, comppot] = complement_pot(clpot{n}, engine.separator{p,n}); + clpot{n} = comppot; + clpot{p} = combine_pots(clpot{p}, margpot); + end +end + +temppot = clpot; +for n=engine.preorder + for c=children(engine.jtree, n) + seppot{n,c} = marginalize_pot(temppot{n}, engine.separator{n,c}); + temppot{c} = direct_combine_pots(temppot{c}, seppot{n,c}); + end +end +engine.clpot = clpot; +engine.seppot = seppot; + +[pot,loglik]=normalize_pot(clpot{engine.root}); + +%%%%%%%%%%%%%%%%%% +function engine = incorporate_dis_evidence(engine, donodes, evidence) +l = length(donodes); +for i=donodes(:)' + node = i; + clqid = engine.clq_ass_to_node(node); + pot = struct(engine.clpot{clqid}); + ns = zeros(1, max(pot.domain)); + ns(pot.ddom) = pot.dsizes; + ns(pot.cheaddom) = pot.cheadsizes; + ns(pot.ctaildom) = pot.ctailsizes; + ddom = pot.ddom; + + potcarray = cell(1, pot.dsize); + for j =1:pot.dsize + tpotc = struct(pot.scgpotc{j}); + potcarray{j} = scgcpot(tpotc.cheadsize, tpotc.ctailsize, 0, tpotc.A, tpotc.B, tpotc.C); + end + + if length(ns(ddom)) == 1 + matrix = pot.scgpotc; + else + matrix = reshape(pot.scgpotc,ns(ddom)); + potcarray = reshape(potcarray, ns(ddom)); + end + + map = find_equiv_posns(node, ddom); + vals = cat(1, evidence{node}); + index = mk_multi_index(length(ddom), map, vals); + potcarray(index{:}) = matrix(index{:}); + potcarray = potcarray(:); + %keyboard; + engine.clpot{clqid} = scgpot(pot.ddom, pot.cheaddom, pot.ctaildom, ns, potcarray); +end + +%%%%%%%%%%%%%%%%%% +function engine = incorporate_singleconts_evidence(engine, node, evidence) +%incorporate_singleconts_evidence incorporate evidence of 1 continuous node +B = engine.cliques_bitv; +clqs_containnode = find(all(B(:,node), 2)); % all selected columns must be 1 +% Every continuous node necessarily apears as head in exactly one clique, +% which is the clique where it appears closest to the strong root. In all other +% clique potentials where it appears, it must be a tail node. +clq_ev_as_head = []; +for i = clqs_containnode(:)' + pot = struct(engine.clpot{i}); + if myismember(node, pot.cheaddom) + clq_ev_as_head = [clq_ev_as_head i]; + break; + end +end + +% If we will incorporate the evidence node which is head of a potential we must rearrange +% the juntion tree by push operation until the tail of the include potential is empty +if ~isempty(clq_ev_as_head) + assert(1 == length(clq_ev_as_head)); + i = clq_ev_as_head; + pot = struct(engine.clpot{i}); + while ~isempty(pot.ctaildom) + [engine, clqtoroot] = push(engine, i, node); + i = clqtoroot; + pot = struct(engine.clpot{i}); + end + B = engine.cliques_bitv; + clqs_containnode = find(all(B(:,node), 2)); +end + +for i = clqs_containnode(:)' + pot = struct(engine.clpot{i}); + if myismember(node, pot.cheaddom) + engine.clpot{i} = incoporate_evidence_headnode(engine.clpot{i}, node, evidence); + else + %assert(myismember(node, pot.ctaildom)); + engine.clpot{i} = incoporate_evidence_tailnode(engine.clpot{i}, node, evidence); + end +end + +%%%%%%%%%%%%%%%%%% +function newscgpot = incoporate_evidence_tailnode(pot, node, evidence) +%ENTER_EVIDENCE_TAILNODE enter the evidence of 1 tailnode of the scgpot +newscgpot = pot; +pot = struct(pot); +%if isempty(pot.ctaildom) +if ~myismember(node, pot.ctaildom) + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + % In this case there is no real dependency of the head nodes % + % on the tail. The potential should be returned unchanged % + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + return; +end +%newscgpot = scgpot([], [], [], []); +assert(myismember(node, pot.ctaildom)); +ni = block(find_equiv_posns(node, pot.ctaildom), pot.ctailsizes); + +ctaildom = mysetdiff(pot.ctaildom, node); +cheaddom = pot.cheaddom; +ddom = pot.ddom; +domain = mysetdiff(pot.domain, node); +dsize = pot.dsize; +ns = zeros(1, max(pot.domain)); +ns(pot.ddom) = pot.dsizes; +ns(pot.cheaddom) = pot.cheadsizes; +ns(pot.ctaildom) = pot.ctailsizes; +cheadsizes = pot.cheadsizes; +cheadsize = pot.cheadsize; +ctailsizes = ns(ctaildom); +ctailsize = sum(ns(ctaildom)); + +potarray = cell(1, dsize); +for i=1:dsize + potc = struct(pot.scgpotc{i}); + B = potc.B; + A = potc.A + B(:, ni)*evidence{node}; + B(:, ni) = []; + potarray{i} = scgcpot(cheadsize, ctailsize, potc.p, A, B, potc.C); +end + +newscgpot = scgpot(ddom, cheaddom, ctaildom, ns, potarray); + +%%%%%%%%%%%%%%%% +function newscgpot = incoporate_evidence_headnode(pot, node, evidence) +%ENTER_EVIDENCE_HEADNODE +pot = struct(pot); +y2 = evidence{node}; +assert(myismember(node, pot.cheaddom)); +assert(isempty(pot.ctaildom)); +ddom = pot.ddom; +cheaddom = mysetdiff(pot.cheaddom, node); +ctaildom = pot.ctaildom; +dsize = pot.dsize; +domain = mysetdiff(pot.domain, node); + +ns = zeros(1, max(pot.domain)); +ns(pot.ddom) = pot.dsizes; +ns(pot.cheaddom) = pot.cheadsizes; +ns(pot.ctaildom) = pot.ctailsizes; +ctailsizes = ns(ctaildom); +ctailsize = sum(ctailsizes); +cheadsizes = ns(cheaddom); +cheadsize = sum(cheadsizes); +onodesize = ns(node); + +p = zeros(1,dsize); +A1 = zeros(cheadsize, dsize); +A2 = zeros(onodesize, dsize); +C11 = zeros(cheadsize, cheadsize, dsize); +C12 = zeros(cheadsize, onodesize, dsize); +C21 = zeros(onodesize, cheadsize, dsize); +C22 = zeros(onodesize, onodesize, dsize); +ZM = zeros(onodesize, onodesize); + +n1i = block(find_equiv_posns(cheaddom, pot.cheaddom), pot.cheadsizes); +n2i = block(find_equiv_posns(node, pot.cheaddom), pot.cheadsizes); + +indic = 0; +for i=1:dsize + potc = struct(pot.scgpotc{i}); + p(i) = potc.p; + if ~isempty(n1i) + A1(:,i) = potc.A(n1i); + end + if ~isempty(n2i) + A2(:,i) = potc.A(n2i); + end + C11(:,:,i) = potc.C(n1i, n1i); + C12(:,:,i) = potc.C(n1i, n2i); + C21(:,:,i) = potc.C(n2i, n1i); + C22(:,:,i) = potc.C(n2i, n2i); + if isequal(0, C22(:,:,i)) & isequal(evidence{node}, A2(:, i)) + indic = i; + end +end + +np = zeros(1,dsize); +nA = zeros(cheadsize, dsize); +nC = zeros(cheadsize, cheadsize, dsize); + +if indic + np(:) = 0; + np(indic) = p(indic); + nA = A1; + nC = C11; +else + for i=1:dsize + if isequal(0, C22(:,:,i)) + p(i) = 0; + nA(:, i) = A1(:, i); + nC(:,:,i) = C11(:,:,i); + else + sq = (y2 - A2(:,i))' * inv(C22(:,:,i)) * (y2 - A2(:,i)); + ex = exp(-0.5*sq); + %np(i) = p(i) * ex / ( (2 * pi)^(-onodesize/2) * sqrt(det(C22(:,:,i))) ); + np(i) = p(i) * ex / ( (2 * pi)^(onodesize/2) * sqrt(det(C22(:,:,i))) ); + nA(:,i) = A1(:,i) + C12(:,:,i) * inv(C22(:,:,i)) * (y2 - A2(:,i)); + tmp1 = C12(:,:,i) * inv(C22(:,:,i)) * C21(:,:,i); + nC(:,:,i) = C11(:,:,i) - tmp1; + end + end +end + +scpot = cell(1, dsize); +W = zeros(cheadsize,ctailsize); +for i=1:dsize + scpot{i} = scgcpot(cheadsize, ctailsize, np(i), nA(:,i), W, nC(:,:,i)); +end +ns(node) = 0; +newscgpot = scgpot(ddom, cheaddom, ctaildom, ns, scpot); diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_difclq_nodes.m b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_difclq_nodes.m new file mode 100644 index 00000000..e1cad6c7 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_difclq_nodes.m @@ -0,0 +1,55 @@ +function marginal = marginal_difclq_nodes(engine, query_nodes) +% MARGINAL_DIFCLQ_NODES get the marginal distribution of nodes which is not in a single clique +% marginal = marginal_difclq_nodes(engine, query_nodes) + +keyboard +num_clique = length(engine.cliques); +B = engine.cliques_bitv; +clqs_containnodes = []; +for i=1:length(query_nodes) + node = query_nodes(i); + tnodes = find(all(B(:, node), 2)); + clqs_containnodes = myunion(clqs_containnodes, tnodes); +end +% get all cliques contains query nodes + +% get the minimal sub tree in junction which contains these cliques and the node closest to the root of jtree +[subtree, nroot_node] = min_subtree_conti_nodes(engine.jtree, engine.root, clqs_containnodes); +if ~mysubset(query_nodes, engine.cliques{nroot_node}); + % if query nodes is not all memers of the clique closest to the root clique performe push operation + engine = push_tree(engine, subtree, query_nodes, nroot_node); +end + +if ~(nroot_node == engine.root) + % if the clique closest to the root clique is not the root clique we must direct combine the + % potential with the potential stored in separator toward to root + p = parents(engine.jtree, nroot_node); + tpot = direct_combine_pots(engine.clpot{nroot_node}, engine.seppot{p, nroot_node}); +else + tpot = engine.clpot{nroot_node}; +end + +pot = marginalize_pot(tpot, query_nodes); +marginal = pot_to_marginal(pot); +marginal.T = normalise(marginal.T); + + + +function engine = push_tree(engine, tree, query_nodes, inode) +% PUSH_TREE recursive perform push opeartion on tree +% engine = push_tree(engine, tree, query_nodes, inode) + +cs = children(tree, inode); +for i = 1:length(cs) + node = cs(i); + push_tree(engine, tree, query_nodes, node); + push_dom = myintersect(engine.cliques{node}, query_nodes); + [engine, clqtoroot] = push(engine, node, push_dom); +end + + + + + + + diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_nodes.m b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_nodes.m new file mode 100644 index 00000000..063c2439 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_nodes.m @@ -0,0 +1,77 @@ +function marginal = marginal_nodes(engine, query, add_ev) +% MARGINAL_NODES Compute the marginal on the specified query nodes (stab_cond_gauss) +% marginal = marginal_nodes(engine, query, add_ev) +% +% 'query' must be a singleton set. +% add_ev is an optional argument; if 1, we will "inflate" the marginal of observed nodes +% to their original size, adding 0s to the positions which contradict the evidence + +if nargin < 3, add_ev = 0; end +if isempty(engine.evidence) + hquery = query; +else + hquery = []; + for i = query + if isempty(engine.evidence{i}) + hquery = [hquery i]; + end + end +end + +bnet = bnet_from_engine(engine); + +nclq = length(engine.cliques); +clique = 0; +for i = 1:nclq + if mysubset(hquery, engine.cliques{i}) + pot = struct(engine.clpot{i}); + %if mysubset(hquery, pot.cheaddom) | mysubset(hquery, pot.ddom) + if mysubset(hquery, pot.domain) + clique = i; + break; + end + end +end + +if isempty(hquery) + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + % If all requested variables are observed, no query is necessary % + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + marginal.mu = []; + marginal.Sigma = []; + marginal.T = 1.0; + marginal.domain = query; +else + if clique == 0 + marginal = marginal_difclq_nodes(engine, hquery); + else + marginal = marginal_singleclq_nodes(engine, clique, hquery); + end + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + % Change the format of output, so that it is identical to the % + % format obtained by the same request for the junction-tree % + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + marginal.domain = query; + bnet = bnet_from_engine(engine); + dquery = myintersect(bnet.dnodes,hquery); + ns = bnet.node_sizes(dquery); + if length(ns) == 0 + marginal.T = 1; + else + if length(ns) == 1 + ns = [1 ns]; + end + marginal.T = reshape(marginal.T,ns); + end +end +if add_ev + bnet = bnet_from_engine(engine); + %marginal = add_ev_to_dmarginal(marginal, engine.evidence, bnet.node_sizes); + marginal = add_evidence_to_gmarginal(marginal, engine.evidence, bnet.node_sizes, bnet.cnodes); +end + + + + + + diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_singleclq_nodes.m b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_singleclq_nodes.m new file mode 100644 index 00000000..d755617f --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_singleclq_nodes.m @@ -0,0 +1,30 @@ +function marginal = marginal_singleclq_nodes(engine, i, query) +% MARGINAL_SINGLECLQ_NODES get the marginal distribution of nodes which is in a single clique +% marginal = marginal_singleclq_nodes(engine, i, query) + +pot = struct(engine.clpot{i}); +if isempty(pot.ctaildom) + if i ~= engine.root + p = parents(engine.jtree, i); + tpot = direct_combine_pots(engine.clpot{i}, engine.seppot{p, i}); + else + tpot = engine.clpot{i}; + end + pot = marginalize_pot(tpot, query); + + marginal = pot_to_marginal(pot); + marginal.T = normalise(marginal.T); +else + [engine, clqtoroot] = push(engine, i, query); + if clqtoroot == engine.root + tpot = engine.clpot{clqtoroot}; + else + p = parents(engine.jtree, clqtoroot); + tpot = direct_combine_pots(engine.clpot{clqtoroot}, engine.seppot{p, clqtoroot}); + end + pot = marginalize_pot(tpot, query); + + marginal = pot_to_marginal(pot); + marginal.T = normalise(marginal.T); +end + diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/problems.txt b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/problems.txt new file mode 100644 index 00000000..fa7c6be8 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/problems.txt @@ -0,0 +1,76 @@ +PROBLEMS WITH STAB_COND_GAUSS_INF_ENGINE + + +- enter_evidence always returns ll=0 + (I set ll=0 since it is not computed) + +- fails on scg_3node, probably because the engine needs to be +re-initialized every time before enter_evidence is called, not just +when the engine is constructed. + +??? Error using ==> assert +assertion violated: + +K>> dbstack +dbstack +> In /home/eecs/murphyk/matlab/BNT/HMM/assert.m at line 9 + In /home/eecs/murphyk/matlab/BNT/examples/static/SCG/scg_3node.m at line 45 + + + +- crashes on scg3 + +Error in ==> /home/eecs/murphyk/matlab/BNT/inference/static/@stab_cond_gauss_inf_engine/stab_cond_gauss_inf_engine.m +On line 77 ==> clpot{cindex} = direct_combine_pots(pot{n}, clpot{cindex}); + +K>> dbstack +dbstack +> In /home/eecs/murphyk/matlab/BNT/inference/static/@stab_cond_gauss_inf_engine/stab_cond_gauss_inf_engine.m at line 77 + In /home/eecs/murphyk/matlab/BNT/examples/static/SCG/scg3.m at line 41 +K>> + + + + + +- fails on scg1 and scg2 + +Warning: One or more output arguments not assigned during call to 'min_subtree_conti_nodes (nearsest_node2)'. +Warning in ==> /home/eecs/murphyk/matlab/BNT/graph/min_subtree_conti_nodes.m (nearsest_node2) +On line 60 ==> nea_node = nearsest_node2(tree, nodes, n); + +K>> dbstack +dbstack +> In /home/eecs/murphyk/matlab/BNT/graph/min_subtree_conti_nodes.m (nearsest_node2) at line 60 + In /home/eecs/murphyk/matlab/BNT/graph/min_subtree_conti_nodes.m (nearest_node) at line 50 + In /home/eecs/murphyk/matlab/BNT/graph/min_subtree_conti_nodes.m at line 11 + In /home/eecs/murphyk/matlab/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_difclq_nodes.m at line 17 + In /home/eecs/murphyk/matlab/BNT/inference/static/@stab_cond_gauss_inf_engine/marginal_nodes.m at line 23 + In /home/eecs/murphyk/matlab/BNT/examples/static/SCG/scg1.m at line 42 + + + + + +- This code fragment, from BNT/graph/min_subtree_conti_nodes, is clearly redundant + +function nea_node = nearest_node(tree, root, nodes) +%get the nearest node to the root in the tree +nea_node = nearsest_node2(tree, nodes, root); + +function nea_node = nearsest_node2(tree, nodes, inode) +if myismember(inode, nodes) + nea_node = inode; + return; +end +cs = children(tree, inode); +for i = 1:length(cs) + n = cs(i); + nea_node = nearsest_node2(tree, nodes, n); +end + + +- Some names are badly chosen. 'nearsest' is a mis-spelling. 'min_subtree_conti_nodes' should be +'min_subtree_containing_nodes' or 'min_subtree_con_nodes'. + +- In general, the code needs some heavy polishing. diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/push.m b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/push.m new file mode 100644 index 00000000..193bf722 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/push.m @@ -0,0 +1,38 @@ +function [engine, clqtoroot] = push(engine, clq, pushdom) +%PUSH_POT push the variables in putshdom which is subset of clq to the clique toword the root and get new engine +%pushdom is pushed variables set +%clq is the index of the clique that pushdom belongs to + +clqdom = engine.cliques{clq}; +assert( mysubset(pushdom, clqdom)); +clqtoroot = parents(engine.jtree, clq); +%sepdom = engine.separator{clq, clqtoroot}; +sepdom = engine.separator{clqtoroot, clq}; + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% Calculate the strong marginal of the union of pushdom and and the separatordomain and % +% the corresponding complement % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +%[margpot, comppot] = complement_pot(engine.clpot{clq}, pushdom); +newsepdom = myunion(pushdom,sepdom); +[margpot,comppot] = complement_pot(engine.clpot{clq}, newsepdom); +engine.clpot{clqtoroot} = direct_combine_pots(engine.clpot{clqtoroot}, margpot); +engine.clpot{clq} = comppot; + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% Calculation of the new separator and separatorpotential of the junction tree % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +engine.seppot{clqtoroot, clq} = direct_combine_pots(engine.seppot{clqtoroot, clq}, margpot); +engine.separator{clqtoroot, clq} = myunion(engine.separator{clqtoroot, clq}, pushdom); + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% Add pushdomain to the clique towards the root % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +engine.cliques{clqtoroot} = myunion(engine.cliques{clqtoroot}, pushdom); + +num_cliques = length(engine.cliques); +B = sparse(num_cliques, 1); +for i=1:num_cliques + B(i, engine.cliques{i}) = 1; +end +engine.cliques_bitv = B; diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/push_pot_toclique.m b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/push_pot_toclique.m new file mode 100644 index 00000000..4bcd0ed0 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/push_pot_toclique.m @@ -0,0 +1,7 @@ +function engine = push_pot_toclique(engine, clqtarget, clq, nodes) +% PUSH_POT push the variables in putshdom which is subset of clq to the target clique toword the root and get new engine +% engine = push_pot_toclique(engine, clqtarget, clq, nodes) +[engine, clqtoroot] = push_pot(engine, clq, nodes) +while clqtoroot ~= clqtarget + [engine, clqtoroot] = push_pot(engine, clqtoroot, nodes) +end \ No newline at end of file diff --git a/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/stab_cond_gauss_inf_engine.m b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/stab_cond_gauss_inf_engine.m new file mode 100644 index 00000000..42c47c6a --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/static/@stab_cond_gauss_inf_engine/stab_cond_gauss_inf_engine.m @@ -0,0 +1,178 @@ +function engine = stab_cond_gauss_inf_engine(bnet) +% STAB_COND_GAUSS_INF_ENGINE Junction tree using stable CG potentials +% engine = cond_gauss_inf_engine(bnet) +% +% This class was written by Shan Huang (shan.huang@intel.com) 2001 +% and fixed by Rainer Deventer deventer@informatik.uni-erlangen.de March 2003 +N = length(bnet.dag); +clusters = {}; +root = N; +stages = { 1:N }; +onodes = []; +engine = init_fields; +engine.evidence = []; +engine = class(engine, 'stab_cond_gauss_inf_engine', inf_engine(bnet)); + +ns = bnet.node_sizes(:); +ns(onodes) = 1; % observed nodes have only 1 possible value + +%[engine.jtree, dummy, engine.cliques, B, w, elim_order, moral_edges, fill_in_edges, strong] = ... +% dag_to_jtree(bnet, onodes, stages, clusters); + + +partial_order = determine_elim_constraints(bnet, onodes); +strong = ~isempty(partial_order); +stages = {}; +clusters = {}; +[engine.jtree, dummy_root, engine.cliques, B, w, elim_order] = + graph_to_jtree(moralize(bnet.dag), ns, partial_order, stages, clusters); + + +engine.cliques_bitv = B; +engine.clique_weight = w; +C = length(engine.cliques); +engine.clpot = cell(1,C); + +% A node can be a member of many cliques, but is assigned to exactly one, to avoid +% double-counting its CPD. We assign node i to clique c if c is the "lightest" clique that +% contains i's family, so it can accomodate its CPD. + +engine.clq_ass_to_node = zeros(1, N); +num_cliques = length(engine.cliques); +for i=1:N + clqs_containing_family = find(all(B(:,family(bnet.dag, i)), 2)); % all selected columns must be 1 + c = clqs_containing_family(argmin(w(clqs_containing_family))); + engine.clq_ass_to_node(i) = c; +end + +% Compute the separators between connected cliques. +[is,js] = find(engine.jtree > 0); +engine.separator = cell(num_cliques, num_cliques); +for k=1:length(is) + i = is(k); j = js(k); + engine.separator{i,j} = find(B(i,:) & B(j,:)); % intersect(cliques{i}, cliques{j}); +end +%keyboard; +engine.seppot = cell(C,C); + +pot_type = 'scg'; +check_for_cd_arcs([], bnet.cnodes, bnet.dag); + +% Make the jtree rooted, so there is a fixed message passing order. +if strong + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + % Start the search for the strong root at the clique with the % + % highest number. % + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + root = length(engine.cliques); + root_found = 0; + + while ((~root_found) & (root >= 1)) + root_found = test_strong_root(engine.jtree,engine.cliques,bnet.dnodes,root); + if ~root_found + root = root - 1; + end + end + assert(root > 0) + engine.root = root; + % the last clique is guaranteed to be a strong root + %engine.root = length(engine.cliques); +else + % jtree_dbn_inf_engine requires the root to contain the interface. + % This may conflict with the strong root requirement! *********** BUG ************* + engine.root = clq_containing_nodes(engine, root); + if engine.root <= 0 + error(['no clique contains ' num2str(root)]); + end +end + +[engine.jtree, engine.preorder, engine.postorder] = mk_rooted_tree(engine.jtree, engine.root); + +% Evaluate CPDs with evidence, and convert to potentials +pot = cell(1, N); +inited = zeros(1, C); +clpot = cell(1, C); +evidence = cell(1, N); +for n=1:N + fam = family(bnet.dag, n); + e = bnet.equiv_class(n); + %pot{n} = CPD_to_scgpot(bnet.CPD{e}, fam, ns, bnet.cnodes, evidence); + pot{n} = convert_to_pot(bnet.CPD{e}, pot_type, fam(:), evidence); + cindex = engine.clq_ass_to_node(n); + if inited(cindex) + clpot{cindex} = direct_combine_pots(pot{n}, clpot{cindex}); + else + clpot{cindex} = pot{n}; + inited(cindex) = 1; + end +end + +for i=1:C + if inited(i) == 0 + clpot{i} = scgpot([], [], [], []); + end +end + +seppot = cell(C, C); +% separators are is not need to initialize + +% collect to root (node to parents) +% Unlike the HUGIN architecture the complements are stored in the cliques during COLLECT +% and the separators are not playing a specific role during this process +for n=engine.postorder(1:end-1) + for p=parents(engine.jtree, n) + if ~isempty(engine.separator{p,n}) + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + % The empty case might happen for unlinked nodes, i.e. the DAG is not % + % a single tree, but a forest % + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + [margpot, comppot] = complement_pot(clpot{n}, engine.separator{p,n}); + clpot{n} = comppot; + clpot{p} = combine_pots(clpot{p}, margpot); + end + end +end + +% distribute message from root +% We have not to store the weak clique marginals and keep the original complement potentials. +% This is a minor variation of HUGIN architecture. +temppot = clpot; +for n=engine.preorder + for c=children(engine.jtree, n) + seppot{n,c} = marginalize_pot(temppot{n}, engine.separator{n,c}); + temppot{c} = direct_combine_pots(temppot{c}, seppot{n,c}); + end +end + +engine.clpot = clpot; +engine.seppot = seppot; + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% init_fields() % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +function engine = init_fields() + +engine.evidence = []; +engine.jtree = []; +engine.cliques = []; +engine.cliques_bitv = []; +engine.clique_weight = []; +engine.preorder = []; +engine.postorder = []; +engine.root = []; +engine.clq_ass_to_node = []; +engine.separator = []; +engine.clpot =[]; +engine.seppot = []; + + + + + + + + + + + + |