From 8070dc963753142bb86c4ed698d91fd623ed28e7 Mon Sep 17 00:00:00 2001 From: ziejd2 Date: Thu, 28 Sep 2017 15:04:40 -0500 Subject: BNW using Octave instead of Matlab. This version of BNW should perform the same as the original version. The only difference is that it uses Octave instead of Matlab when running BayesNet Toolbox during parameter learning. I am calling this BNW_1.02. It can be accessed at: compbio.uthsc.edu/BNW_1.02 --- .../dynamic/@frontier_inf_engine/CVS/Entries | 7 + .../dynamic/@frontier_inf_engine/CVS/Repository | 1 + .../dynamic/@frontier_inf_engine/CVS/Root | 1 + .../dynamic/@frontier_inf_engine/enter_evidence.m | 44 +++++++ .../@frontier_inf_engine/enter_soft_evidence.m | 142 +++++++++++++++++++++ .../@frontier_inf_engine/frontier_inf_engine.m | 121 ++++++++++++++++++ .../dynamic/@frontier_inf_engine/marginal_family.m | 7 + .../dynamic/@frontier_inf_engine/marginal_nodes.m | 21 +++ .../dynamic/@frontier_inf_engine/set_fwdback.m | 8 ++ 9 files changed, 352 insertions(+) create mode 100644 sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Entries create mode 100644 sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Repository create mode 100644 sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Root create mode 100644 sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/enter_evidence.m create mode 100644 sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/enter_soft_evidence.m create mode 100644 sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/frontier_inf_engine.m create mode 100644 sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/marginal_family.m create mode 100644 sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/marginal_nodes.m create mode 100644 sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/set_fwdback.m (limited to 'sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine') diff --git a/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Entries b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Entries new file mode 100644 index 00000000..79297e05 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Entries @@ -0,0 +1,7 @@ +/enter_evidence.m/1.1.1.1/Wed May 29 15:59:56 2002// +/enter_soft_evidence.m/1.1.1.1/Wed May 29 15:59:56 2002// +/frontier_inf_engine.m/1.1.1.1/Wed May 29 15:59:56 2002// +/marginal_family.m/1.1.1.1/Wed May 29 15:59:56 2002// +/marginal_nodes.m/1.1.1.1/Wed May 29 15:59:56 2002// +/set_fwdback.m/1.1.1.1/Wed May 29 15:59:56 2002// +D diff --git a/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Repository b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Repository new file mode 100644 index 00000000..0e85f66f --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/inference/dynamic/@frontier_inf_engine diff --git a/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Root b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/enter_evidence.m b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/enter_evidence.m new file mode 100644 index 00000000..bd30a57c --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/enter_evidence.m @@ -0,0 +1,44 @@ +function [engine, loglik] = enter_evidence(engine, evidence, varargin) +% ENTER_EVIDENCE Add the specified evidence to the network (frontier) +% [engine, loglik] = enter_evidence(engine, evidence, ...) +% +% evidence{i,t} = [] if if X(i,t) is hidden, and otherwise contains its observed value (scalar or column vector) +% +% The following optional arguments can be specified in the form of name/value pairs: +% [default value in brackets] +% +% maximize - if 1, does max-product (not yet supported), else sum-product [0] +% filter - if 1, do filtering, else smoothing [0] +% +% e.g., engine = enter_evidence(engine, ev, 'maximize', 1) + +maximize = 0; +filter = 0; + +% parse optional params +args = varargin; +nargs = length(args); +if nargs > 0 + for i=1:2:nargs + switch args{i}, + case 'maximize', maximize = args{i+1}; + case 'filter', filter = args{i+1}; + otherwise, + error(['invalid argument name ' args{i}]); + end + end +end + +assert(~maximize); + +[ss T] = size(evidence); +bnet = bnet_from_engine(engine); +onodes = find(~isemptycell(evidence)); +cnodes = unroll_set(bnet.cnodes(:), ss, T); +pot_type = determine_pot_type(bnet, onodes); + +CPDpot = convert_dbn_CPDs_to_pots(bnet, evidence, pot_type); + +[engine.fwdback, loglik, engine.fwd_frontier, engine.back_frontier] = ... + enter_soft_evidence(engine, CPDpot, onodes, pot_type, filter); + diff --git a/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/enter_soft_evidence.m b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/enter_soft_evidence.m new file mode 100644 index 00000000..8da339c7 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/enter_soft_evidence.m @@ -0,0 +1,142 @@ +function [fwdback, loglik, fwd_frontier, back_frontier] = enter_soft_evidence(engine, CPD, onodes, pot_type, filter) +% ENTER_SOFT_EVIDENCE Add soft evidence to network (frontier) +% [fwdback, loglik] = enter_soft_evidence(engine, CPDpot, onodes, filter) + +if nargin < 3, filter = 0; end + +[ss T] = size(CPD); +bnet = bnet_from_engine(engine); +ns = repmat(bnet.node_sizes_slice(:), 1, T); +cnodes = unroll_set(bnet.cnodes(:), ss, T); + +% FORWARDS +fwd = cell(ss,T); +ll = zeros(1,T); +S = 2*ss; % num. intermediate frontiers to get from t to t+1 +frontier = cell(S,T); + +% Start with empty frontier, and add each node in slice 1 +init = mk_initial_pot(pot_type, [], ns, cnodes, onodes); +t = 1; +s = 1; +j = 1; +frontier{s,t} = update(init, j, 1, CPD{j}, engine.fdom1{s}, pot_type, ns, cnodes, onodes); +fwd{j} = frontier{s,t}; +for s=2:ss + j = s; % add node j at step s + frontier{s,t} = update(frontier{s-1,t}, j, 1, CPD{j}, engine.fdom1{s}, pot_type, ns, cnodes, onodes); + fwd{j} = frontier{s,t}; +end +frontier{S,t} = frontier{ss,t}; +[frontier{S,t}, ll(1)] = normalize_pot(frontier{S,t}); + +% Now move frontier from slice to slice +OPS = engine.ops; +add = OPS>0; +nodes = [zeros(S,1) unroll_set(abs(OPS(:)), ss, T-1)]; +for t=2:T + offset = (t-2)*ss; + for s=1:S + if s==1 + prev_ndx = (t-2)*S + S; % S,t-1 + else + prev_ndx = (t-1)*S + s-1; % s-1,t + end + j = nodes(s,t); + frontier{s,t} = update(frontier{prev_ndx}, j, add(s), CPD{j}, engine.fdom{s}+offset, pot_type, ns, cnodes, onodes); + if add(s) + fwd{j} = frontier{s,t}; + end + end + [frontier{S,t}, ll(t)] = normalize_pot(frontier{S,t}); +end +loglik = sum(ll); + + +fwd_frontier = frontier; + +if filter + fwdback = fwd; + return; +end + + +% BACKWARDS +back = cell(ss,T); +add = ~add; % forwards add = backwards remove +frontier = cell(S,T+1); +t = T; +dom = (1:ss) + (t-1)*ss; +frontier{1,T+1} = mk_initial_pot(pot_type, dom, ns, cnodes, onodes); % all 1s for last slice +for t=T:-1:2 + offset = (t-2)*ss; + for s=S:-1:1 % reverse order + if s==S + prev_ndx = t*S + 1; % 1,t+1 + else + prev_ndx = (t-1)*S + (s+1); % s+1,t + end + j = nodes(s,t); + if ~add(s) + back{j} = frontier{prev_ndx}; % save frontier before removing + end + frontier{s,t} = rev_update(frontier{prev_ndx}, t, s, j, add(s), CPD{j}, engine.fdom{s}+offset, pot_type, ns, cnodes, onodes); + end + frontier{1,t} = normalize_pot(frontier{1,t}); +end +% Remove each node in first slice until left with empty set +t = 1; +frontier{ss+1,t} = frontier{1,2}; +add = 0; +for s=ss:-1:1 + j = s; % remove node j at step s + back{j} = frontier{s+1,t}; + frontier{s,t} = rev_update(frontier{s+1,t}, t, s, j, add, CPD{j}, 1:s, pot_type, ns, cnodes, onodes); +end + +% COMBINE +for t=1:T + for i=1:ss + %fwd{i,t} = multiply_by_pot(fwd{i,t}, back{i,t}); + %fwdback{i,t} = normalize_pot(fwd{i,t}); + fwdback{i,t} = normalize_pot(multiply_pots(fwd{i,t}, back{i,t})); + end +end + +back_frontier = frontier; + +%%%%%%%%%% +function new_frontier = update(old_frontier, j, add, CPD, newdom, pot_type, ns, cnodes, onodes) + +if add + new_frontier = mk_initial_pot(pot_type, newdom, ns, cnodes, onodes); + new_frontier = multiply_by_pot(new_frontier, old_frontier); + new_frontier = multiply_by_pot(new_frontier, CPD); +else + new_frontier = marginalize_pot(old_frontier, mysetdiff(domain_pot(old_frontier), j)); +end + + +%%%%%% +function new_frontier = rev_update(old_frontier, t, s, j, add, CPD, junk, pot_type, ns, cnodes, onodes) + +olddom = domain_pot(old_frontier); +assert(isequal(junk, olddom)); + +if add + % add: extend domain to include j by multiplying by 1 + newdom = myunion(olddom, j); + new_frontier = mk_initial_pot(pot_type, newdom, ns, cnodes, onodes); + new_frontier = multiply_by_pot(new_frontier, old_frontier); + %fprintf('t=%d, s=%d, add %d to %s to make %s\n', t, s, j, num2str(olddom), num2str(newdom)); +else + % remove: multiply in CPT and then marginalize out j + % parents of j are guaranteed to be in old_frontier, else couldn't have added j on fwds pass + old_frontier = multiply_by_pot(old_frontier, CPD); + newdom = mysetdiff(olddom, j); + new_frontier = marginalize_pot(old_frontier, newdom); + %newdom2 = domain_pot(new_frontier); + %fprintf('t=%d, s=%d, rem %d from %s to make %s\n', t, s, j, num2str(olddom), num2str(newdom2)); +end + + diff --git a/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/frontier_inf_engine.m b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/frontier_inf_engine.m new file mode 100644 index 00000000..fd550579 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/frontier_inf_engine.m @@ -0,0 +1,121 @@ +function engine = frontier_inf_engine(bnet) +% FRONTIER_INF_ENGINE Inference engine for DBNs which which uses the frontier algorithm. +% engine = frontier_inf_engine(bnet) +% +% The frontier algorithm extends the forwards-backwards algorithm to DBNs in the obvious way, +% maintaining a joint distribution (frontier) over all the nodes in a time slice. +% When all the hidden nodes in the DBN are persistent (have children in the next time slice), +% its theoretical running time is often similar to that of the junction tree algorithm, +% although in practice, this algorithm seems to very slow (at least in matlab). +% However, it is extremely simple to describe and implement. +% +% Suppose there are n binary nodes per slice, so the frontier takes O(2^n) space. +% Each time step takes between O(n 2^{n+1}) and O(n 2^{2n}) operations, depending on the graph structure. +% The lower bound is achieved by a set of n independent chains, as in a factorial HMM. +% The upper bound is achieved by a set of n fully interconnected chains, as in an HMM. +% +% The factor of n arises because we need to multiply in each CPD from slice t+1. +% The second factor depends on the size of the frontier to which we add the new node. +% In an FHMM, once we have added X(i,t+1), we can marginalize out X(i,t) from the frontier, since +% no other nodes depend on it; hence the frontier never contains more than n+1 nodes. +% In a fully coupled HMM, we must leave X(i,t) in the frontier until all X(j,t+1) have been +% added; hence the frontier will contain 2*n nodes at its peak. +% +% For details, see +% "The Factored Frontier Algorithm for Approximate Inference in DBNs", +% Kevin Murphy and Yair Weiss, UAI 01. + +ns = bnet.node_sizes_slice; +onodes = bnet.observed; +ns(onodes) = 1; +ss = length(bnet.intra); + +[engine.ops, engine.fdom] = best_first_frontier_seq(ns, bnet.dag); +engine.ops1 = 1:ss; + +engine.fwdback = []; +engine.fwd_frontier = []; +engine.back_frontier = []; + +engine.fdom1 = cell(1,ss); +for s=1:ss + engine.fdom1{s} = 1:s; +end + +engine = class(engine, 'frontier_inf_engine', inf_engine(bnet)); + + +%%%%%%%%% + +function [ops, frontier_set] = best_first_frontier_seq(ns, dag) +% BEST_FIRST_FRONTIER_SEQ Do a greedy search for the sequence of additions/removals to the frontier. +% [ops, frontier_set] = best_first_frontier_seq(ns, dag) +% +% We maintain 3 sets: the frontier (F), the right set (R), and the left set (L). +% The invariant is that the nodes in R are d-separated from L given F. +% We start with slice 1 in F and slice 2 in R. +% The goal is to move slice 1 from F to L, and slice 2 from R to F, so as to minimize the size +% of the frontier at each step, where the size(F) = product of the node-sizes of nodes in F. +% A node may be removed (from F to L) if it has no children in R. +% A node may be added (from R to F) if its parents are in F. +% +% ns(i) = num. discrete values node i can take on (i=1..ss, where ss = slice size) +% dag is the (2*ss) x (2*ss) adjacency matrix for the 2-slice DBN. + +% Example: +% +% 4 9 +% ^ ^ +% | | +% 2 -> 7 +% ^ ^ +% | | +% 1 -> 6 +% | | +% v v +% 3 -> 8 +% | | +% v V +% 5 10 +% +% ops = -4, -5, 6, -1, 7, -2, 8, -3, 9, 10 + +ss = length(ns); +ns = [ns(:)' ns(:)']; +ops = zeros(1,ss); +L = []; F = 1:ss; R = (1:ss)+ss; +frontier_set = cell(1,2*ss); +for s=1:2*ss + remcost = inf*ones(1,2*ss); + %disp(['L: ' num2str(L) ', F: ' num2str(F) ', R: ' num2str(R)]); + maybe_removable = myintersect(F, 1:ss); + for n=maybe_removable(:)' + cs = children(dag, n); + if isempty(myintersect(cs, R)) + remcost(n) = prod(ns(mysetdiff(F, n))); + end + end + %remcost + if any(remcost < inf) + n = argmin(remcost); + ops(s) = -n; + L = myunion(L, n); + F = mysetdiff(F, n); + else + addcost = inf*ones(1,2*ss); + for n=R(:)' + ps = parents(dag, n); + if mysubset(ps, F) + addcost(n) = prod(ns(myunion(F, [ps n]))); + end + end + %addcost + assert(any(addcost < inf)); + n = argmin(addcost); + ops(s) = n; + R = mysetdiff(R, n); + F = myunion(F, n); + end + %fprintf('op at step %d = %d\n\n', s, ops(s)); + frontier_set{s} = F; +end diff --git a/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/marginal_family.m b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/marginal_family.m new file mode 100644 index 00000000..4d28263b --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/marginal_family.m @@ -0,0 +1,7 @@ +function marginal = marginal_family(engine, i, t) +% MARGINAL_FAMILY Compute the marginal on node i in slice t and its parents (frontier) +% marginal = marginal_family(engine, i, t) + +bnet = bnet_from_engine(engine); +fam = family(bnet.dag, i, t); +marginal = pot_to_marginal(normalize_pot(marginalize_pot(engine.fwdback{i,t}, fam))); diff --git a/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/marginal_nodes.m b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/marginal_nodes.m new file mode 100644 index 00000000..898d1130 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/marginal_nodes.m @@ -0,0 +1,21 @@ +function marginal = marginal_nodes(engine, nodes, t) +% MARGINAL_NODES Compute the marginal on the specified query nodes (frontier) +% marginal = marginal_nodes(engine, nodes, t) +% +% 't' specifies the time slice of the earliest node in 'nodes'. +% 'nodes' cannot span more than 2 time slices. +% +% Example: +% Consider a DBN with 2 nodes per slice. +% Then t=2, nodes=[1 3] refers to node 1 in slice 2 and node 1 in slice 3, +% i.e., nodes 3 and 5 in the unrolled network, + +if nargin < 3, t = 1; end +assert(length(nodes)==1); +i = nodes(1); +bigpot = engine.fwdback{i,t}; +bnet = bnet_from_engine(engine); +ss = length(bnet.intra); +nodes = nodes + (t-1)*ss; +%if t > 1, nodes = nodes + ss; end +marginal = pot_to_marginal(marginalize_pot(bigpot, nodes)); diff --git a/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/set_fwdback.m b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/set_fwdback.m new file mode 100644 index 00000000..6752d827 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/inference/dynamic/@frontier_inf_engine/set_fwdback.m @@ -0,0 +1,8 @@ +function engine = set_fwdback(engine, fb) +% SET_FWDBACK Set the field 'fwdback', which contains the frontiers after propagation +% engine = set_fwdback(engine, fb) +% +% This is used by frontier_fast_inf_engine/enter_evidence +% as a workaround for Matlab's annoying privacy control + +engine.fwdback = fb; -- cgit 1.4.1