/* 

  Bugs book example 2.7.1 in Picat.

  Example 2.7.1 Surgery (continued): prediction

  Cf my Gamble model gamble_bugs_book_2_7_1.rkt

  This program was created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Picat page: http://www.hakank.org/picat/

*/

import ppl_distributions, ppl_utils.
import util.
% import ordset.

main => go.

/*

  var : p6
  Probabilities:
  false: 0.9540000000000000
  true: 0.0460000000000000
  mean = [false = 0.954,true = 0.046]
  HPD intervals:
  show_hpd_intervals: data is not numeric

  var : y
  Probabilities (truncated):
  1: 0.2920000000000000
  2: 0.2000000000000000
  0: 0.1910000000000000
  3: 0.1370000000000000
  .........
  7: 0.0080000000000000
  8: 0.0060000000000000
  9: 0.0030000000000000
  11: 0.0010000000000000
  mean = 1.99
  HPD intervals:
  HPD interval (0.84): 0.00000000000000..4.00000000000000

  var : theta
  Probabilities (truncated):
  0.38222921182783: 0.0010000000000000
  0.31779901837293: 0.0010000000000000
  0.31377295882435: 0.0010000000000000
  0.310500683865886: 0.0010000000000000
  .........
  0.011423218919007: 0.0010000000000000
  0.010191017942465: 0.0010000000000000
  0.004757395106247: 0.0010000000000000
  0.004503938265627: 0.0010000000000000
  mean = 0.0984157
  HPD intervals:
  HPD interval (0.84): 0.02393567223937..0.15888419988151

*/
go ?=>
  reset_store,
  run_model(1_000,$model,[show_probs_trunc,mean,show_hpd_intervals,hpd_intervals=[0.84]]),
  nl,
  % show_store_lengths,nl,
  % fail,
  nl.
go => true.
  
  
model() =>
  Theta = beta_dist(3,27),     % Prior distribution
  Y = binomial_dist(20,Theta), % Sampling distribution
  P6 = check(Y >= 6),
  
  add("p6",P6),
  add("y",Y),
  add("theta",Theta).