/* nag_chi_sq_2_way_table (g11aac) Example Program.

 *

 * Copyright 1996 Numerical Algorithms Group.

 *

 * Mark 4, 1996.

 * Mark 8 revised, 2004.

 *

 */



#include <nag.h>

#include <stdio.h>

#include <nag_stdlib.h>

#include <nagg11.h>



#define EXPT(I,J) expt[(I)*tdt + J]

#define CHIST(I,J) chist[(I)*tdt + J]

#define NOBST(I,J) nobst[(I)*tdt + J]

int main(void)

{



  Integer exit_status=0, i, j, ncol, *nobst=0, nrow, tdt;

  NagError fail;

  double chi, *chist=0, df, *expt=0, g, prob;



  INIT_FAIL(fail);

  Vprintf("nag_chi_sq_2_way_table (g11aac) Example Program Results\n\n");



  /* Skip heading in data file */

  Vscanf("%*[^\n] ");



  Vscanf("%ld %ld %*[^\n] ", &nrow, &ncol);

  if (nrow>=2  && ncol>=2)

    {

      if ( !( expt = NAG_ALLOC((nrow)*(ncol), double)) ||

           !( chist = NAG_ALLOC((nrow)*(ncol), double))  ||

           !( nobst = NAG_ALLOC((nrow)*(ncol), Integer)) )

        {

          Vprintf("Allocation failure\n");

          exit_status = -1;

          goto END;

        }

      tdt = ncol;

    }

  else

    {

      Vprintf("Invalid nrow or ncol.\n");

      exit_status = 1;

      return exit_status;

    }

  for (i = 0; i < nrow; ++i)

    {

      for (j = 0; j < ncol; ++j)

        Vscanf("%ld", &NOBST(i,j));

      Vscanf("%*[^\n]");

    }



  /* nag_chi_sq_2_way_table (g11aac).

   * chi^2 statistics for two-way contingency table

   */

  nag_chi_sq_2_way_table(nrow, ncol, nobst, tdt, expt,

                         chist, &prob, &chi, &g, &df, &fail);

  if (fail.code != NE_NOERROR)

    {

      Vprintf("Error from nag_chi_sq_2_way_table (g11aac).\n%s\n",

              fail.message);

      exit_status = 1;

      goto END;

    }



  Vprintf("Probability = %6.4f\n", prob);

  Vprintf("Pearson Chi-square statistic = %8.3f\n", chi);

  Vprintf("Likelihood ratio test statistic = %8.3f\n", g);

  Vprintf("Degrees of freedom = %4.0f\n", df);

 END:

  if (expt) NAG_FREE(expt);

  if (chist) NAG_FREE(chist);

  if (nobst) NAG_FREE(nobst);

  return exit_status;

}