| [pdag ,GG, maxScore,ischanged]=MGraph_BGscore_simpleSearch(d,pdag,G)
|
function [pdag ,GG, maxScore,ischanged]=MGraph_BGscore_simpleSearch(d,pdag,G)
%Input: d is the data, G is the initial graph
%Output: G best graph
[nr nc]=size(d);
u0=median(d);
n=nc;
L=nr;
sigma=nc+2;
alpha=nc+2;
v=ones(1,n);
%test oct 12
all_edges=wang_kSubset(1:nc,2);
total_number_of_edges=(nc^2-nc)/2;
number_of_repeated=20;
for repeated=1:number_of_repeated
repeated
%end test
b=gnt_graph_to_coeffiecent_b(pdag);
[current_score, T0, TL]=gnt_scoring_completeG(u0,v,b,sigma,alpha,n,L,d);
max_score=current_score
done=0;
ischanged=0;
%check this line
%while repeated<10
while max_score<=current_score
%add path to graph
A=[];
B=[];
temp_edges={};
idx=1;
pscore=[];
done=0;
temp_graphs={};
%clear Z1 Z2 ZZ1 ZZ2 ZZZ
[A B]=find(abs(setdiag(G,1))==0); %potential add edges
if isempty(A) break; end
for i=1:length(A);
if (G(A(i),B(i))==0 & (G(B(i),A(i))==0 & A(i)<B(i)))
tempG=G;
%try to find unshielded at Y
x=A(i);
y=B(i);
ZZ=[];
Z1=[];
Z12=[];
Z1 = find(tempG(:,y))';
Z12=find(tempG(y,:));
ZZ=unique([Z1,Z12]);
ZZ1 = mysetdiff(ZZ, x);
if ~isempty(ZZ1)
len_z=length(ZZ1);
for ii=1:len_z
z=ZZ1(ii);
tempG=G;
if ( (tempG(x,y)==0 & tempG(y,x)==0) )
%fprintf('%d -> %d <- %d\n', x, y, z);
if (x<y & z<y )
tempG(x,y) = -1; tempG(y,x) = 0;
tempG(z,y) = -1; tempG(y,z) = 0;
else
tempG(x,y) = 1; tempG(y,x) = 1;
tempG(z,y) = 1; tempG(y,z) = 1;
end
%if ~isempty(find(tempG==-1))
%tempG,idx
temp_graphs{idx}=tempG;
%disp('at y');
%pause
pscore(idx) = (gnt_scoring_uncompleteG(tempG,sigma,alpha,L,T0,TL));
idx=idx+1;
% end
end
end
end
%try to find unshielded at X
ZZ=[];
tempG=G;
Z2 = find(tempG(:,x))';
Z22=find(tempG(x,:));
ZZ=unique([Z2,Z22]);
ZZ2 = mysetdiff(ZZ, y);
if ~isempty(ZZ2)
len_z=length(ZZ2);
for ii=1:len_z
z=ZZ2(ii);
tempG=G;
if ( (tempG(y,x)==0 & tempG(x,y)==0) )
%fprintf('%d -> %d <- %d\n', x, y, z);
if (y<x & z<x)
tempG(y,x) = -1; tempG(x,y) = 0;
tempG(z,x) = -1; tempG(x,z) = 0;
else
tempG(y,x) = 1; tempG(x,y) = 1;
tempG(z,x) = 1; tempG(x,z) = 1;
end
% if ~isempty(find(tempG==-1))
%tempG , idx
temp_graphs{idx}=tempG;
%disp('at x');
%pause
pscore(idx) = (gnt_scoring_uncompleteG(tempG,sigma,alpha,L,T0,TL));
idx=idx+1;
% end
end
end
end
%add undirected edge
if (tempG(x,y)==0 & tempG(y,x)==0)
tempG=G;
tempG(x,y)=1;
tempG(y,x)=1;
%if ~isempty(find(tempG==-1))
%tempG,idx
%disp('not x y')
%pause
temp_graphs{idx}=tempG;
pscore(idx) = (gnt_scoring_uncompleteG(tempG,sigma,alpha,L,T0,TL));
idx=idx+1;
%end
end
end %end if
end %end for
if ~isempty(temp_graphs)
[best_pscore, best_p] = max(pscore);
best_graph=temp_graphs{best_p};
current_score=best_pscore;
%[pdag2, G2] = wang_learn_struct_undirect2pdag(1,'cond_indep_fisher_z', nc, 5,best_graph,cov(d), nr,0.05);
%current_score=(gnt_scoring_uncompleteG(pdag2,sigma,alpha,L,T0,TL));
current_G=best_graph;
done=1;
else
done=1;
end
if current_score>max_score
disp('add edge')
max_score=current_score
G=current_G;
ischanged=1;
elseif current_score==max_score & done
break;
end
end
current_score=max_score;
done=0;
while max_score<=current_score
rA=[];
rB=[];
temp_r_edges={};
rpscore=[];
rev_score=[];
temp_rev_edges={};
both_p=[];
idx=1;
idx2=1;
temp_graphs={};
need_r_edge=[];
done=0;
%star test
[rA rB]=find(abs(G)==1); %potential remove edge/ another is revers edge ??
%[rA,rB]
%disp('Start remove edge ...')
%pause
for i=1:length(rA)
if rA(i)<rB(i)
temp_rG=G;
%remove edge
%try to find unshielded between a b
x=rA(i);y=rB(i);
Z1=[];
Z12=[];
ZZ=[];
Z1=[];
Z1 = find(temp_rG(x,:));
Z12=find(temp_rG(:,x))';
ZZ=unique([Z1,Z12]);
ZZ1 = mysetdiff(ZZ, y);
Z2=[];
Z22=[];
ZZ=[];
ZZ2=[];
Z2 = find(temp_rG(y,:));
Z22=find(temp_rG(:,y))';
ZZ=unique([Z2,Z22]);
ZZ2 = mysetdiff(ZZ, x);
ZZZ=intersect(ZZ1,ZZ2);
if ~isempty(ZZZ)
%have unshielder
len_zzz=length(ZZZ);
for ii=1:len_zzz
z=ZZZ(ii);
temp_rG=G;
temp_rG(x,y)=0;
temp_rG(y,x)=0;
%fprintf('%d -> %d <- %d\n', x, y, z);
if (x<z & y<z)
temp_rG(x,z) = -1; temp_rG(z,x) = 0;
temp_rG(y,z) = -1; temp_rG(z,y) = 0;
else
temp_rG(x,z) = 1; temp_rG(z,x) = 1;
temp_rG(y,z) = 1; temp_rG(z,y) = 1;
end
% if ~isempty(find(tempG==-1))
% temp_rG
% disp('Remove unshilder');
% pause
temp_graphs{idx}=temp_rG;
rpscore(idx) = (gnt_scoring_uncompleteG(temp_rG,sigma,alpha,L,T0,TL));
idx=idx+1;
% end
% temp_rG=G;
% temp_rG(x,z) = 1; temp_rG(z,x) = 1;
% temp_rG(y,z) = 1; temp_rG(z,y) = 1;
% if ~isempty(find(tempG==-1))
% temp_graphs{idx}=temp_rG;
% rpscore(idx) = (gnt_scoring_uncompleteG(temp_rG,sigma,alpha,L,T0,TL));
% idx=idx+1;
% end
%temp_rG(x,z) = 1; temp_rG(z,x) = 1;
%temp_rG(y,z) = 1; temp_rG(z,y) = 1;
%temp_graphs{idx}=temp_rG;
%rpscore(idx) = (gnt_scoring_uncompleteG(temp_rG,sigma,alpha,L,T0,TL));
%idx=idx+1;
end
else
%non unshielder
temp_rG=G;
temp_rG(x,y)=0;
temp_rG(y,x)=0;
% if ~MGraph_isundirected(temp_rG)
%if ~isempty(find(tempG==-1))
% temp_rG
% disp('remove edg')
temp_graphs{idx}=temp_rG;
rpscore(idx) = (gnt_scoring_uncompleteG(temp_rG,sigma,alpha,L,T0,TL));
idx=idx+1;
%end
end
end %end if
end %end for
if ~isempty(temp_graphs)
[best_rpscore best_rp]=max(rpscore);
best_graph=temp_graphs{best_rp};
%[pdag2, G2] = wang_learn_struct_undirect2pdag(1,'cond_indep_fisher_z', nc, 5,best_graph,cov(d), nr,0.05);
%current_score=(gnt_scoring_uncompleteG(pdag2,sigma,alpha,L,T0,TL));
current_score=best_rpscore;
current_G=best_graph;
done=1;
else
done=1;
end
if current_score>max_score;
disp('Remove edge')
max_score=current_score
G=current_G;
ischanged=1;
%need_r_edge
%best_r_edge
elseif current_score==max_score & done
break;
end
end
% repeated=repeated+1
% if repeated<10
% rand('seed',repeated);
% idx_n=randperm(n);
% if mod(idx_n,2)
% G(idx_n(1),idx_n(2))=0;
% G(idx_n(2),idx_n(1))=0;
% G(idx_n(2),idx_n(3))=0;
% G(idx_n(3),idx_n(2))=0;
% if length(idx_n)>3
% G(idx_n(3),idx_n(4))=0;
% G(idx_n(4),idx_n(3))=0;
% end
% else
% G(idx_n(1),idx_n(2))=1;
% G(idx_n(2),idx_n(1))=1;
% G(idx_n(2),idx_n(3))=1;
% G(idx_n(3),idx_n(2))=1;
% if length(idx_n)>3
% G(idx_n(3),idx_n(4))=1;
% G(idx_n(4),idx_n(3))=1;
% end
% end
% %[pdag2, G2] = wang_learn_struct_undirect2pdag(1,'cond_indep_fisher_z', nc, 5,G,cov(d), nr,0.05);
% %current_score=(gnt_scoring_uncompleteG(pdag2,sigma,alpha,L,T0,TL));
% current_score = (gnt_scoring_uncompleteG(G,sigma,alpha,L,T0,TL));
% max_score=current_score;
% end
%end %end all
%change pattern to DA
%test oct 12
%repeated
all_maxScore(repeated)=max_score
all_G{repeated}=G;
if repeated <number_of_repeated
rand('seed',repeated);
idx_n=randperm(total_number_of_edges);
number_of_changes=ceil(total_number_of_edges*0.5);
%if start with zero G or very sparse graph, the number of changes should
%increase otherwise it should enoutgh
%
for jj=1:number_of_changes
temp_edge=all_edges(idx_n(jj),:);
if mod(idx_n(jj),2)
G(temp_edge(1),temp_edge(2))=0;
G(temp_edge(2),temp_edge(1))=0;
else
% if temp_edge(1)<temp_edge(2)
G(temp_edge(1),temp_edge(2))=1;
G(temp_edge(2),temp_edge(1))=1;
% else
% G(temp_edge(2),temp_edge(1))=-1;
% G(temp_edge(1),temp_edge(2))=0;
% end
end
end
%try this line
[pdag, G2] = wang_learn_struct_undirect2pdag(1,'cond_indep_fisher_z', nc, 5,G,cov(d), nr,0.05);
%check whether have edges go from high order to lower order and remove it
[pa pb]=find(pdag);
isreverse=pa>pb;
idx_remove=find(isreverse==1);
for i=1:length(idx_remove)
pdag(pa(idx_remove(i)),pb(idx_remove(i))) =0;
pdag(pb(idx_remove(i)),pa(idx_remove(i))) =-1;
end
end
end %end all
%end test
GG=[];
[maxScore idx_maxScore]=max(all_maxScore);
GG=all_G{idx_maxScore}
%GG=G;
[pdag, tG] = wang_learn_struct_undirect2pdag(1,'cond_indep_fisher_z', n, 5,GG,cov(d), L, 0.05);
|
|
