diff --git a/packages/CLPBN/clpbn/bdd.yap b/packages/CLPBN/clpbn/bdd.yap index ded653947..1f37c5c21 100644 --- a/packages/CLPBN/clpbn/bdd.yap +++ b/packages/CLPBN/clpbn/bdd.yap @@ -43,6 +43,8 @@ Va <- P*X1*Y1 + Q*X2*Y2 + ... :- use_module(library(atts)). +:- use_module(library(hacks)). + :- use_module(library(lists)). :- use_module(library(dgraphs)). @@ -65,19 +67,23 @@ bdd([QueryVars], AllVars, AllDiffs) :- finalize_bdd_solver(BayesNet), clpbn_bind_vals([QueryVars], [LPs], AllDiffs). -init_bdd_solver(_, AllVars0, _, bdd(Term, Leaves)) :- +init_bdd_solver(_, AllVars0, _, bdd(Term, Leaves, Tops)) :- check_for_agg_vars(AllVars0, AllVars1), sort_vars(AllVars1, AllVars, Leaves), rb_new(Vars0), rb_new(Pars0), - get_vars_info(AllVars, Vars0, _Vars, Pars0, _Pars, Term, []). + init_tops(Leaves,Tops), + get_vars_info(AllVars, Vars0, _Vars, Pars0, _Pars, Leaves, Tops, Term, []). + +init_tops([],[]). +init_tops(_.Leaves,_.Tops) :- + init_tops(Leaves,Tops). sort_vars(AllVars0, AllVars, Leaves) :- dgraph_new(Graph0), build_graph(AllVars0, Graph0, Graph), dgraph_leaves(Graph, Leaves), - dgraph_top_sort(Graph, RAllVars), - reverse(RAllVars, AllVars). + dgraph_top_sort(Graph, AllVars). build_graph([], Graph, Graph). build_graph(V.AllVars0, Graph0, Graph) :- @@ -93,10 +99,10 @@ add_parents(V0.Parents, V, Graph0, GraphF) :- dgraph_add_edge(Graph0, V0, V, GraphI), add_parents(Parents, V, GraphI, GraphF). - -get_vars_info([], Vs, Vs, Ps, Ps) --> []. -get_vars_info([V|MoreVs], Vs, VsF, Ps, PsF) --> +get_vars_info([], Vs, Vs, Ps, Ps, _, _) --> []. +get_vars_info([V|MoreVs], Vs, VsF, Ps, PsF, Lvs, Outs) --> { clpbn:get_atts(V, [dist(DistId, Parents)]) }, !, +%{ clpbn:get_atts(V, [key(K)]), writeln(V:K:DistId:Parents) }, [DIST], { check_p(DistId, Parms, _ParmVars, Ps, Ps1), unbound_parms(Parms, ParmVars), @@ -104,12 +110,13 @@ get_vars_info([V|MoreVs], Vs, VsF, Ps, PsF) --> DIST = info(V, Tree, Ev, Values, Formula, ParmVars, Parms), get_parents(Parents, PVars, Vs1, Vs2), cross_product(Values, Ev, PVars, ParmVars, Formula0), - get_evidence(V, Tree, Ev, Formula0, Formula) -% (numbervars(Formula,0,_),writeln(formula:Formula), fail ; true) +% (numbervars(Formula0,0,_),writeln(formula0:Ev:Formula0), fail ; true), + get_evidence(V, Tree, Ev, Formula0, Formula, Lvs, Outs) +%, (numbervars(Formula,0,_),writeln(formula:Formula), fail ; true) }, - get_vars_info(MoreVs, Vs2, VsF, Ps1, PsF). -get_vars_info([_|MoreVs], Vs0, VsF, Ps0, PsF, VarsInfo) :- - get_vars_info(MoreVs, Vs0, VsF, Ps0, PsF, VarsInfo). + get_vars_info(MoreVs, Vs2, VsF, Ps1, PsF, Lvs, Outs). +get_vars_info([_|MoreVs], Vs0, VsF, Ps0, PsF, VarsInfo, Lvs, Outs) :- + get_vars_info(MoreVs, Vs0, VsF, Ps0, PsF, VarsInfo, Lvs, Outs). % % look for parameters in the rb-tree, or add a new. @@ -125,7 +132,8 @@ check_p(DistId, Parms, ParmVars, Ps, PsF) :- L is L0-L1, initial_maxes(L1, Multipliers), copy(L, Multipliers, NextMults, NextMults, Parms0, Parms, ParmVars), - rb_insert(Ps, DistId, theta(DistId, Parms, ParmVars), PsF). +%writeln(t:Size:Parms0:Parms:ParmVars), + rb_insert(Ps, DistId, theta(Parms, ParmVars), PsF). % % we are using switches by two @@ -138,10 +146,29 @@ initial_maxes(Size, [1.0|Multipliers]) :- !, copy(0, [], [], _, _Parms0, [], []) :- !. copy(N, [], [], Ms, Parms0, Parms, ParmVars) :-!, copy(N, Ms, NewMs, NewMs, Parms0, Parms, ParmVars). -copy(N, D.Ds, ND.NDs, New, El.Parms0, NEl.Parms, _.ParmVars) :- +copy(N, D.Ds, ND.NDs, New, El.Parms0, NEl.Parms, V.ParmVars) :- N1 is N-1, - NEl is El/D, - ND is D-El, + (El == 0.0 -> + NEl = 0, + ND = D, + V = NEl + ;El == 1.0 -> + NEl = 1, + ND = 0.0, + V = NEl + ;El == 0 -> + NEl = 0, + ND = D, + V = NEl + ;El =:= 1 -> + NEl = 1, + ND = 0.0, + V = NEl + ; + NEl is El/D, + ND is D-El, + V = NEl + ), copy(N1, Ds, NDs, New, Parms0, Parms, ParmVars). unbound_parms([], []). @@ -225,6 +252,14 @@ apply_first_parent(Parents.PVars, Disj+Conj, Theta.TheseParents) :- parents_to_conj(Parents,Theta,Conj), apply_first_parent(PVars, Disj, TheseParents). +apply_middle_parent([Parents], Other, Conj, [ThetaPar]) :- !, + skim_for_theta(Other, Theta, _, ThetaPar), + parents_to_conj(Parents,Theta,Conj). +apply_middle_parent(Parents.PVars, Other, Conj+Disj, ThetaPar.TheseParents) :- + skim_for_theta(Other, Theta, Remaining, ThetaPar), + parents_to_conj(Parents,(Theta),Conj), + apply_middle_parent(PVars, Remaining, Disj, TheseParents). + apply_last_parent([Parents], Other, Conj) :- !, parents_to_conj(Parents,(Theta),Conj), skim_for_theta(Other, Theta, _, _). @@ -233,16 +268,6 @@ apply_last_parent(Parents.PVars, Other, Conj+Disj) :- skim_for_theta(Other, Theta, Remaining, _), apply_last_parent(PVars, Remaining, Disj). -apply_middle_parent([Parents], Other, Conj, [ThetaPar]) :- !, - skim_for_theta(Other, Theta, _, ThetaPar), - parents_to_conj(Parents,Theta,Conj), - commit_deterministic(ThetaPar, Theta). -apply_middle_parent(Parents.PVars, Other, Conj+Disj, ThetaPar.TheseParents) :- - skim_for_theta(Other, Theta, Remaining, ThetaPar), - parents_to_conj(Parents,(Theta),Conj), - commit_deterministic(ThetaPar, Theta), - apply_middle_parent(PVars, Remaining, Disj, TheseParents). - % % % simplify stuff, removing process that is cancelled by 0s @@ -255,15 +280,6 @@ parents_to_conj2([P],P) :- !. parents_to_conj2(P.Ps,Conj*P) :- parents_to_conj2(Ps,Conj). -% -% don't need variables for deterministic parameters of CPTs. -% -commit_deterministic(1, 1) :- !. -commit_deterministic(0, 0) :- !. -commit_deterministic(1.0, 1) :- !. -commit_deterministic(0.0, 0) :- !. -commit_deterministic(_ThetaPar, _Theta). - % % first case we haven't reached the end of the list so we need % to create a new parameter variable @@ -279,13 +295,23 @@ skim_for_theta([[P|Other]], not(P), [Other], _) :- !. skim_for_theta([[P|Other]|More], Ps*not(P), [Other|Left], New ) :- skim_for_theta(More, Ps, Left, New ). -get_evidence(V, Tree, Values, F0, (Tree=Outs).F0) :- +get_evidence(V, Tree, Ev, F0, F, Leaves, Finals) :- clpbn:get_atts(V, [evidence(Pos)]), !, - zero_pos(0, Pos, Values), - get_outs(F0, Outs). + zero_pos(0, Pos, Ev), + insert_output(Leaves, V, Finals, Tree, Outs, SendOut), + get_outs(F0, F, SendOut, Outs). +% hidden deterministic node, can be removed. +get_evidence(_V, _Tree, Ev, F0, [], _Leaves, _Finals) :- + clpbn:get_atts(V, [key(K)]), + functor(K, Name, 2), + ( Name = 'AVG' ; Name = 'MAX' ; Name = 'MIN' ), + !, + one_list(Ev), + eval_outs(F0). %% no evidence !!! -get_evidence(_V, Tree, _Values, F0, (Tree=Outs).F0) :- - get_outs(F0, Outs). +get_evidence(V, Tree, _Values, F0, F1, Leaves, Finals) :- + insert_output(Leaves, V, Finals, Tree, Outs, SendOut), + get_outs(F0, F1, SendOut, Outs). zero_pos(_, _Pos, []). zero_pos(Pos, Pos, 1.Values) :- !, @@ -295,31 +321,95 @@ zero_pos(I0, Pos, 0.Values) :- I is I0+1, zero_pos(I, Pos, Values). -get_outs([V=_F], V) :- !. -get_outs((V=_F).Outs, (V + F0)) :- - get_outs(Outs, F0). +one_list([]). +one_list(1.Ev) :- + one_list(Ev). -run_bdd_solver([[V]], LPs, bdd(Term, Leaves)) :- - build_out_node(Term, Leaves, Node), +% +% insert a node with the disj of all alternatives, this is only done if node ends up to be in the output +% +insert_output([], _V, [], _Out, _Outs, []). +insert_output(V.Lvs, V0, [Top|_], Top, Outs, [Top = Outs]) :- V == V0, !. +insert_output(_.Leaves, V, _.Finals, Top, Outs, SendOut) :- + insert_output(Leaves, V, Finals, Top, Outs, SendOut). + + +get_outs([V=F], [V=NF|End], End, V) :- !, + simplify_exp(F,NF). +get_outs((V=F).Outs, (V=NF).NOuts, End, (F0 + V)) :- + simplify_exp(F,NF), + get_outs(Outs, NOuts, End, F0). + +eval_outs([]). +eval_outs((V=F).Outs) :- + simplify_exp(F,NF), + V = NF, + get_outs(Outs). + +simplify_exp(V,V) :- var(V), !. +simplify_exp(S1+S2,NS) :- !, + simplify_exp(S1, SS1), + simplify_exp(S2, SS2), + simplify_sum(SS1, SS2, NS). +simplify_exp(S1*S2,NS) :- !, + simplify_exp(S1, SS1), + simplify_exp(S2, SS2), + simplify_prod(SS1, SS2, NS). +simplify_exp(not(S),NS) :- !, + simplify_exp(S, SS), + simplify_not(SS, NS). +simplify_exp(S,S). + +simplify_sum(V1, V2, O) :- + ( var(V1) -> + ( var(V2) -> + ( V1 == V2 -> O = V1 ; O = V1+V2 ) ; /* var(V1) , var(V2) */ + ( V2 == 0 -> O = V1 ; V2 == 1 -> O = 1 ; O = V1+V2 ) /* var(V1) , nonvar(V2) */ + ) ; + ( var(V2) -> + ( V1 == 0 -> O = V2 ; V1 == 1 -> O = 1 ; O = V2+V1 ) ; /* nonvar(V1) , var(V2) */ + ( V2 == 0 -> O = V1 ; V2 == 1 -> O = 1 ; V1 == 0 -> O = V2 ; V1 == 1 -> O = 1; O = V1+V2 ) /* nonvar(V1) , nonvar(V2) */ + ) + ). + +simplify_prod(V1, V2, O) :- + ( var(V1) -> + ( var(V2) -> + ( V1 == V2 -> O = V1 ; O = V1*V2 ) ; /* var(V1) , var(V2) */ + ( V2 == 0 -> O = 0 ; V2 == 1 -> O = V1 ; O = V1*V2 ) /* var(V1) , nonvar(V2) */ + ) ; + ( var(V2) -> + ( V1 == 0 -> O = 0 ; V1 == 1 -> O = V2 ; O = V2*V1 ) ; /* nonvar(V1) , var(V2) */ + ( V2 == 0 -> O = 0 ; V2 == 1 -> O = V1 ; V1 == 0 -> O = 0 ; V1 == 1 -> O = V2; V1 == V2 -> O = V1 ; O = V1*V2 ) /* nonvar(V1) , nonvar(V2) */ + ) + ). + + +simplify_not(V, not(V)) :- var(V), !. +simplify_not(0, 1) :- !. +simplify_not(1, 0) :- !. +simplify_not(SS, not(SS)). + + +run_bdd_solver([[V]], LPs, bdd(Term, Leaves, Nodes)) :- + build_out_node(Nodes, Node), findall(Prob, get_prob(Term, Node, V, Prob),TermProbs), sumlist(TermProbs, Sum), normalise(TermProbs, Sum, LPs). -build_out_node(Term, [Leaf], Top) :- !, - find_exp(Leaf, Term, Top). -build_out_node(Term, [Leaf|Leaves], Tops*Top) :- - find_exp(Leaf, Term, Top), - build_out_node(Term, Leaves, Tops). +build_out_node([Top], []). +build_out_node([T,T1|Tops], [Top = T*Top]) :- + build_out_node2(T1.Tops, Top). -find_exp(Leaf, info(V, Top, _Ev, _Values, _Formula, _ParmVars, _Parms)._, Top) :- - V == Leaf, !. -find_exp(Leaf, _.Term, Top) :- - find_exp(Leaf, Term, Top). +build_out_node2([Top], Top). +build_out_node2([T,T1|Tops], T*Top) :- + build_out_node2(T1.Tops, Top). -get_prob(Term, Top, V, SP) :- - bind_all(Term, [], Bindings, V, AllParms, AllParmValues), - term_variables(AllParms, NVs0), - reverse(NVs0, NVs), + +get_prob(Term, Node, V, SP) :- + bind_all(Term, Node, Bindings, V, AllParms, AllParmValues), +% reverse(AllParms, RAllParms), + term_variables(AllParms, NVs), build_bdd(Bindings, NVs, AllParms, AllParmValues, Bdd), bdd_to_probability_sum_product(Bdd, SP), bdd_close(Bdd). @@ -327,27 +417,27 @@ get_prob(Term, Top, V, SP) :- build_bdd(Bindings, NVs, VTheta, Theta, Bdd) :- bdd_from_list(Bindings, NVs, Bdd), bdd_tree(Bdd, bdd(_F,Tree,_Vs)), length(Tree, Len), - VTheta = Theta, - writeln(length=Len). - -bind_all([], Binds, Binds, _V, [], []). -bind_all(info(V, _Tree, Ev, _Values, Formula, ParmVars, Parms).Term, Binds, BindsF, V0, ParmVars.AllParms, Parms.AllTheta) :- + writeln(length=Len), + VTheta = Theta. + +bind_all([], End, End, _V, [], []). +bind_all(info(V, _Tree, Ev, _Values, Formula, ParmVars, Parms).Term, End, BindsF, V0, ParmVars.AllParms, Parms.AllTheta) :- V0 == V, !, set_to_one_zeros(Ev), - bind_formula(Formula, Binds, BindsI), - bind_all(Term, BindsI, BindsF, V0, AllParms, AllTheta). -bind_all(info(_V, _Tree, Ev, _Values, Formula, ParmVars, Parms).Term, Binds, BindsF, V0, ParmVars.AllParms, Parms.AllTheta) :- + bind_formula(Formula, BindsF, BindsI), + bind_all(Term, End, BindsI, V0, AllParms, AllTheta). +bind_all(info(_V, _Tree, Ev, _Values, Formula, ParmVars, Parms).Term, End, BindsF, V0, ParmVars.AllParms, Parms.AllTheta) :- set_to_ones(Ev),!, - bind_formula(Formula, Binds, BindsI), - bind_all(Term, BindsI, BindsF, V0, AllParms, AllTheta). + bind_formula(Formula, BindsF, BindsI), + bind_all(Term, End, BindsI, V0, AllParms, AllTheta). % evidence: no need to add any stuff. -bind_all(info(_V, _Tree, _Ev, _Values, Formula, ParmVars, Parms).Term, Binds, BindsF, V0, ParmVars.AllParms, Parms.AllTheta) :- - bind_formula(Formula, Binds, BindsI), - bind_all(Term, BindsI, BindsF, V0, AllParms, AllTheta). +bind_all(info(_V, _Tree, _Ev, _Values, Formula, ParmVars, Parms).Term, End, BindsF, V0, ParmVars.AllParms, Parms.AllTheta) :- + bind_formula(Formula, BindsF, BindsI), + bind_all(Term, End, BindsI, V0, AllParms, AllTheta). bind_formula([], L, L). -bind_formula(B.Formula, Bs0, BsF) :- - bind_formula(Formula, B.Bs0, BsF). +bind_formula(B.Formula, B.BsF, Bs0) :- + bind_formula(Formula, BsF, Bs0). set_to_one_zeros([1|Values]) :- set_to_zeros(Values).