enter_evidence2.m

来自「Bayes网络工具箱」· M 代码 · 共 152 行

function [engine, loglik] = enter_evidence(engine, evidence, filter)
% ENTER_EVIDENCE Add the specified evidence to the network (frontier)
% [engine, loglik] = enter_evidence(engine, evidence, filter)
%
% evidence{i,t} = [] if if X(i,t) is hidden, and otherwise contains its observed value (scalar or column vector)
% If filter = 1, we do filtering, otherwise smoothing (default).

if nargin < 3, filter = 0; end

[ss T] = size(evidence);
bnet = bnet_from_engine(engine);
ns = repmat(bnet.node_sizes_slice(:), 1, T);
eclass = [bnet.equiv_class(:,1) repmat(bnet.equiv_class(:,2), 1, T-1)];
onodes = find(~isemptycell(evidence));
cnodes = unroll_set(bnet.cnodes(:), ss, T);
big_dag = unroll_dbn_topology(bnet.intra1, bnet.intra, bnet.inter, T);
pot_type = determine_pot_type(onodes, cnodes, big_dag);

ns = ns(:)';
cnodes = cnodes(:)';

% Convert evidence-specific CPDs to potentials
CPD = cell(ss,T);
for t=1:T
  for i=1:ss
    fam = family(bnet.dag, i, t);
    CPD{i,t} = CPD_to_pot(pot_type, bnet.CPD{eclass(i,t)}, fam, ns, cnodes, evidence);   
  end
end


% FORWARDS
OPS = engine.ops;
fwd = cell(ss,T);
ll = zeros(1,T);
S = 2*ss; 


frontier = cell(S,T);
% start with empty frontier, and add each node one at a time
t = 1;
s = 1;
frontier{s,t} = mk_initial_pot(pot_type, [], ns, cnodes, onodes);  
for s=1:ss
  j = s; % add node j at step s
  frontier{s+1,t} = update(frontier{s,t}, j, 1, CPD{j});
  fwd{j} = frontier{s+1,t};
end
frontier{S,t} = frontier{ss+1,t};
[frontier{S,t}, ll(1)] = normalize_pot(frontier{S,t});

add = OPS>0;
nodes = [zeros(S,1) unroll_set(abs(OPS(:)), ss, T-1)];
for t=2:T
  for s=1:S
    j = nodes(s,t);
    if s==1
      %frontier{s,t} = update(frontier{S,t-1}, j, add(s), CPD{j});        
      prev_ndx = (t-2)*S + S; % S,t-1
    else
      %frontier{s,t} = update(frontier{s-1,t}, j, add(s), CPD{j});          
      prev_ndx = (t-1)*S + s-1; % s-1,t
    end
    frontier{s,t} = update(frontier{prev_ndx}, j, add(s), CPD{j});
    if add(s)
      fwd{j} = frontier{s,t};
    end
  end
  [frontier{S,t}, ll(t)] = normalize_pot(frontier{S,t});
end


loglik = sum(ll);


if filter
  engine.fwdback = fwd;
  return;
end


% BACKWARDS
back = cell(ss,T);
OPS = -OPS(end:-1:1);
OPS1 = 1:ss;   % add all nodes in topological order
OPS1 = -OPS1(end:-1:1);
clear frontier;

if 0
add = OPS>0;
nodes = [zeros(S,1) unroll_set(abs(OPS(:)), ss, T-1)];
frontier = cell(T+1,S);
% start with frontier of all 1s
t = T;
dom = (1:ss) + (t-1)*ss;
frontier{T+1,S} = mk_initial_pot(pot_type, dom, ns, cnodes, onodes);
end

t = T;
% start with frontier of all 1s
fc = 1;
dom = (1:ss) + (t-1)*ss;
frontier{fc} = mk_initial_pot(pot_type, dom, ns, cnodes, onodes);  

for t=T:-1:1
  offset = max(0,t-2);
  if t==1
    ops = OPS1;
  else
    ops = OPS;
  end           
  for s=1:length(ops)
    i = ops(s);
    if i > 0
      j = i + offset*ss;
      % add: extend domain to include i by multiplying by 1
      pot = mk_initial_pot(pot_type, j, ns, cnodes, onodes);      
      frontier{fc+1} = multiply_pots(frontier{fc}, pot);
    else 
      % remove: multiply in CPT and then marginalize it out
      j = -i + offset*ss;
      back{j} = frontier{fc};
      frontier{fc+1} = multiply_pots(frontier{fc}, CPD{j});
      frontier{fc+1} = marginalize_pot(frontier{fc+1}, mysetdiff(domain_pot(frontier{fc+1}), j));
    end
    fc = fc + 1;
  end
  frontier{fc} = normalize_pot(frontier{fc});
end
% end with empty frontier


% COMBINE
for t=1:T
  for i=1:ss
   engine.fwdback{i,t} = normalize_pot(multiply_pots(fwd{i,t}, back{i,t}));
  end
end




%%%%%%%%%%
function new_frontier = update(old_frontier, j, add, CPD)

if add
  new_frontier = multiply_pots(old_frontier, CPD);
else
  new_frontier = marginalize_pot(old_frontier, mysetdiff(domain_pot(old_frontier), j));    
end