/* 

  Poisson mean inference in Picat.

  From the SPPL model poisson_mean_inference.pynb

  Here is the PDF of Y = 0..20 (exact probabilities) from the SPPL model:
  """
  (0, 0.19572329209779615),
  (1, 0.24541518485806976),
  (2, 0.15633363576277434),
  (3, 0.07275764007128807),
  (4, 0.03749784880330323),
  (5, 0.03225414602603766),
  (6, 0.036272539178041456),
  (7, 0.03919655164760693),
  (8, 0.03782130413444809),
  (9, 0.032547494609049904),
  (10, 0.025221805385485633),
  (11, 0.017769927037907036),
  (12, 0.011477227802063514),
  (13, 0.006844292955361095),
  (14, 0.00379364785941744),
  (15, 0.0019704640279916506),
  (16, 0.000975431521715664),
  (17, 0.0004844975064651556),
  (18, 0.0002799391781190319),
  (19, 0.00023526933546341676)
  (20, 0.0002864142016352986)
  """


  Cf my Gamble model gamble_poisson_mean_inference.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.

main => go.

/*

  var : x
  ...
  mean = 0.0125491

  var : y
  Probabilities:
  1: 0.2451000000000000
  0: 0.1852000000000000
  2: 0.1676000000000000
  3: 0.1033000000000000
  4: 0.0622000000000000
  5: 0.0398000000000000
  6: 0.0320000000000000
  7: 0.0244000000000000
  8: 0.0190000000000000
  9: 0.0158000000000000
  10: 0.0117000000000000
  11: 0.0106000000000000
  12: 0.0088000000000000
  13: 0.0079000000000000
  14: 0.0059000000000000
  16: 0.0056000000000000
  17: 0.0052000000000000
  15: 0.0052000000000000
  19: 0.0045000000000000
  18: 0.0040000000000000
  22: 0.0036000000000000
  21: 0.0032000000000000
  23: 0.0028000000000000
  24: 0.0025000000000000
  27: 0.0022000000000000
  20: 0.0021000000000000
  32: 0.0014000000000000
  25: 0.0014000000000000
  28: 0.0013000000000000
  26: 0.0013000000000000
  30: 0.0012000000000000
  29: 0.0011000000000000
  31: 0.0010000000000000
  37: 0.0009000000000000
  33: 0.0007000000000000
  41: 0.0006000000000000
  39: 0.0006000000000000
  49: 0.0005000000000000
  44: 0.0005000000000000
  35: 0.0005000000000000
  61: 0.0004000000000000
  47: 0.0004000000000000
  43: 0.0004000000000000
  38: 0.0004000000000000
  36: 0.0004000000000000
  34: 0.0004000000000000
  53: 0.0003000000000000
  51: 0.0003000000000000
  45: 0.0003000000000000
  40: 0.0003000000000000
  76: 0.0002000000000000
  65: 0.0002000000000000
  64: 0.0002000000000000
  60: 0.0002000000000000
  50: 0.0002000000000000
  48: 0.0002000000000000
  46: 0.0002000000000000
  42: 0.0002000000000000
  116: 0.0001000000000000
  97: 0.0001000000000000
  80: 0.0001000000000000
  79: 0.0001000000000000
  75: 0.0001000000000000
  73: 0.0001000000000000
  72: 0.0001000000000000
  69: 0.0001000000000000
  66: 0.0001000000000000
  63: 0.0001000000000000
  62: 0.0001000000000000
  59: 0.0001000000000000
  57: 0.0001000000000000
  56: 0.0001000000000000
  55: 0.0001000000000000
  54: 0.0001000000000000
  mean = 4.0863

*/
go ?=>
  reset_store,
  run_model(10_000,$model,[show_probs,mean]),
  nl,
  % show_store_lengths,
  % fail,
  nl.
go => true.
  
  
model() =>
  X = normal_dist(0,1),
  observe(X >= -4, X <=4 ),
  
  Y = poisson_dist(1 + abs(pow(X,3))*2),
  
  add("x",X),
  add("y",Y).