```    1/*
2Factory producing widgets.
3Consider a factory with two machines a and b. Each machine produces a widget
4with a continuous feature. A widget is produced by machine a with probability
50.7 and by machine b with probability b.
6If the widget is produced by machine a, the feature is distributed as a
7Gaussian with mean 2.0 and variance 1.0.
8If the widget is produced by machine b, the feature is distributed as a
9Gaussian with mean 3.0 and variance 1.0.
10The widget then is processed by a third machine that adds a random quantity to
11the feature distributed as a Gaussian with mean 0.5 and variance 0.1.
12What is the distribution of the feature?
13What is the distribution of the feature given that the widget was procuded
14by machine a?
15What is the distribution of the feature given that the third machine added a
16quantity greater than 0.2?
17What is the distribution of the feature given that the third machine added
18a quantity of 2.0?
20Islam, Muhammad Asiful, C. R. Ramakrishnan, and I. V. Ramakrishnan.
21"Inference in probabilistic logic programs with continuous random variables."
22Theory and Practice of Logic Programming 12.4-5 (2012): 505-523.
23http://arxiv.org/pdf/1112.2681v3.pdf
24*/
25:- use_module(library(mcintyre)).   26:- use_module(library(clpr)).   27
28:- if(current_predicate(use_rendering/1)).   29:- use_rendering(c3).   30:- endif.   31:- mc.   32:- begin_lpad.   33
34widget(X) :-
35  machine(M),
36  st(M,Z),
37  pt(Y),
38  {X =:= Y + Z}.
39
40machine(a):0.3;machine(b):0.7.
41
42st(a,X): gaussian(X,2.0, 1.0).
43
44st(b,X): gaussian(X,3.0, 1.0).
45
46pt(X): gaussian(X,0.5, 0.1).
47
48
51hist_uncond(Samples,NBins,Chart):-
52  mc_sample_arg(widget(X),Samples,X,L0),
53  histogram(L0,Chart,[nbins(NBins)]).
54% What is the distribution of the feature?
55
56hist_rej_macha(Samples,NBins,Chart):-
57  mc_rejection_sample_arg(widget(X),machine(a),Samples,X,L0),
58  histogram(L0,Chart,[nbins(NBins)]).
59% What is the distribution of the feature given that the widget was procuded
60% by machine a, computed by taking Samples samples with rejection sampling and
61% drawing a histogram with NBins bins?
62
63hist_mh_macha(Samples,Lag,NBins,Chart):-
64  mc_mh_sample_arg(widget(X),machine(a),Samples,X,L0,[lag(Lag)]),
65  histogram(L0,Chart,[nbins(NBins)]).
66% What is the distribution of the feature given that the widget was procuded
67% by machine a, computed by taking Samples samples with Metropolis-Hastings
68% (lag=Lag) and drawing a histogram with NBins bins?
69
70hist_rej_dis(Samples,NBins,Chart):-
71  mc_rejection_sample_arg(widget(X),(pt(Y),Y>0.2),Samples,X,L0),
72  histogram(L0,Chart,[nbins(NBins)]).
73% What is the distribution of the feature given that the third machine added a
74% quantity greater than 0.2, computed by taking Samples samples with rejection
75% sampling and
76% drawing a histogram with NBins bins?
77
78hist_mh_dis(Samples,Lag,NBins,Chart):-
79  mc_mh_sample_arg(widget(X),(pt(Y),Y>0.2),Samples,X,L0,[lag(Lag)]),
80  histogram(L0,Chart,[nbins(NBins)]).
81% What is the distribution of the feature given that the third machine added a
82% quantity greater than 0.2, computed by taking Samples samples with
83% Metropolis-Hastings and
84% drawing a histogram with NBins bins?
85
86hist_lw(Samples,NBins,Chart):-
87  mc_sample_arg(widget(Y),Samples,Y,L0),
88  mc_lw_sample_arg(widget(X),pt(2.0),Samples,X,L),
89  densities(L0,L,Chart,[nbins(NBins)]).
90% What is the distribution of the feature given that the third machine added
91% a quantity of 2.0, computed by taking Samples samples with likelihood weighting
92% and drawing a density with NBins bins?
```

?- `hist_lw(1000,40,G)`. % What is the distribution of the feature given that the third machine added % a quantity of 2.0, computed by taking 1000 samples with likelihood weighting % and drawing a density with 40 bins?

?- `hist_uncond(10000,40,G)`. % What is the distribution of the feature?

?- `hist_rej_macha(10000,40,G)`. % What is the distribution of the feature given that the widget was procuded % by machine a, computed by taking 10000 samples with rejection sampling and % drawing a histogram with 40 bins?

?- `hist_mh_macha(10000,2,40,G)`. % What is the distribution of the feature given that the widget was procuded % by machine a, computed by taking 10000 samples with Metropolis-Hastings % (lag=Lag) and drawing a histogram with 40 bins?

?- `hist_rej_dis(10000,40,G)`. % What is the distribution of the feature given that the third machine added a % quantity greater than 0.2, computed by taking 10000 samples with rejection % sampling and % drawing a histogram with 40 bins?

?- `hist_mh_dis(10000,2,40,G)`. % What is the distribution of the feature given that the third machine added a % quantity greater than 0.2, computed by taking 10000 samples with % Metropolis-Hastings and % drawing a histogram with 40 bins?

*/