/* 

  Chi square test in Picat.

  From Mathematica ChiSquareDistribution
  """
  The weight in grams of a particular boxed cereal product is known to 
  follow a normal distribution. A quality assurance team samples 15 boxes at 
  random and records their weights. Test the hypothesis that the standard 
  deviation of the product weight is less than 36:

  weights = 
  Quantity[{367.9, 384.7, 353.8, 334.7, 450.9, 390.6, 422.6, 352.2, 
    330.9, 342.0, 388.9, 386.7, 374.2, 388.5, 382.2}, "Grams"];

  {Mean[weights], StandardDeviation[weights]}
  -> {Quantity[376.72, "Grams"], Quantity[32.2513, "Grams"]}

  Under the null hypothesis of Sigma^2>=36^2, the following statistic 
  follows ChiSquareDistribution
 
  Sigma0 = Quantity[36, "Grams"];
  Chi2stat = (15 - 1) Variance[weights]/Sigma0^2
  ->11.2362

  The null hypothesis cannot be rejected at the 5% level:
  InverseCDF[ChiSquareDistribution[14], 0.05]
  -> 6.57063

  Chi2stat < InverseCDF[ChiSquareDistribution[14], 0.05]
  -> False
  """

  Cf my Gamble model gamble_chi_squared_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.
% import ordset.

main => go.

/*
  weights = [367.9,384.7,353.8,334.7,450.9,390.6,422.6,352.2,330.9,342.0,388.9,386.7,374.2,388.5,382.2]
  [len = 15,min = 330.9,mean = 376.72,median = 382.2,max = 450.9,variance = 970.804,stdev = 31.1577]
  [sigma_theta = 36,chi2_stat = 11.2362]
  quantile0_05 = 6.57063
  null_hypotheses_tested = false


*/
go ?=>
  Weights = [367.9,384.7,353.8,334.7,450.9,390.6,422.6,352.2,
             330.9,342.0,388.9,386.7,374.2,388.5,382.2],
  println(weights=Weights),
  show_simple_stats(Weights),
  SigmaTheta = 36,
  Chi2Stat = Weights.len * Weights.variance / SigmaTheta**2,
  println(sigma_theta=SigmaTheta),
  println(chi2_stat=Chi2Stat),
  Q0_05 = chi_square_dist_quantile(Weights.len - 1, 0.05),
  println(quantile0_05=Q0_05),
  println(null_hypotheses_tested=(cond(Chi2Stat < Q0_05,true,false))),
  nl.