src/examine.q

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