Daniel@0: function [bnet, vals] = mk_minimal_qmr_bnet(G, inhibit, leak, prior, pos, neg, pos_only)
Daniel@0: % MK_MINIMAL_QMR_BNET Make a QMR model which only contains the observed findings
Daniel@0: % [bnet, vals] = mk_minimal_qmr_bnet(G, inhibit, prior, leak, pos, neg)
Daniel@0: %
Daniel@0: % Input:
Daniel@0: % G(i,j) = 1 iff there is an arc from disease i to finding j
Daniel@0: % inhibit(i,j) = inhibition probability on i->j arc
Daniel@0: % leak(j) = inhibition prob. on leak->j arc
Daniel@0: % prior(i) = prob. disease i is on
Daniel@0: % pos = list of leaves that have positive observations
Daniel@0: % neg = list of leaves that have negative observations
Daniel@0: % pos_only = 1 means only include positively observed leaves in the model - the negative
Daniel@0: %   ones are absorbed into the prior terms
Daniel@0: %
Daniel@0: % Output:
Daniel@0: % bnet
Daniel@0: % vals is their value
Daniel@0: 
Daniel@0: if pos_only
Daniel@0:   obs = pos;
Daniel@0: else
Daniel@0:   obs = myunion(pos, neg);
Daniel@0: end
Daniel@0: Nfindings = length(obs);
Daniel@0: [Ndiseases maxNfindings] = size(inhibit);
Daniel@0: N = Ndiseases + Nfindings;
Daniel@0: finding_node = Ndiseases+1:N;
Daniel@0: 
Daniel@0: % j = finding_node(i) means the i'th finding node is the j'th node in the bnet
Daniel@0: % k = obs(i) means the i'th observed (positive) finding is the k'th finding overall
Daniel@0: % If all findings are observed, and posonly = 0, we have i = obs(i) for all i.
Daniel@0: 
Daniel@0: %dag = sparse(N, N);
Daniel@0: dag = zeros(N, N);
Daniel@0: dag(1:Ndiseases, Ndiseases+1:N) = G(:,obs);
Daniel@0: 
Daniel@0: ns = 2*ones(1,N);
Daniel@0: bnet = mk_bnet(dag, ns, 'observed', finding_node);
Daniel@0: 
Daniel@0: CPT = cell(1, Ndiseases);
Daniel@0: for d=1:Ndiseases
Daniel@0:   CPT{d} = [1-prior(d) prior(d)];
Daniel@0: end
Daniel@0: 
Daniel@0: if pos_only
Daniel@0:   % Fold in the negative evidence into the prior
Daniel@0:   for i=1:length(neg)
Daniel@0:     n = neg(i);
Daniel@0:     ps = parents(G,n);
Daniel@0:     for pi=1:length(ps)
Daniel@0:       p = ps(pi);
Daniel@0:       q = inhibit(p,n);
Daniel@0:       CPT{p} = CPT{p} .* [1 q];
Daniel@0:     end
Daniel@0:     % Arbitrarily attach the leak term to the first parent
Daniel@0:     p = ps(1);
Daniel@0:     q = leak(n);
Daniel@0:     CPT{p} = CPT{p} .* [q q];
Daniel@0:   end
Daniel@0: end
Daniel@0: 
Daniel@0: for d=1:Ndiseases
Daniel@0:   bnet.CPD{d} = tabular_CPD(bnet, d, CPT{d}');
Daniel@0: end
Daniel@0: 
Daniel@0: for i=1:Nfindings
Daniel@0:   fnode = finding_node(i);
Daniel@0:   fid = obs(i);
Daniel@0:   ps = parents(G, fid);
Daniel@0:   bnet.CPD{fnode} = noisyor_CPD(bnet, fnode, leak(fid), inhibit(ps, fid));
Daniel@0: end
Daniel@0: 
Daniel@0: obs_nodes = finding_node;
Daniel@0: vals = sparse(1, maxNfindings);
Daniel@0: vals(pos) = 2;
Daniel@0: vals(neg) = 1;
Daniel@0: vals = full(vals(obs));
Daniel@0: 
Daniel@0: 
Daniel@0: 
Daniel@0: 
Daniel@0: