/* 

  Bayesian null hypothesis test in Picat.

  From https://mhtess.github.io/bdappl/chapters/04-hypothesisTesting.html
  """
  '''
  Let’s see an example: Consider our children helping study from the previous study. 
  Rather than estimate the propensity_to_help, we might want to perform a 
  "Null Hypothesis" test: What is the probability that the behavior observed in 
  our experiment (that 15 out of 20 kids helped) could be explained by true 
  randomness (i.e., propensity_to_help = 0.5)?  
 
  ...
 
  Bayesian Null Hypothesis Testing requires you to specify both hypotheses: 
  It is not enough to define what the null hypothesis is. The alternative 
  hypothesis is intuitively M1:θ≠0.5 but we must be more specific. One standard 
  way of defining the null hypothesis is to put say that the parameter is 
  drawn from a uniform distribution.(Note that this was the model from the 
  previous chapter.)
  """

  Cf my Gamble model gamble_bayesian_null_hypothesis_test.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 : the better model
  Probabilities:
  alternative: 0.7694704049844237
  null: 0.2305295950155763
  mean = [alternative = 0.76947,null = 0.23053]

*/
go ?=>
  reset_store,
  run_model(10_000,$model,[show_probs_trunc,mean]),
  nl,
  % fail,
  nl.
go => true.
  
  
model() =>
  K = 15,
  N = 20,
  NullModel = 0.5,
  AlternativeModel = uniform(0,1),
  
  TheBetterModel = uniform_draw([null,alternative]),

  P = cond(TheBetterModel == null,
           NullModel,
           AlternativeModel),

  observe(binomial_dist(N, P) == K),
  
  if observed_ok then
    add("the better model",TheBetterModel),
  end.