    1:- module(prism, [
    2  % session management
    3		prism_start/0,	prism_start/1,	prism_start/2
    4	,	prism_restart/0
    5	,	prism_recover/0
    6	,	prism_close/0
    7   ,  prism_is_running/0
    8
    9	% state
   10	,	prism_state_get/1
   11	,	prism_state_set/1
   12
   13	% loading source files
   14	,	load_bprolog/1
   15	,	load_prism/1,	load_prism/2
   16
   17	% sampling
   18	,	prism_sample/1
   19	,	prism_sample_n/3
   20
   21	% inference
   22	,	prob/2
   23	,	explain/2
   24   ,  explain_probs/3
   25	,	hindsight/3
   26	,	hindsight_agg/3
   27	,	chindsight/3
   28	,	chindsight_agg/3
   29	,	viterbi/4, viterbi_expl/5, viterbi_tree/5
   30	,	viterbi_graph_to_tree/2
   31
   32	% learning
   33	,	prism_learn/4
   34
   35	% pretty printing
   36	,	print_graph/1
   37	,	print_tree/1
   38
   39	% switch management
   40	,	sw_set/3
   41	,	sw_get/3 
   42	,	sw_fix/2
   43	,	sw_unfix/2
   44	,	sw_values/2
   45	,	sw_set_to_sample/1, sw_set_to_sample/2
   46	,	sw_set_to_mean/1, sw_set_to_mean/2
   47
   48	% flags
   49	,	prism_flag/3
   50	,	prism_flag_affects/2
   51	,	prism_flag_set/2
   52	,	prism_flag_get/2
   53
   54	% information
   55	,	prism_statistics/3,	prism_statistics/2
   56	,	prism_show/1
   57
   58	% low level interface
   59	,	prism/1
   60	,	prism_nd/1
   61	,	(prism_dynamic)/1
   62
   63	,	(#)/1
   64	,	(##)/1
   65
   66	% operators
   67	,	op(1150,fx,prism_dynamic)
   68	,	op(950,fx,#)
   69	,	op(950,fx,##)
   70	]).

Using PRISM as a child process

This module provides access PRISM, a Prolog-based probabilistic programming system. PRISM runs under B-Prolog, so this module provides ways to manage and communicate with a PRISM/B-Prolog child process.

Most of the flags that affect learning and inference are managed explicitly and statelessly by the relevant procedures. The states of switches (distributions and pseudocounts) must still be managed by the user, but they can easily be maniplated using sw_get/3 and sw_set/3.

Flags which are still effective statefully are:

flags to do with information and progress display
default parameters for switches
explanation graph housekeeping
sort_hindsight
log_scale affects explanation, viterbi and learning

Starting PRISM

?- use_module(library(prism)).
?- prism_start(path(prism),'prism.log').
?- #X is 2+2. % as PRISM - sanity check
?- load_prism('some_test_program.psm').
?- prism_show(_). % check prism.log for output

You may find it useful to run tail -f prism.log in a terminal to keep an eye on PRISMs output messages.

Types used in this module:

prism_goal   % a callable term
prism_state  % term describing state of PRISM
filename     % a literal absolute or relative filename
filespec     % path to file using search path mechanism.

To be done: - Add some sample PRISM programs and utilities. */

  112:- use_module(library(plrand), [with_rnd_state/1]).  113:- use_module(library(typedef)).  114
  115:- dynamic current_prism/3.  116:- dynamic current_prism_file/1.  117:- dynamic (prism_dynamic)/2.  118:- dynamic saved_state/2.  119
  120:- prolog_load_context(directory, Dir),
  121   directory_file_path(Dir,'psm',PSM),
  122   assert(user:file_search_path(prism,PSM)).  123
  124% ------------------------------- TYPES ----------------------------------
  125
  126:- type pair(A,B) ---> A-B.
  127:- type prism_goal == callable.
  128:- type flag(_)    == atom. % not specified any further 
  129:- type switch(_)  == term. % not specified any further 
  130:- type switch     == switch(_). 
  131:- type param_type ---> probs; counts.
  132:- type fixness    ---> fixed; unfixed.
  133:- type switch_params == pair(fixness, list(float)).
  134
  135:- type switch_distribution ---> default
  136                               ; uniform
  137                               ; f_geometric
  138                               ; f_geometric(natural)
  139                               ; f_geometric(natural, oneof([asc,desc])).
  140
  141:- type prism_state ---> ps( list(switch_values) % ordinary switches with values
  142                           , list(switch_spec)   % switch probabilities
  143                           , list(switch_spec)   % switch pseudocounts
  144                           , list(flag_value)    % flags and values
  145                           , list(switch_values) % switches with dynamic values
  146                           , list(dynamic_pred)  % dynamic predicates
  147                           ).
  148
  149:- type switch_values ---> sw(switch(A), list(A)).
  150:- type switch_state  ---> unset; set(fixness, list(float)).
  151:- type switch_spec   ---> sw(switch(_), oneof([p,a]), switch_state).
  152:- type flag_value    ---> flag(flag(A), A).
  153:- type dynamic_pred  ---> dyn(atom, natural, list(term)).
  154
  155
  156% ------------------------------------------------------------------------
  157%                             MAIN INTERFACE 
  158% ------------------------------------------------------------------------
  159
  160% ---------------------- SESSION CONTROL----------------------------------

prism_start(+Exec:filespec, +LogFile:filename) is det
prism_start(+LogFile:filename) is det
prism_start is det: Start PRISM as a child process. Exec should be a standard filespec pointing to the PRISM executable. % Output from PRISM is recorded to LogFile. PRISM immediately runs bprepl.pl, which receives queries on stdin and replies on stdout. Any currently running PRISM % is closed first. If not supplied, executable is set to path(prism). If not supplied, log file is 'prism.log'.

  171prism_start :- prism_start('prism.log').
  172prism_start(Log) :- 
  173   prism_start(path(prism), Log).
  174prism_start(ExecSpec,Log) :-
  175	absolute_file_name(ExecSpec, [access(execute)],Exec),
  176	absolute_file_name(prism('bprepl.pl'),[access(read)],S),
  177	format(atom(Args),"cl('~w'),main(file('~w'))",[S,Log]),
  178   start_prism_with_args(Exec, ["-g", Args]).
  179
  180start_prism_with_args(Exec,Args) :-
  181	process_create(Exec,Args,[stdin(pipe(In)),stdout(pipe(Out)),process(PID)]),
  182	(current_prism(_,_,_) -> prism_close; true),
  183	set_streams(In,Out,child(Exec,Args,PID)),
  184	nl, once((repeat, wait(infinite), recv(term(repl(ready),_)))),
  185	foreach(prism_dynamic(FF,AA),prism(dynamic(FF/AA))).
  186
  187set_streams(In,Out,Descriptor) :-
  188	%	set_stream(Out,timeout(0)),
  189	set_stream(Out,close_on_abort(false)),
  190	set_stream(In,close_on_abort(false)),
  191	assert(current_prism(Descriptor,In,Out)),
  192	set_prolog_flag(float_format,'%.15g'). % for accurate exchange of floats

prism_is_running is semidet: Succeeds if an instance of prism is running.

  198prism_is_running :- current_prism(_,_,_).

prism_restart is det: Restart the current PRISM process, which must be a child as created by prism_start/0 or prism_start/2. The current state is saved and restored to the new PRISM process.

  205prism_restart :-
  206	prism_full_state(S),
  207	prism_start_with(S).

prism_restore_state(+ID) is det: Retrieves PRISM state previously stored under ID and installs it into PRISM. See prism_save_state/1.

  212prism_restore_state(ID) :-
  213	saved_state(ID,S),
  214	prism_start_with(S).

prism_save_state(+ID) is det: Get current PRISM state and save it under the name ID. Can later be restored used prism_restore_state(ID).

  220prism_save_state(ID) :-
  221	prism_full_state(S),
  222	retractall(saved_state(ID,_)),
  223	assert(saved_state(ID,S)).
  224
  225prism_full_state(state(Exec,Args,Files,State)) :-
  226	prism_state_get(State),
  227	findall(F,current_prism_file(F),Files),
  228	current_prism(child(Exec,Args,_),_,_).

prism_recover is det: Try to restart from saved state if there is one. Recovery state can be saved from prism_restart/0 and prism_restore_state/1.

  237prism_recover :- 
  238	recovery_state(nothing,nothing), !,
  239	format('Nothing to recover from\n').
  240
  241prism_recover :-
  242	recovery_state(just(State),nothing),
  243	prism_start_with(State).
  244
  245recovery_state(S1,S2) :-
  246	( saved_state(recovery,S) -> S1=just(S); S1=nothing),
  247	retractall(saved_state(recovery,_)),
  248	(	S2=just(SS2) 
  249	->	assert(saved_state(recovery,SS2))
  250	;	true).

prism_start_with(+S:prism_state) is det: Tries to start PRISM with the given state, saving the state as recovery_state if it fails.

  257prism_start_with(S) :-
  258	S=state(Exec,Args,Files,State),
  259	catch((
  260			start_prism_with_args(Exec,Args),
  261			maplist(prism,Files),
  262			prism_state_set(State)
  263		), Ex, (
  264			recovery_state(_,just(S)),
  265			format('PRISM exception: ~w.\n',[Ex]),
  266			format('State was saved. Use prism_recover when the problem has been fixed.\n'),
  267			throw(Ex)
  268		)
  269	).

prism_close is det: Close the current PRISM process.

  275prism_close :-
  276	(	current_prism(Desc,In,Out) 
  277	-> writeln('Closing PRISM'),
  278		catch((close(In), close(Out)), Ex, print_message(informational,Ex)),
  279      (Desc=child(_,_,PID) -> process_wait(PID,_); true),
  280		retractall(current_prism(_,_,_))
  281	;	writeln('No active PRISM.')
  282	).
  283
  284% ---------------------- Low level interface ----------------------------

prism(Goal:prism_goal) is semidet: Submit query to PRISM process and get reply. Goal is a normal Prolog goal which submitted to PRISM. It is run as by once/1 and any bindings returned.

  290prism(Goal) :- 
  291	gensym(q,Id),
  292	term_variables(Goal,Vars),
  293	send(query(Id,Goal,Vars)), wait(infinite), 
  294	once((repeat,recv(term(query(Id,Vars,Reply),_)))),
  295	interp(Reply).
  296
  297interp(true).
  298interp(fail) :- !, fail.
  299interp(throw(Ex)) :- throw(Ex). %format('PRISM exception: ~w.\n',[Ex]).

prism_nd(Goal:prism_goal) is nondet: Submit nondeterministic query to PRISM. The goal in run in PRISM via findall/3. All the results are returned and retrieved on backtracking in this predicate.

  307prism_nd(Q) :- term_variables(Q,V), prism(findall(V,Q,VX)), member(V,VX).

# Goal:prism_goal is semidet: Unary prefix operator, equivalent to prism/1.

  312#(G) :- prism(G).

##(Goal:prism_goal) is semidet: Unary prefix operator, equivalent to prism_nd/1.

  316##(G) :- prism_nd(G).
  317
  318
  319% ------------------------------------------------------------------------
  320%                   Higher level prism operations 
  321% ------------------------------------------------------------------------
  322
  323% ---------------------- PRISM state -------------------------------------

prism_state(-S1:prism_state, +S2:prism_state) is det: Perform a transition of the global PRISM state, retreiving the current state into S1 and then setting it to S2.

  329prism_state(S1,S2) :- prism_state_get(S1), prism_state_set(S2).

prism_state_get(-S:prism_state) is det: Get the current state of the PRISM system into a Prolog term S. Can be restored using prism_state_set/1.

  334prism_state_get(ps(SX,PX,CX,FX,VX,DX)) :-
  335	prism(findall(sw(I,V), (get_reg_sw(I), get_values1(I,V)),SX)),           % get all switch values
  336	prism(findall(sw(I,V), dynamic_values(I,V),VX)),   % just the dynamic ones
  337	prism(findall(sw(I,p,set(F,P)),(get_reg_sw(I),get_sw(I,[F,_,P])),PX)),   % all switch probabilities
  338   % !! NB I changed this from h to a as h is a synonym for d, not a
  339	prism(findall(sw(I,a,set(F,C)),(get_reg_sw(I),get_sw_a(I,[F,_,C])),CX)), % all switch pseudocounts
  340	prism(findall(flag(F,V), get_prism_flag(F,V),FX)), % all prism flags
  341	findall( dyn(FF,AA,CCX), (
  342			prism_dynamic(FF,AA), functor(PP,FF,AA),
  343			prism(findall(PP,PP,CCX))
  344		), DX).

prism_state_set(+S:prism_state) is det: Set the current state of the PRISM system to that described in S. State can be obtained using prism_state_get/1 or prism_state/2.

  350prism_state_set(ps(SX,PX,CX,FX,VX,DX)) :-
  351	foreach((member(flag(F,V),FX), V\='$disabled'), prism_flag_set(F,V)),       % restore flags
  352	foreach(member(sw(I,V),VX), prism(set_dynamic_values(I,V))), % assert dynamic values if necessary
  353	foreach(member(sw(I,V),SX), check(prism(get_values(I,V)))),      % check all switch values
  354	forall( member(dyn(_,_,CCX),DX), 
  355		prism(forall(member(CC,CCX),assert(CC)))),
  356	maplist(set_switch,PX), % set switch probabilities
  357	maplist(set_switch,CX). % set switch pseudocounts
  358
  359set_switch(sw(I,p,set(F,P))) :- !, sw_set(probs,I,F-P).
  360set_switch(sw(I,a,set(F,C))) :- !, sw_set(counts,I,F-C).
  361set_switch(sw(I,h,set(F,C))) :- !, 
  362   maplist(add(1),C,A), % incr and use set_sw_a in case pseudocounts are negative
  363   set_switch(sw(I,a,set(F,A))).
  364set_switch(I) :- format('*** unrecognised switch setting: ~q.\n',[I]).
  365
  366check(G) :- call(G) -> true; format('*** check failed: ~w.\n',[G]), fail.

prism_dynamic(+Predicates) is det: Registers predicates as dynamic predicates in PRISM. These are then declared as dynamic automatically whenever PRISM is started.

  372prism_dynamic(F/A) :- prism_dynamic(F,A), !.
  373prism_dynamic(F/A) :- !, assert(prism_dynamic(F,A)).
  374prism_dynamic((P,PX)) :- prism_dynamic(P), prism_dynamic(PX).
  375
  376stoch(H,P,N) :- sumlist(H,N), N>0, maplist(divby(N),H,P).
  377divby(N,X,Y) :- Y is X/N.
  378
  379
  380% ---------------------- PROGRAM LOADING ----------------------------

load_bprolog(+Name:filename) is det: Compile and link a B-Prolog source file.

  384load_bprolog(Name)  :- prism(cl(Name)), assert(current_prism_file(cl(Name))).

load_prism(+Name:filespec, +Opts:list(oneof([neg]))) is det
load_prism(+Name:filespec) is det: Compile and link a PRISM source file. If file extension is ommitted, then '.psm' is assumed. There is only one option currently: if neg is included in the option list, the program is assumed to include negation a failure predicate, and is loaded via prismn/1.

  393load_prism(Name) :- load_prism(Name, []). 
  394load_prism(Name,Opts) :- 
  395	absolute_file_name(Name,[access(read),extensions([psm])],Filename),
  396	file_directory_name(Filename,Directory),
  397	prism(cd(Directory)),
  398	(  member(neg,Opts) 
  399	-> prism(prismn(Filename))
  400	;  prism(prism(Filename))
  401	),
  402	retractall(current_prism_file(prism(_))), % only one prism file active at a time
  403	assert(current_prism_file(prism(Name))).
  404
  405% ---------------------- MANAGING SWITCHES ----------------------------

sw_values(+S:switch(A), -V:list(A)) is semidet
sw_values(-S:switch(A), -V:list(A)) is nondet: True when S is a registered switch and V is its list of values.

  410sw_values(S,V) :- prism_nd(get_reg_sw(S)), prism(get_values1(S,V)).

sw_get(+PType:param_type, +S:switch, -I:switch_params) is semidet
sw_get(-PType:param_type, -S:switch, -I:switch_params) is nondet: Gets information on probability distribution or pseudocounts for a switch. Fails if the switch parameters are unset.

  417sw_get(T,S,F-Vals)  :- getter(T,S,Vals,F,G), prism_nd(get_reg_sw(S)), prism(G).

sw_set(+PType:param_type, +S:switch, +H:switch_params) is det: Set a switch's parameters depending on PType.

  422sw_set(Type,S,F-V)  :- setter(Type,S,V,Setter), fixer(Type,S,F,Fixer), prism((Setter,Fixer)).

sw_fix(+PType:param_type, S:switch) is det: Fixes the distribution or pseudocounts for the named switch. This prevents the setting from changing during learning.

  427sw_fix(Type,S)   :- fixer(Type, S, fixed, G), prism(G).

sw_unfix(+PType:param_type, S:switch) is det: Unfixes the distribution or pseudocounts for the named switch. This allows the setting to change during learning.

  432sw_unfix(Type,S)   :- fixer(Type, S, unfixed, G), prism(G).
  433
  434
  435% tables of PRISM commands
  436setter(probs,  S, V, set_sw(S,V)).
  437setter(counts, S, V, set_sw_a(S,V)).
  438
  439getter(probs, S, V, F, get_sw(S,[F,_,V])).
  440getter(counts, S, V, F, get_sw_a(S,[F,_,V])).
  441
  442fixer(probs,  S, fixed,   fix_sw(S)).
  443fixer(probs,  S, unfixed, unfix_sw(S)).
  444fixer(counts, S, fixed,   fix_sw_a(S)).
  445fixer(counts, S, unfixed, unfix_sw_a(S)).

sw_set_to_mean(+S:switch) is det
sw_set_to_mean(+S:switch, +Offset:float) is det: Sets a switch's distribution to the mean of the Dirichlet distribution whose parameters come from the current 'counts' setting of the switch. Off, if present, is added to all the pseudocounts. Also sets the probabilities' 'fixed' states to match that of the pseudocounts.

  455sw_set_to_mean(S) :- sw_set_to_mean(S,0).
  456sw_set_to_mean(S,Off) :-
  457	sw_get(counts,S,F-D),
  458	maplist(add(Off+1),D,A),
  459   stoch(A,Probs,_),
  460	sw_set(probs,S,F-Probs).

sw_set_to_sample(+S:switch) is det
sw_set_to_sample(+S:switch, Offset:float) is det: Sample a new discrete distribution for a switch's parameters from a Dirichlet distribution whose parameters come from the current 'counts' setting of the switch. Off, if present, is added to all the pseudocounts. This uses the undocumented sample_Dirichlet/5 predicate from library(plrand).

  470sw_set_to_sample(S) :- sw_set_to_sample(S,0).
  471sw_set_to_sample(S,Off) :-
  472	sw_get(counts,S,_-D),
  473	maplist(add(Off+1),D,A),
  474   length(A,N),
  475	with_rnd_state(plrand:sample_Dirichlet(N,A,P)),
  476	sw_set(probs,S,P).
  477
  478add(Z,X,Y) :- Y is X+Z. 
  479
  480
  481% ---------------------- SAMPLING EXECUTION ----------------------------

prism_sample(Goal:prism_goal) is semidet: Sample from distribution specified by PRISM goal.

  485prism_sample(Q) :- prism(sample(Q)).

prism_sample_n(+N:nat, Goal:cond_goal, Results:list(prism_goal)) is det: Multiple samples from distribution specified by PRISM goal. Goal can be an ordinary PRISM goal or a conditional goal of the form (Goal|Cond), where Cond is an ordinary Prolog goal. Returns a list of copies of the called goal with return values filled in.

  494prism_sample_n(N,Q,X) :- Q\=(_|_), !, prism(get_samples(N,Q,X)).
  495prism_sample_n(N,(Q|C),X) :- prism(get_samples_c(N,Q,C,X)).
  496
  497% ---------------------- INFERENCE ----------------------------

prob(Goal:prism_goal, -P:float) is det: Compute probability of given goal.

  501prob(Q,P) :- prism(prob(Q,P)).

explain(Goal:prism_goal, -Graph) is det: Gets explanation graph for given goal.

  505explain(Q,G)          :- prism(probf(Q,G)).

explain_probs(+Type:oneof([inside,outside,viterbi,in_out]), Goal:prism_goal, -Graph) is det

Gets explanation graph annoted with probabilities for given goal. The explanation types map to a PRISM commands as follows:

inside -> probfi/2
outside -> probfo/2
in_out -> probfio/2
viterbi -> probvf/2

  514explain_probs(inside,Q,G)   :- prism(probfi(Q,G)).
  515explain_probs(outside,Q,G)  :- prism(probfo(Q,G)).
  516explain_probs(viterbi,Q,G)  :- prism(probfv(Q,G)).
  517explain_probs(in_out,Q,G)   :- prism(probfio(Q,G)).

hindsight(TopGoal:prism_goal, SubGoal:prism_goal, -Probs:list(pair(prism_goal,float))) is nondet: Computes probabilities of subgoals unifying with SubGoal for a given top goal. Matches different subgoals on backtracking. hindsight(Q,R,P) :- prism(hindsight(Q,R,P1)), maplist(list_pair,P1,P).

  524hindsight(Q,R,P) :-      prism(hindsight(Q,R,P1)), member([R,P],P1).
  525chindsight(Q,R,P) :-     prism(chindsight(Q,R,P1)), member([R,P],P1).

hindsight_agg(TopGoal:prism_goal, SubGoal:prism_goal, -Probs:list(list(pair(prism_goal,float)))) is det: Computes total probability of all subgoals unifying with SubGoal for a given top goal. NB. this functionality can be replicated in a more idiomatic way using hindsight/3 and SWI Prolog's library(aggregate).

  532hindsight_agg(Q,R,P) :- prism(hindsight_agg(Q,R,P1)), maplist((maplist(list_pair)),P1,P).

chindsight(TopGoal:prism_goal, SubGoal:prism_goal, -Probs:list(pair(prism_goal,float))) is nondet: Computes conditional probabilities of matching subgoals given that the top goal is true. Same as hindsight/3 but with probabilities divided by the probabality of the top goal. chindsight(Q,R,P) :- prism(chindsight(Q,R,P1)), maplist(list_pair,P1,P).

chindsight_agg(TopGoal:prism_goal, SubGoal:prism_goal, -Probs:list(list(pair(prism_goal,float)))) is det: Computes aggregate conditional probability of matching subgoals given that the top goal is true. Same as chindsight_agg/3 but with probabilities divided by the probabality of the top goal. NB. this functionality can be replicated in a more idiomatic way using hindsight/3 and SWI Prolog's library(aggregate).

  546chindsight_agg(Q,R,P) :- prism(chindsight_agg(Q,R,P1)), maplist((maplist(list_pair)),P1,P).
  547
  548goal_probs(P1,P2) :- maplist(list_pair,P1,P2).
  549list_pair([X,Y],X-Y).
  550
  551% expectation
  552expect(X,G,EX) :-
  553	chindsight(G,G,PX),
  554	accum(X,G,PX,EX).
  555
  556accum(_,_,[], 0).
  557accum(X,G,[G1-P1 | PX], Tot) :-
  558	accum(X,G,PX,TotX),
  559	copy_term(X/G,X1/G1),
  560	Tot is P1*X1 + TotX.

viterbi(+N:nat, Goal:prism_goal, -P:float, +Opts:list(viterbi_opt)) is nondet

Compute probabilities of N most probable explanation for Goal. Options can be as follows:

viterbi_opt ---> ground(bool)        % [false] whether or not to ground variable in Goal
               ; mode(oneof([ml,vb]) % [ml]  use probs (ml) or counts (vb)
               ; rerank(nat)         % [10]  number of explanations for reranking
               .

  573viterbi(N,Q,P,Opts) :- 
  574	process_viterbi_opts(Opts,Ground),
  575	viterbi_prob(Ground,N,Q,P).
  576
  577viterbi_prob(false,1,Q,P) :- !, prism(viterbi(Q,P)).
  578viterbi_prob(false,N,Q,P) :- prism(n_viterbi(N,Q,PX)), member(P,PX).
  579viterbi_prob(true,1,Q,P)  :- !, prism(viterbig(Q,P)).
  580viterbi_prob(true,N,Q,P)  :- prism(n_viterbig(N,Q,P)).

viterbi_expl(+N:nat, Goal:prism_goal, -P:float, +Opts:list(viterbi_opt), -G:graph) is nondet: Compute probabilities of N most probable explanation for Goal and also produce graphs for each. See viterbi/4 for options.

  586viterbi_expl(N,Q,P,Opts,G) :- 
  587	process_viterbi_opts(Opts,Ground),
  588	viterbi_expl1(Ground,N,Q,P,G).
  589
  590viterbi_expl1(false,1,Q,P,G) :- !, prism(viterbif(Q,P,G)).
  591viterbi_expl1(false,N,Q,P,G) :- prism(n_viterbif(N,Q,VX)), member(v_expl(_,P,G),VX).
  592viterbi_expl1(true,1,Q,P,G)  :- !, prism(viterbig(Q,P,G)).
  593viterbi_expl1(true,N,Q,P,G)  :- prism(n_viterbig(N,Q,P,G)).

viterbi_tree(+N:nat, Goal:prism_goal, -P:float, +Opts:list(viterbi_opt), -T:tree) is nondet: Compute probabilities of N most probable explanation for Goal, as viterbi_expl/5, but produce a tree instead of a graph by using viterbi_graph_to_tree/2. also produce tree for each. See viterbi/4 for options.

  599viterbi_tree(N,Q,P,Opts,T) :- 
  600   viterbi_expl(N,Q,P,Opts,G),
  601   viterbi_graph_to_tree(G,T).
  602
  603process_viterbi_opts(Opts,Ground) :-
  604	select_option(ground(Ground),Opts,O1,false),
  605	select_option(mode(Mode),O1,O2,ml),
  606	select_option(rerank(Rerank),O2,O3,10),
  607	(	O3=[] -> true; throw(error(unrecognised_options(O2)))),
  608	prism_flag_set(viterbi_mode,Mode),
  609	prism_flag_set(rerank,Rerank).

viterbi_graph_to_tree(+G:graph, -T:tree) is det: Computes Viterbi tree for a given Viterbi explanation graph. Graph comes from one of the viterbi or n_viterbi predicates.

  616viterbi_graph_to_tree(G,T) :- prism(viterbi_tree(G,T)).

print_graph(+G:graph) is det: Prints an explanation graph on PRISM's user output stream.

  620print_graph(G) :- prism(print_graph(user_output,G,[])).

print_tree(+G:tree) is det: Prints a Viterbi tree PRISM's user output stream.

  624print_tree(T)  :- prism(print_tree(user_output,T,[])).
  625
  626
  627
  628% ---------------------- LEARNING ------------------------

prism_learn(+Method:learn_method, +GG:list(goal), +Opts:list(learn_opt), -Scores:list(learn_score)) is det

Learn model parameters from a list of goals.

 learn_method ---> map(Init:oneof([none,random,noisy(Std)]))
                 ; vb(Init:oneof([none,perturb(Std),reset,noisy(Std),Pt:oneof([on,off]))
                 ; vb_pm(Init:oneof([none,perturb(Std),reset,noisy(Std))
                 .

For MAP learning, the Init option determines how switch probabilities are initialised; the values have the following meanings:

none: Switch probabilities begin with their current values.
random: Probabilities are set to random values
noisy(Std:nonneg): Probabilities are drawn from rectified Gaussian with given standard deviation.

For VB learning, initialisation methods are as follows

none: Learning continues from current pseudocounts (ie evidence accumulates).
perturb(Std:nonneg): Current values have rectified Gaussian noise added
reset: Hyperparameters are set to value of default_sw_a flag.
noisy(Std:nonneg): Hyperparameters are set to value of default_sw_a plus some noise.

Method vb_pm is like vb except that switch parameters are set to their posterior means after learning.

Valid options are (with defaults in comments):

 learn_opt ---> daem(oneof([off,sched(Init,Rate)])) % on
              ; epsilon(nonneg)                     % 0.0001
              ; max_iterate(nonneg)                 % default
              ; restart(natural)                    % 1
              .

 learn_score ---> log_prior(float)    % only for map learning
                ; log_post(float)     % only for map learning
                ; log_lik(float)
                ; free_energy(float)  % only for VB learning
                .

  678prism_learn(map(Init),GG,Opts,[log_prior(LPr),log_lik(LL),log_post(LPo),bic(BIC)]) :- 
  679	set_prism_option(learn_mode,ml),
  680	set_prism_option(init,Init),
  681	process_learn_options(Opts),
  682	prism_flush(learn(GG)),
  683   prism_statistics(log_likelihood,LL),
  684   (  maplist(prism_statistics,[log_prior,log_post,bic],[LPr,LPo,BIC]) -> true
  685   ;  prism_statistics(num_parameters,K),
  686      length(GG,N),
  687      LPr=0, LPo=LL, 
  688      BIC is LL - (K/2)*log(N)
  689   ).
  690
  691prism_learn(vb(Init),GG,Opts,[free_energy(FE)]) :- 
  692	set_prism_option(vb_init,Init),
  693	set_prism_option(learn_mode,vb),
  694	process_learn_options(Opts),
  695	prism_flush(learn(GG)),
  696	prism_statistics(free_energy,FE).
  697
  698
  699prism_learn(vb_pm(Init),GG,Opts,[free_energy(FE)]) :-
  700	set_prism_option(vb_init,Init),
  701	set_prism_option(learn_mode,both),
  702	process_learn_options(Opts),
  703	prism_flush(learn(GG)),
  704	prism_statistics(free_energy,FE).
  705
  706
  707process_learn_options(Options) :-
  708	foldl(process_option, 
  709		[	fix_init_order	/on
  710		,	epsilon  		    /0.0001
  711		,	max_iterate 	  /default
  712		, daem            /off
  713		,	restart 			  /1
  714		],Options,O2),
  715		(	O2=[] -> true; throw(error(unrecognised_options(O2)))).
  716
  717process_option(Name/Default,O1,O2) :- 
  718	functor(Opt,Name,1), arg(1,Opt,Val),
  719	select_option(Opt,O1,O2,Default),
  720	set_prism_option(Name,Val).
  721
  722
  723set_prism_option(vb_init,Init) :- !,
  724	member(Init/(Reset,Std),[none/(off,0.000000001), perturb(U)/(off,U), 
  725									 reset/(on,0.000000001), noisy(U)/(on,U)]),
  726	prism_flag_set(reset_hparams,Reset),
  727	prism_flag_set(std_ratio,Std).
  728
  729set_prism_option(init,noisy(Std)) :- !,
  730	prism_flag_set(init,noisy_u), 
  731	prism_flag_set(std_ratio,Std).
  732
  733set_prism_option(daem,off) :- !, prism_flag_set(daem,off).
  734set_prism_option(daem,sched(Bstart,Brate)) :- !, 
  735	prism_flag_set(daem,on),
  736	prism_flag_set(itemp_init,Bstart),
  737	prism_flag_set(itemp_rate,Brate).
  738	
  739set_prism_option(Name,Val) :- prism_flag_set(Name,Val).
  740
  741% --------------------- PRISM FLAGS ---------------------------

prism_flag_set(+Name:flag(A), +Value:A) is det: Set value of the named PRISM flag.

  746prism_flag_set(F,V) :- prism(set_prism_flag(F,V)).

prism_flag_get(+Name:flag(A), -Value:A) is det
prism_flag_get(-Name:flag(A), -Value:A) is nondet: Get value of the named PRISM flag or of all flags in turn.

  752prism_flag_get(F,V) :- prism_nd(get_prism_flag(F,V)).

prism_flag(+Name:flag(A), -A:type, -Description) is det
prism_flag(?Name:flag(A), ?A:type, ?Description) is nondet: Contains information about all the PRISM flags.

  759% cosmetic - affect display only
  760prism_flag( warn,       	   oneof([on,off]), 'controls display of warning messages').
  761prism_flag( verb,      		   oneof([none,graph,em,full]), 'control verbosity in learning').
  762prism_flag( show_itemp,       oneof([on,off]), 'print inverse temperature rate in DAEM algorithm').
  763prism_flag( search_progress, 	natural,         'frequency of message printing in explanation finding').
  764prism_flag( em_progress,      natural,         'frequence of message printing in EM algorirhm').
  765prism_flag( mcmc_progress, 	natural,         'frequency of message printing in MCMC learning').
  766prism_flag( learn_message, 	oneof([search,em,stats,misc,none,all]),      'messages during learning (can combine with +)').
  767prism_flag( mcmc_message,   	oneof([search,em,mcmc,stats,misc,none,all]), 'messages during MCMC learning (can combine with +)').
  768prism_flag( write_call_events,  oneof([none, off, all, _]), 'no idea').
  769
  770% model initialisation
  771prism_flag( default_sw_a,  oneof([none,uniform,uniform(nonneg)]), 'default values for switch pseudo-counts').
  772prism_flag( default_sw_d,  oneof([none,uniform,uniform(nonneg)]), 'default values for switch pseudo-counts').
  773prism_flag( default_sw, 	switch_distribution,                   'default distribution for new switches').
  774
  775% affect inference, explanation generation
  776prism_flag( clean_table, 	oneof([on,off]), 'dispose of explanation table after use').
  777prism_flag( force_gc,   	oneof([on,off]), 'trigger garbage collection after inference operations').
  778prism_flag( reduce_copy, 	oneof([on,off]), 'affects copying of terms in explanation generation').
  779prism_flag( error_on_cycle,oneof([on,off]), 'checks for cycles in derivations').
  780prism_flag( log_scale,  	oneof([on,off]), 'log scaling in inference and learning').
  781prism_flag( explicit_empty_expls, oneof([on,off]), 'affects structure of explanation graphs').
  782prism_flag( sort_hindsight,oneof([by_goal,by_prob]),   'how to sort subgoals in hindsight').
  783
  784% affect viterbi algorithm
  785prism_flag( viterbi_mode,  oneof([ml,vb]),  'which parameters to use in Viterbi (probs or counts)').
  786prism_flag( log_viterbi,	oneof([on,off]), 'use log probabilities in Viterbi algorithm').
  787prism_flag( rerank, 			natural,         'affects Viterbi algorithm with hyperparameters').
  788
  789% affect learning
  790prism_flag( learn_mode, 	oneof([ml, vb, both, ml_vt, vb_vt, both_vt]), 'switch between ML/EM and VB/EM'). 
  791prism_flag( data_source, 	oneof([data/1, file(_), none]), 'data source for learning').
  792prism_flag( init,     		oneof([none,random,noisy_u]), 'initialisation method for EM algorithm').
  793prism_flag( fix_init_order,oneof([on,off]), 'fixes switch initialisation order in EM algorithm').
  794prism_flag( std_ratio,  	nonneg,          'variance of noise when EM initialisation method is noisy_u').
  795prism_flag( restart,    	natural,         'number of restarts in EM algorithm').
  796prism_flag( epsilon,       nonneg,          'convergence threshold for EM algorirhm').
  797prism_flag( max_iterate,   natural,         'maximum iterations in EM algorithm (can be default)').
  798prism_flag( params_after_vbem, oneof([none,mean,max]),  'method of parameter estimation after VB/EM').
  799prism_flag( reset_hparams, oneof([on,off]), 'reset pseudocounts to default before VB/EM').
  800
  801% DAEM learning
  802prism_flag( daem,     		oneof([on,off]), 'enable deterministic annealing EM algorithm').
  803prism_flag( itemp_init,    range(0,1),      'initial inverse temperature for DAEM algorithm').
  804prism_flag( itemp_rate,    range(1,inf),    'inverse temperature rate for DAEM algorithm').
  805
  806prism_flag( mcmc_b,    	   natural,         'MCMC burn in').
  807prism_flag( mcmc_e,    	   natural,         'MCMC chain length').
  808prism_flag( mcmc_s,    	   natural,         'MCMC cycle length').

prism_flag_affects(-Name:flag(_), -Subsystem) is nondet: Relation between flag names and the subsytem affected, which can be one of display, initialisation (of switches), explanation, vitberbi, hindsight or learning.

  815prism_flag_affects( warn,	             display).
  816prism_flag_affects( verb,	             display).
  817prism_flag_affects( show_itemp,	       display).
  818prism_flag_affects( em_progress,        display).
  819prism_flag_affects( search_progress,	 display).
  820prism_flag_affects( learn_progress,	    display).
  821prism_flag_affects( mcmc_progress,	    display).
  822prism_flag_affects( write_call_events,  display).
  823
  824prism_flag_affects( default_sw,           initialisation).
  825prism_flag_affects( default_sw_a,         initialisation).
  826prism_flag_affects( default_sw_d,         initialisation).
  827
  828prism_flag_affects( clean_table,    explanation).
  829prism_flag_affects( reduce_copy,    explanation).
  830prism_flag_affects( log_scale,      explanation).
  831prism_flag_affects( error_on_cycle, explanation).
  832prism_flag_affects( sort_hindsight, hindsight).
  833
  834prism_flag_affects( viterbi_mode,  viterbi).
  835prism_flag_affects( log_viterbi,   viterbi).
  836prism_flag_affects( rerank,        viterbi).
  837prism_flag_affects( log_scale,     viterbi).
  838
  839prism_flag_affects( daem,           learning).
  840prism_flag_affects( itemp_init,     learning).
  841prism_flag_affects( itemp_rate,     learning).
  842prism_flag_affects( learn_mode, 	   learning).	
  843prism_flag_affects( data_source, 	learning).
  844prism_flag_affects( init,     		learning).
  845prism_flag_affects( fix_init_order, learning).
  846prism_flag_affects( std_ratio,     	learning).
  847prism_flag_affects( restart,    	   learning).
  848prism_flag_affects( epsilon,        learning).
  849prism_flag_affects( max_iterate,    learning).
  850prism_flag_affects( reset_hparams,  learning). 
  851prism_flag_affects( log_scale,      learning).
  852prism_flag_affects( mcmc_b,         learning).
  853prism_flag_affects( mcmc_e,         learning).
  854prism_flag_affects( mcmc_s,         learning).
  855
  856% --------------------- INFORMATION ---------------------------

prism_show(+Subsystem:oneof([values,probs,counts,goals,flags])) is det
prism_show(-Subsystem:oneof([values,probs,counts,goals,flags,stats])) is nondet: Causes PRISM to print information to its output stream. By default, this appears in a file 'prism.log' in the directory where PRISM was started.

  864prism_show(values) :- prism_flush(show_values).
  865prism_show(probs)  :- prism_flush(show_sw).
  866prism_show(counts) :- prism_flush(show_sw_a).
  867prism_show(goals)  :- prism_flush(show_goals).
  868prism_show(flags)  :- prism_flush(show_prism_flags).
  869prism_show(stats)  :- prism_flush(prism_statistics).
  870
  871prism_flush(G) :- prism(G), prism(flush_output).

prism_statistics(+StatName, -Value) is semidet
prism_statistics(-StatName, -Value) is nondet: Get values of various statistics on inference, learning, and the explanation graph.

  878prism_statistics(S,V) :- prism_statistics(_,S,V).

prism_statistics(+Subsytem, +StatName, -Value) is semidet
prism_statistics(-Subsystem, -StatName, -Value) is nondet: Get values of various statistics on a PRISM subsytem, which can be one of infer, learn or graph.

  885prism_statistics(infer,S,V) :- prism_nd(infer_statistics(S,V)).
  886prism_statistics(learn,S,V) :- prism_nd(learn_statistics(S,V)).
  887prism_statistics(graph,S,V) :- prism_nd(graph_statistics(S,V)).
  888
  889
  890% ------------------ SUPPORT PREDICATES -----------------------
  891
  892
  893send(Term) :-
  894	current_prism(_,Stream,_),
  895	write_canon(Stream,Term),
  896	write(Stream,'.\n'),
  897	flush_output(Stream).
  898
  899recv(Term) :-
  900	current_prism(_,_,Stream),
  901	read_line_to_string(Stream, Input),
  902	parse(Term,Input).
  903
  904parse(line(Input),Input).
  905parse(term(T,V),Input) :-
  906	(	catch(term_string(T,Input,[variable_names(V)]),_,fail) -> true
  907	;	(string_concat("query",_,Input) -> parse_bug(Input,T); true),
  908		format('PRISM> ~s\n',[Input]),fail).
  909
  910parse_bug(Input,Term) :-
  911	term_string(Parse, Input, []),
  912	assert(parse_bug(Input,Parse,Term)).
  913
  914compare(A,B) :- atomic(A), !, (A=B -> true; writeln(mismatch(A,B)), fail).
  915compare(A,B) :- 
  916	A =.. [FA|AA], B=.. [FB|BB], 
  917	compare(FA,FB),
  918	maplist(compare,AA,BB).
  919	
  920read_line(S,Line) :-
  921	wait_for_input([S],[_],0.001),
  922	read_line_to_codes(S,Line).
  923
  924wait(TO) :-
  925	current_prism(_,_,Stream),
  926	wait_for_input([Stream],[_],TO).
  927
  928% PRISM can't read 1e7. It wants 1.0e7.
  929% Therefore we have to reimplement write_canonical
  930write_canon(S,T) :- 
  931	(	numbervars(T,0,_,[singletons(true)]), 
  932		(canon(T,Codes,[]) -> true; throw(error(write_canon(T)))), 
  933		% format('~s',[Codes]), nl,
  934		format(S,'~s',[Codes]), fail
  935	;	true).
  936
  937wrq(A,C1,C2) :- format(codes(C1,C2),'~q',[A]).
  938
  939canon('$VAR'(N)) --> !, wrq('$VAR'(N)). 
  940canon([])        --> !, "[]".
  941canon([H|T])     --> !, "[", canon(H), comma_list(T), "]".
  942canon(X)         --> {float(X)}, !, {format(codes(C,T),'~15g',[X])}, float_codes(C,T).
  943canon(A)         --> {atomic(A)}, !, wrq(A).
  944canon(A)         --> {A=..[F,H|T]}, wrq(F), "(", canon(H), comma_list(T), ")".  
  945
  946-(A,B,A,B).
  947float_codes(H,U)       --> {H==U}, !.
  948float_codes([46|T],U)  --> !, ".", T-U.
  949float_codes([101|T],U) --> !, ".0e", T-U.
  950float_codes([C|T],U)   --> [C], float_codes(T,U).
  951
  952comma_list('$VAR'(N)) --> !, " | ", wrq('$VAR'(N)).
  953comma_list([])        --> !, {true}.
  954comma_list([H|T])     --> ", ", canon(H), comma_list(T)