1/*
    2Posterior estimation in Bayesian models.
    3We are trying to estimate the true value of a Gaussian distributed random
    4variable, given some observed data. The variance is known (2) and we
    5suppose that the mean has a Gaussian distribution with mean 1 and variance
    65. We take different measurement (e.g. at different times), indexed
    7with an integer.
    8Given that we observe 9 and 8 at indexes 1 and 2, how does the distribution
    9of the random variable (value at index 0) changes with respect to the case of
   10no observations?
   11From
   12http://www.robots.ox.ac.uk/~fwood/anglican/examples/viewer/?worksheet=gaussian-posteriors
   13*/
   14:- use_module(library(mcintyre)).   15
   16:- if(current_predicate(use_rendering/1)).   17:- use_rendering(c3).   18:- endif.   19
   20
   21
   22
   23:- mc.   24:- begin_lpad.   25
   26val(I,X) :-
   27  mean(M),
   28  val(I,M,X).
   29% at time I we see X sampled from a Gaussian with mean M and variamce 2.0
   30
   31mean(M): user(M,gauss(1.0, 5.0)).
   32% Gaussian distribution of the mean of the Gaussian of the variable
   33
   34val(_,M,X): user(X,gauss(M, 2.0)).
   35% Gaussian distribution of the variable
   36
   37
   38:- end_lpad.   39
   40gauss(Mean,Variance,S):-
   41  number(Mean),!,
   42  random(U1),
   43  random(U2),
   44  R is sqrt(-2*log(U1)),
   45  Theta is 2*pi*U2,
   46  S0 is R*cos(Theta),
   47  StdDev is sqrt(Variance),
   48  S is StdDev*S0+Mean.
   49
   50gauss(Mean,Variance,S,D):-
   51  StdDev is sqrt(Variance),
   52  D is 1/(StdDev*sqrt(2*pi))*exp(-(S-Mean)*(S-Mean)/(2*Variance)).
   53
   54hist_uncond(Samples,NBins,Chart):-
   55  mc_sample_arg(val(0,X),Samples,X,L0),
   56  histogram(L0,Chart,[nbins(NBins)]).
   57% plot an histogram of the density of the random variable before any
   58% observations by taking Samples samples and by dividing the domain
   59% in NBins bins
   60
   61dens_lw(Samples,NBins,Chart):-
   62  mc_sample_arg(val(0,Y),Samples,Y,L0),
   63  mc_lw_sample_arg(val(0,X),(val(1,9),val(2,8)),Samples,X,L),
   64  densities(L0,L,Chart,[nbins(NBins)]).
   65% plot the densities of the random variable before and after
   66% observing 9 and 8 by taking Samples samples using likelihood weighting
   67% and by dividing the domain
   68% in NBins bins
   69
   70dens_part(Samples,NBins,Chart):-
   71  mc_sample_arg(val(0,Y),Samples,Y,L0),
   72  mc_particle_sample_arg(val(0,X),[val(1,9),val(2,8)],Samples,X,L),
   73  densities(L0,L,Chart,[nbins(NBins)]).
   74% plot the densities of the random variable before and after
   75% observing 9 and 8 by taking Samples samples using particle filtering
   76% and by dividing the domain
   77% in NBins bins