1/* 2Mixture of two Gaussians. A biased coin is thrown, if it lands heads X in mix(X) 3is sampled from a Gaussian with mean 0 and variance 1. if it lands tails X is 4sampled from a Gaussian with mean 5 and variance 2. 5The example illustrates the use of continuous random variables and 6the use of sampling, including 7rejection sampling and Metropolis/Hastings. Moreover the example 8illustrates the use of the predicate histogram/3 for graphing the 9probability density function of continuous random variables. 10*/ 11:- use_module(library(mcintyre)). 12 13:- if(current_predicate(use_rendering/1)). 14:- use_rendering(c3). 15:- endif. 16:- mc. 17:- begin_lpad. 18 19heads:0.6;tails:0.4. 20% a coin is thrown. The coin is biased: with probability 0.6 it lands heads, 21% with probabiity 0.4 it lands tails 22 23g ~ uniform_dens(0, 1). 24% X in g(X) follows a Gaussian distribution with mean 0 and variance 1 25h ~ uniform_dens(5, 2). 26% X in h(X) follows a Gaussian distribution with mean 5 and variance 2 27 28mix(X) := heads, g ~= X. 29% if the coin lands heads, X in mix(X) is given by g(X) 30mix(X) := tails, h ~= X. 31% if the coin lands tails, X in mix(X) is given by h(X) 32 33:- end_lpad. 34 35hist_uncond(Samples,NBins,Chart):- 36 mc_sample_arg(mix(X),Samples,X,L0), 37 histogram(L0,Chart,[nbins(NBins)]). 38% take SAmples samples of X in mix(X) and draw a histogram with NBins bins representing 39% the probability density of X 40 41hist_rej_heads(Samples,NBins,Chart):- 42 mc_rejection_sample_arg(mix(X),heads,Samples,X,L0), 43 histogram(L0,Chart,[nbins(NBins)]). 44% take Samples samples of X in mix(X) given that heads was true using 45% rejection sampling and draw an 46% histogram with NBins bins representing the probability density of X 47 48hist_mh_heads(Samples,Lag,NBins,Chart):- 49 mc_mh_sample_arg(mix(X),heads,Samples,Lag,X,L0), 50 histogram(L0,Chart,[nbins(NBins)]). 51% take Samples samples of X in mix(X) given that heads was true using 52% Metropolis-Hastings and draw an 53% histogram with NBins bins representing the probability density of X 54 55hist_rej_dis(Samples,NBins,Chart):- 56 mc_rejection_sample_arg(mix(X),(mix(Y),Y>2),Samples,X,L0), 57 histogram(L0,Chart,[nbins(NBins)]). 58% take Samples samples of X in mix(X) given that X>2 was true using 59% rejection sampling and draw an 60% histogram with NBins bins representing the probability density of X 61 62hist_mh_dis(Samples,Lag,NBins,Chart):- 63 mc_mh_sample_arg(mix(X),(mix(Y),Y>2),Samples,Lag,X,L0), 64 histogram(L0,Chart,[nbins(NBins)]). 65% take Samples samples of X in mix(X) given that X>2 was true using 66% Metropolis-Hastings and draw an 67% histogram with NBins bins representing the probability density of X
?-
hist_uncond(1000,40,G)
. % take 1000 samples of X inmix(X)
and draw a histogram with 40 bins representing % the probability density of X ?-mc_sample_arg(mix(X),1000,X,L)
,histogram(L,40,Chart)
. % take 1000 samples of X inmix(X)
and draw a histogram with 40 bins representing % the probability density of X ?-mc_expectation(mix(X),1000,X,E)
. % E=2.017964749114414 ?-hist_rej_heads(1000,40,G)
. % take 1000 samples of X inmix(X)
given that heads was true using % rejection sampling and draw an % histogram with 40 bins representing the probability density of X ?-hist_mh_heads(1000,2,40,G)
. % take 1000 samples of X inmix(X)
given that heads was true using % Metropolis-Hastings and draw an % histogram with 40 bins representing the probability density of X ?-mc_mh_expectation(mix(X),heads,1000,2,X,E)
. % E=-0.018433307290594284 ?-hist_rej_dis(1000,40,G)
. % take 1000 samples of X inmix(X)
given that X>2 was true using % rejection sampling and draw an % histogram with 40 bins representing the probability density of X ?-hist_mh_dis(1000,2,40,G)
. % take 1000 samples of X inmix(X)
given that X>2 was true using % Metropolis-Hastings and draw an % histogram with 40 bins representing the probability density of X ?-mc_mh_expectation(mix(X),(mix(Y),Y>2),1000,2,X,E)
.*/