cs_score.m.svn-base

bayesian network structrue learning matlab program

function [PDG, PD2G, PD2Gtheta, PDGtheta] = cs(data, dag, ns, data_mat)
%  [PDG, PD2G, PD2Gtheta, PDGtheta] = cs(data, dag, node_sizes, data_mat)
%
% Compute ln P(D|G) = ln P(D'|G) - ln P(D'|G,^Theta) + ln P(D|G,^Theta)
%
% where D' is a completion of the cell array data
%  that respect the essential statistics in data
%
% INPUTS :
%   - data{i,m},    node i in case m (cell array)
%   - dag,          adjacent graph to score
%   - node_sizes,   vector containing sizes of nodes
%   - data_mat,     a completion of data that respect the essential statistics in data (optionnal)
%               --> not to recompute this all of time
%
% OUTPUTS :
%   - PDG       = ln P(D|G)         --> MAP approximation of Cheesement and Stuntz formula
%   - PD2G      = ln P(D'|G)        --> compute using Heckerman
%   - PD2Gtheta = ln P(D'|G,^theta) --> sum of logs of thetas
%   - PDGtheta  = ln P(D|G,^theta)  --> sum of logs of sum of marginals
% 
% for discretes variables
%

[N, m]=size(data);
verbose = 1;
max_iter = 8;
thresh = 2e-3;
%test = 'LRT'; % 'pearson'

mm=0;
PM=1;

if nargin < 3 | mm==1,
  misv = -9999;
  data_m = bnt_to_mat(data,misv);
  ns=max(data_m');
  %clear data_m misv
  long = size(data_m,2);
end

% -> apprentissage des params

bnet = mk_bnet(dag, ns);
for i = 1:N
  bnet.CPD{i}=tabular_CPD(bnet, i, 'prior_type', 'dirichlet', 'dirichlet_type', 'unif');
end
engine = jtree_inf_engine(bnet);
[bnet, LLtrace, engine] = learn_params_em(engine, data, max_iter, thresh);

% -> faire data2
if nargin < 4,
  [data_mat, LLik, logP] = fill_in_by_inference(data, engine);
end

bnet2 = mk_bnet(dag, ns);
for i = 1:N
  bnet2.CPD{i}=tabular_CPD(bnet, i, 'prior_type', 'dirichlet', 'dirichlet_type', 'unif');
end
bnet2 = learn_params(bnet, data_mat);

[CPT, counts, nsamples] = CPT_from_bnet(bnet, verbose);
[CPT2, counts2, nsamples] = CPT_from_bnet(bnet, verbose);

% POUR V镽IFIER QUE L'ON A BIEN LES MM STATISTIQUES ESSENTIELLES
if 0,
  res=1;
  for i =1:N
    res=res*(prod(counts{i}==counts2{i}));
  end
  if res == 0, error('ca deconne ... ??'); end
end

% -> CALCUL DU SCORE
% ln P(D|G) = ln P(D'|G) - ln P(D'|G,^Theta) + ln P(D|G,^Theta)

%ln P(D'|G)
PD2G = 0;
% Utilisation de Heckerman
PD2G = score_dags(data_mat, ns, {dag}); % see log_marg_prob_node.m

% ln P(D'|G,^Theta)
PD2Gtheta = 0;
% incorpor