20 :- use_module(library(mcintyre)). 21:- if(current_predicate(use_rendering/1)). 22:- use_rendering(c3). 23:- endif. 24:- mc. 25:- begin_lpad. 26
27dp_n_values(N,N,_Alpha,[]):-!.
28
29dp_n_values(N0,N,Alpha,[[V]-1|Vs]):-
30 N0<N,
31 dp_value(N0,Alpha,V),
32 N1 is N0+1,
33 dp_n_values(N1,N,Alpha,Vs).
34
35dp_value(NV,Alpha,V):-
36 dp_stick_index(NV,Alpha,I),
37 dp_pick_value(I,NV,V).
38
39dp_pick_value(I,NV,V):-
40 ivar(I,IV),
41 Var is 1.0/IV,
42 mean(I,Var,M),
43 value(NV,M,Var,V).
44
45ivar(_,IV):gamma(IV,1,0.1).
46
47mean(_,V0,M):gaussian(M,0,V):-V is V0*30.
48
49value(_,M,V,Val):gaussian(Val,M,V).
50
51dp_stick_index(NV,Alpha,I):-
52 dp_stick_index(1,NV,Alpha,I).
53
54dp_stick_index(N,NV,Alpha,V):-
55 stick_proportion(N,Alpha,P),
56 choose_prop(N,NV,Alpha,P,V).
57
58choose_prop(N,NV,_Alpha,P,N):-
59 pick_portion(N,NV,P).
60
61choose_prop(N,NV,Alpha,P,V):-
62 neg_pick_portion(N,NV,P),
63 N1 is N+1,
64 dp_stick_index(N1,NV,Alpha,V).
65
66
67
68stick_proportion(_,Alpha,P):beta(P,1,Alpha).
69
70pick_portion(N,NV,P):P;neg_pick_portion(N,NV,P):1-P.
71
72:- end_lpad. 73
74
75obs([-1,7,3]).
76
77prior(Samples,NBins,Chart):-
78 mc_sample_arg_first(dp_n_values(0,Samples,10.0,V),1,V,L),
79 L=[Vs-_],
80 histogram(Vs,Chart,[nbins(NBins)]).
81
82post(Samples,NBins,Chart):-
83 obs(O),
84 maplist(to_val,O,O1),
85 length(O1,N),
86 mc_lw_sample_arg_log(dp_value(0,10.0,T),dp_n_values(0,N,10.0,O1),Samples,T,L),
87 maplist(keys,L,LW),
88 min_list(LW,Min),
89 maplist(exp(Min),L,L1),
90 density(L1,Chart,[nbins(NBins),min(-8),max(15)]).
91
92keys(_-W,W).
93
94exp(Min,L-W,L-W1):- W1 is exp(W-Min).
95
96to_val(V,[V]-1)
?-
prior(200,100,G)
. % draw the prior density ?-post(200,100,G)
. % draw the posterior density*/