From 1bd96722de5e0a5ad183b3e275b3baccfb8458fe Mon Sep 17 00:00:00 2001 From: vsc Date: Wed, 28 Nov 2007 23:52:14 +0000 Subject: [PATCH] junction tree algorithm git-svn-id: https://yap.svn.sf.net/svnroot/yap/trunk@2031 b08c6af1-5177-4d33-ba66-4b1c6b8b522a --- C/absmi.c | 10 +- C/cdmgr.c | 9 +- CLPBN/Makefile.in | 1 + CLPBN/clpbn.yap | 10 +- CLPBN/clpbn/jt.yap | 469 +++++++++++++++++++++++++++++++ CLPBN/clpbn/matrix_cpt_utils.yap | 104 +++++-- CLPBN/clpbn/vel.yap | 12 +- changes-5.1.html | 1 + library/dgraphs.yap | 22 ++ library/matrix.yap | 4 + library/matrix/matrix.c | 247 +++++++++++++++- library/ordsets.yap | 8 +- library/rbtrees.yap | 86 +++--- library/wdgraphs.yap | 14 +- library/wundgraphs.yap | 16 ++ 15 files changed, 925 insertions(+), 88 deletions(-) create mode 100644 CLPBN/clpbn/jt.yap diff --git a/C/absmi.c b/C/absmi.c index 335f61b9b..dd689c197 100644 --- a/C/absmi.c +++ b/C/absmi.c @@ -10,8 +10,11 @@ * * * File: absmi.c * * comments: Portable abstract machine interpreter * -* Last rev: $Date: 2007-11-26 23:43:07 $,$Author: vsc $ * +* Last rev: $Date: 2007-11-28 23:52:14 $,$Author: vsc $ * * $Log: not supported by cvs2svn $ +* Revision 1.231 2007/11/26 23:43:07 vsc +* fixes to support threads and assert correctly, even if inefficiently. +* * Revision 1.230 2007/11/08 15:52:15 vsc * fix some bugs in new dbterm code. * @@ -2688,11 +2691,6 @@ Yap_absmi(int inp) #ifndef NO_CHECKING check_stack(NoStackCall, H); #endif - if (pt->PredFlags & LogUpdatePredFlag) { - if (pt->OpcodeOfPred != LOCKPRED_OPCODE && - pt->ModuleOfPred != IDB_MODULE && pt->OpcodeOfPred != UNDEF_OPCODE) - fprintf(stderr,"OOPS\n"); - } ENV = ENV_YREG; /* Try to preserve the environment */ ENV_YREG = (CELL *) (((char *) ENV_YREG) + PREG->u.sla.s); diff --git a/C/cdmgr.c b/C/cdmgr.c index 9237da6bb..1a6ec093d 100644 --- a/C/cdmgr.c +++ b/C/cdmgr.c @@ -11,8 +11,11 @@ * File: cdmgr.c * * comments: Code manager * * * -* Last rev: $Date: 2007-11-26 23:43:07 $,$Author: vsc $ * +* Last rev: $Date: 2007-11-28 23:52:14 $,$Author: vsc $ * * $Log: not supported by cvs2svn $ +* Revision 1.213 2007/11/26 23:43:07 vsc +* fixes to support threads and assert correctly, even if inefficiently. +* * Revision 1.212 2007/11/16 14:58:40 vsc * implement sophisticated operations with matrices. * @@ -2966,7 +2969,7 @@ p_undefined(void) { /* '$undefined'(P,Mod) */ PredEntry *pe; - pe = get_pred(Deref(ARG1), CurrentModule, "undefined/1"); + pe = get_pred(Deref(ARG1), Deref(ARG2), "undefined/1"); if (EndOfPAEntr(pe)) return TRUE; LOCK(pe->PELock); @@ -2992,7 +2995,7 @@ p_kill_dynamic(void) { /* '$kill_dynamic'(P,M) */ PredEntry *pe; - pe = get_pred(Deref(ARG1), CurrentModule, "kill_dynamic/1"); + pe = get_pred(Deref(ARG1), Deref(ARG2), "kill_dynamic/1"); if (EndOfPAEntr(pe)) return TRUE; LOCK(pe->PELock); diff --git a/CLPBN/Makefile.in b/CLPBN/Makefile.in index 5d44adf2a..c5bbb9e41 100644 --- a/CLPBN/Makefile.in +++ b/CLPBN/Makefile.in @@ -40,6 +40,7 @@ CLPBN_PROGRAMS= \ $(CLPBN_SRCDIR)/graphs.yap \ $(CLPBN_SRCDIR)/graphviz.yap \ $(CLPBN_SRCDIR)/hmm.yap \ + $(CLPBN_SRCDIR)/jt.yap \ $(CLPBN_SRCDIR)/matrix_cpt_utils.yap \ $(CLPBN_SRCDIR)/topsort.yap \ $(CLPBN_SRCDIR)/utils.yap \ diff --git a/CLPBN/clpbn.yap b/CLPBN/clpbn.yap index 7bea2badf..e0a39ba5f 100644 --- a/CLPBN/clpbn.yap +++ b/CLPBN/clpbn.yap @@ -29,6 +29,9 @@ check_if_vel_done/1 ]). +:- use_module('clpbn/jt', [jt/3 + ]). + :- use_module('clpbn/bnt', [do_bnt/3, check_if_bnt_done/1 ]). @@ -57,7 +60,7 @@ :- dynamic solver/1,output/1,use/1. -solver(vel). +solver(jt). %output(xbif(user_error)). %output(gviz(user_error)). @@ -142,6 +145,8 @@ get_clpbn_vars([_|GVars],CLPBNGVars) :- write_out(vel, GVars, AVars, DiffVars) :- vel(GVars, AVars, DiffVars). +write_out(jt, GVars, AVars, DiffVars) :- + jt(GVars, AVars, DiffVars). write_out(gibbs, GVars, AVars, DiffVars) :- gibbs(GVars, AVars, DiffVars). write_out(bnt, GVars, AVars, DiffVars) :- @@ -224,6 +229,9 @@ bind_clpbn(_, Var, _, _, _, _) :- bind_clpbn(_, Var, _, _, _, _) :- use(vel), check_if_vel_done(Var), !. +bind_clpbn(_, Var, _, _, _, _) :- + use(jt), + check_if_vel_done(Var), !. bind_clpbn(T, Var, Key0, _, _, _) :- get_atts(Var, [key(Key)]), !, ( diff --git a/CLPBN/clpbn/jt.yap b/CLPBN/clpbn/jt.yap new file mode 100644 index 000000000..2c1ca2c09 --- /dev/null +++ b/CLPBN/clpbn/jt.yap @@ -0,0 +1,469 @@ + + +:- use_module(library(dgraphs), + [dgraph_new/1, + dgraph_add_edges/3, + dgraph_add_vertex/3, + dgraph_add_vertices/3, + dgraph_edges/2, + dgraph_vertices/2, + dgraph_transpose/2, + dgraph_to_ugraph/2, + ugraph_to_dgraph/2, + dgraph_neighbors/3 + ]). + +:- use_module(library(undgraphs), + [undgraph_new/1, + undgraph_add_edge/4, + undgraph_add_edges/3, + undgraph_del_vertex/3, + undgraph_del_vertices/3, + undgraph_vertices/2, + undgraph_edges/2, + undgraph_neighbors/3, + undgraph_edge/3, + dgraph_to_undgraph/2 + ]). + +:- use_module(library(wundgraphs), + [wundgraph_new/1, + wundgraph_max_tree/3, + wundgraph_add_edges/3, + wundgraph_add_vertices/3, + wundgraph_to_undgraph/2 + ]). + +:- use_module(library(rbtrees), + [rb_new/1, + rb_insert/4, + rb_lookup/3]). + +:- use_module(library(ordsets), + [ord_subset/2, + ord_insert/3, + ord_intersection/3, + ord_del_element/3, + ord_memberchk/2]). + +:- use_module(library(lists), + [reverse/2]). + +:- use_module(library('clpbn/dists'), + [get_dist_domain_size/2, + get_dist_domain/2, + get_dist_matrix/5]). + +:- use_module(library('clpbn/matrix_cpt_utils'), + [project_from_CPT/3, + reorder_CPT/5, + unit_CPT/2, + multiply_CPTs/4, + divide_CPTs/3, + normalise_CPT/2, + expand_CPT/4, + get_CPT_sizes/2, + reset_CPT_that_disagrees/5, + sum_out_from_CPT/4, + list_from_CPT/2]). + +:- use_module(library('clpbn/display'), [ + clpbn_bind_vals/3]). + +jt(LVs,Vs0,AllDiffs) :- + get_graph(Vs0, BayesNet, CPTs, Evidence), + build_jt(BayesNet, CPTs, JTree), + % JTree is a dgraph + % now our tree has cpts + fill_with_cpts(JTree, NewTree), + propagate_evidence(Evidence, NewTree, EvTree), + message_passing(EvTree, MTree), + get_margin(MTree, LVs, LPs), + clpbn_bind_vals(LVs,LPs,AllDiffs). + + +get_graph(LVs, BayesNet, CPTs, Evidence) :- + run_vars(LVs, Edges, Vertices, CPTs, Evidence), + dgraph_new(V0), + dgraph_add_edges(Edges, V0, V1), + dgraph_add_vertices(Vertices, V1, V2), + dgraph_to_ugraph(V2, BayesNet). + +run_vars([], [], [], [], []). +run_vars([V|LVs], Edges, [V|Vs], [CPTVars-dist([V|Parents],Id)|CPTs], Ev) :- + clpbn:get_atts(V, [dist(Id,Parents)]), + add_evidence(V, Id, Ev, Ev0), + sort([V|Parents],CPTVars), + add_edges(Parents, V, Edges, Edges0), + run_vars(LVs, Edges0, Vs, CPTs, Ev0). + +add_evidence(V, Id, [e(V,P)|Evs], Evs) :- + clpbn:get_atts(V, [evidence(Ev)]), !, + get_dist_domain(Id, D), + find_nth0(D, Ev, 0, P). +add_evidence(_, _, Evs, Evs). + +find_nth0([Id|_], Id, P, P) :- !. +find_nth0([_|D], Id, P0, P) :- + P1 is P0+1, + find_nth0(D, Id, P1, P). + +add_edges([], _, Edges, Edges). +add_edges([P|Parents], V, [V-P|Edges], Edges0) :- + add_edges(Parents, V, Edges, Edges0). + +build_jt(BayesNet, CPTs, Tree) :- + init_undgraph(BayesNet, Moral0), + moralised(BayesNet, Moral0, Markov), + undgraph_vertices(Markov, Vertices), + triangulate(Vertices, Markov, Markov, _, Cliques0), + cliques(Cliques0, EndCliques), + wundgraph_max_tree(EndCliques, J0Tree, _), + root(J0Tree, JTree), + populate(CPTs, JTree, Tree). + +initial_graph(_,Parents, CPTs) :- + test_graph(0, Graph0, CPTs), + dgraph_new(V0), + dgraph_add_edges(Graph0, V0, V1), + % OK, this is a bit silly, I could have written the transposed graph + % from the very beginning. + dgraph_transpose(V1, V2), + dgraph_to_ugraph(V2, Parents). + + +problem_graph([], []). +problem_graph([V|BNet], GraphF) :- + clpbn:get_atts(V, [dist(_,_,Parents)]), + add_parents(Parents, V, Graph0, GraphF), + problem_graph(BNet, Graph0). + +add_parents([], _, Graph, Graph). +add_parents([P|Parents], V, Graph0, [P-V|GraphF]) :- + add_parents(Parents, V, Graph0, GraphF). + + +% From David Page's lectures +test_graph(0, + [1-3,2-3,2-4,5-4,5-7,10-7,10-9,11-9,3-6,4-6,7-8,9-8,6-12,8-12], + [[1]-a, + [2]-b, + [1,2,3]-c, + [2,4,5]-d, + [5]-e, + [3,4,6]-f, + [5,7,10]-g, + [7,8,9]-h, + [9]-i, + [10]-j, + [11]-k, + [6,8,12]-l + ]). +test_graph(1,[a-b,a-c,b-d,c-e,d-f,e-f], + []). + + +init_undgraph(Parents, UndGraph) :- + ugraph_to_dgraph(Parents, DGraph), + dgraph_to_undgraph(DGraph, UndGraph). + +get_par_keys([], []). +get_par_keys([P|Parents],[K|KPars]) :- + clpbn:get_atts(P, [key(K)]), + get_par_kets(Parents,KPars). + +moralised([],Moral,Moral). +moralised([_-KPars|Ks],Moral0,MoralF) :- + add_moral_edges(KPars, Moral0, MoralI), + moralised(Ks,MoralI,MoralF). + +add_moral_edges([], Moral, Moral). +add_moral_edges([_], Moral, Moral). +add_moral_edges([K1,K2|KPars], Moral0, MoralF) :- + undgraph_add_edge(K1,K2,Moral0,MoralI), + add_moral_edges([K1|KPars], MoralI, MoralJ), + add_moral_edges([K2|KPars],MoralJ,MoralF). + +triangulate([], _, Triangulated, Triangulated, []) :- !. +triangulate(Vertices, S0, T0, Tf, Cliques) :- + choose(Vertices, S0, +inf, [], -1, BestVertex, _, Cliques0, Cliques, Edges), + ord_del_element(Vertices, BestVertex, NextVertices), + undgraph_add_edges(Edges, T0, T1), + undgraph_del_vertex(BestVertex, S0, Si), + undgraph_add_edges(Edges, Si, Si2), + triangulate(NextVertices, Si2, T1, Tf, Cliques0). + +choose([], _, _, NewEdges, Best, Best, Clique, Cliques0, [Clique|Cliques0], NewEdges). +choose([V|Vertices], Graph, Score0, _, _, Best, _, Cliques0, Cliques, EdgesF) :- + undgraph_neighbors(V, Graph, Neighbors), + ord_insert(Neighbors, V, PossibleClique), + new_edges(Neighbors, Graph, NewEdges), + ( + % simplicial edge + NewEdges == [] + -> + !, + Best = V, + NewEdges = EdgesF, + length(PossibleClique,L), + Cliques = [L-PossibleClique|Cliques0] + ; + length(PossibleClique,CL), + CL < Score0, !, + choose(Vertices,Graph,CL,NewEdges, V, Best, CL-PossibleClique, Cliques0,Cliques,EdgesF) + ). +choose([_|Vertices], Graph, Score0, Edges0, BestSoFar, Best, Clique, Cliques0, Cliques, EdgesF) :- + choose(Vertices,Graph,Score0,Edges0, BestSoFar, Best, Clique, Cliques0,Cliques,EdgesF). + +new_edges([], _, []). +new_edges([N|Neighbors], Graph, NewEdgesF) :- + new_edges(Neighbors,N,Graph,NewEdges0, NewEdgesF), + new_edges(Neighbors, Graph, NewEdges0). + +new_edges([],_,_,NewEdges, NewEdges). +new_edges([N1|Neighbors],N,Graph,NewEdges0, NewEdgesF) :- + undgraph_edge(N, N1, Graph), !, + new_edges(Neighbors,N,Graph,NewEdges0, NewEdgesF). +new_edges([N1|Neighbors],N,Graph,NewEdges0, [N-N1|NewEdgesF]) :- + new_edges(Neighbors,N,Graph,NewEdges0, NewEdgesF). + + + +% +% This is simple stuff, I just have to remove cliques that +% are subset of the others. +% +cliques(CliqueList, CliquesF) :- + wundgraph_new(Cliques0), + % first step, order by size, + keysort(CliqueList,Sort), + reverse(Sort, Rev), + get_links(Rev, [], Vertices, [], Edges), + wundgraph_add_vertices(Vertices, Cliques0, CliquesI), + wundgraph_add_edges(Edges, CliquesI, CliquesF). + +% stupid quadratic algorithm, needs to be improved. +get_links([], Vertices, Vertices, Edges, Edges). +get_links([Sz-Clique|Cliques], SoFar, Vertices, Edges0, Edges) :- + add_clique_edges(SoFar, Clique, Sz, Edges0, EdgesI), !, + get_links(Cliques, [Clique|SoFar], Vertices, EdgesI, Edges). +get_links([_|Cliques], SoFar, Vertices, Edges0, Edges) :- + get_links(Cliques, SoFar, Vertices, Edges0, Edges). + +add_clique_edges([], _, _, Edges, Edges). +add_clique_edges([Clique1|Cliques], Clique, Sz, Edges0, EdgesF) :- + ord_intersection(Clique1, Clique, Int), + Int \== Clique, + ( + Int = [] -> + add_clique_edges(Cliques, Clique, Sz, Edges0, EdgesF) + ; + % we connect + length(Int, LSz), + add_clique_edges(Cliques, Clique, Sz, [Clique-(Clique1-LSz)|Edges0], EdgesF) + ). + +root(WTree, JTree) :- + wundgraph_to_undgraph(WTree, Tree), + remove_leaves(Tree, SmallerTree), + undgraph_vertices(SmallerTree, InnerVs), + pick_root(InnerVs, Root), + rb_new(M0), + build_tree(Root, M0, Tree, JTree, _). + +remove_leaves(Tree, SmallerTree) :- + undgraph_vertices(Tree, Vertices), + Vertices = [_,_,_|_], + get_leaves(Vertices, Tree, Leaves), + Leaves = [_|_], !, + undgraph_del_vertices(Leaves, Tree, NTree), + remove_leaves(NTree, SmallerTree). +remove_leaves(Tree, Tree). + +get_leaves([], _, []). +get_leaves([V|Vertices], Tree, [V|Leaves]) :- + undgraph_neighbors(V, Tree, [_]), !, + get_leaves(Vertices, Tree, Leaves). +get_leaves([_|Vertices], Tree, Leaves) :- + get_leaves(Vertices, Tree, Leaves). + +pick_root([V|_],V). + +direct_edges([], _, [], []) :- !. +direct_edges([], NewVs, RemEdges, Directed) :- + direct_edges(RemEdges, NewVs, [], Directed). +direct_edges([V1-V2|Edges], NewVs0, RemEdges, [V1-V2|Directed]) :- + ord_memberchk(V1, NewVs0), !, + ord_insert(NewVs0, V2, NewVs), + direct_edges(Edges, NewVs, RemEdges, Directed). +direct_edges([V1-V2|Edges], NewVs0, RemEdges, [V2-V1|Directed]) :- + ord_memberchk(V2, NewVs0), !, + ord_insert(NewVs0, V1, NewVs), + direct_edges(Edges, NewVs, RemEdges, Directed). +direct_edges([Edge|Edges], NewVs, RemEdges, Directed) :- + direct_edges(Edges, NewVs, [Edge|RemEdges], Directed). + + +populate(CPTs, JTree, NewJTree) :- + keysort(CPTs, KCPTs), + populate_cliques(JTree, KCPTs, NewJTree, []). + +populate_cliques(tree(Clique,Kids), CPTs, tree(Clique-MyCPTs,NewKids), RemCPTs) :- + get_cpts(CPTs, Clique, MyCPTs, MoreCPTs), + populate_trees_with_cliques(Kids, MoreCPTs, NewKids, RemCPTs). + +populate_trees_with_cliques([], MoreCPTs, [], MoreCPTs). +populate_trees_with_cliques([Node|Kids], MoreCPTs, [NewNode|NewKids], RemCPts) :- + populate_cliques(Node, MoreCPTs, NewNode, ExtraCPTs), + populate_trees_with_cliques(Kids, ExtraCPTs, NewKids, RemCPts). + + +get_cpts([], _, [], []). +get_cpts([CPT|CPts], [], [], [CPT|CPts]) :- !. +get_cpts([[I|MCPT]-Info|CPTs], [J|Clique], MyCPTs, MoreCPTs) :- + compare(C,I,J), + ( C == < -> + % our CPT cannot be a part of the clique. + MoreCPTs = [[I|MCPT]-Info|LeftoverCPTs], + get_cpts(CPTs, [J|Clique], MyCPTs, LeftoverCPTs) + ; + C == = -> + % our CPT cannot be a part of the clique. + get_cpt(MCPT, Clique, I, Info, MyCPTs, MyCPTs0, MoreCPTs, MoreCPTs0), + get_cpts(CPTs, [J|Clique], MyCPTs0, MoreCPTs0) + ; + % the first element in our CPT may not be in a clique + get_cpts([[I|MCPT]-Info|CPTs], Clique, MyCPTs, MoreCPTs) + ). + +get_cpt(MCPT, Clique, I, Info, [[I|MCPT]-Info|MyCPTs], MyCPTs, MoreCPTs, MoreCPTs) :- + ord_subset(MCPT, Clique), !. +get_cpt(MCPT, _, I, Info, MyCPTs, MyCPTs, [[I|MCPT]-Info|MoreCPTs], MoreCPTs). + + +translate_edges([], [], []). +translate_edges([E1-E2|Edges], [(E1-A)-(E2-B)|NEdges], [E1-A,E2-B|Vs]) :- + translate_edges(Edges, NEdges, Vs). + +match_vs(_,[]). +match_vs([K-A|Cls],[K1-B|KVs]) :- + compare(C, K, K1), + (C == = -> + A = B, + match_vs([K-A|Cls], KVs) + ; + C = < -> + match_vs(Cls,[K1-B|KVs]) + ; + match_vs([K-A|Cls],KVs) + ). + +fill_with_cpts(tree(Clique-Dists,Leafs), tree(Clique-NewDists,NewLeafs)) :- + compile_cpts(Dists, Clique, NewDists), + fill_tree_with_cpts(Leafs, NewLeafs). + + +fill_tree_with_cpts([], []). +fill_tree_with_cpts([L|Leafs], [NL|NewLeafs]) :- + fill_with_cpts(L, NL), + fill_tree_with_cpts(Leafs, NewLeafs). + +transform([], []). +transform([Clique-Dists|Nodes],[Clique-NewDist|NewNodes]) :- + compile_cpts(Dists, Clique, NewDist), + transform(Nodes, NewNodes). + +compile_cpts([Vs-dist(OVs,Id)|Dists], Clique, TAB) :- + OVs = [_|Ps], !, + get_dist_matrix(Id, Ps, _, _, TAB0), + reorder_CPT(OVs, TAB0, Vs, TAB1, Sz1), + multiply_dists(Dists,Vs,TAB1,Sz1,Vars2,ITAB), + expand_CPT(ITAB,Vars2,Clique,TAB). +compile_cpts([], [V|Clique], TAB) :- + unit_CPT(V, CPT0), + expand_CPT(CPT0, [V], [V|Clique], TAB). + +multiply_dists([],Vs,TAB,_,Vs,TAB). +multiply_dists([Vs-dist(OVs,Id)|Dists],MVs,TAB2,Sz2,FVars,FTAB) :- + OVs = [_|Ps], + get_dist_matrix(Id, Ps, _, _, TAB0), + reorder_CPT(OVs, TAB0, Vs, TAB1, Sz1), + multiply_CPTs(tab(TAB1,Vs,Sz1),tab(TAB2,MVs,Sz2),tab(TAB3,NVs,Sz),_), + multiply_dists(Dists,NVs,TAB3,Sz,FVars,FTAB). + +build_tree(Root, Leafs, WTree, tree(Root,Leaves), NewLeafs) :- + rb_insert(Leafs, Root, [], Leafs0), + undgraph_neighbors(Root, WTree, Children), + build_trees(Children, Leafs0, WTree, Leaves, NewLeafs). + +build_trees( [], Leafs, _, [], Leafs). +build_trees([V|Children], Leafs, WTree, NLeaves, NewLeafs) :- + % back pointer + rb_lookup(V, _, Leafs), !, + build_trees(Children, Leafs, WTree, NLeaves, NewLeafs). +build_trees([V|Children], Leafs, WTree, [VT|NLeaves], NewLeafs) :- + build_tree(V, Leafs, WTree, VT, Leafs1), + build_trees(Children, Leafs1, WTree, NLeaves, NewLeafs). + + +propagate_evidence([], NewTree, NewTree). +propagate_evidence([e(V,P)|Evs], Tree0, NewTree) :- + add_evidence_to_matrix(Tree0, V, P, Tree1), !, + propagate_evidence(Evs, Tree1, NewTree). + +add_evidence_to_matrix(tree(Clique-Dist,Kids), V, P, tree(Clique-NDist,Kids)) :- + ord_memberchk(V, Clique), !, + reset_CPT_that_disagrees(Dist, Clique, V, P, NDist). +add_evidence_to_matrix(tree(C,Kids), V, P, tree(C,NKids)) :- + add_evidence_to_kids(Kids, V, P, NKids). + +add_evidence_to_kids([K|Kids], V, P, [NK|Kids]) :- + add_evidence_to_matrix(K, V, P, NK), !. +add_evidence_to_kids([K|Kids], V, P, [K|NNKids]) :- + add_evidence_to_kids(Kids, V, P, NNKids). + +message_passing(tree(Clique-Dist,Kids), tree(Clique-NDist,NKids)) :- + get_CPT_sizes(Dist, Sizes), + upward(Kids, Clique, tab(Dist, Clique, Sizes), IKids, ITab), + ITab = tab(NDist, _, _), + downward(IKids, Clique, ITab, NKids). + +upward([], _, Dist, [], Dist). +upward([tree(Clique1-Dist1,DistKids)|Kids], Clique, Tab, [tree(Clique1-(NewDist1,EDist1),NDistKids)|Kids], NewTab) :- + get_CPT_sizes(Dist1, Sizes1), + upward(DistKids, Clique1, tab(Dist1,Clique1,Sizes1), NDistKids, NewTab1), + NewTab1 = tab(NewDist1,_,_), + ord_intersection(Clique1, Clique, Int), + sum_out_from_CPT(Int, NewDist1, Clique1, Tab1), + multiply_CPTs(Tab, Tab1, NewTab, EDist1). + +downward([], _, _, []). +downward([tree(Clique1-(Dist1,Msg1),DistKids)|Kids], Clique, Tab, [tree(Clique1-NDist1,NDistKids)|Kids]) :- + get_CPT_sizes(Dist1, Sizes1), + ord_intersection(Clique1, Clique, Int), + Tab = tab(Dist,_,_), + divide_CPTs(Dist, Msg1, Div), + sum_out_from_CPT(Int, Div, Clique, STab), + multiply_CPTs(STab, tab(Dist1, Clique1, Sizes1), NewTab, _), + NewTab = tab(NDist1,_,_), + downward(DistKids, Clique1, NewTab, NDistKids). + + +get_margin(NewTree, LVs0, LPs) :- + sort(LVs0, LVs), + find_clique(NewTree, LVs, Clique, Dist), + sum_out_from_CPT(LVs, Dist, Clique, tab(TAB,_,_)), + reorder_CPT(LVs, TAB, LVs0, NTAB, _), + normalise_CPT(NTAB, Ps), + list_from_CPT(Ps, LPs). + +find_clique(tree(Clique-Dist,_), LVs, Clique, Dist) :- + ord_subset(LVs, Clique), !. +find_clique(tree(_,Kids), LVs, Clique, Dist) :- + find_clique_from_kids(Kids, LVs, Clique, Dist). + +find_clique_from_kids([K|_], LVs, Clique, Dist) :- + find_clique(K, LVs, Clique, Dist), !. +find_clique_from_kids([_|Kids], LVs, Clique, Dist) :- + find_clique_from_kids(Kids, LVs, Clique, Dist). + diff --git a/CLPBN/clpbn/matrix_cpt_utils.yap b/CLPBN/clpbn/matrix_cpt_utils.yap index e21d6264c..7c503928c 100644 --- a/CLPBN/clpbn/matrix_cpt_utils.yap +++ b/CLPBN/clpbn/matrix_cpt_utils.yap @@ -1,25 +1,36 @@ -:- module(clpbn_matrix_utils, [init_CPT/2, - project_from_CPT/3, - reorder_CPT/5, - get_dist_size/2, - normalise_CPT/2, - multiply_CPTs/3, - list_from_CPT/2]). +:- module(clpbn_matrix_utils, + [init_CPT/2, + project_from_CPT/3, + reorder_CPT/5, + get_CPT_sizes/2, + normalise_CPT/2, + multiply_CPTs/4, + divide_CPTs/3, + expand_CPT/4, + reset_CPT_that_disagrees/5, + unit_CPT/2, + sum_out_from_CPT/4, + list_from_CPT/2]). -:- use_module(dists, [get_dist_domain_size/2, - get_dist_domain/2]). +:- use_module(dists, + [get_dist_domain_size/2, + get_dist_domain/2]). -:- use_module(library(matrix), [matrix_new/4, - matrix_select/4, - matrix_dims/2, - matrix_shuffle/3, - matrix_expand/3, - matrix_op/4, - matrix_dims/2, - matrix_sum/2, - matrix_sum_out/3, - matrix_op_to_all/4, - matrix_to_list/2]). +:- use_module(library(matrix), + [matrix_new/4, + matrix_new_set/4, + matrix_select/4, + matrix_dims/2, + matrix_shuffle/3, + matrix_expand/3, + matrix_op/4, + matrix_dims/2, + matrix_sum/2, + matrix_sum_out/3, + matrix_sum_out_several/3, + matrix_op_to_all/4, + matrix_set_all_that_disagree/5, + matrix_to_list/2]). :- use_module(library(lists), [nth0/3]). @@ -51,21 +62,21 @@ reorder_CPT(Vs0,T0,Vs,TF,Sizes) :- var(Vs), !, order_vec(Vs0,Vs,Map), ( - Vs == V0 + Vs == Vs0 -> - matrix_shuffle(T0,Map,TF) - ; TF = T0 + ; + matrix_shuffle(T0,Map,TF) ), matrix_dims(TF, Sizes). reorder_CPT(Vs0,T0,Vs,TF,Sizes) :- mapping(Vs0,Vs,Map), ( - Vs == V0 + Vs == Vs0 -> - matrix_shuffle(T0,Map,TF) - ; TF = T0 + ; + matrix_shuffle(T0,Map,TF) ), matrix_dims(TF, Sizes). @@ -98,7 +109,11 @@ split_map([], []). split_map([_-M|Is], [M|Map]) :- split_map(Is, Map). -multiply_CPTs(tab(Tab1, Deps1, Sz1), tab(Tab2, Deps2, Sz2), tab(OT, NDeps, NSz)) :- +divide_CPTs(Tab1, Tab2, OT) :- + matrix_op(Tab1,Tab2,/,OT). + + +multiply_CPTs(tab(Tab1, Deps1, Sz1), tab(Tab2, Deps2, Sz2), tab(OT, NDeps, NSz), NTab2) :- expand_tabs(Deps1, Sz1, Deps2, Sz2, Map1, Map2, NDeps), matrix_expand(Tab1, Map1, NTab1), matrix_expand(Tab2, Map2, NTab2), @@ -140,4 +155,39 @@ normalise_CPT(MAT,NMAT) :- list_from_CPT(MAT, List) :- matrix_to_list(MAT, List). +expand_CPT(MAT0, Dims0, DimsNew, MAT) :- + generate_map(DimsNew, Dims0, Map), + matrix_expand(MAT0, Map, MAT). +generate_map([], [], []). +generate_map([V|DimsNew], [V0|Dims0], [0|Map]) :- V == V0, !, + generate_map(DimsNew, Dims0, Map). +generate_map([V|DimsNew], Dims0, [Sz|Map]) :- + clpbn:get_atts(V, [dist(Id,_)]), + get_dist_domain_size(Id, Sz), + generate_map(DimsNew, Dims0, Map). + +unit_CPT(V,CPT) :- + clpbn:get_atts(V, [dist(Id,_)]), + get_dist_domain_size(Id, Sz), + matrix_new_set(floats,[Sz],1.0,CPT). + +reset_CPT_that_disagrees(CPT, Vars, V, Pos, NCPT) :- + vnth(Vars, 0, V, Dim, _), + matrix_set_all_that_disagree(CPT, Dim, Pos, 0.0, NCPT). + +sum_out_from_CPT(Vs,Table,Deps,tab(NewTable,Vs,Sz)) :- + conversion_matrix(Vs, Deps, Conv), + matrix_sum_out_several(Table, Conv, NewTable), + matrix_dims(NewTable, Sz). + +conversion_matrix([], [], []). +conversion_matrix([], [_|Deps], [1|Conv]) :- + conversion_matrix([], Deps, Conv). +conversion_matrix([V|Vs], [V1|Deps], [0|Conv]) :- V==V1, !, + conversion_matrix(Vs, Deps, Conv). +conversion_matrix([V|Vs], [_|Deps], [1|Conv]) :- + conversion_matrix([V|Vs], Deps, Conv). + +get_CPT_sizes(CPT, Sizes) :- + matrix_dims(CPT, Sizes). diff --git a/CLPBN/clpbn/vel.yap b/CLPBN/clpbn/vel.yap index 3d4e0306b..cc4cb347d 100644 --- a/CLPBN/clpbn/vel.yap +++ b/CLPBN/clpbn/vel.yap @@ -25,9 +25,10 @@ :- use_module(library('clpbn/graphviz'), [clpbn2gviz/4]). -:- use_module(library('clpbn/dists'), [ - get_dist_domain_size/2, - get_dist_matrix/5]). +:- use_module(library('clpbn/dists'), + [ + get_dist_domain_size/2, + get_dist_matrix/5]). :- use_module(library('clpbn/utils'), [ clpbn_not_var_member/2, @@ -39,8 +40,7 @@ :- use_module(library('clpbn/matrix_cpt_utils'), [project_from_CPT/3, reorder_CPT/5, - get_dist_size/2, - multiply_CPTs/3, + multiply_CPTs/4, normalise_CPT/2, list_from_CPT/2]). @@ -178,7 +178,7 @@ find_best([V|LV], V0, Threshold, VF, WorkTables, [V|LVF], Inputs) :- multiply_tables([Table], Table) :- !. multiply_tables([TAB1, TAB2| Tables], Out) :- - multiply_CPTs(TAB1, TAB2, TAB), + multiply_CPTs(TAB1, TAB2, TAB, _), multiply_tables([TAB| Tables], Out). diff --git a/changes-5.1.html b/changes-5.1.html index 48b24b3d3..cb1f4756b 100644 --- a/changes-5.1.html +++ b/changes-5.1.html @@ -17,6 +17,7 @@

Yap-5.1.3: