src/language/stats/examine.q

   1 /* PSPP - a program for statistical analysis.
   2    Copyright (C) 2004 Free Software Foundation, Inc.
   3
   4    This program is free software: you can redistribute it and/or modify
   5    it under the terms of the GNU General Public License as published by
   6    the Free Software Foundation, either version 3 of the License, or
   7    (at your option) any later version.
   8
   9    This program is distributed in the hope that it will be useful,
  10    but WITHOUT ANY WARRANTY; without even the implied warranty of
  11    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  12    GNU General Public License for more details.
  13
  14    You should have received a copy of the GNU General Public License
  15    along with this program.  If not, see <http://www.gnu.org/licenses/>. */
  16
  17 #include <config.h>
  18
  19 #include <gsl/gsl_cdf.h>
  20 #include <libpspp/message.h>
  21 #include <math.h>
  22 #include <stdio.h>
  23 #include <stdlib.h>
  24
  25 #include <data/case.h>
  26 #include <data/casegrouper.h>
  27 #include <data/casereader.h>
  28 #include <data/dictionary.h>
  29 #include <data/procedure.h>
  30 #include <data/value-labels.h>
  31 #include <data/variable.h>
  32 #include <language/command.h>
  33 #include <language/dictionary/split-file.h>
  34 #include <language/lexer/lexer.h>
  35 #include <libpspp/compiler.h>
  36 #include <libpspp/hash.h>
  37 #include <libpspp/message.h>
  38 #include <libpspp/misc.h>
  39 #include <libpspp/str.h>
  40 #include <math/factor-stats.h>
  41 #include <math/moments.h>
  42 #include <math/percentiles.h>
  43 #include <output/charts/box-whisker.h>
  44 #include <output/charts/cartesian.h>
  45 #include <output/manager.h>
  46 #include <output/table.h>
  47
  48 #include "minmax.h"
  49 #include "xalloc.h"
  50
  51 #include "gettext.h"
  52 #define _(msgid) gettext (msgid)
  53 #define N_(msgid) msgid
  54
  55 /* (headers) */
  56 #include <output/chart.h>
  57 #include <output/charts/plot-hist.h>
  58 #include <output/charts/plot-chart.h>
  59
  60 /* (specification)
  61    "EXAMINE" (xmn_):
  62    *^variables=custom;
  63    +total=custom;
  64    +nototal=custom;
  65    missing=miss:pairwise/!listwise,
  66            rep:report/!noreport,
  67            incl:include/!exclude;
  68    +compare=cmp:variables/!groups;
  69    +percentiles=custom;
  70    +id=var;
  71    +plot[plt_]=stemleaf,boxplot,npplot,:spreadlevel(*d:n),histogram,all,none;
  72    +cinterval=double;
  73    +statistics[st_]=descriptives,:extreme(*d:n),all,none.
  74 */
  75
  76 /* (declarations) */
  77
  78 /* (functions) */
  79
  80
  81
  82 static struct cmd_examine cmd;
  83
  84 static const struct variable **dependent_vars;
  85
  86 static size_t n_dependent_vars;
  87
  88
  89 struct factor
  90 {
  91   /* The independent variable */
  92   struct variable *indep_var[2];
  93
  94
  95   /* Hash table of factor stats indexed by 2 values */
  96   struct hsh_table *fstats;
  97
  98   /* The hash table after it has been crunched */
  99   struct factor_statistics **fs;
 100
 101   struct factor *next;
 102
 103 };
 104
 105 /* Linked list of factors */
 106 static struct factor *factors = 0;
 107
 108 static struct metrics *totals = 0;
 109
 110 /* Parse the clause specifying the factors */
 111 static int examine_parse_independent_vars (struct lexer *lexer, const struct dictionary *dict, struct cmd_examine *cmd);
 112
 113
 114
 115 /* Output functions */
 116 static void show_summary (const struct variable **dependent_var, int n_dep_var,
 117                          const struct factor *f);
 118
 119 static void show_extremes (const struct variable **dependent_var,
 120                           int n_dep_var,
 121                           const struct factor *factor,
 122                           int n_extremities);
 123
 124 static void show_descriptives (const struct variable **dependent_var,
 125                               int n_dep_var,
 126                               struct factor *factor);
 127
 128 static void show_percentiles (const struct variable **dependent_var,
 129                              int n_dep_var,
 130                              struct factor *factor);
 131
 132
 133
 134
 135 void np_plot (const struct metrics *m, const char *factorname);
 136
 137
 138 void box_plot_group (const struct factor *fctr,
 139                     const struct variable **vars, int n_vars,
 140                     const struct variable *id
 141                     ) ;
 142
 143
 144 void box_plot_variables (const struct factor *fctr,
 145                         const struct variable **vars, int n_vars,
 146                         const struct variable *id
 147                         );
 148
 149
 150
 151 /* Per Split function */
 152 static void run_examine (struct cmd_examine *, struct casereader *,
 153                          struct dataset *);
 154
 155 static void output_examine (void);
 156
 157
 158 void factor_calc (const struct ccase *c, int case_no,
 159                   double weight, int case_missing);
 160
 161
 162 /* Represent a factor as a string, so it can be
 163    printed in a human readable fashion */
 164 static void factor_to_string (const struct factor *fctr,
 165                                const struct factor_statistics *fs,
 166                               const struct variable *var,
 167                               struct string *str
 168                               );
 169
 170 /* Represent a factor as a string, so it can be
 171    printed in a human readable fashion,
 172    but sacrificing some readablility for the sake of brevity */
 173 static void factor_to_string_concise (const struct factor *fctr,
 174                                       const struct factor_statistics *fs,
 175                                       struct string *);
 176
 177
 178
 179
 180 /* Categories of missing values to exclude. */
 181 static enum mv_class exclude_values;
 182
 183 /* PERCENTILES */
 184
 185 static subc_list_double percentile_list;
 186
 187 static enum pc_alg percentile_algorithm;
 188
 189 static short sbc_percentile;
 190
 191
 192 int
 193 cmd_examine (struct lexer *lexer, struct dataset *ds)
 194 {
 195   struct casegrouper *grouper;
 196   struct casereader *group;
 197   bool ok;
 198
 199   subc_list_double_create (&percentile_list);
 200   percentile_algorithm = PC_HAVERAGE;
 201
 202   if ( !parse_examine (lexer, ds, &cmd, NULL) )
 203     {
 204       subc_list_double_destroy (&percentile_list);
 205       return CMD_FAILURE;
 206     }
 207
 208   /* If /MISSING=INCLUDE is set, then user missing values are ignored */
 209   exclude_values = cmd.incl == XMN_INCLUDE ? MV_SYSTEM : MV_ANY;
 210
 211   if ( cmd.st_n == SYSMIS )
 212     cmd.st_n = 5;
 213
 214   if ( ! cmd.sbc_cinterval)
 215     cmd.n_cinterval[0] = 95.0;
 216
 217   /* If descriptives have been requested, make sure the
 218      quartiles are calculated */
 219   if ( cmd.a_statistics[XMN_ST_DESCRIPTIVES] )
 220     {
 221       subc_list_double_push (&percentile_list, 25);
 222       subc_list_double_push (&percentile_list, 50);
 223       subc_list_double_push (&percentile_list, 75);
 224     }
 225
 226   grouper = casegrouper_create_splits (proc_open (ds), dataset_dict (ds));
 227   while (casegrouper_get_next_group (grouper, &group))
 228     run_examine (&cmd, group, ds);
 229   ok = casegrouper_destroy (grouper);
 230   ok = proc_commit (ds) && ok;
 231
 232   if ( totals )
 233     {
 234       free ( totals );
 235     }
 236
 237   if ( dependent_vars )
 238     free (dependent_vars);
 239
 240   {
 241     struct factor *f = factors ;
 242     while ( f )
 243       {
 244         struct factor *ff = f;
 245
 246         f = f->next;
 247         free ( ff->fs );
 248         hsh_destroy ( ff->fstats ) ;
 249         free ( ff ) ;
 250       }
 251     factors = 0;
 252   }
 253
 254   subc_list_double_destroy (&percentile_list);
 255
 256   return ok ? CMD_SUCCESS : CMD_CASCADING_FAILURE;
 257 };
 258
 259
 260
 261 /* Show all the appropriate tables */
 262 static void
 263 output_examine (void)
 264 {
 265   struct factor *fctr;
 266
 267   /* Show totals if appropriate */
 268   if ( ! cmd.sbc_nototal || factors == 0 )
 269     {
 270       show_summary (dependent_vars, n_dependent_vars, 0);
 271
 272       if ( cmd.sbc_statistics )
 273         {
 274           if ( cmd.a_statistics[XMN_ST_EXTREME])
 275             show_extremes (dependent_vars, n_dependent_vars, 0, cmd.st_n);
 276
 277           if ( cmd.a_statistics[XMN_ST_DESCRIPTIVES])
 278             show_descriptives (dependent_vars, n_dependent_vars, 0);
 279
 280         }
 281       if ( sbc_percentile )
 282         show_percentiles (dependent_vars, n_dependent_vars, 0);
 283
 284       if ( cmd.sbc_plot)
 285         {
 286           int v;
 287           if ( cmd.a_plot[XMN_PLT_STEMLEAF] )
 288             msg (SW, _ ("%s is not currently supported."), "STEMLEAF");
 289
 290           if ( cmd.a_plot[XMN_PLT_SPREADLEVEL] )
 291             msg (SW, _ ("%s is not currently supported."), "SPREADLEVEL");
 292
 293           if ( cmd.a_plot[XMN_PLT_NPPLOT] )
 294             {
 295               for ( v = 0 ; v < n_dependent_vars; ++v )
 296                 np_plot (&totals[v], var_to_string (dependent_vars[v]));
 297             }
 298
 299           if ( cmd.a_plot[XMN_PLT_BOXPLOT] )
 300             {
 301               if ( cmd.cmp == XMN_GROUPS )
 302                 {
 303                   box_plot_group (0, (const struct variable **) dependent_vars,
 304                                   n_dependent_vars, cmd.v_id);
 305                 }
 306               else
 307                 box_plot_variables (0,
 308                                     (const struct variable **) dependent_vars,
 309                                     n_dependent_vars, cmd.v_id);
 310             }
 311
 312           if ( cmd.a_plot[XMN_PLT_HISTOGRAM] )
 313             {
 314               for ( v = 0 ; v < n_dependent_vars; ++v )
 315                 {
 316                   struct normal_curve normal;
 317
 318                   normal.N      = totals[v].n;
 319                   normal.mean   = totals[v].mean;
 320                   normal.stddev = totals[v].stddev;
 321
 322                   histogram_plot (totals[v].histogram,
 323                                  var_to_string (dependent_vars[v]),
 324                                  &normal, 0);
 325                 }
 326             }
 327
 328         }
 329
 330     }
 331
 332
 333   /* Show grouped statistics  as appropriate */
 334   fctr = factors;
 335   while ( fctr )
 336     {
 337       show_summary (dependent_vars, n_dependent_vars, fctr);
 338
 339       if ( cmd.sbc_statistics )
 340         {
 341           if ( cmd.a_statistics[XMN_ST_EXTREME])
 342             show_extremes (dependent_vars, n_dependent_vars, fctr, cmd.st_n);
 343
 344           if ( cmd.a_statistics[XMN_ST_DESCRIPTIVES])
 345             show_descriptives (dependent_vars, n_dependent_vars, fctr);
 346         }
 347
 348       if ( sbc_percentile )
 349         show_percentiles (dependent_vars, n_dependent_vars, fctr);
 350
 351
 352       if ( cmd.sbc_plot)
 353         {
 354           size_t v;
 355
 356           struct factor_statistics **fs = fctr->fs ;
 357
 358           if ( cmd.a_plot[XMN_PLT_BOXPLOT] )
 359             {
 360               if ( cmd.cmp == XMN_VARIABLES )
 361                 box_plot_variables (fctr,
 362                                     (const struct variable **) dependent_vars,
 363                                     n_dependent_vars, cmd.v_id);
 364               else
 365                 box_plot_group (fctr,
 366                                 (const struct variable **) dependent_vars,
 367                                 n_dependent_vars, cmd.v_id);
 368             }
 369
 370           for ( v = 0 ; v < n_dependent_vars; ++v )
 371             {
 372
 373               for ( fs = fctr->fs ; *fs ; ++fs )
 374                 {
 375                   struct string str;
 376                   ds_init_empty (&str);
 377                   factor_to_string (fctr, *fs, dependent_vars[v], &str);
 378
 379                   if ( cmd.a_plot[XMN_PLT_NPPLOT] )
 380                     np_plot (& (*fs)->m[v], ds_cstr (&str));
 381
 382                   if ( cmd.a_plot[XMN_PLT_HISTOGRAM] )
 383                     {
 384                       struct normal_curve normal;
 385
 386                       normal.N      = (*fs)->m[v].n;
 387                       normal.mean   = (*fs)->m[v].mean;
 388                       normal.stddev = (*fs)->m[v].stddev;
 389
 390                       histogram_plot ((*fs)->m[v].histogram,
 391                                      ds_cstr (&str) ,  &normal, 0);
 392                     }
 393
 394                   ds_destroy (&str);
 395
 396                 } /* for ( fs .... */
 397
 398             } /* for ( v = 0 ..... */
 399
 400         }
 401
 402       fctr = fctr->next;
 403     }
 404
 405 }
 406
 407
 408 /* Create a hash table of percentiles and their values from the list of
 409    percentiles */
 410 static struct hsh_table *
 411 list_to_ptile_hash (const subc_list_double *l)
 412 {
 413   int i;
 414
 415   struct hsh_table *h ;
 416
 417   h = hsh_create (subc_list_double_count (l),
 418                  (hsh_compare_func *) ptile_compare,
 419                  (hsh_hash_func *) ptile_hash,
 420                  (hsh_free_func *) free,
 421                  0);
 422
 423
 424   for ( i = 0 ; i < subc_list_double_count (l) ; ++i )
 425     {
 426       struct percentile *p = xmalloc (sizeof *p);
 427
 428       p->p = subc_list_double_at (l,i);
 429       p->v = SYSMIS;
 430
 431       hsh_insert (h, p);
 432
 433     }
 434
 435   return h;
 436
 437 }
 438
 439 /* Parse the PERCENTILES subcommand */
 440 static int
 441 xmn_custom_percentiles (struct lexer *lexer, struct dataset *ds UNUSED,
 442                        struct cmd_examine *p UNUSED, void *aux UNUSED)
 443 {
 444   sbc_percentile = 1;
 445
 446   lex_match (lexer, '=');
 447
 448   lex_match (lexer, '(');
 449
 450   while ( lex_is_number (lexer) )
 451     {
 452       subc_list_double_push (&percentile_list, lex_number (lexer));
 453
 454       lex_get (lexer);
 455
 456       lex_match (lexer, ',') ;
 457     }
 458   lex_match (lexer, ')');
 459
 460   lex_match (lexer, '=');
 461
 462   if ( lex_match_id (lexer, "HAVERAGE"))
 463     percentile_algorithm = PC_HAVERAGE;
 464
 465   else if ( lex_match_id (lexer, "WAVERAGE"))
 466     percentile_algorithm = PC_WAVERAGE;
 467
 468   else if ( lex_match_id (lexer, "ROUND"))
 469     percentile_algorithm = PC_ROUND;
 470
 471   else if ( lex_match_id (lexer, "EMPIRICAL"))
 472     percentile_algorithm = PC_EMPIRICAL;
 473
 474   else if ( lex_match_id (lexer, "AEMPIRICAL"))
 475     percentile_algorithm = PC_AEMPIRICAL;
 476
 477   else if ( lex_match_id (lexer, "NONE"))
 478     percentile_algorithm = PC_NONE;
 479
 480
 481   if ( 0 == subc_list_double_count (&percentile_list))
 482     {
 483       subc_list_double_push (&percentile_list, 5);
 484       subc_list_double_push (&percentile_list, 10);
 485       subc_list_double_push (&percentile_list, 25);
 486       subc_list_double_push (&percentile_list, 50);
 487       subc_list_double_push (&percentile_list, 75);
 488       subc_list_double_push (&percentile_list, 90);
 489       subc_list_double_push (&percentile_list, 95);
 490     }
 491
 492   return 1;
 493 }
 494
 495 /* TOTAL and NOTOTAL are simple, mutually exclusive flags */
 496 static int
 497 xmn_custom_total (struct lexer *lexer UNUSED, struct dataset *ds UNUSED, struct cmd_examine *p, void *aux UNUSED)
 498 {
 499   if ( p->sbc_nototal )
 500     {
 501       msg (SE, _ ("%s and %s are mutually exclusive"),"TOTAL","NOTOTAL");
 502       return 0;
 503     }
 504
 505   return 1;
 506 }
 507
 508 static int
 509 xmn_custom_nototal (struct lexer *lexer UNUSED, struct dataset *ds UNUSED,
 510                     struct cmd_examine *p, void *aux UNUSED)
 511 {
 512   if ( p->sbc_total )
 513     {
 514       msg (SE, _ ("%s and %s are mutually exclusive"),"TOTAL","NOTOTAL");
 515       return 0;
 516     }
 517
 518   return 1;
 519 }
 520
 521
 522
 523 /* Parser for the variables sub command
 524    Returns 1 on success */
 525 static int
 526 xmn_custom_variables (struct lexer *lexer, struct dataset *ds, struct cmd_examine *cmd, void *aux UNUSED)
 527 {
 528   const struct dictionary *dict = dataset_dict (ds);
 529   lex_match (lexer, '=');
 530
 531   if ( (lex_token (lexer) != T_ID || dict_lookup_var (dict, lex_tokid (lexer)) == NULL)
 532       && lex_token (lexer) != T_ALL)
 533     {
 534       return 2;
 535     }
 536
 537   if (!parse_variables_const (lexer, dict, &dependent_vars, &n_dependent_vars,
 538                         PV_NO_DUPLICATE | PV_NUMERIC | PV_NO_SCRATCH) )
 539     {
 540       free (dependent_vars);
 541       return 0;
 542     }
 543
 544   assert (n_dependent_vars);
 545
 546   totals = xnmalloc (n_dependent_vars, sizeof *totals);
 547
 548   if ( lex_match (lexer, T_BY))
 549     {
 550       int success ;
 551       success =  examine_parse_independent_vars (lexer, dict, cmd);
 552       if ( success != 1 ) {
 553         free (dependent_vars);
 554         free (totals) ;
 555       }
 556       return success;
 557     }
 558
 559   return 1;
 560 }
 561
 562
 563
 564 /* Parse the clause specifying the factors */
 565 static int
 566 examine_parse_independent_vars (struct lexer *lexer, const struct dictionary *dict, struct cmd_examine *cmd)
 567 {
 568   int success;
 569   struct factor *sf = xmalloc (sizeof *sf);
 570
 571   if ( (lex_token (lexer) != T_ID || dict_lookup_var (dict, lex_tokid (lexer)) == NULL)
 572       && lex_token (lexer) != T_ALL)
 573     {
 574       free ( sf ) ;
 575       return 2;
 576     }
 577
 578
 579   sf->indep_var[0] = parse_variable (lexer, dict);
 580   sf->indep_var[1] = 0;
 581
 582   if ( lex_token (lexer) == T_BY )
 583     {
 584
 585       lex_match (lexer, T_BY);
 586
 587       if ( (lex_token (lexer) != T_ID || dict_lookup_var (dict, lex_tokid (lexer)) == NULL)
 588           && lex_token (lexer) != T_ALL)
 589         {
 590           free ( sf ) ;
 591           return 2;
 592         }
 593
 594       sf->indep_var[1] = parse_variable (lexer, dict);
 595
 596     }
 597
 598
 599   sf->fstats = hsh_create (4,
 600                           (hsh_compare_func *) factor_statistics_compare,
 601                           (hsh_hash_func *) factor_statistics_hash,
 602                           (hsh_free_func *) factor_statistics_free,
 603                           0);
 604
 605   sf->next = factors;
 606   factors = sf;
 607
 608   lex_match (lexer, ',');
 609
 610   if ( lex_token (lexer) == '.' || lex_token (lexer) == '/' )
 611     return 1;
 612
 613   success =  examine_parse_independent_vars (lexer, dict, cmd);
 614
 615   if ( success != 1 )
 616     free ( sf ) ;
 617
 618   return success;
 619 }
 620
 621
 622
 623
 624 void populate_percentiles (struct tab_table *tbl, int col, int row,
 625                           const struct metrics *m);
 626
 627 void populate_descriptives (struct tab_table *t, int col, int row,
 628                            const struct metrics *fs);
 629
 630 void populate_extremes (struct tab_table *t, int col, int row, int n,
 631                        const struct metrics *m);
 632
 633 void populate_summary (struct tab_table *t, int col, int row,
 634                       const struct metrics *m);
 635
 636
 637
 638
 639 /* Perform calculations for the sub factors */
 640 void
 641 factor_calc (const struct ccase *c, int case_no, double weight,
 642              int case_missing)
 643 {
 644   size_t v;
 645   struct factor *fctr = factors;
 646
 647   while ( fctr)
 648     {
 649       struct factor_statistics **foo ;
 650       union value *indep_vals[2] ;
 651
 652       indep_vals[0] = value_dup (
 653                                  case_data (c, fctr->indep_var[0]),
 654                                  var_get_width (fctr->indep_var[0])
 655                                  );
 656
 657       if ( fctr->indep_var[1] )
 658         indep_vals[1] = value_dup (
 659                                    case_data (c, fctr->indep_var[1]),
 660                                    var_get_width (fctr->indep_var[1])
 661                                    );
 662       else
 663         {
 664           const union value sm = {SYSMIS};
 665           indep_vals[1] = value_dup (&sm, 0);
 666         }
 667
 668       assert (fctr->fstats);
 669
 670       foo = ( struct factor_statistics ** )
 671         hsh_probe (fctr->fstats, (void *) &indep_vals);
 672
 673       if ( !*foo )
 674         {
 675
 676           *foo = create_factor_statistics (n_dependent_vars,
 677                                           indep_vals[0],
 678                                           indep_vals[1]);
 679
 680           for ( v =  0 ; v  < n_dependent_vars ; ++v )
 681             {
 682               metrics_precalc ( & (*foo)->m[v] );
 683             }
 684
 685         }
 686       else
 687         {
 688           free (indep_vals[0]);
 689           free (indep_vals[1]);
 690         }
 691
 692       for ( v =  0 ; v  < n_dependent_vars ; ++v )
 693         {
 694           const struct variable *var = dependent_vars[v];
 695           union value *val = value_dup (
 696                                         case_data (c, var),
 697                                         var_get_width (var)
 698                                         );
 699
 700           if (case_missing || var_is_value_missing (var, val, exclude_values))
 701             {
 702               free (val);
 703               continue;
 704             }
 705
 706           metrics_calc ( & (*foo)->m[v], val, weight, case_no);
 707
 708           free (val);
 709         }
 710
 711       fctr = fctr->next;
 712     }
 713 }
 714
 715 static void
 716 run_examine (struct cmd_examine *cmd, struct casereader *input,
 717              struct dataset *ds)
 718 {
 719   struct dictionary *dict = dataset_dict (ds);
 720   casenumber case_no;
 721   struct ccase c;
 722   int v;
 723   bool ok;
 724
 725   struct factor *fctr;
 726
 727   if (!casereader_peek (input, 0, &c))
 728     {
 729       casereader_destroy (input);
 730       return;
 731     }
 732   output_split_file_values (ds, &c);
 733   case_destroy (&c);
 734
 735   input = casereader_create_filter_weight (input, dict, NULL, NULL);
 736   input = casereader_create_counter (input, &case_no, 0);
 737
 738   /* Make sure we haven't got rubbish left over from a
 739      previous split. */
 740   fctr = factors;
 741   while (fctr)
 742     {
 743       struct factor *next = fctr->next;
 744
 745       hsh_clear (fctr->fstats);
 746
 747       fctr->fs = 0;
 748
 749       fctr = next;
 750     }
 751
 752   for ( v = 0 ; v < n_dependent_vars ; ++v )
 753     metrics_precalc (&totals[v]);
 754
 755   for (; casereader_read (input, &c); case_destroy (&c))
 756     {
 757       int case_missing = 0;
 758       const double weight = dict_get_case_weight (dict, &c, NULL);
 759
 760       if ( cmd->miss == XMN_LISTWISE )
 761         {
 762           for ( v = 0 ; v < n_dependent_vars ; ++v )
 763             {
 764               const struct variable *var = dependent_vars[v];
 765               union value *val = value_dup (
 766                                                   case_data (&c, var),
 767                                                   var_get_width (var)
 768                                                   );
 769
 770               if ( var_is_value_missing (var, val, exclude_values))
 771                 case_missing = 1;
 772
 773               free (val);
 774             }
 775         }
 776
 777       for ( v = 0 ; v < n_dependent_vars ; ++v )
 778         {
 779           const struct variable *var = dependent_vars[v];
 780           union value *val = value_dup (
 781                                         case_data (&c, var),
 782                                         var_get_width (var)
 783                                         );
 784
 785           if ( var_is_value_missing (var, val, exclude_values)
 786                || case_missing )
 787             {
 788               free (val) ;
 789               continue ;
 790             }
 791
 792           metrics_calc (&totals[v], val, weight, case_no);
 793
 794           free (val);
 795         }
 796
 797       factor_calc (&c, case_no, weight, case_missing);
 798     }
 799   ok = casereader_destroy (input);
 800
 801   for ( v = 0 ; v < n_dependent_vars ; ++v)
 802     {
 803       fctr = factors;
 804       while ( fctr )
 805         {
 806           struct hsh_iterator hi;
 807           struct factor_statistics *fs;
 808
 809           for ( fs = hsh_first (fctr->fstats, &hi);
 810                 fs != 0 ;
 811                 fs = hsh_next (fctr->fstats, &hi))
 812             {
 813
 814               fs->m[v].ptile_hash = list_to_ptile_hash (&percentile_list);
 815               fs->m[v].ptile_alg = percentile_algorithm;
 816               metrics_postcalc (&fs->m[v]);
 817             }
 818
 819           fctr = fctr->next;
 820         }
 821
 822       totals[v].ptile_hash = list_to_ptile_hash (&percentile_list);
 823       totals[v].ptile_alg = percentile_algorithm;
 824       metrics_postcalc (&totals[v]);
 825     }
 826
 827
 828   /* Make sure that the combination of factors are complete */
 829
 830   fctr = factors;
 831   while ( fctr )
 832     {
 833       struct hsh_iterator hi;
 834       struct hsh_iterator hi0;
 835       struct hsh_iterator hi1;
 836       struct factor_statistics *fs;
 837
 838       struct hsh_table *idh0=0;
 839       struct hsh_table *idh1=0;
 840       union value *val0;
 841       union value *val1;
 842
 843       idh0 = hsh_create (4, (hsh_compare_func *) compare_values,
 844                          (hsh_hash_func *) hash_value,
 845                         0,0);
 846
 847       idh1 = hsh_create (4, (hsh_compare_func *) compare_values,
 848                          (hsh_hash_func *) hash_value,
 849                         0,0);
 850
 851
 852       for ( fs = hsh_first (fctr->fstats, &hi);
 853             fs != 0 ;
 854             fs = hsh_next (fctr->fstats, &hi))
 855         {
 856           hsh_insert (idh0, (void *) &fs->id[0]);
 857           hsh_insert (idh1, (void *) &fs->id[1]);
 858         }
 859
 860       /* Ensure that the factors combination is complete */
 861       for ( val0 = hsh_first (idh0, &hi0);
 862             val0 != 0 ;
 863             val0 = hsh_next (idh0, &hi0))
 864         {
 865           for ( val1 = hsh_first (idh1, &hi1);
 866                 val1 != 0 ;
 867                 val1 = hsh_next (idh1, &hi1))
 868             {
 869               struct factor_statistics **ffs;
 870               union value key[2];
 871               key[0] = *val0;
 872               key[1] = *val1;
 873
 874               ffs = (struct factor_statistics **)
 875                 hsh_probe (fctr->fstats, (void *) &key );
 876
 877               if ( !*ffs ) {
 878                 size_t i;
 879                  (*ffs) = create_factor_statistics (n_dependent_vars,
 880                                                    &key[0], &key[1]);
 881                 for ( i = 0 ; i < n_dependent_vars ; ++i )
 882                   metrics_precalc ( & (*ffs)->m[i]);
 883               }
 884             }
 885         }
 886
 887       hsh_destroy (idh0);
 888       hsh_destroy (idh1);
 889
 890       fctr->fs = (struct factor_statistics **) hsh_sort_copy (fctr->fstats);
 891
 892       fctr = fctr->next;
 893     }
 894
 895   if (ok)
 896     output_examine ();
 897
 898
 899   if ( totals )
 900     {
 901       size_t i;
 902       for ( i = 0 ; i < n_dependent_vars ; ++i )
 903         {
 904           metrics_destroy (&totals[i]);
 905         }
 906     }
 907 }
 908
 909
 910 static void
 911 show_summary (const struct variable **dependent_var, int n_dep_var,
 912              const struct factor *fctr)
 913 {
 914   static const char *subtitle[]=
 915     {
 916       N_ ("Valid"),
 917       N_ ("Missing"),
 918       N_ ("Total")
 919     };
 920
 921   int i;
 922   int heading_columns ;
 923   int n_cols;
 924   const int heading_rows = 3;
 925   struct tab_table *tbl;
 926
 927   int n_rows ;
 928   int n_factors = 1;
 929
 930   if ( fctr )
 931     {
 932       heading_columns = 2;
 933       n_factors = hsh_count (fctr->fstats);
 934       n_rows = n_dep_var * n_factors ;
 935
 936       if ( fctr->indep_var[1] )
 937         heading_columns = 3;
 938     }
 939   else
 940     {
 941       heading_columns = 1;
 942       n_rows = n_dep_var;
 943     }
 944
 945   n_rows += heading_rows;
 946
 947   n_cols = heading_columns + 6;
 948
 949   tbl = tab_create (n_cols,n_rows,0);
 950   tab_headers (tbl, heading_columns, 0, heading_rows, 0);
 951
 952   tab_dim (tbl, tab_natural_dimensions);
 953
 954   /* Outline the box */
 955   tab_box (tbl,
 956            TAL_2, TAL_2,
 957            -1, -1,
 958            0, 0,
 959            n_cols - 1, n_rows - 1);
 960
 961   /* Vertical lines for the data only */
 962   tab_box (tbl,
 963            -1, -1,
 964            -1, TAL_1,
 965            heading_columns, 0,
 966            n_cols - 1, n_rows - 1);
 967
 968
 969   tab_hline (tbl, TAL_2, 0, n_cols - 1, heading_rows );
 970   tab_hline (tbl, TAL_1, heading_columns, n_cols - 1, 1 );
 971   tab_hline (tbl, TAL_1, heading_columns, n_cols - 1, heading_rows -1 );
 972
 973   tab_vline (tbl, TAL_2, heading_columns, 0, n_rows - 1);
 974
 975
 976   tab_title (tbl, _ ("Case Processing Summary"));
 977
 978
 979   tab_joint_text (tbl, heading_columns, 0,
 980                  n_cols -1, 0,
 981                  TAB_CENTER | TAT_TITLE,
 982                  _ ("Cases"));
 983
 984   /* Remove lines ... */
 985   tab_box (tbl,
 986            -1, -1,
 987            TAL_0, TAL_0,
 988            heading_columns, 0,
 989            n_cols - 1, 0);
 990
 991   for ( i = 0 ; i < 3 ; ++i )
 992     {
 993       tab_text (tbl, heading_columns + i*2 , 2, TAB_CENTER | TAT_TITLE,
 994                 _ ("N"));
 995
 996       tab_text (tbl, heading_columns + i*2 + 1, 2, TAB_CENTER | TAT_TITLE,
 997                 _ ("Percent"));
 998
 999       tab_joint_text (tbl, heading_columns + i*2 , 1,
1000                      heading_columns + i*2 + 1, 1,
1001                      TAB_CENTER | TAT_TITLE,
1002                      subtitle[i]);
1003
1004       tab_box (tbl, -1, -1,
1005                TAL_0, TAL_0,
1006                heading_columns + i*2, 1,
1007                heading_columns + i*2 + 1, 1);
1008
1009     }
1010
1011
1012   /* Titles for the independent variables */
1013   if ( fctr )
1014     {
1015       tab_text (tbl, 1, heading_rows - 1, TAB_CENTER | TAT_TITLE,
1016                 var_to_string (fctr->indep_var[0]));
1017
1018       if ( fctr->indep_var[1] )
1019         {
1020           tab_text (tbl, 2, heading_rows - 1, TAB_CENTER | TAT_TITLE,
1021                     var_to_string (fctr->indep_var[1]));
1022         }
1023
1024     }
1025
1026
1027   for ( i = 0 ; i < n_dep_var ; ++i )
1028     {
1029       int n_factors = 1;
1030       if ( fctr )
1031         n_factors = hsh_count (fctr->fstats);
1032
1033
1034       if ( i > 0 )
1035         tab_hline (tbl, TAL_1, 0, n_cols -1 , i * n_factors + heading_rows);
1036
1037       tab_text (tbl,
1038                 0, i * n_factors + heading_rows,
1039                 TAB_LEFT | TAT_TITLE,
1040                 var_to_string (dependent_var[i])
1041                 );
1042
1043
1044       if ( !fctr )
1045         populate_summary (tbl, heading_columns,
1046                          (i * n_factors) + heading_rows,
1047                          &totals[i]);
1048
1049
1050       else
1051         {
1052           struct factor_statistics **fs = fctr->fs;
1053           int count = 0 ;
1054           const union value *prev = NULL;
1055
1056           while (*fs)
1057             {
1058               if ( !prev ||
1059                    0 != compare_values (prev, (*fs)->id[0],
1060                                    var_get_width (fctr->indep_var[0])))
1061                 {
1062                   struct string vstr;
1063                   ds_init_empty (&vstr);
1064                   var_append_value_name (fctr->indep_var[0],
1065                                       (*fs)->id[0], &vstr);
1066
1067                   tab_text (tbl,
1068                             1,
1069                             (i * n_factors ) + count +
1070                             heading_rows,
1071                             TAB_LEFT | TAT_TITLE,
1072                             ds_cstr (&vstr)
1073                             );
1074
1075                   ds_destroy (&vstr);
1076
1077                   if (fctr->indep_var[1] && count > 0 )
1078                     tab_hline (tbl, TAL_1, 1, n_cols - 1,
1079                               (i * n_factors ) + count + heading_rows);
1080
1081                 }
1082
1083               prev = (*fs)->id[0];
1084
1085               if ( fctr->indep_var[1])
1086                 {
1087                   struct string vstr;
1088                   ds_init_empty (&vstr);
1089                   var_append_value_name (fctr->indep_var[1],
1090                                          (*fs)->id[1], &vstr);
1091                 tab_text (tbl,
1092                           2,
1093                           (i * n_factors ) + count +
1094                           heading_rows,
1095                           TAB_LEFT | TAT_TITLE,
1096                             ds_cstr (&vstr)
1097                           );
1098                   ds_destroy (&vstr);
1099             }
1100
1101               populate_summary (tbl, heading_columns,
1102                                (i * n_factors) + count
1103                                + heading_rows,
1104                                & (*fs)->m[i]);
1105
1106               count++ ;
1107               fs++;
1108             }
1109         }
1110     }
1111
1112   tab_submit (tbl);
1113 }
1114
1115
1116 void
1117 populate_summary (struct tab_table *t, int col, int row,
1118                  const struct metrics *m)
1119
1120 {
1121   const double total = m->n + m->n_missing ;
1122
1123   tab_float (t, col + 0, row + 0, TAB_RIGHT, m->n, 8, 0);
1124   tab_float (t, col + 2, row + 0, TAB_RIGHT, m->n_missing, 8, 0);
1125   tab_float (t, col + 4, row + 0, TAB_RIGHT, total, 8, 0);
1126
1127
1128   if ( total > 0 ) {
1129     tab_text (t, col + 1, row + 0, TAB_RIGHT | TAT_PRINTF, "%2.0f%%",
1130               100.0 * m->n / total );
1131
1132     tab_text (t, col + 3, row + 0, TAB_RIGHT | TAT_PRINTF, "%2.0f%%",
1133               100.0 * m->n_missing / total );
1134
1135     /* This seems a bit pointless !!! */
1136     tab_text (t, col + 5, row + 0, TAB_RIGHT | TAT_PRINTF, "%2.0f%%",
1137               100.0 * total / total );
1138
1139
1140   }
1141
1142
1143 }
1144
1145
1146
1147 static void
1148 show_extremes (const struct variable **dependent_var, int n_dep_var,
1149               const struct factor *fctr, int n_extremities)
1150 {
1151   int i;
1152   int heading_columns ;
1153   int n_cols;
1154   const int heading_rows = 1;
1155   struct tab_table *tbl;
1156
1157   int n_factors = 1;
1158   int n_rows ;
1159
1160   if ( fctr )
1161     {
1162       heading_columns = 2;
1163       n_factors = hsh_count (fctr->fstats);
1164
1165       n_rows = n_dep_var * 2 * n_extremities * n_factors;
1166
1167       if ( fctr->indep_var[1] )
1168         heading_columns = 3;
1169     }
1170   else
1171     {
1172       heading_columns = 1;
1173       n_rows = n_dep_var * 2 * n_extremities;
1174     }
1175
1176   n_rows += heading_rows;
1177
1178   heading_columns += 2;
1179   n_cols = heading_columns + 2;
1180
1181   tbl = tab_create (n_cols,n_rows,0);
1182   tab_headers (tbl, heading_columns, 0, heading_rows, 0);
1183
1184   tab_dim (tbl, tab_natural_dimensions);
1185
1186   /* Outline the box, No internal lines*/
1187   tab_box (tbl,
1188            TAL_2, TAL_2,
1189            -1, -1,
1190            0, 0,
1191            n_cols - 1, n_rows - 1);
1192
1193   tab_hline (tbl, TAL_2, 0, n_cols - 1, heading_rows );
1194
1195   tab_title (tbl, _ ("Extreme Values"));
1196
1197   tab_vline (tbl, TAL_2, n_cols - 2, 0, n_rows -1);
1198   tab_vline (tbl, TAL_1, n_cols - 1, 0, n_rows -1);
1199
1200   if ( fctr )
1201     {
1202       tab_text (tbl, 1, heading_rows - 1, TAB_CENTER | TAT_TITLE,
1203                 var_to_string (fctr->indep_var[0]));
1204
1205       if ( fctr->indep_var[1] )
1206         tab_text (tbl, 2, heading_rows - 1, TAB_CENTER | TAT_TITLE,
1207                   var_to_string (fctr->indep_var[1]));
1208     }
1209
1210   tab_text (tbl, n_cols - 1, 0, TAB_CENTER | TAT_TITLE, _ ("Value"));
1211   tab_text (tbl, n_cols - 2, 0, TAB_CENTER | TAT_TITLE, _ ("Case Number"));
1212
1213   for ( i = 0 ; i < n_dep_var ; ++i )
1214     {
1215
1216       if ( i > 0 )
1217         tab_hline (tbl, TAL_1, 0, n_cols -1 ,
1218                   i * 2 * n_extremities * n_factors + heading_rows);
1219
1220       tab_text (tbl, 0,
1221                 i * 2 * n_extremities * n_factors  + heading_rows,
1222                 TAB_LEFT | TAT_TITLE,
1223                 var_to_string (dependent_var[i])
1224                 );
1225
1226
1227       if ( !fctr )
1228         populate_extremes (tbl, heading_columns - 2,
1229                           i * 2 * n_extremities * n_factors  + heading_rows,
1230                           n_extremities, &totals[i]);
1231
1232       else
1233         {
1234           struct factor_statistics **fs = fctr->fs;
1235           int count = 0 ;
1236           const union value *prev  = NULL;
1237
1238           while (*fs)
1239             {
1240               const int row = heading_rows + ( 2 * n_extremities )  *
1241                  ( ( i  * n_factors  ) +  count );
1242
1243
1244               if ( !prev || 0 != compare_values (prev, (*fs)->id[0],
1245                                         var_get_width (fctr->indep_var[0])))
1246                 {
1247                   struct string vstr;
1248                   ds_init_empty (&vstr);
1249                   var_append_value_name (fctr->indep_var[0],
1250                                       (*fs)->id[0], &vstr);
1251
1252                   if ( count > 0 )
1253                     tab_hline (tbl, TAL_1, 1, n_cols - 1, row);
1254
1255                   tab_text (tbl,
1256                             1, row,
1257                             TAB_LEFT | TAT_TITLE,
1258                             ds_cstr (&vstr)
1259                             );
1260
1261                   ds_destroy (&vstr);
1262                 }
1263
1264               prev = (*fs)->id[0];
1265
1266               if (fctr->indep_var[1] && count > 0 )
1267                 tab_hline (tbl, TAL_1, 2, n_cols - 1, row);
1268
1269               if ( fctr->indep_var[1])
1270                 {
1271                   struct string vstr;
1272                   ds_init_empty (&vstr);
1273                   var_append_value_name (fctr->indep_var[1], (*fs)->id[1], &vstr);
1274
1275                 tab_text (tbl, 2, row,
1276                           TAB_LEFT | TAT_TITLE,
1277                             ds_cstr (&vstr)
1278                           );
1279
1280                   ds_destroy (&vstr);
1281                 }
1282
1283               populate_extremes (tbl, heading_columns - 2,
1284                                 row, n_extremities,
1285                                 & (*fs)->m[i]);
1286
1287               count++ ;
1288               fs++;
1289             }
1290         }
1291     }
1292
1293   tab_submit (tbl);
1294 }
1295
1296
1297
1298 /* Fill in the extremities table */
1299 void
1300 populate_extremes (struct tab_table *t,
1301                   int col, int row, int n, const struct metrics *m)
1302 {
1303   int extremity;
1304   int idx=0;
1305
1306
1307   tab_text (t, col, row,
1308            TAB_RIGHT | TAT_TITLE ,
1309            _ ("Highest")
1310            );
1311
1312   tab_text (t, col, row + n ,
1313            TAB_RIGHT | TAT_TITLE ,
1314            _ ("Lowest")
1315            );
1316
1317
1318   tab_hline (t, TAL_1, col, col + 3, row + n );
1319
1320   for (extremity = 0; extremity < n ; ++extremity )
1321     {
1322       /* Highest */
1323       tab_float (t, col + 1, row + extremity,
1324                 TAB_RIGHT,
1325                 extremity + 1, 8, 0);
1326
1327
1328       /* Lowest */
1329       tab_float (t, col + 1, row + extremity + n,
1330                 TAB_RIGHT,
1331                 extremity + 1, 8, 0);
1332
1333     }
1334
1335
1336   /* Lowest */
1337   for (idx = 0, extremity = 0; extremity < n && idx < m->n_data ; ++idx )
1338     {
1339       int j;
1340       const struct weighted_value *wv = m->wvp[idx];
1341       struct case_node *cn = wv->case_nos;
1342
1343
1344       for (j = 0 ; j < wv->w ; ++j  )
1345         {
1346           if ( extremity + j >= n )
1347             break ;
1348
1349           tab_float (t, col + 3, row + extremity + j  + n,
1350                     TAB_RIGHT,
1351                     wv->v.f, 8, 2);
1352
1353           tab_float (t, col + 2, row + extremity + j  + n,
1354                     TAB_RIGHT,
1355                     cn->num, 8, 0);
1356
1357           if ( cn->next )
1358             cn = cn->next;
1359
1360         }
1361
1362       extremity +=  wv->w ;
1363     }
1364
1365
1366   /* Highest */
1367   for (idx = m->n_data - 1, extremity = 0; extremity < n && idx >= 0; --idx )
1368     {
1369       int j;
1370       const struct weighted_value *wv = m->wvp[idx];
1371       struct case_node *cn = wv->case_nos;
1372
1373       for (j = 0 ; j < wv->w ; ++j  )
1374         {
1375           if ( extremity + j >= n )
1376             break ;
1377
1378           tab_float (t, col + 3, row + extremity + j,
1379                     TAB_RIGHT,
1380                     wv->v.f, 8, 2);
1381
1382           tab_float (t, col + 2, row + extremity + j,
1383                     TAB_RIGHT,
1384                     cn->num, 8, 0);
1385
1386           if ( cn->next )
1387             cn = cn->next;
1388
1389         }
1390
1391       extremity +=  wv->w ;
1392     }
1393 }
1394
1395
1396 /* Show the descriptives table */
1397 void
1398 show_descriptives (const struct variable **dependent_var,
1399                   int n_dep_var,
1400                   struct factor *fctr)
1401 {
1402   int i;
1403   int heading_columns ;
1404   int n_cols;
1405   const int n_stat_rows = 13;
1406
1407   const int heading_rows = 1;
1408
1409   struct tab_table *tbl;
1410
1411   int n_factors = 1;
1412   int n_rows ;
1413
1414   if ( fctr )
1415     {
1416       heading_columns = 4;
1417       n_factors = hsh_count (fctr->fstats);
1418
1419       n_rows = n_dep_var * n_stat_rows * n_factors;
1420
1421       if ( fctr->indep_var[1] )
1422         heading_columns = 5;
1423     }
1424   else
1425     {
1426       heading_columns = 3;
1427       n_rows = n_dep_var * n_stat_rows;
1428     }
1429
1430   n_rows += heading_rows;
1431
1432   n_cols = heading_columns + 2;
1433
1434
1435   tbl = tab_create (n_cols, n_rows, 0);
1436
1437   tab_headers (tbl, heading_columns + 1, 0, heading_rows, 0);
1438
1439   tab_dim (tbl, tab_natural_dimensions);
1440
1441   /* Outline the box and have no internal lines*/
1442   tab_box (tbl,
1443            TAL_2, TAL_2,
1444            -1, -1,
1445            0, 0,
1446            n_cols - 1, n_rows - 1);
1447
1448   tab_hline (tbl, TAL_2, 0, n_cols - 1, heading_rows );
1449
1450   tab_vline (tbl, TAL_1, heading_columns, 0, n_rows - 1);
1451   tab_vline (tbl, TAL_2, n_cols - 2, 0, n_rows - 1);
1452   tab_vline (tbl, TAL_1, n_cols - 1, 0, n_rows - 1);
1453
1454   tab_text (tbl, n_cols - 2, 0, TAB_CENTER | TAT_TITLE, _ ("Statistic"));
1455   tab_text (tbl, n_cols - 1, 0, TAB_CENTER | TAT_TITLE, _ ("Std. Error"));
1456
1457   tab_title (tbl, _ ("Descriptives"));
1458
1459
1460   for ( i = 0 ; i < n_dep_var ; ++i )
1461     {
1462       const int row = heading_rows + i * n_stat_rows * n_factors ;
1463
1464       if ( i > 0 )
1465         tab_hline (tbl, TAL_1, 0, n_cols - 1, row );
1466
1467       tab_text (tbl, 0,
1468                 i * n_stat_rows * n_factors  + heading_rows,
1469                 TAB_LEFT | TAT_TITLE,
1470                 var_to_string (dependent_var[i])
1471                 );
1472
1473
1474       if ( fctr  )
1475         {
1476           const union value *prev = NULL;
1477
1478           struct factor_statistics **fs = fctr->fs;
1479           int count = 0;
1480
1481           tab_text (tbl, 1, heading_rows - 1, TAB_CENTER | TAT_TITLE,
1482                     var_to_string (fctr->indep_var[0]));
1483
1484
1485           if ( fctr->indep_var[1])
1486             tab_text (tbl, 2, heading_rows - 1, TAB_CENTER | TAT_TITLE,
1487                       var_to_string (fctr->indep_var[1]));
1488
1489           while ( *fs )
1490             {
1491               const int row = heading_rows + n_stat_rows  *
1492                  ( ( i  * n_factors  ) +  count );
1493
1494
1495               if ( !prev || 0 != compare_values (prev, (*fs)->id[0],
1496                                         var_get_width (fctr->indep_var[0])))
1497                 {
1498                   struct string vstr;
1499                   ds_init_empty (&vstr);
1500                   var_append_value_name (fctr->indep_var[0],
1501                                       (*fs)->id[0], &vstr);
1502
1503                   if ( count > 0 )
1504                     tab_hline (tbl, TAL_1, 1, n_cols - 1, row);
1505
1506                   tab_text (tbl,
1507                             1, row,
1508                             TAB_LEFT | TAT_TITLE,
1509                             ds_cstr (&vstr)
1510                             );
1511
1512                   ds_destroy (&vstr);
1513                 }
1514
1515               prev = (*fs)->id[0];
1516
1517               if (fctr->indep_var[1] && count > 0 )
1518                 tab_hline (tbl, TAL_1, 2, n_cols - 1, row);
1519
1520               if ( fctr->indep_var[1])
1521                 {
1522                   struct string vstr;
1523                   ds_init_empty (&vstr);
1524                   var_append_value_name (fctr->indep_var[1], (*fs)->id[1], &vstr);
1525
1526                 tab_text (tbl, 2, row,
1527                           TAB_LEFT | TAT_TITLE,
1528                             ds_cstr (&vstr)
1529                           );
1530
1531                   ds_destroy (&vstr);
1532                 }
1533
1534               populate_descriptives (tbl, heading_columns - 2,
1535                                     row, & (*fs)->m[i]);
1536
1537               count++ ;
1538               fs++;
1539             }
1540
1541         }
1542
1543       else
1544         {
1545
1546           populate_descriptives (tbl, heading_columns - 2,
1547                                 i * n_stat_rows * n_factors  + heading_rows,
1548                                 &totals[i]);
1549         }
1550     }
1551
1552   tab_submit (tbl);
1553
1554 }
1555
1556
1557 /* Fill in the descriptives data */
1558 void
1559 populate_descriptives (struct tab_table *tbl, int col, int row,
1560                       const struct metrics *m)
1561 {
1562   const double t = gsl_cdf_tdist_Qinv ((1 - cmd.n_cinterval[0] / 100.0)/2.0,
1563                                       m->n -1);
1564
1565   tab_text (tbl, col,
1566             row,
1567             TAB_LEFT | TAT_TITLE,
1568             _ ("Mean"));
1569
1570   tab_float (tbl, col + 2,
1571              row,
1572              TAB_CENTER,
1573              m->mean,
1574              8,2);
1575
1576   tab_float (tbl, col + 3,
1577              row,
1578              TAB_CENTER,
1579              m->se_mean,
1580              8,3);
1581
1582
1583   tab_text (tbl, col,
1584             row + 1,
1585             TAB_LEFT | TAT_TITLE | TAT_PRINTF,
1586             _ ("%g%% Confidence Interval for Mean"), cmd.n_cinterval[0]);
1587
1588
1589   tab_text (tbl, col + 1,
1590             row  + 1,
1591             TAB_LEFT | TAT_TITLE,
1592             _ ("Lower Bound"));
1593
1594   tab_float (tbl, col + 2,
1595              row + 1,
1596              TAB_CENTER,
1597              m->mean - t * m->se_mean,
1598              8,3);
1599
1600   tab_text (tbl, col + 1,
1601             row + 2,
1602             TAB_LEFT | TAT_TITLE,
1603             _ ("Upper Bound"));
1604
1605
1606   tab_float (tbl, col + 2,
1607              row + 2,
1608              TAB_CENTER,
1609              m->mean + t * m->se_mean,
1610              8,3);
1611
1612   tab_text (tbl, col,
1613             row + 3,
1614             TAB_LEFT | TAT_TITLE | TAT_PRINTF,
1615             _ ("5%% Trimmed Mean"));
1616
1617   tab_float (tbl, col + 2,
1618              row + 3,
1619              TAB_CENTER,
1620              m->trimmed_mean,
1621              8,2);
1622
1623   tab_text (tbl, col,
1624             row + 4,
1625             TAB_LEFT | TAT_TITLE,
1626             _ ("Median"));
1627
1628   {
1629     struct percentile *p;
1630     double d = 50;
1631
1632     p = hsh_find (m->ptile_hash, &d);
1633
1634     assert (p);
1635
1636
1637     tab_float (tbl, col + 2,
1638                row + 4,
1639                TAB_CENTER,
1640                p->v,
1641                8, 2);
1642   }
1643
1644
1645   tab_text (tbl, col,
1646             row + 5,
1647             TAB_LEFT | TAT_TITLE,
1648             _ ("Variance"));
1649
1650   tab_float (tbl, col + 2,
1651              row + 5,
1652              TAB_CENTER,
1653              m->var,
1654              8,3);
1655
1656
1657   tab_text (tbl, col,
1658             row + 6,
1659             TAB_LEFT | TAT_TITLE,
1660             _ ("Std. Deviation"));
1661
1662
1663   tab_float (tbl, col + 2,
1664              row + 6,
1665              TAB_CENTER,
1666              m->stddev,
1667              8,3);
1668
1669
1670   tab_text (tbl, col,
1671             row + 7,
1672             TAB_LEFT | TAT_TITLE,
1673             _ ("Minimum"));
1674
1675   tab_float (tbl, col + 2,
1676              row + 7,
1677              TAB_CENTER,
1678              m->min,
1679              8,3);
1680
1681   tab_text (tbl, col,
1682             row + 8,
1683             TAB_LEFT | TAT_TITLE,
1684             _ ("Maximum"));
1685
1686   tab_float (tbl, col + 2,
1687              row + 8,
1688              TAB_CENTER,
1689              m->max,
1690              8,3);
1691
1692
1693   tab_text (tbl, col,
1694             row + 9,
1695             TAB_LEFT | TAT_TITLE,
1696             _ ("Range"));
1697
1698
1699   tab_float (tbl, col + 2,
1700              row + 9,
1701              TAB_CENTER,
1702              m->max - m->min,
1703              8,3);
1704
1705   tab_text (tbl, col,
1706             row + 10,
1707             TAB_LEFT | TAT_TITLE,
1708             _ ("Interquartile Range"));
1709
1710   {
1711     struct percentile *p1;
1712     struct percentile *p2;
1713
1714     double d = 75;
1715     p1 = hsh_find (m->ptile_hash, &d);
1716
1717     d = 25;
1718     p2 = hsh_find (m->ptile_hash, &d);
1719
1720     assert (p1);
1721     assert (p2);
1722
1723     tab_float (tbl, col + 2,
1724                row + 10,
1725                TAB_CENTER,
1726                p1->v - p2->v,
1727                8, 2);
1728   }
1729
1730
1731
1732   tab_text (tbl, col,
1733             row + 11,
1734             TAB_LEFT | TAT_TITLE,
1735             _ ("Skewness"));
1736
1737
1738   tab_float (tbl, col + 2,
1739              row + 11,
1740              TAB_CENTER,
1741              m->skewness,
1742              8,3);
1743
1744   /* stderr of skewness */
1745   tab_float (tbl, col + 3,
1746              row + 11,
1747              TAB_CENTER,
1748              calc_seskew (m->n),
1749              8,3);
1750
1751
1752   tab_text (tbl, col,
1753             row + 12,
1754             TAB_LEFT | TAT_TITLE,
1755             _ ("Kurtosis"));
1756
1757
1758   tab_float (tbl, col + 2,
1759              row + 12,
1760              TAB_CENTER,
1761              m->kurtosis,
1762              8,3);
1763
1764   /* stderr of kurtosis */
1765   tab_float (tbl, col + 3,
1766              row + 12,
1767              TAB_CENTER,
1768              calc_sekurt (m->n),
1769              8,3);
1770
1771
1772 }
1773
1774
1775
1776 void
1777 box_plot_variables (const struct factor *fctr,
1778                    const struct variable **vars, int n_vars,
1779                    const struct variable *id)
1780 {
1781
1782   int i;
1783   struct factor_statistics **fs ;
1784
1785   if ( ! fctr )
1786     {
1787       box_plot_group (fctr, vars, n_vars, id);
1788       return;
1789     }
1790
1791   for ( fs = fctr->fs ; *fs ; ++fs )
1792     {
1793       struct string str;
1794       double y_min = DBL_MAX;
1795       double y_max = -DBL_MAX;
1796       struct chart *ch = chart_create ();
1797       ds_init_empty (&str);
1798       factor_to_string (fctr, *fs, 0, &str );
1799
1800       chart_write_title (ch, ds_cstr (&str));
1801
1802       for ( i = 0 ; i < n_vars ; ++i )
1803         {
1804           y_max = MAX (y_max, (*fs)->m[i].max);
1805           y_min = MIN (y_min, (*fs)->m[i].min);
1806         }
1807
1808       boxplot_draw_yscale (ch, y_max, y_min);
1809
1810       for ( i = 0 ; i < n_vars ; ++i )
1811         {
1812
1813           const double box_width = (ch->data_right - ch->data_left)
1814             / (n_vars * 2.0 ) ;
1815
1816           const double box_centre = ( i * 2 + 1) * box_width
1817             + ch->data_left;
1818
1819           boxplot_draw_boxplot (ch,
1820                                box_centre, box_width,
1821                                & (*fs)->m[i],
1822                                var_to_string (vars[i]));
1823
1824
1825         }
1826
1827       chart_submit (ch);
1828       ds_destroy (&str);
1829     }
1830 }
1831
1832
1833
1834 /* Do a box plot, grouping all factors into one plot ;
1835    each dependent variable has its own plot.
1836 */
1837 void
1838 box_plot_group (const struct factor *fctr,
1839                const struct variable **vars,
1840                int n_vars,
1841                const struct variable *id UNUSED)
1842 {
1843
1844   int i;
1845
1846   for ( i = 0 ; i < n_vars ; ++i )
1847     {
1848       struct factor_statistics **fs ;
1849       struct chart *ch;
1850
1851       ch = chart_create ();
1852
1853       boxplot_draw_yscale (ch, totals[i].max, totals[i].min);
1854
1855       if ( fctr )
1856         {
1857           int n_factors = 0;
1858           int f=0;
1859           for ( fs = fctr->fs ; *fs ; ++fs )
1860             ++n_factors;
1861
1862           chart_write_title (ch, _ ("Boxplot of %s vs. %s"),
1863                             var_to_string (vars[i]), var_to_string (fctr->indep_var[0]) );
1864
1865           for ( fs = fctr->fs ; *fs ; ++fs )
1866             {
1867               struct string str;
1868               const double box_width = (ch->data_right - ch->data_left)
1869                 / (n_factors * 2.0 ) ;
1870
1871               const double box_centre = ( f++ * 2 + 1) * box_width
1872                 + ch->data_left;
1873
1874               ds_init_empty (&str);
1875               factor_to_string_concise (fctr, *fs, &str);
1876
1877               boxplot_draw_boxplot (ch,
1878                                    box_centre, box_width,
1879                                    & (*fs)->m[i],
1880                                    ds_cstr (&str));
1881               ds_destroy (&str);
1882             }
1883         }
1884       else if ( ch )
1885         {
1886           const double box_width = (ch->data_right - ch->data_left) / 3.0;
1887           const double box_centre = (ch->data_right + ch->data_left) / 2.0;
1888
1889           chart_write_title (ch, _ ("Boxplot"));
1890
1891           boxplot_draw_boxplot (ch,
1892                                box_centre,    box_width,
1893                                &totals[i],
1894                                var_to_string (vars[i]) );
1895
1896         }
1897
1898       chart_submit (ch);
1899     }
1900 }
1901
1902
1903 /* Plot the normal and detrended normal plots for m
1904    Label the plots with factorname */
1905 void
1906 np_plot (const struct metrics *m, const char *factorname)
1907 {
1908   int i;
1909   double yfirst=0, ylast=0;
1910
1911   /* Normal Plot */
1912   struct chart *np_chart;
1913
1914   /* Detrended Normal Plot */
1915   struct chart *dnp_chart;
1916
1917   /* The slope and intercept of the ideal normal probability line */
1918   const double slope = 1.0 / m->stddev;
1919   const double intercept = - m->mean / m->stddev;
1920
1921   /* Cowardly refuse to plot an empty data set */
1922   if ( m->n_data == 0 )
1923     return ;
1924
1925   np_chart = chart_create ();
1926   dnp_chart = chart_create ();
1927
1928   if ( !np_chart || ! dnp_chart )
1929     return ;
1930
1931   chart_write_title (np_chart, _ ("Normal Q-Q Plot of %s"), factorname);
1932   chart_write_xlabel (np_chart, _ ("Observed Value"));
1933   chart_write_ylabel (np_chart, _ ("Expected Normal"));
1934
1935
1936   chart_write_title (dnp_chart, _ ("Detrended Normal Q-Q Plot of %s"),
1937                     factorname);
1938   chart_write_xlabel (dnp_chart, _ ("Observed Value"));
1939   chart_write_ylabel (dnp_chart, _ ("Dev from Normal"));
1940
1941   yfirst = gsl_cdf_ugaussian_Pinv (m->wvp[0]->rank / ( m->n + 1));
1942   ylast =  gsl_cdf_ugaussian_Pinv (m->wvp[m->n_data-1]->rank / ( m->n + 1));
1943
1944
1945   {
1946     /* Need to make sure that both the scatter plot and the ideal fit into the
1947        plot */
1948     double x_lower = MIN (m->min, (yfirst - intercept) / slope) ;
1949     double x_upper = MAX (m->max, (ylast  - intercept) / slope) ;
1950     double slack = (x_upper - x_lower)  * 0.05 ;
1951
1952     chart_write_xscale (np_chart, x_lower - slack, x_upper + slack, 5);
1953
1954     chart_write_xscale (dnp_chart, m->min, m->max, 5);
1955
1956   }
1957
1958   chart_write_yscale (np_chart, yfirst, ylast, 5);
1959
1960   {
1961     /* We have to cache the detrended data, beacause we need to
1962        find its limits before we can plot it */
1963     double *d_data = xnmalloc (m->n_data, sizeof *d_data);
1964     double d_max = -DBL_MAX;
1965     double d_min = DBL_MAX;
1966     for ( i = 0 ; i < m->n_data; ++i )
1967       {
1968         const double ns = gsl_cdf_ugaussian_Pinv (m->wvp[i]->rank / ( m->n + 1));
1969
1970         chart_datum (np_chart, 0, m->wvp[i]->v.f, ns);
1971
1972         d_data[i] = (m->wvp[i]->v.f - m->mean) / m->stddev  - ns;
1973
1974         if ( d_data[i] < d_min ) d_min = d_data[i];
1975         if ( d_data[i] > d_max ) d_max = d_data[i];
1976       }
1977     chart_write_yscale (dnp_chart, d_min, d_max, 5);
1978
1979     for ( i = 0 ; i < m->n_data; ++i )
1980       chart_datum (dnp_chart, 0, m->wvp[i]->v.f, d_data[i]);
1981
1982     free (d_data);
1983   }
1984
1985   chart_line (np_chart, slope, intercept, yfirst, ylast , CHART_DIM_Y);
1986   chart_line (dnp_chart, 0, 0, m->min, m->max , CHART_DIM_X);
1987
1988   chart_submit (np_chart);
1989   chart_submit (dnp_chart);
1990 }
1991
1992
1993
1994
1995 /* Show the percentiles */
1996 void
1997 show_percentiles (const struct variable **dependent_var,
1998                  int n_dep_var,
1999                  struct factor *fctr)
2000 {
2001   struct tab_table *tbl;
2002   int i;
2003
2004   int n_cols, n_rows;
2005   int n_factors;
2006
2007   struct hsh_table *ptiles ;
2008
2009   int n_heading_columns;
2010   const int n_heading_rows = 2;
2011   const int n_stat_rows = 2;
2012
2013   int n_ptiles ;
2014
2015   if ( fctr )
2016     {
2017       struct factor_statistics **fs = fctr->fs ;
2018       n_heading_columns = 3;
2019       n_factors = hsh_count (fctr->fstats);
2020
2021       ptiles = (*fs)->m[0].ptile_hash;
2022
2023       if ( fctr->indep_var[1] )
2024         n_heading_columns = 4;
2025     }
2026   else
2027     {
2028       n_factors = 1;
2029       n_heading_columns = 2;
2030
2031       ptiles = totals[0].ptile_hash;
2032     }
2033
2034   n_ptiles = hsh_count (ptiles);
2035
2036   n_rows = n_heading_rows + n_dep_var * n_stat_rows * n_factors;
2037
2038   n_cols = n_heading_columns + n_ptiles ;
2039
2040   tbl = tab_create (n_cols, n_rows, 0);
2041
2042   tab_headers (tbl, n_heading_columns + 1, 0, n_heading_rows, 0);
2043
2044   tab_dim (tbl, tab_natural_dimensions);
2045
2046   /* Outline the box and have no internal lines*/
2047   tab_box (tbl,
2048            TAL_2, TAL_2,
2049            -1, -1,
2050            0, 0,
2051            n_cols - 1, n_rows - 1);
2052
2053   tab_hline (tbl, TAL_2, 0, n_cols - 1, n_heading_rows );
2054
2055   tab_vline (tbl, TAL_2, n_heading_columns, 0, n_rows - 1);
2056
2057
2058   tab_title (tbl, _ ("Percentiles"));
2059
2060
2061   tab_hline (tbl, TAL_1, n_heading_columns, n_cols - 1, 1 );
2062
2063
2064   tab_box (tbl,
2065            -1, -1,
2066            -1, TAL_1,
2067            0, n_heading_rows,
2068            n_heading_columns - 1, n_rows - 1);
2069
2070
2071   tab_box (tbl,
2072            -1, -1,
2073            -1, TAL_1,
2074            n_heading_columns, n_heading_rows - 1,
2075            n_cols - 1, n_rows - 1);
2076
2077   tab_joint_text (tbl, n_heading_columns + 1, 0,
2078                  n_cols - 1 , 0,
2079                  TAB_CENTER | TAT_TITLE ,
2080                  _ ("Percentiles"));
2081
2082
2083   {
2084     /* Put in the percentile break points as headings */
2085
2086     struct percentile **p = (struct percentile **) hsh_sort (ptiles);
2087
2088     i = 0;
2089     while ( (*p)  )
2090       {
2091         tab_float (tbl, n_heading_columns + i++ , 1,
2092                   TAB_CENTER,
2093                   (*p)->p, 8, 0);
2094
2095         p++;
2096       }
2097
2098   }
2099
2100   for ( i = 0 ; i < n_dep_var ; ++i )
2101     {
2102       const int n_stat_rows = 2;
2103       const int row = n_heading_rows + i * n_stat_rows * n_factors ;
2104
2105       if ( i > 0 )
2106         tab_hline (tbl, TAL_1, 0, n_cols - 1, row );
2107
2108       tab_text (tbl, 0,
2109                 i * n_stat_rows * n_factors  + n_heading_rows,
2110                 TAB_LEFT | TAT_TITLE,
2111                 var_to_string (dependent_var[i])
2112                 );
2113
2114       if ( fctr  )
2115         {
2116           const union value *prev  = NULL ;
2117           struct factor_statistics **fs = fctr->fs;
2118           int count = 0;
2119
2120           tab_text (tbl, 1, n_heading_rows - 1,
2121                     TAB_CENTER | TAT_TITLE,
2122                     var_to_string (fctr->indep_var[0]));
2123
2124
2125           if ( fctr->indep_var[1])
2126             tab_text (tbl, 2, n_heading_rows - 1, TAB_CENTER | TAT_TITLE,
2127                       var_to_string (fctr->indep_var[1]));
2128
2129           while ( *fs )
2130             {
2131               const int row = n_heading_rows + n_stat_rows  *
2132                  ( ( i  * n_factors  ) +  count );
2133
2134
2135               if ( !prev || 0 != compare_values (prev, (*fs)->id[0],
2136                                         var_get_width (fctr->indep_var[0])))
2137                 {
2138                   struct string vstr;
2139                   ds_init_empty (&vstr);
2140                   var_append_value_name (fctr->indep_var[0],
2141                                       (*fs)->id[0], &vstr);
2142
2143
2144                   if ( count > 0 )
2145                     tab_hline (tbl, TAL_1, 1, n_cols - 1, row);
2146
2147                   tab_text (tbl,
2148                             1, row,
2149                             TAB_LEFT | TAT_TITLE,
2150                             ds_cstr (&vstr)
2151                             );
2152
2153                   ds_destroy (&vstr);
2154                 }
2155
2156               prev = (*fs)->id[0];
2157
2158               if (fctr->indep_var[1] && count > 0 )
2159                 tab_hline (tbl, TAL_1, 2, n_cols - 1, row);
2160
2161               if ( fctr->indep_var[1])
2162                 {
2163                   struct string vstr;
2164                   ds_init_empty (&vstr);
2165                   var_append_value_name (fctr->indep_var[1], (*fs)->id[1], &vstr);
2166
2167                 tab_text (tbl, 2, row,
2168                           TAB_LEFT | TAT_TITLE,
2169                             ds_cstr (&vstr)
2170                           );
2171
2172                   ds_destroy (&vstr);
2173                 }
2174
2175
2176               populate_percentiles (tbl, n_heading_columns - 1,
2177                                    row, & (*fs)->m[i]);
2178
2179
2180               count++ ;
2181               fs++;
2182             }
2183
2184
2185         }
2186       else
2187         {
2188           populate_percentiles (tbl, n_heading_columns - 1,
2189                                i * n_stat_rows * n_factors  + n_heading_rows,
2190                                &totals[i]);
2191         }
2192
2193
2194     }
2195
2196
2197   tab_submit (tbl);
2198
2199
2200 }
2201
2202
2203
2204
2205 void
2206 populate_percentiles (struct tab_table *tbl, int col, int row,
2207                      const struct metrics *m)
2208 {
2209   int i;
2210
2211   struct percentile **p = (struct percentile **) hsh_sort (m->ptile_hash);
2212
2213   tab_text (tbl,
2214             col, row + 1,
2215             TAB_LEFT | TAT_TITLE,
2216             _ ("Tukey\'s Hinges")
2217             );
2218
2219   tab_text (tbl,
2220             col, row,
2221             TAB_LEFT | TAT_TITLE,
2222             ptile_alg_desc[m->ptile_alg]
2223             );
2224
2225
2226   i = 0;
2227   while ( (*p)  )
2228     {
2229       tab_float (tbl, col + i + 1 , row,
2230                 TAB_CENTER,
2231                  (*p)->v, 8, 2);
2232       if ( (*p)->p == 25 )
2233         tab_float (tbl, col + i + 1 , row + 1,
2234                   TAB_CENTER,
2235                   m->hinge[0], 8, 2);
2236
2237       if ( (*p)->p == 50 )
2238         tab_float (tbl, col + i + 1 , row + 1,
2239                   TAB_CENTER,
2240                   m->hinge[1], 8, 2);
2241
2242       if ( (*p)->p == 75 )
2243         tab_float (tbl, col + i + 1 , row + 1,
2244                   TAB_CENTER,
2245                   m->hinge[2], 8, 2);
2246
2247
2248       i++;
2249
2250       p++;
2251     }
2252
2253 }
2254
2255 static void
2256 factor_to_string (const struct factor *fctr,
2257                   const struct factor_statistics *fs,
2258                   const struct variable *var,
2259                   struct string *str
2260                   )
2261 {
2262   if (var)
2263     ds_put_format (str, "%s (",var_to_string (var) );
2264
2265
2266   ds_put_format (str,  "%s = ",
2267                  var_to_string (fctr->indep_var[0]));
2268
2269   var_append_value_name (fctr->indep_var[0], fs->id[0], str);
2270
2271   if ( fctr->indep_var[1] )
2272     {
2273       ds_put_format (str, "; %s = )",
2274                      var_to_string (fctr->indep_var[1]));
2275
2276       var_append_value_name (fctr->indep_var[1], fs->id[1], str);
2277     }
2278   else
2279     {
2280       if ( var )
2281         ds_put_cstr (str, ")");
2282     }
2283 }
2284
2285
2286 static void
2287 factor_to_string_concise (const struct factor *fctr,
2288                           const struct factor_statistics *fs,
2289                           struct string *str
2290                           )
2291
2292 {
2293   var_append_value_name (fctr->indep_var[0], fs->id[0], str);
2294
2295   if ( fctr->indep_var[1] )
2296     {
2297       ds_put_cstr (str, ",");
2298
2299       var_append_value_name (fctr->indep_var[1],fs->id[1], str);
2300
2301       ds_put_cstr (str, ")");
2302     }
2303 }
2304
2305 /*
2306   Local Variables:
2307   mode: c
2308   End:
2309 */