1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
|
function [cpdag, best_score, cache] = learn_struct_ges(data, nodesizes, varargin)
%
% LEARN_STRUCT_GES learns a structure of Bayesian net by Greedy Equivalence Search.
% cpdag = learn_struct_ges(Data, Nodesizes, 'cache', cache, 'scoring_fn', 'bic', 'verbose', 'yes')
%
% cpdag: the final cpdag
% Data : training data, data(i,m) is the m obsevation of node i
% Nodesizes: the size array of different nodes
% cache : data structure used to memorize local score computations
% (cf. SCORE_INIT_CACHE function)
%
% V1.1 : 28 july 2003 (Ph. Leray - philippe.leray@univ-nantes.fr, O. francois - francois.olivier.c.h@gmail.com)
%
% Ref:
% Optimal Structure Identification with Greedy Search, Chickering 2002
%
[N ncases] = size(data);
seeddag = zeros(N,N);
% set default params
scoring_fn = 'bayesian';
verbose = 0;
cache=[];
% get params
args = varargin;
nargs = length(args);
if length(args) > 0
if isstr(args{1})
for i = 1:2:nargs
switch args{i}
case 'scoring_fn', scoring_fn = args{i+1};
case 'verbose', verbose = strcmp(args{i+1},'yes');
case 'cache', cache=args{i+1} ;
end;
end;
end;
end;
if verbose
names=cellstr(int2str((1:N)'));
carre=zeros(N,1);
end
done = 0;
[best_score cache] = score_dags(data,nodesizes, {seeddag},'scoring_fn',scoring_fn,'cache',cache);
cptt=0;
% First step : INSERT
while ~done
cptt=cptt+1;
[pdags,nodes] = mk_nbrs_of_pdag_add(seeddag);
seedold=seeddag;
sold=best_score;
nbrs = length(pdags);
dags=pdag_to_dag(pdags);
[scores cache] = score_dags(data, nodesizes, dags,'scoring_fn',scoring_fn,'cache',cache);
max_score = max(scores);
new = find(scores == max_score );
if ~isempty(new) & (max_score > best_score)
p = sample_discrete(normalise(ones(1, length(new))));
best_score = max_score;
seeddag = dag_to_cpdag(dags{new(p)});
new=new(p);
if verbose
figure;
subplot(1,2,1), [xx yy]=draw_graph(seedold,names,carre);
set(gca,'color',[1 1 0]);
title(sprintf('current CPDAG (Smax=%5.2f)',sold));
subplot(1,2,2), draw_graph(seeddag,names,carre,xx,yy);
s=sprintf(' %d',nodes{new,3});
title([sprintf('Best in N+ = INSERT(%d, %d,',nodes{new,1},nodes{new,2}) s ')' sprintf(' S=%5.2f',max_score)]);
drawnow;
end
else
done = 1;
end
end;
done = 0;
%[best_score cache] = score_dags(data,nodesizes, {seeddag},'scoring_fn',scoring_fn,'cache',cache);
cptt=0;
if sum(sum(seeddag))==0, done=1;end
% Second step : DELETE
while ~done
cptt=cptt+1;
[pdags,nodes] = mk_nbrs_of_pdag_del(seeddag);
seedold=seeddag; sold=best_score;
nbrs = length(pdags);
dags=pdag_to_dag(pdags);
[scores cache] = score_dags(data, nodesizes, dags,'scoring_fn',scoring_fn,'cache',cache);
max_score = max(scores);
new = find(scores == max_score );
if ~isempty(new) & (max_score > best_score)
p = sample_discrete(normalise(ones(1, length(new))));
best_score = max_score;
seeddag = dag_to_cpdag(dags{new(p)});
new=new(p);
if verbose
cpdags=dag_to_cpdag(dags);
figure;
subplot(1,2,1), [xx yy]=draw_graph(seedold,names,carre);
set(gca,'color',[1 1 0]);
title(sprintf('current CPDAG (Smax=%5.2f)',best_score));
subplot(1,2,2), draw_graph(seeddag,names,carre,xx,yy);
s=sprintf('%d',nodes{new,3});
title([sprintf('Best in N- = DELETE(%d, %d,',nodes{new,1},nodes{new,2}) s ')' sprintf(' S=%5.2f',max_score)]);
drawnow;
end
else
done = 1;
end
end
cpdag = seeddag;
|
