393 lines
9.9 KiB
Perl
393 lines
9.9 KiB
Perl
|
/*
|
||
|
|
LPAD and CP-Logic reasoning suite
|
||
|
|
File mcintyre.pl
|
||
|
|
Monte Carlo inference for LPADs
|
||
|
|
Copyright (c) 2011, Fabrizio Riguzzi
|
||
|
|
*/
|
||
|
|
|
||
|
|
|
||
|
|
:-dynamic rule_n/1,setting/2.
|
||
|
|
:-use_module(library(random)).
|
||
|
|
:-use_module(library(lists)).
|
||
|
|
:- set_prolog_flag(unknown,fail).
|
||
|
|
:-source.
|
||
|
|
|
||
|
|
rule_n(1).
|
||
|
|
|
||
|
|
setting(epsilon_parsing, 0.0001).
|
||
|
|
|
||
|
|
setting(compiling,true).
|
||
|
|
/* values: true, failse */
|
||
|
|
|
||
|
|
/* k
|
||
|
|
*
|
||
|
|
* numnber of samples
|
||
|
|
*
|
||
|
|
* Default value: 1000
|
||
|
|
*/
|
||
|
|
setting(k, 1000).
|
||
|
|
/* min_error
|
||
|
|
* ---------
|
||
|
|
* This parameter shows the threshold for the probability interval.
|
||
|
|
*
|
||
|
|
* Default value: 0.01
|
||
|
|
*/
|
||
|
|
setting(min_error, 0.01).
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
s(Goals, Samples, CPUTime, WallTime, Lower, Prob, Upper):-
|
||
|
|
solve(Goals, Samples, CPUTime, WallTime, Lower, Prob, Upper).
|
||
|
|
|
||
|
|
solve(Goals, Samples, CPUTime, WallTime, Lower, Prob, Upper) :-
|
||
|
|
% Retrieving functional parameters...
|
||
|
|
setting(k, K),
|
||
|
|
setting(min_error, MinError),
|
||
|
|
% Resetting the clocks...
|
||
|
|
statistics(cputime,[_,_]),
|
||
|
|
statistics(walltime,[_,_]),
|
||
|
|
% Performing resolution...
|
||
|
|
montecarlo_cycle(0, 0, Goals, K, MinError, Samples, Lower, Prob, Upper),
|
||
|
|
% Taking elapsed times...
|
||
|
|
statistics(cputime,[_,CT1]),
|
||
|
|
CPUTime is CT1/1000,
|
||
|
|
statistics(walltime,[_,WT1]),
|
||
|
|
WallTime is WT1/1000.
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
/* montecarlo(Count, Success, Goals, K, MinError, Samples, Lower, Prob, Upper)
|
||
|
|
* ---------------------------------------------------------------------------
|
||
|
|
* This tail recursive predicate solves the problem currently in memory with a
|
||
|
|
* Monte Carlo approach.
|
||
|
|
* It requires the number of samples and successes so far (Count and Success),
|
||
|
|
* the desired list of Goals to fulfil, the number K of samples to consider at
|
||
|
|
* once and the threshold MinError for the binomial proportion confidence
|
||
|
|
* interval.
|
||
|
|
* It returns the total number of Samples considered, the Lower and Upper ends
|
||
|
|
* of the the binomial proportion confidence interval and the extimated Prob.
|
||
|
|
*
|
||
|
|
* INPUT
|
||
|
|
* - Count: number of samples considered so far.
|
||
|
|
* - Success: number of successfull samples considered so far.
|
||
|
|
* - Goals: list of goals to fulfil.
|
||
|
|
* - K: number of samples to consider at once.
|
||
|
|
* - MinError: threshold for the binomial proportion confidence interval.
|
||
|
|
*
|
||
|
|
* OUTPUT
|
||
|
|
* - Samples: total number of samples considered.
|
||
|
|
* - Lower: lower end of the the binomial proportion confidence interval.
|
||
|
|
* - Prob: extimated probability.
|
||
|
|
* - Upper: upper end of the the binomial proportion confidence interval.
|
||
|
|
*
|
||
|
|
* NB: This method is based on the binomial proportion confidence interval and
|
||
|
|
* the Brown's rule of thumb to avoid the case the sample proportion is
|
||
|
|
* exactly 0.0 or 1.0 and doesn't make use of BDDs.
|
||
|
|
*/
|
||
|
|
montecarlo_cycle(N0, S0, Goals, K, MinError, Samples, Lower, Prob, Upper):-!,
|
||
|
|
montecarlo(K,N0, S0, Goals, N, S),
|
||
|
|
P is S / N,
|
||
|
|
D is N - S,
|
||
|
|
Semi is 1.95996 * sqrt(P * (1 - P) / N),
|
||
|
|
Int is 2 * Semi,
|
||
|
|
/* N * P > 5; N * S / N > 5; S > 5
|
||
|
|
* N (1 - P) > 5; N (1 - S / N) > 5; N (N - S) / N > 5; N - S > 5
|
||
|
|
*/
|
||
|
|
%format("Batch: samples ~d positive ~d interval ~f~n",[N,S,Int]),
|
||
|
|
%flush_output,
|
||
|
|
%((S > 5, D > 5,
|
||
|
|
((Int < MinError; Int =:= 0) ->
|
||
|
|
Samples is N,
|
||
|
|
Lower is P - Semi,
|
||
|
|
Prob is P,
|
||
|
|
Upper is P + Semi
|
||
|
|
;
|
||
|
|
montecarlo_cycle(N, S, Goals, K, MinError, Samples, Lower, Prob, Upper)
|
||
|
|
).
|
||
|
|
|
||
|
|
montecarlo(0,N,S , _Goals,N,S):-!.
|
||
|
|
|
||
|
|
montecarlo(K1,Count, Success, Goals,N1,S1):-
|
||
|
|
abolish_all_tables,
|
||
|
|
eraseall(exp),
|
||
|
|
(Goals->
|
||
|
|
Valid=1
|
||
|
|
;
|
||
|
|
Valid=0
|
||
|
|
),
|
||
|
|
N is Count + 1,
|
||
|
|
S is Success + Valid,
|
||
|
|
K2 is K1-1,
|
||
|
|
montecarlo(K2,N, S, Goals, N1,S1).
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
member_eq(Item, [Head|_Tail]) :-
|
||
|
|
Item==Head, !.
|
||
|
|
|
||
|
|
member_eq(Item, [_Head|Tail]) :-
|
||
|
|
member_eq(Item, Tail).
|
||
|
|
|
||
|
|
|
||
|
|
list2and([X], X) :-
|
||
|
|
X\= (_, _) , !.
|
||
|
|
|
||
|
|
list2and([H|T], (H, Ta)) :- !,
|
||
|
|
list2and(T, Ta).
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
list2or([X], X) :-
|
||
|
|
X \= ( _ ; _ ) , !.
|
||
|
|
|
||
|
|
list2or([H|T], (H ; Ta)) :- !,
|
||
|
|
list2or(T, Ta).
|
||
|
|
|
||
|
|
get_var_n(R,S,Probs,V):-
|
||
|
|
(v(R,S,V)->
|
||
|
|
true
|
||
|
|
;
|
||
|
|
length(Probs,L),
|
||
|
|
add_var(L,Probs,V),
|
||
|
|
assert(v(R,S,V))
|
||
|
|
).
|
||
|
|
sample_head(_HeadList,R,VC,NH):-
|
||
|
|
recorded(exp,(R,VC,NH),_),!.
|
||
|
|
|
||
|
|
sample_head(HeadList,R,VC,NH):-
|
||
|
|
sample(HeadList,NH),
|
||
|
|
recorda(exp,(R,VC,NH),_).
|
||
|
|
|
||
|
|
generate_rules_fact([],_HeadList,_VC,_R,_Probs,_N,[]).
|
||
|
|
|
||
|
|
generate_rules_fact([(Head:_P1),(null:_P2)],_HeadList,VC,R,Probs,N,[Clause]):-!,
|
||
|
|
(Clause=(Head:-(sample_head(Probs,R,VC,NH),NH=N))).
|
||
|
|
|
||
|
|
|
||
|
|
generate_rules_fact([(Head:_P)|T],HeadList,VC,R,Probs,N,[Clause|Clauses]):-
|
||
|
|
(Clause=(Head:-(sample_head(Probs,R,VC,NH),NH=N))),
|
||
|
|
N1 is N+1,
|
||
|
|
generate_rules_fact(T,HeadList,VC,R,Probs,N1,Clauses).
|
||
|
|
|
||
|
|
generate_clause(Head,Body,_HeadList,VC,R,Probs,_BDDAnd,N,_Builtin,Clause):-
|
||
|
|
(ground(VC)->
|
||
|
|
Clause=(Head:-(sample_head(Probs,R,VC,NH),NH=N,Body))
|
||
|
|
;
|
||
|
|
Clause=(Head:-(Body,sample_head(Probs,R,VC,NH),NH=N))
|
||
|
|
).
|
||
|
|
|
||
|
|
generate_rules([],_Body,_HeadList,_VC,_R,_Probs,_BDDAnd,_N,_Builtin,[]).
|
||
|
|
|
||
|
|
generate_rules([(Head:_P1),(null:_P2)],Body,HeadList,VC,R,Probs,BDDAnd,N,Builtin,[Clause]):-!,
|
||
|
|
generate_clause(Head,Body,HeadList,VC,R,Probs,BDDAnd,N,Builtin,Clause).
|
||
|
|
|
||
|
|
generate_rules([(Head:_P)|T],Body,HeadList,VC,R,Probs,BDDAnd,N,Builtin,[Clause|Clauses]):-
|
||
|
|
generate_clause(Head,Body,HeadList,VC,R,Probs,BDDAnd,N,Builtin,Clause),
|
||
|
|
N1 is N+1,
|
||
|
|
generate_rules(T,Body,HeadList,VC,R,Probs,BDDAnd,N1,Builtin,Clauses).
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
/* set(Par,Value) can be used to set the value of a parameter */
|
||
|
|
set(Parameter,Value):-
|
||
|
|
retract(setting(Parameter,_)),
|
||
|
|
assert(setting(Parameter,Value)).
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
extract_vars(Variable, Var0, Var1) :-
|
||
|
|
var(Variable), !,
|
||
|
|
(member_eq(Variable, Var0) ->
|
||
|
|
Var1 = Var0;
|
||
|
|
append(Var0, [Variable], Var1)).
|
||
|
|
|
||
|
|
extract_vars(Term, Var0, Var1) :-
|
||
|
|
Term=..[_F|Args],
|
||
|
|
extract_vars_list(Args, Var0, Var1).
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
extract_vars_list([], Var, Var).
|
||
|
|
|
||
|
|
extract_vars_list([Term|Tail], Var0, Var1) :-
|
||
|
|
extract_vars(Term, Var0, Var),
|
||
|
|
extract_vars_list(Tail, Var, Var1).
|
||
|
|
|
||
|
|
|
||
|
|
difference([],_,[]).
|
||
|
|
|
||
|
|
difference([H|T],L2,L3):-
|
||
|
|
member_eq(H,L2),!,
|
||
|
|
difference(T,L2,L3).
|
||
|
|
|
||
|
|
difference([H|T],L2,[H|L3]):-
|
||
|
|
difference(T,L2,L3).
|
||
|
|
|
||
|
|
|
||
|
|
process_head(HeadList, GroundHeadList) :-
|
||
|
|
ground_prob(HeadList), !,
|
||
|
|
process_head_ground(HeadList, 0, GroundHeadList).
|
||
|
|
|
||
|
|
process_head(HeadList, HeadList).
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
/* process_head_ground([Head:ProbHead], Prob, [Head:ProbHead|Null])
|
||
|
|
* ----------------------------------------------------------------
|
||
|
|
*/
|
||
|
|
process_head_ground([Head:ProbHead], Prob, [Head:ProbHead|Null]) :-!,
|
||
|
|
ProbLast is 1 - Prob - ProbHead,
|
||
|
|
setting(epsilon_parsing, Eps),
|
||
|
|
EpsNeg is - Eps,
|
||
|
|
ProbLast > EpsNeg,
|
||
|
|
(ProbLast > Eps ->
|
||
|
|
Null = [null:ProbLast];
|
||
|
|
Null = []).
|
||
|
|
|
||
|
|
process_head_ground([Head:ProbHead|Tail], Prob, [Head:ProbHead|Next]) :-
|
||
|
|
ProbNext is Prob + ProbHead,
|
||
|
|
process_head_ground(Tail, ProbNext, Next).
|
||
|
|
|
||
|
|
ground_prob([]).
|
||
|
|
|
||
|
|
ground_prob([_Head:ProbHead|Tail]) :-
|
||
|
|
ground(ProbHead), % Succeeds if there are no free variables in the term ProbHead.
|
||
|
|
ground_prob(Tail).
|
||
|
|
|
||
|
|
get_probs([], []).
|
||
|
|
|
||
|
|
get_probs([_H:P|T], [P1|T1]) :-
|
||
|
|
P1 is P,
|
||
|
|
get_probs(T, T1).
|
||
|
|
|
||
|
|
builtin(_A is _B).
|
||
|
|
builtin(_A > _B).
|
||
|
|
builtin(_A < _B).
|
||
|
|
builtin(_A >= _B).
|
||
|
|
builtin(_A =< _B).
|
||
|
|
builtin(_A =:= _B).
|
||
|
|
builtin(_A =\= _B).
|
||
|
|
builtin(true).
|
||
|
|
builtin(false).
|
||
|
|
builtin(_A = _B).
|
||
|
|
builtin(_A==_B).
|
||
|
|
builtin(_A\=_B).
|
||
|
|
builtin(_A\==_B).
|
||
|
|
builtin(length(_L, _N)).
|
||
|
|
builtin(member(_El, _L)).
|
||
|
|
builtin(average(_L, _Av)).
|
||
|
|
builtin(max_list(_L, _Max)).
|
||
|
|
builtin(min_list(_L, _Max)).
|
||
|
|
builtin(nth0(_, _, _)).
|
||
|
|
builtin(nth(_, _, _)).
|
||
|
|
builtin(eraseall(_Id)).
|
||
|
|
builtin(recordzifnot(_Id, _Item, _)).
|
||
|
|
|
||
|
|
parse(FileIn,FileOut):-
|
||
|
|
open(FileIn,read,SI),
|
||
|
|
read_clauses(SI,C),
|
||
|
|
close(SI),
|
||
|
|
process_clauses(C,[],C1),
|
||
|
|
open(FileOut,write,SO),
|
||
|
|
write_clauses(C1,SO),
|
||
|
|
close(SO).
|
||
|
|
|
||
|
|
process_clauses([end_of_file],C,C).
|
||
|
|
|
||
|
|
process_clauses([H|T],C0,C1):-
|
||
|
|
(expand_term(H,H1)->
|
||
|
|
true
|
||
|
|
;
|
||
|
|
H1=H
|
||
|
|
),
|
||
|
|
(H1=[_|_]->
|
||
|
|
append(C0,H1,C2)
|
||
|
|
;
|
||
|
|
append(C0,[H1],C2)
|
||
|
|
),
|
||
|
|
process_clauses(T,C2,C1).
|
||
|
|
|
||
|
|
read_clauses(S,[Cl|Out]):-
|
||
|
|
read_term(S,Cl,[]),
|
||
|
|
(Cl=end_of_file->
|
||
|
|
Out=[]
|
||
|
|
;
|
||
|
|
read_clauses(S,Out)
|
||
|
|
).
|
||
|
|
|
||
|
|
write_clauses([],_).
|
||
|
|
|
||
|
|
write_clauses([H|T],S):-
|
||
|
|
write(S,H),
|
||
|
|
write(S,'.'),
|
||
|
|
nl(S),
|
||
|
|
write_clauses(T,S).
|
||
|
|
|
||
|
|
sample(HeadList, HeadId) :-
|
||
|
|
random(Prob),
|
||
|
|
sample(HeadList, 0, 0, Prob, HeadId), !.
|
||
|
|
|
||
|
|
sample([HeadProb|Tail], Index, Prev, Prob, HeadId) :-
|
||
|
|
Succ is Index + 1,
|
||
|
|
Next is Prev + HeadProb,
|
||
|
|
(Prob =< Next ->
|
||
|
|
HeadId = Index;
|
||
|
|
sample(Tail, Succ, Next, Prob, HeadId)).
|
||
|
|
|
||
|
|
get_next_rule_number(R):-
|
||
|
|
retract(rule_n(R)),
|
||
|
|
R1 is R+1,
|
||
|
|
assert(rule_n(R1)).
|
||
|
|
|
||
|
|
|
||
|
|
user:term_expansion((Head :- Body), Clauses):-
|
||
|
|
% disjunctive clause with more than one head atom
|
||
|
|
setting(compiling,true),
|
||
|
|
Head = (_;_), !,
|
||
|
|
list2or(HeadListOr, Head),
|
||
|
|
process_head(HeadListOr, HeadList),
|
||
|
|
get_next_rule_number(R),
|
||
|
|
get_probs(HeadList,Probs),
|
||
|
|
extract_vars((Head:-Body),[],VC),
|
||
|
|
generate_rules(HeadList,Body,HeadList,VC,R,Probs,_BDDAnd,0,_Builtin,Clauses).
|
||
|
|
|
||
|
|
user:term_expansion((Head :- Body), Clauses) :-
|
||
|
|
% disjunctive clause with a single head atom
|
||
|
|
setting(compiling,true),
|
||
|
|
((Head:-Body) \= ((user:term_expansion(_,_) ):- _ )),
|
||
|
|
Head = (H:_), !,
|
||
|
|
list2or(HeadListOr, Head),
|
||
|
|
process_head(HeadListOr, HeadList),
|
||
|
|
get_next_rule_number(R),
|
||
|
|
get_probs(HeadList,Probs),
|
||
|
|
extract_vars((Head:-Body),[],VC),
|
||
|
|
generate_clause(H,Body,HeadList,VC,R,Probs,_BDDAnd,0,_Builtin,Clauses).
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
user:term_expansion(Head,Clauses) :-
|
||
|
|
% disjunctive fact with more than one head atom
|
||
|
|
setting(compiling,true),
|
||
|
|
Head = (_;_),!,
|
||
|
|
list2or(HeadListOr, Head),
|
||
|
|
process_head(HeadListOr, HeadList),
|
||
|
|
extract_vars((Head),[],VC),
|
||
|
|
get_next_rule_number(R),
|
||
|
|
get_probs(HeadList,Probs),
|
||
|
|
generate_rules_fact(HeadList,HeadList,VC,R,Probs,0,Clauses).
|
||
|
|
|
||
|
|
user:term_expansion(Head,Clause) :-
|
||
|
|
% disjunctive fact with a single head atom
|
||
|
|
setting(compiling,true),
|
||
|
|
(Head \= ((term_expansion(_,_)) :- _ )),
|
||
|
|
Head = (Head1:_), !,
|
||
|
|
list2or(HeadListOr, Head),
|
||
|
|
process_head(HeadListOr, HeadList),
|
||
|
|
get_probs(HeadList,Probs),
|
||
|
|
extract_vars((Head),[],VC),
|
||
|
|
get_next_rule_number(R),
|
||
|
|
(Clause=(Head1:-(sample_head(Probs,R,VC,NH),NH=0))).
|
||
|
|
|