Mercurial > hg > camir-aes2014
diff toolboxes/FullBNT-1.0.7/bnt/inference/static/@belprop_inf_engine/belprop_inf_engine.m @ 0:e9a9cd732c1e tip
first hg version after svn
author | wolffd |
---|---|
date | Tue, 10 Feb 2015 15:05:51 +0000 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/toolboxes/FullBNT-1.0.7/bnt/inference/static/@belprop_inf_engine/belprop_inf_engine.m Tue Feb 10 15:05:51 2015 +0000 @@ -0,0 +1,90 @@ +function engine = belprop_inf_engine(bnet, varargin) +% BELPROP_INF_ENGINE Make a loopy belief propagation inference engine +% engine = belprop_inf_engine(bnet, ...) +% +% This is like pearl_inf_engine, except it uses potential objects, +% instead of lambda/pi structs. Hence it is slower. +% +% The following optional arguments can be specified in the form of name/value pairs: +% [default in brackets] +% +% protocol - 'tree' means send messages up then down the tree, +% 'parallel' means use synchronous updates ['parallel'] +% max_iter - max. num. iterations [ 2*num_nodes ] +% momentum - weight assigned to old message in convex combination (useful for damping oscillations) [0] +% tol - tolerance used to assess convergence [1e-3] +% maximize - 1 means use max-product, 0 means use sum-product [0] +% filename - name of file to write beliefs to after each iteration within enter_evidence [ [] ] +% +% e.g., engine = belprop_inf_engine(bnet, 'maximize', 1, 'max_iter', 10) + +% gdl = general distributive law +engine.gdl = bnet_to_gdl(bnet); + +% set default params +N = length(engine.gdl.G); +engine.protocol = 'parallel'; +engine.max_iter = 2*N; +engine.momentum = 0; +engine.tol = 1e-3; +engine.maximize = 0; +engine.filename = []; +engine.fid = []; + +args = varargin; +nargs = length(args); +for i=1:2:nargs + switch args{i}, + case 'max_iter', engine.max_iter = args{i+1}; + case 'momentum', engine.momentum = args{i+1}; + case 'tol', engine.tol = args{i+1}; + case 'protocol', engine.protocol = args{i+1}; + case 'filename', engine.filename = args{i+1}; + otherwise, + error(['invalid argument name ' args{i}]); + end +end + + +if strcmp(engine.protocol, 'tree') + % Make a rooted tree, so there is a fixed message passing order. + root = N; + [engine.tree, engine.preorder, engine.postorder, height, cyclic] = mk_rooted_tree(engine.gdl.G, root); + assert(~cyclic); +end + +% store results computed by enter_evidence here +engine.marginal_domains = cell(1, N); + +engine.niter = []; + +engine = class(engine, 'belprop_inf_engine', inf_engine(bnet)); + +%%%%%%%%% + +function gdl = bnet_to_gdl(bnet) + +gdl.G = mk_undirected(bnet.dag); +N = length(bnet.dag); +gdl.doms = cell(1,N); +for i=1:N + gdl.doms{i} = family(bnet.dag, i); +end + +% Compute a bit vector representation of the set of domains +% dom_bitv(i,j) = 1 iff variable j occurs in domain i +gdl.dom_bitv = zeros(N, N); +for i=1:N + gdl.dom_bitv(i, gdl.doms{i}) = 1; +end + +% compute the interesection of the domains on either side of each edge (separating set) +gdl.sepset = cell(N, N); +gdl.nbrs = cell(1,N); +for i=1:N + nbrs = neighbors(gdl.G, i); + gdl.nbrs{i} = nbrs; + for j = nbrs(:)' + gdl.sepset{i,j} = myintersect(gdl.doms{i}, gdl.doms{j}); + end +end