function bnode = maximize_params_noisyfn(bnode)
% MAXIMIZE_PARAMS_NOISYFN Find the ML parameters of a noisy-fn node.
% bnode = maximize_params_noisyfn(bnode)
%
% We follow the technique discussed in
% "Learning causal interaction models", Meek and Heckerman, UAI 97.
% See the README file for details.

% The CPT of node Xi looks like this, where p(j) = pfail(i,j),
% and . means 1 - the sum of the other cols (child vals)
% C\X  1     2   3
%  1   .   p(2) p(3)
%  2  p(1) .    p(3)
%  3  p(1) p(2) .
% Hence we can estimate pfail(i,j) by finding a p(j) in this table.
% For non-binary nodes, this is not optimal, since there is more than one such p(j):
% we should have pooled the stats of each copy, since they are tied.
% If this is statistical inefficiency is a problem, you should create a noisy_node type to model
% a noisy channel, and use it instead of tabular nodes for the Xi's.
% (A noisy_node behaves just like a noisyfn_node where F is the identity fn.)

% We treat clamped (non-adjustable) parameters in a rather cheesy way:
% Inside the EM routine, they are unconstrained, but afterwards, we set them to
% their original values.

bnet = mk_bnet_from_noisyfn(bnode.pfail, bnode.fn, bnode.self_size);
np = length(bnode.parents); % num. parents
Cs = 1:np;
Xs = np+1:2*np;
E = 2*np+1;

soft_nodes = num2cell(repmat([Cs E]', 1, bnode.num_samples), 1);   % each row gets its own cell
%soft_evidence = num2cell([bnode.parent_distrib bnode.self_distrib], 1);
N = bnode.num_samples;
soft_evidence = cell(1, N);
for l=1:N
  ev = cell(1,np+1);
  for i=1:np
    dist = bnode.parent_distrib(l, i, :);
    ev{i} = dist(:)';
  end
  dist = bnode.self_distrib(l, :);
  ev{end} = dist(:)';
  soft_evidence{l} = ev;
end
ev = struct('soft_nodes', soft_nodes, 'soft_vals', soft_evidence);

%fprintf('about to do M step for noisyfn node %d\n', bnode.self);
bnet = learn_params(bnet, ev, bnode.options);
%fprintf('done\n');

for i=1:np
  pf = zeros(1, bnode.parent_size);
  pf(1) = bnet.bnodes{Xs(i)}.CPT(2,1);
  pf(2:end) = bnet.bnodes{Xs(i)}.CPT(1,2:end);
  mask = logical(bnode.adjustable_mask(i,:));
  bnode.pfail(i,mask) = pf(mask);
end

