/* 

  Covid probabilities in Picat.

  From Howie Hua https://twitter.com/howie_hua/status/1421502809862664197
  """
  New TikTok video: Doing my part in helping people understand the difference 
  between P(vacc|infected) and P(infected|vacc)
  """

  Population of 100 people where 90 are vaccinated and 10 are not vaccinated.
  There are 4 infections:
  - 3 that are vaccinated
  - 1 that are not vaccinated

  Cf my Gamble model gamble_covid_prob.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 : p(infected|not_vaccinated)
  Probabilities:
  false: 0.8983000000000000
  true: 0.1017000000000000
  mean = [false = 0.8983,true = 0.1017]

  var : p(infected|vaccinated)
  Probabilities:
  false: 0.9686000000000000
  true: 0.0314000000000000
  mean = [false = 0.9686,true = 0.0314]

  var : p(vaccinated|infected)
  Probabilities:
  true: 0.7500000000000000
  false: 0.2500000000000000
  mean = [true = 0.75,false = 0.25]


*/
go ?=>
  reset_store,
  run_model(10_000,$model,[show_probs_trunc,mean]),
  nl,
  % show_store_lengths,
  % fail,
  nl.
go => true.
  
  
model() =>
  N = 100,
  NumVaccinated = 90,
  TotalInfected = 4,
  InfectedInVaccinated = 3,
 
  %  P(vaccinated|infected): 3/4
  PVaccinatedGivenInfected = flip(InfectedInVaccinated / TotalInfected),
   
  %  P(infected|vaccinated): 3/90
  PInfectedGivenVaccinated = flip(InfectedInVaccinated / NumVaccinated),
   
  %  P(infected|not vaccinated): 1/10
  PInfectedGivenNotVaccinated = flip((TotalInfected - InfectedInVaccinated) / (N - NumVaccinated)),

  add("p(vaccinated|infected)",PVaccinatedGivenInfected),
  add("p(infected|vaccinated)",PInfectedGivenVaccinated),
  add("p(infected|not_vaccinated)",PInfectedGivenNotVaccinated).