pintos-os.org Git - pspp/blob - src/language/stats/examine.c

   1 /*
   2   PSPP - a program for statistical analysis.
   3   Copyright (C) 2012 Free Software Foundation, Inc.
   4
   5   This program is free software: you can redistribute it and/or modify
   6   it under the terms of the GNU General Public License as published by
   7   the Free Software Foundation, either version 3 of the License, or
   8   (at your option) any later version.
   9
  10   This program is distributed in the hope that it will be useful,
  11   but WITHOUT ANY WARRANTY; without even the implied warranty of
  12   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  13   GNU General Public License for more details.
  14
  15   You should have received a copy of the GNU General Public License
  16   along with this program.  If not, see <http://www.gnu.org/licenses/>.
  17 */
  18
  19 #include <config.h>
  20
  21 #include <math.h>
  22 #include <gsl/gsl_cdf.h>
  23
  24 #include "libpspp/assertion.h"
  25 #include "libpspp/message.h"
  26 #include "libpspp/pool.h"
  27
  28
  29 #include "data/dataset.h"
  30 #include "data/dictionary.h"
  31 #include "data/casegrouper.h"
  32 #include "data/casereader.h"
  33 #include "data/casewriter.h"
  34 #include "data/caseproto.h"
  35 #include "data/subcase.h"
  36
  37
  38 #include "data/format.h"
  39
  40 #include "math/interaction.h"
  41 #include "math/box-whisker.h"
  42 #include "math/categoricals.h"
  43 #include "math/chart-geometry.h"
  44 #include "math/histogram.h"
  45 #include "math/moments.h"
  46 #include "math/np.h"
  47 #include "math/sort.h"
  48 #include "math/order-stats.h"
  49 #include "math/percentiles.h"
  50 #include "math/tukey-hinges.h"
  51 #include "math/trimmed-mean.h"
  52
  53 #include "output/charts/boxplot.h"
  54 #include "output/charts/np-plot.h"
  55 #include "output/charts/plot-hist.h"
  56
  57 #include "language/command.h"
  58 #include "language/lexer/lexer.h"
  59 #include "language/lexer/value-parser.h"
  60 #include "language/lexer/variable-parser.h"
  61
  62 #include "output/tab.h"
  63
  64 #include "gettext.h"
  65 #define _(msgid) gettext (msgid)
  66 #define N_(msgid) msgid
  67
  68 enum bp_mode
  69   {
  70     BP_GROUPS,
  71     BP_VARIABLES
  72   };
  73
  74
  75 /* Indices for the ex_proto member (below) */
  76 enum
  77   {
  78     EX_VAL,  /* value */
  79     EX_ID,   /* identity */
  80     EX_WT    /* weight */
  81   };
  82
  83
  84 struct examine
  85 {
  86   struct pool *pool;
  87
  88   /* A caseproto used to contain the data subsets under examination,
  89      see (enum above)   */
  90   struct caseproto *ex_proto;
  91
  92   size_t n_dep_vars;
  93   const struct variable **dep_vars;
  94
  95   size_t n_iacts;
  96   struct interaction **iacts;
  97
  98   enum mv_class exclude;
  99
 100   const struct dictionary *dict;
 101
 102   struct categoricals *cats;
 103
 104   /* how many extremities to display */
 105   int disp_extremes;
 106   int calc_extremes;
 107   bool descriptives;
 108
 109   double conf;
 110
 111   bool missing_pw;
 112
 113   /* The case index of the ID value (or -1) if not applicable */
 114   size_t id_idx;
 115   int id_width;
 116
 117   enum pc_alg pc_alg;
 118   double *ptiles;
 119   size_t n_percentiles;
 120
 121   bool npplot;
 122   bool histogram;
 123   bool boxplot;
 124
 125   enum bp_mode boxplot_mode;
 126
 127   const struct variable *id_var;
 128
 129   const struct variable *wv;
 130 };
 131
 132 struct extremity
 133 {
 134   /* The value of this extremity */
 135   double val;
 136
 137   /* Either the casenumber or the value of the variable specified
 138      by the /ID subcommand which corresponds to this extremity */
 139   union value identity;
 140 };
 141
 142 struct exploratory_stats
 143 {
 144   double missing;
 145   double non_missing;
 146
 147   struct moments *mom;
 148
 149   /* Most operations need a sorted reader/writer */
 150   struct casewriter *sorted_writer;
 151   struct casereader *sorted_reader;
 152
 153   struct extremity *minima;
 154   struct extremity *maxima;
 155
 156   /*
 157      Minimum should alway equal mimima[0].val.
 158      Likewise, maximum should alway equal maxima[0].val.
 159      This redundancy exists as an optimisation effort.
 160      Some statistics (eg histogram) require early calculation
 161      of the min and max
 162   */
 163   double minimum;
 164   double maximum;
 165
 166   struct trimmed_mean *trimmed_mean;
 167   struct percentile *quartiles[3];
 168   struct percentile **percentiles;
 169
 170   struct tukey_hinges *hinges;
 171
 172   /* The data for the NP Plots */
 173   struct np *np;
 174
 175   struct histogram *histogram;
 176
 177   /* The data for the box plots */
 178   struct box_whisker *box_whisker;
 179
 180   /* Total weight */
 181   double cc;
 182
 183   /* The minimum weight */
 184   double cmin;
 185 };
 186
 187
 188 /* Returns an array of (iact->n_vars) pointers to union value initialised to NULL.
 189    The caller must free this array when no longer required. */
 190 static const union value **
 191 previous_value_alloc (const struct interaction *iact)
 192 {
 193   int ivar_idx;
 194
 195   const union value **prev_val = xcalloc (iact->n_vars, sizeof (*prev_val));
 196
 197   for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
 198     prev_val[ivar_idx] = NULL;
 199
 200   return prev_val;
 201 }
 202
 203 /* Set the contents of PREV_VAL to the values of C indexed by the variables of IACT */
 204 static int
 205 previous_value_record (const struct interaction *iact, const struct ccase *c, const union value **prev_val)
 206 {
 207   int ivar_idx;
 208   int diff_idx = -1;
 209
 210   for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
 211     {
 212       const struct variable *ivar = iact->vars[ivar_idx];
 213       const int width = var_get_width (ivar);
 214       const union value *val = case_data (c, ivar);
 215
 216       if (prev_val[ivar_idx])
 217         if (! value_equal (prev_val[ivar_idx], val, width))
 218           {
 219             diff_idx = ivar_idx;
 220             break;
 221           }
 222     }
 223
 224   for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
 225     {
 226       const struct variable *ivar = iact->vars[ivar_idx];
 227       const union value *val = case_data (c, ivar);
 228
 229       prev_val[ivar_idx] = val;
 230     }
 231   return diff_idx;
 232 }
 233
 234
 235 static void
 236 show_boxplot_grouped (const struct examine *cmd, int iact_idx)
 237 {
 238   int v;
 239
 240   const struct interaction *iact = cmd->iacts[iact_idx];
 241   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
 242
 243   for (v = 0; v < cmd->n_dep_vars; ++v)
 244     {
 245       double y_min = DBL_MAX;
 246       double y_max = -DBL_MAX;
 247       int grp;
 248       struct boxplot *boxplot;
 249       struct string title;
 250       ds_init_empty (&title);
 251
 252       if (iact->n_vars > 0)
 253         {
 254           struct string istr;
 255           ds_init_empty (&istr);
 256           interaction_to_string (iact, &istr);
 257           ds_put_format (&title, _("Boxplot of %s vs. %s"),
 258                          var_to_string (cmd->dep_vars[v]),
 259                          ds_cstr (&istr));
 260           ds_destroy (&istr);
 261         }
 262       else
 263         ds_put_format (&title, _("Boxplot of %s"), var_to_string (cmd->dep_vars[v]));
 264
 265       for (grp = 0; grp < n_cats; ++grp)
 266         {
 267           const struct exploratory_stats *es =
 268             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
 269
 270           if ( y_min > es[v].minimum)
 271             y_min = es[v].minimum;
 272
 273           if ( y_max < es[v].maximum)
 274             y_max = es[v].maximum;
 275         }
 276
 277       boxplot = boxplot_create (y_min, y_max, ds_cstr (&title));
 278
 279       ds_destroy (&title);
 280
 281       for (grp = 0; grp < n_cats; ++grp)
 282         {
 283           int ivar_idx;
 284           struct string label;
 285
 286           const struct ccase *c =
 287             categoricals_get_case_by_category_real (cmd->cats,  iact_idx, grp);
 288
 289           const struct exploratory_stats *es =
 290             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
 291
 292           ds_init_empty (&label);
 293           for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
 294             {
 295               const struct variable *ivar = iact->vars[ivar_idx];
 296               const union value *val = case_data (c, ivar);
 297
 298               ds_put_cstr (&label, var_to_string (ivar));
 299               ds_put_cstr (&label, " = ");
 300               var_append_value_name (ivar, val, &label);
 301               ds_put_cstr (&label, "; ");
 302             }
 303
 304           boxplot_add_box (boxplot, es[v].box_whisker, ds_cstr (&label));
 305
 306           ds_destroy (&label);
 307         }
 308
 309       boxplot_submit (boxplot);
 310     }
 311 }
 312
 313 static void
 314 show_boxplot_variabled (const struct examine *cmd, int iact_idx)
 315 {
 316   int grp;
 317   const struct interaction *iact = cmd->iacts[iact_idx];
 318   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
 319
 320   for (grp = 0; grp < n_cats; ++grp)
 321     {
 322       struct boxplot *boxplot;
 323       int v;
 324       double y_min = DBL_MAX;
 325       double y_max = -DBL_MAX;
 326
 327       const struct ccase *c =
 328         categoricals_get_case_by_category_real (cmd->cats,  iact_idx, grp);
 329
 330       struct string title;
 331       ds_init_empty (&title);
 332
 333       for (v = 0; v < cmd->n_dep_vars; ++v)
 334         {
 335           const struct exploratory_stats *es =
 336             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
 337
 338           if ( y_min > es[v].minimum)
 339             y_min = es[v].minimum;
 340
 341           if ( y_max < es[v].maximum)
 342             y_max = es[v].maximum;
 343         }
 344
 345       if ( iact->n_vars == 0)
 346         ds_put_format (&title, _("Boxplot"));
 347       else
 348         {
 349           int ivar_idx;
 350           struct string label;
 351           ds_init_empty (&label);
 352           for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
 353             {
 354               const struct variable *ivar = iact->vars[ivar_idx];
 355               const union value *val = case_data (c, ivar);
 356
 357               ds_put_cstr (&label, var_to_string (ivar));
 358               ds_put_cstr (&label, " = ");
 359               var_append_value_name (ivar, val, &label);
 360               ds_put_cstr (&label, "; ");
 361             }
 362
 363           ds_put_format (&title, _("Boxplot of %s"),
 364                          ds_cstr (&label));
 365
 366           ds_destroy (&label);
 367         }
 368
 369       boxplot = boxplot_create (y_min, y_max, ds_cstr (&title));
 370
 371       ds_destroy (&title);
 372
 373       for (v = 0; v < cmd->n_dep_vars; ++v)
 374         {
 375           const struct exploratory_stats *es =
 376             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
 377
 378           boxplot_add_box (boxplot, es[v].box_whisker,
 379                            var_to_string (cmd->dep_vars[v]));
 380         }
 381
 382       boxplot_submit (boxplot);
 383     }
 384 }
 385
 386
 387 static void
 388 show_npplot (const struct examine *cmd, int iact_idx)
 389 {
 390   const struct interaction *iact = cmd->iacts[iact_idx];
 391   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
 392
 393   int v;
 394
 395   for (v = 0; v < cmd->n_dep_vars; ++v)
 396     {
 397       int grp;
 398       for (grp = 0; grp < n_cats; ++grp)
 399         {
 400           struct chart_item *npp, *dnpp;
 401           struct casereader *reader;
 402           struct np *np;
 403
 404           int ivar_idx;
 405           const struct ccase *c =
 406             categoricals_get_case_by_category_real (cmd->cats,
 407                                                     iact_idx, grp);
 408
 409           const struct exploratory_stats *es =
 410             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
 411
 412           struct string label;
 413           ds_init_cstr (&label,
 414                         var_to_string (cmd->dep_vars[v]));
 415
 416           if ( iact->n_vars > 0)
 417             {
 418               ds_put_cstr (&label, " (");
 419               for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
 420                 {
 421                   const struct variable *ivar = iact->vars[ivar_idx];
 422                   const union value *val = case_data (c, ivar);
 423
 424                   ds_put_cstr (&label, var_to_string (ivar));
 425                   ds_put_cstr (&label, " = ");
 426                   var_append_value_name (ivar, val, &label);
 427                   ds_put_cstr (&label, "; ");
 428
 429                 }
 430               ds_put_cstr (&label, ")");
 431             }
 432
 433           np = es[v].np;
 434           reader = casewriter_make_reader (np->writer);
 435           np->writer = NULL;
 436
 437           npp = np_plot_create (np, reader, ds_cstr (&label));
 438           dnpp = dnp_plot_create (np, reader, ds_cstr (&label));
 439
 440           if (npp == NULL || dnpp == NULL)
 441             {
 442               msg (MW, _("Not creating NP plot because data set is empty."));
 443               chart_item_unref (npp);
 444               chart_item_unref (dnpp);
 445             }
 446           else
 447             {
 448               chart_item_submit (npp);
 449               chart_item_submit (dnpp);
 450             }
 451           casereader_destroy (reader);
 452
 453           ds_destroy (&label);
 454         }
 455     }
 456 }
 457
 458
 459 static void
 460 show_histogram (const struct examine *cmd, int iact_idx)
 461 {
 462   const struct interaction *iact = cmd->iacts[iact_idx];
 463   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
 464
 465   int v;
 466
 467   for (v = 0; v < cmd->n_dep_vars; ++v)
 468     {
 469       int grp;
 470       for (grp = 0; grp < n_cats; ++grp)
 471         {
 472           double n, mean, var;
 473           int ivar_idx;
 474           const struct ccase *c =
 475             categoricals_get_case_by_category_real (cmd->cats,
 476                                                     iact_idx, grp);
 477
 478           const struct exploratory_stats *es =
 479             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
 480
 481           struct string label;
 482
 483           if (es[v].histogram == NULL)
 484             continue;
 485
 486           ds_init_cstr (&label,
 487                         var_to_string (cmd->dep_vars[v]));
 488
 489           if ( iact->n_vars > 0)
 490             {
 491               ds_put_cstr (&label, " (");
 492               for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
 493                 {
 494                   const struct variable *ivar = iact->vars[ivar_idx];
 495                   const union value *val = case_data (c, ivar);
 496
 497                   ds_put_cstr (&label, var_to_string (ivar));
 498                   ds_put_cstr (&label, " = ");
 499                   var_append_value_name (ivar, val, &label);
 500                   ds_put_cstr (&label, "; ");
 501
 502                 }
 503               ds_put_cstr (&label, ")");
 504             }
 505
 506
 507           moments_calculate (es[v].mom, &n, &mean, &var, NULL, NULL);
 508
 509           chart_item_submit
 510             ( histogram_chart_create (es[v].histogram->gsl_hist,
 511                                       ds_cstr (&label), n, mean,
 512                                       sqrt (var), false));
 513
 514
 515           ds_destroy (&label);
 516         }
 517     }
 518 }
 519
 520 static void
 521 percentiles_report (const struct examine *cmd, int iact_idx)
 522 {
 523   const struct interaction *iact = cmd->iacts[iact_idx];
 524   int i, v;
 525   const int heading_columns = 1 + iact->n_vars + 1;
 526   const int heading_rows = 2;
 527   struct tab_table *t;
 528
 529   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
 530
 531   const int rows_per_cat = 2;
 532   const int rows_per_var = n_cats * rows_per_cat;
 533
 534   const int nr = heading_rows + cmd->n_dep_vars * rows_per_var;
 535   const int nc = heading_columns + cmd->n_percentiles;
 536
 537   t = tab_create (nc, nr);
 538   tab_title (t, _("Percentiles"));
 539
 540   tab_headers (t, heading_columns, 0, heading_rows, 0);
 541
 542   /* Internal Vertical lines */
 543   tab_box (t, -1, -1, -1, TAL_1,
 544            heading_columns, 0, nc - 1, nr - 1);
 545
 546   /* External Frame */
 547   tab_box (t, TAL_2, TAL_2, -1, -1,
 548            0, 0, nc - 1, nr - 1);
 549
 550   tab_hline (t, TAL_2, 0, nc - 1, heading_rows);
 551   tab_vline (t, TAL_2, heading_columns, 0, nr - 1);
 552
 553   tab_joint_text (t, heading_columns, 0,
 554                   nc - 1, 0,
 555                   TAT_TITLE | TAB_CENTER,
 556                   _("Percentiles")
 557                   );
 558
 559   tab_hline (t, TAL_1, heading_columns, nc - 1, 1);
 560
 561
 562   for (i = 0; i < cmd->n_percentiles; ++i)
 563     {
 564       tab_text_format (t, heading_columns + i, 1,
 565                        TAT_TITLE | TAB_CENTER,
 566                        _("%g"), cmd->ptiles[i]);
 567     }
 568
 569   for (i = 0; i < iact->n_vars; ++i)
 570     {
 571       tab_text (t,
 572                 1 + i, 1,
 573                 TAT_TITLE,
 574                 var_to_string (iact->vars[i])
 575                 );
 576     }
 577
 578
 579
 580   if (n_cats > 0)
 581     {
 582       tab_vline (t, TAL_1, heading_columns - 1, heading_rows, nr - 1);
 583
 584       for (v = 0; v < cmd->n_dep_vars; ++v)
 585         {
 586           const union value **prev_vals = previous_value_alloc (iact);
 587
 588           int ivar_idx;
 589           if ( v > 0 )
 590             tab_hline (t, TAL_1, 0, nc - 1, heading_rows + v * rows_per_var);
 591
 592           tab_text (t,
 593                     0, heading_rows + v * rows_per_var,
 594                     TAT_TITLE | TAB_LEFT,
 595                     var_to_string (cmd->dep_vars[v])
 596                     );
 597
 598           for (i = 0; i < n_cats; ++i)
 599             {
 600               const struct ccase *c =
 601                 categoricals_get_case_by_category_real (cmd->cats,
 602                                                         iact_idx, i);
 603
 604               const struct exploratory_stats *ess =
 605                 categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, i);
 606
 607               const struct exploratory_stats *es = ess + v;
 608
 609               int diff_idx = previous_value_record (iact, c, prev_vals);
 610
 611               double hinges[3];
 612               int p;
 613
 614               for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
 615                 {
 616                   const struct variable *ivar = iact->vars[ivar_idx];
 617                   const union value *val = case_data (c, ivar);
 618
 619                   if (( diff_idx != -1 && diff_idx <= ivar_idx)
 620                       || i == 0)
 621                     {
 622                       struct string str;
 623                       ds_init_empty (&str);
 624                       var_append_value_name (ivar, val, &str);
 625
 626                       tab_text (t,
 627                                 1 + ivar_idx,
 628                                 heading_rows + v * rows_per_var + i * rows_per_cat,
 629                                 TAT_TITLE | TAB_LEFT,
 630                                 ds_cstr (&str)
 631                                 );
 632
 633                       ds_destroy (&str);
 634                     }
 635                 }
 636
 637               if ( diff_idx != -1 && diff_idx < iact->n_vars)
 638                 {
 639                   tab_hline (t, TAL_1, 1 + diff_idx, nc - 1,
 640                              heading_rows + v * rows_per_var + i * rows_per_cat
 641                              );
 642                 }
 643
 644               tab_text (t, heading_columns - 1,
 645                         heading_rows + v * rows_per_var + i * rows_per_cat,
 646                         TAT_TITLE | TAB_LEFT,
 647                         gettext (ptile_alg_desc [cmd->pc_alg]));
 648
 649               tukey_hinges_calculate (es->hinges, hinges);
 650
 651               for (p = 0; p < cmd->n_percentiles; ++p)
 652                 {
 653                   tab_double (t, heading_columns + p,
 654                               heading_rows + v * rows_per_var + i * rows_per_cat,
 655                               0,
 656                               percentile_calculate (es->percentiles[p], cmd->pc_alg),
 657                               0);
 658
 659                   if (cmd->ptiles[p] == 25.0)
 660                     {
 661                       tab_double (t, heading_columns + p,
 662                                   heading_rows + v * rows_per_var + i * rows_per_cat + 1,
 663                                   0,
 664                                   hinges[0],
 665                                   0);
 666                     }
 667                   else if (cmd->ptiles[p] == 50.0)
 668                     {
 669                       tab_double (t, heading_columns + p,
 670                                   heading_rows + v * rows_per_var + i * rows_per_cat + 1,
 671                                   0,
 672                                   hinges[1],
 673                                   0);
 674                     }
 675                   else if (cmd->ptiles[p] == 75.0)
 676                     {
 677                       tab_double (t, heading_columns + p,
 678                                   heading_rows + v * rows_per_var + i * rows_per_cat + 1,
 679                                   0,
 680                                   hinges[2],
 681                                   0);
 682                     }
 683                 }
 684
 685
 686               tab_text (t, heading_columns - 1,
 687                         heading_rows + v * rows_per_var + i * rows_per_cat + 1,
 688                         TAT_TITLE | TAB_LEFT,
 689                         _("Tukey's Hinges"));
 690
 691             }
 692
 693           free (prev_vals);
 694         }
 695     }
 696   tab_submit (t);
 697 }
 698
 699 static void
 700 descriptives_report (const struct examine *cmd, int iact_idx)
 701 {
 702   const struct interaction *iact = cmd->iacts[iact_idx];
 703   int i, v;
 704   const int heading_columns = 1 + iact->n_vars + 2;
 705   const int heading_rows = 1;
 706   struct tab_table *t;
 707
 708   size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
 709
 710   const int rows_per_cat = 13;
 711   const int rows_per_var = n_cats * rows_per_cat;
 712
 713   const int nr = heading_rows + cmd->n_dep_vars * rows_per_var;
 714   const int nc = 2 + heading_columns;
 715
 716   t = tab_create (nc, nr);
 717   tab_title (t, _("Descriptives"));
 718
 719   tab_headers (t, heading_columns, 0, heading_rows, 0);
 720
 721   /* Internal Vertical lines */
 722   tab_box (t, -1, -1, -1, TAL_1,
 723            heading_columns, 0, nc - 1, nr - 1);
 724
 725   /* External Frame */
 726   tab_box (t, TAL_2, TAL_2, -1, -1,
 727            0, 0, nc - 1, nr - 1);
 728
 729   tab_hline (t, TAL_2, 0, nc - 1, heading_rows);
 730   tab_vline (t, TAL_2, heading_columns, 0, nr - 1);
 731
 732
 733   tab_text (t, heading_columns, 0, TAB_CENTER | TAT_TITLE,
 734             _("Statistic"));
 735
 736   tab_text (t, heading_columns + 1, 0, TAB_CENTER | TAT_TITLE,
 737             _("Std. Error"));
 738
 739   for (i = 0; i < iact->n_vars; ++i)
 740     {
 741       tab_text (t,
 742                 1 + i, 0,
 743                 TAT_TITLE,
 744                 var_to_string (iact->vars[i])
 745                 );
 746     }
 747
 748   for (v = 0; v < cmd->n_dep_vars; ++v)
 749     {
 750       const union value **prev_val = previous_value_alloc (iact);
 751
 752       int ivar_idx;
 753       if ( v > 0 )
 754         tab_hline (t, TAL_1, 0, nc - 1, heading_rows + v * rows_per_var);
 755
 756       tab_text (t,
 757                 0, heading_rows + v * rows_per_var,
 758                 TAT_TITLE | TAB_LEFT,
 759                 var_to_string (cmd->dep_vars[v])
 760                 );
 761
 762       for (i = 0; i < n_cats; ++i)
 763         {
 764           const struct ccase *c =
 765             categoricals_get_case_by_category_real (cmd->cats,
 766                                                     iact_idx, i);
 767
 768           const struct exploratory_stats *ess =
 769             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, i);
 770
 771           const struct exploratory_stats *es = ess + v;
 772
 773           const int diff_idx = previous_value_record (iact, c, prev_val);
 774
 775           double m0, m1, m2, m3, m4;
 776           double tval;
 777
 778           moments_calculate (es->mom, &m0, &m1, &m2, &m3, &m4);
 779
 780           tval = gsl_cdf_tdist_Qinv ((1.0 - cmd->conf) / 2.0, m0 - 1.0);
 781
 782           for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
 783             {
 784               const struct variable *ivar = iact->vars[ivar_idx];
 785               const union value *val = case_data (c, ivar);
 786
 787               if (( diff_idx != -1 && diff_idx <= ivar_idx)
 788                   || i == 0)
 789                 {
 790                   struct string str;
 791                   ds_init_empty (&str);
 792                   var_append_value_name (ivar, val, &str);
 793
 794                   tab_text (t,
 795                             1 + ivar_idx,
 796                             heading_rows + v * rows_per_var + i * rows_per_cat,
 797                             TAT_TITLE | TAB_LEFT,
 798                             ds_cstr (&str)
 799                             );
 800
 801                   ds_destroy (&str);
 802                 }
 803             }
 804
 805           if ( diff_idx != -1 && diff_idx < iact->n_vars)
 806             {
 807               tab_hline (t, TAL_1, 1 + diff_idx, nc - 1,
 808                          heading_rows + v * rows_per_var + i * rows_per_cat
 809                          );
 810             }
 811
 812           tab_text (t,
 813                     1 + iact->n_vars,
 814                     heading_rows + v * rows_per_var + i * rows_per_cat,
 815                     TAB_LEFT,
 816                     _("Mean")
 817                     );
 818
 819           tab_double (t,
 820                       1 + iact->n_vars + 2,
 821                       heading_rows + v * rows_per_var + i * rows_per_cat,
 822                       0, m1, 0);
 823
 824           tab_double (t,
 825                       1 + iact->n_vars + 3,
 826                       heading_rows + v * rows_per_var + i * rows_per_cat,
 827                       0, calc_semean (m2, m0), 0);
 828
 829           tab_text_format (t,
 830                            1 + iact->n_vars,
 831                            heading_rows + v * rows_per_var + i * rows_per_cat + 1,
 832                            TAB_LEFT,
 833                            _("%g%% Confidence Interval for Mean"),
 834                            cmd->conf * 100.0
 835                            );
 836
 837           tab_text (t,
 838                     1 + iact->n_vars + 1,
 839                     heading_rows + v * rows_per_var + i * rows_per_cat + 1,
 840                     TAB_LEFT,
 841                     _("Lower Bound")
 842                     );
 843
 844           tab_double (t,
 845                       1 + iact->n_vars + 2,
 846                       heading_rows + v * rows_per_var + i * rows_per_cat + 1,
 847                       0, m1 - tval * calc_semean (m2, m0), 0);
 848
 849
 850           tab_text (t,
 851                     1 + iact->n_vars + 1,
 852                     heading_rows + v * rows_per_var + i * rows_per_cat + 2,
 853                     TAB_LEFT,
 854                     _("Upper Bound")
 855                     );
 856
 857           tab_double (t,
 858                       1 + iact->n_vars + 2,
 859                       heading_rows + v * rows_per_var + i * rows_per_cat + 2,
 860                       0, m1 + tval * calc_semean (m2, m0), 0);
 861
 862
 863           tab_text (t,
 864                     1 + iact->n_vars,
 865                     heading_rows + v * rows_per_var + i * rows_per_cat + 3,
 866                     TAB_LEFT,
 867                     _("5% Trimmed Mean")
 868                     );
 869
 870           tab_double (t,
 871                       1 + iact->n_vars + 2,
 872                       heading_rows + v * rows_per_var + i * rows_per_cat + 3,
 873                       0,
 874                       trimmed_mean_calculate (es->trimmed_mean),
 875                       0);
 876
 877           tab_text (t,
 878                     1 + iact->n_vars,
 879                     heading_rows + v * rows_per_var + i * rows_per_cat + 4,
 880                     TAB_LEFT,
 881                     _("Median")
 882                     );
 883
 884           tab_double (t,
 885                       1 + iact->n_vars + 2,
 886                       heading_rows + v * rows_per_var + i * rows_per_cat + 4,
 887                       0,
 888                       percentile_calculate (es->quartiles[1], cmd->pc_alg),
 889                       0);
 890
 891
 892           tab_text (t,
 893                     1 + iact->n_vars,
 894                     heading_rows + v * rows_per_var + i * rows_per_cat + 5,
 895                     TAB_LEFT,
 896                     _("Variance")
 897                     );
 898
 899           tab_double (t,
 900                       1 + iact->n_vars + 2,
 901                       heading_rows + v * rows_per_var + i * rows_per_cat + 5,
 902                       0, m2, 0);
 903
 904           tab_text (t,
 905                     1 + iact->n_vars,
 906                     heading_rows + v * rows_per_var + i * rows_per_cat + 6,
 907                     TAB_LEFT,
 908                     _("Std. Deviation")
 909                     );
 910
 911           tab_double (t,
 912                       1 + iact->n_vars + 2,
 913                       heading_rows + v * rows_per_var + i * rows_per_cat + 6,
 914                       0, sqrt (m2), 0);
 915
 916           tab_text (t,
 917                     1 + iact->n_vars,
 918                     heading_rows + v * rows_per_var + i * rows_per_cat + 7,
 919                     TAB_LEFT,
 920                     _("Minimum")
 921                     );
 922
 923           tab_double (t,
 924                       1 + iact->n_vars + 2,
 925                       heading_rows + v * rows_per_var + i * rows_per_cat + 7,
 926                       0,
 927                       es->minima[0].val,
 928                       0);
 929
 930           tab_text (t,
 931                     1 + iact->n_vars,
 932                     heading_rows + v * rows_per_var + i * rows_per_cat + 8,
 933                     TAB_LEFT,
 934                     _("Maximum")
 935                     );
 936
 937           tab_double (t,
 938                       1 + iact->n_vars + 2,
 939                       heading_rows + v * rows_per_var + i * rows_per_cat + 8,
 940                       0,
 941                       es->maxima[0].val,
 942                       0);
 943
 944           tab_text (t,
 945                     1 + iact->n_vars,
 946                     heading_rows + v * rows_per_var + i * rows_per_cat + 9,
 947                     TAB_LEFT,
 948                     _("Range")
 949                     );
 950
 951           tab_double (t,
 952                       1 + iact->n_vars + 2,
 953                       heading_rows + v * rows_per_var + i * rows_per_cat + 9,
 954                       0,
 955                       es->maxima[0].val - es->minima[0].val,
 956                       0);
 957
 958           tab_text (t,
 959                     1 + iact->n_vars,
 960                     heading_rows + v * rows_per_var + i * rows_per_cat + 10,
 961                     TAB_LEFT,
 962                     _("Interquartile Range")
 963                     );
 964
 965
 966           tab_double (t,
 967                       1 + iact->n_vars + 2,
 968                       heading_rows + v * rows_per_var + i * rows_per_cat + 10,
 969                       0,
 970                       percentile_calculate (es->quartiles[2], cmd->pc_alg) -
 971                       percentile_calculate (es->quartiles[0], cmd->pc_alg),
 972                       0);
 973
 974
 975
 976
 977           tab_text (t,
 978                     1 + iact->n_vars,
 979                     heading_rows + v * rows_per_var + i * rows_per_cat + 11,
 980                     TAB_LEFT,
 981                     _("Skewness")
 982                     );
 983
 984           tab_double (t,
 985                       1 + iact->n_vars + 2,
 986                       heading_rows + v * rows_per_var + i * rows_per_cat + 11,
 987                       0, m3, 0);
 988
 989           tab_double (t,
 990                       1 + iact->n_vars + 3,
 991                       heading_rows + v * rows_per_var + i * rows_per_cat + 11,
 992                       0, calc_seskew (m0), 0);
 993
 994           tab_text (t,
 995                     1 + iact->n_vars,
 996                     heading_rows + v * rows_per_var + i * rows_per_cat + 12,
 997                     TAB_LEFT,
 998                     _("Kurtosis")
 999                     );
1000
1001           tab_double (t,
1002                       1 + iact->n_vars + 2,
1003                       heading_rows + v * rows_per_var + i * rows_per_cat + 12,
1004                       0, m4, 0);
1005
1006           tab_double (t,
1007                       1 + iact->n_vars + 3,
1008                       heading_rows + v * rows_per_var + i * rows_per_cat + 12,
1009                       0, calc_sekurt (m0), 0);
1010         }
1011
1012       free (prev_val);
1013     }
1014   tab_submit (t);
1015 }
1016
1017
1018 static void
1019 extremes_report (const struct examine *cmd, int iact_idx)
1020 {
1021   const struct interaction *iact = cmd->iacts[iact_idx];
1022   int i, v;
1023   const int heading_columns = 1 + iact->n_vars + 2;
1024   const int heading_rows = 1;
1025   struct tab_table *t;
1026
1027   size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
1028
1029   const int rows_per_cat = 2 * cmd->disp_extremes;
1030   const int rows_per_var = n_cats * rows_per_cat;
1031
1032   const int nr = heading_rows + cmd->n_dep_vars * rows_per_var;
1033   const int nc = 2 + heading_columns;
1034
1035   t = tab_create (nc, nr);
1036   tab_title (t, _("Extreme Values"));
1037
1038   tab_headers (t, heading_columns, 0, heading_rows, 0);
1039
1040   /* Internal Vertical lines */
1041   tab_box (t, -1, -1, -1, TAL_1,
1042            heading_columns, 0, nc - 1, nr - 1);
1043
1044   /* External Frame */
1045   tab_box (t, TAL_2, TAL_2, -1, -1,
1046            0, 0, nc - 1, nr - 1);
1047
1048   tab_hline (t, TAL_2, 0, nc - 1, heading_rows);
1049   tab_vline (t, TAL_2, heading_columns, 0, nr - 1);
1050
1051
1052   if ( cmd->id_var )
1053     tab_text (t, heading_columns, 0, TAB_CENTER | TAT_TITLE,
1054               var_to_string (cmd->id_var));
1055   else
1056     tab_text (t, heading_columns, 0, TAB_CENTER | TAT_TITLE,
1057               _("Case Number"));
1058
1059   tab_text (t, heading_columns + 1, 0, TAB_CENTER | TAT_TITLE,
1060             _("Value"));
1061
1062   for (i = 0; i < iact->n_vars; ++i)
1063     {
1064       tab_text (t,
1065                 1 + i, 0,
1066                 TAT_TITLE,
1067                 var_to_string (iact->vars[i])
1068                 );
1069     }
1070
1071   for (v = 0; v < cmd->n_dep_vars; ++v)
1072     {
1073       const union value **prev_val = previous_value_alloc (iact);
1074
1075       int ivar_idx;
1076       if ( v > 0 )
1077         tab_hline (t, TAL_1, 0, nc - 1, heading_rows + v * rows_per_var);
1078
1079       tab_text (t,
1080                 0, heading_rows + v * rows_per_var,
1081                 TAT_TITLE,
1082                 var_to_string (cmd->dep_vars[v])
1083                 );
1084
1085       for (i = 0; i < n_cats; ++i)
1086         {
1087           int e;
1088           const struct ccase *c =
1089             categoricals_get_case_by_category_real (cmd->cats, iact_idx, i);
1090
1091           const struct exploratory_stats *ess =
1092             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, i);
1093
1094           const struct exploratory_stats *es = ess + v;
1095
1096           int diff_idx = previous_value_record (iact, c, prev_val);
1097
1098           for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
1099             {
1100               const struct variable *ivar = iact->vars[ivar_idx];
1101               const union value *val = case_data (c, ivar);
1102
1103               if (( diff_idx != -1 && diff_idx <= ivar_idx)
1104                   || i == 0)
1105                 {
1106                   struct string str;
1107                   ds_init_empty (&str);
1108                   var_append_value_name (ivar, val, &str);
1109
1110                   tab_text (t,
1111                             1 + ivar_idx,
1112                             heading_rows + v * rows_per_var + i * rows_per_cat,
1113                             TAT_TITLE | TAB_LEFT,
1114                             ds_cstr (&str)
1115                             );
1116
1117                   ds_destroy (&str);
1118                 }
1119             }
1120
1121           if ( diff_idx != -1 && diff_idx < iact->n_vars)
1122             {
1123               tab_hline (t, TAL_1, 1 + diff_idx, nc - 1,
1124                          heading_rows + v * rows_per_var + i * rows_per_cat
1125                          );
1126             }
1127
1128           tab_text (t,
1129                     heading_columns - 2,
1130                     heading_rows + v * rows_per_var + i * rows_per_cat,
1131                     TAB_RIGHT,
1132                     _("Highest"));
1133
1134
1135           tab_hline (t, TAL_1, heading_columns - 2, nc - 1,
1136                      heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes
1137                      );
1138
1139           tab_text (t,
1140                     heading_columns - 2,
1141                     heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes,
1142                     TAB_RIGHT,
1143                     _("Lowest"));
1144
1145           for (e = 0 ; e < cmd->disp_extremes; ++e)
1146             {
1147               tab_double (t,
1148                           heading_columns - 1,
1149                           heading_rows + v * rows_per_var + i * rows_per_cat + e,
1150                           TAB_RIGHT,
1151                           e + 1,
1152                           &F_8_0);
1153
1154               /* The casenumber */
1155               if (cmd->id_var)
1156                 tab_value (t,
1157                            heading_columns,
1158                            heading_rows + v * rows_per_var + i * rows_per_cat + e,
1159                            TAB_RIGHT,
1160                            &es->maxima[e].identity,
1161                            cmd->id_var,
1162                            NULL);
1163               else
1164                 tab_double (t,
1165                           heading_columns,
1166                             heading_rows + v * rows_per_var + i * rows_per_cat + e,
1167                             TAB_RIGHT,
1168                             es->maxima[e].identity.f,
1169                             &F_8_0);
1170
1171               tab_double (t,
1172                           heading_columns + 1,
1173                           heading_rows + v * rows_per_var + i * rows_per_cat + e,
1174                           0,
1175                           es->maxima[e].val,
1176                           0);
1177
1178
1179
1180               tab_double (t,
1181                           heading_columns - 1,
1182                           heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes + e,
1183                           TAB_RIGHT,
1184                           e + 1,
1185                           &F_8_0);
1186
1187               /* The casenumber */
1188               if (cmd->id_var)
1189                 tab_value (t,
1190                            heading_columns,
1191                            heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes + e,
1192                            TAB_RIGHT,
1193                            &es->minima[e].identity,
1194                            cmd->id_var,
1195                            NULL);
1196               else
1197                 tab_double (t,
1198                             heading_columns,
1199                             heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes + e,
1200                             TAB_RIGHT,
1201                             es->minima[e].identity.f,
1202                             &F_8_0);
1203
1204               tab_double (t,
1205                           heading_columns + 1,
1206                           heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes + e,
1207                           0,
1208                           es->minima[e].val,
1209                           0);
1210             }
1211         }
1212       free (prev_val);
1213     }
1214
1215   tab_submit (t);
1216 }
1217
1218
1219 static void
1220 summary_report (const struct examine *cmd, int iact_idx)
1221 {
1222   const struct interaction *iact = cmd->iacts[iact_idx];
1223   int i, v;
1224   const int heading_columns = 1 + iact->n_vars;
1225   const int heading_rows = 3;
1226   struct tab_table *t;
1227
1228   const struct fmt_spec *wfmt = cmd->wv ? var_get_print_format (cmd->wv) : &F_8_0;
1229
1230   size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
1231
1232   const int nr = heading_rows + n_cats * cmd->n_dep_vars;
1233   const int nc = 6 + heading_columns;
1234
1235   t = tab_create (nc, nr);
1236   tab_title (t, _("Case Processing Summary"));
1237
1238   tab_headers (t, heading_columns, 0, heading_rows, 0);
1239
1240   /* Internal Vertical lines */
1241   tab_box (t, -1, -1, -1, TAL_1,
1242            heading_columns, 0, nc - 1, nr - 1);
1243
1244   /* External Frame */
1245   tab_box (t, TAL_2, TAL_2, -1, -1,
1246            0, 0, nc - 1, nr - 1);
1247
1248   tab_hline (t, TAL_2, 0, nc - 1, heading_rows);
1249   tab_vline (t, TAL_2, heading_columns, 0, nr - 1);
1250
1251   tab_joint_text (t, heading_columns, 0,
1252                   nc - 1, 0, TAB_CENTER | TAT_TITLE, _("Cases"));
1253   tab_joint_text (t,
1254                   heading_columns, 1,
1255                   heading_columns + 1, 1,
1256                   TAB_CENTER | TAT_TITLE, _("Valid"));
1257
1258   tab_joint_text (t,
1259                   heading_columns + 2, 1,
1260                   heading_columns + 3, 1,
1261                   TAB_CENTER | TAT_TITLE, _("Missing"));
1262
1263   tab_joint_text (t,
1264                   heading_columns + 4, 1,
1265                   heading_columns + 5, 1,
1266                   TAB_CENTER | TAT_TITLE, _("Total"));
1267
1268   for (i = 0; i < 3; ++i)
1269     {
1270       tab_text (t, heading_columns + i * 2, 2, TAB_CENTER | TAT_TITLE,
1271                 _("N"));
1272       tab_text (t, heading_columns + i * 2 + 1, 2, TAB_CENTER | TAT_TITLE,
1273                 _("Percent"));
1274     }
1275
1276   for (i = 0; i < iact->n_vars; ++i)
1277     {
1278       tab_text (t,
1279                 1 + i, 2,
1280                 TAT_TITLE,
1281                 var_to_string (iact->vars[i])
1282                 );
1283     }
1284
1285   if (n_cats > 0)
1286     for (v = 0; v < cmd->n_dep_vars; ++v)
1287       {
1288         int ivar_idx;
1289         const union value **prev_values = previous_value_alloc (iact);
1290
1291         if ( v > 0 )
1292           tab_hline (t, TAL_1, 0, nc - 1, heading_rows + v * n_cats);
1293
1294         tab_text (t,
1295                   0, heading_rows + n_cats * v,
1296                   TAT_TITLE,
1297                   var_to_string (cmd->dep_vars[v])
1298                   );
1299
1300
1301         for (i = 0; i < n_cats; ++i)
1302           {
1303             double total;
1304             const struct exploratory_stats *es;
1305
1306             const struct ccase *c =
1307               categoricals_get_case_by_category_real (cmd->cats,
1308                                                       iact_idx, i);
1309             if (c)
1310               {
1311                 int diff_idx = previous_value_record (iact, c, prev_values);
1312
1313                 if ( diff_idx != -1 && diff_idx < iact->n_vars - 1)
1314                   tab_hline (t, TAL_1, 1 + diff_idx, nc - 1,
1315                              heading_rows + n_cats * v + i );
1316
1317                 for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
1318                   {
1319                     const struct variable *ivar = iact->vars[ivar_idx];
1320                     const union value *val = case_data (c, ivar);
1321
1322                     if (( diff_idx != -1 && diff_idx <= ivar_idx)
1323                         || i == 0)
1324                       {
1325                         struct string str;
1326                         ds_init_empty (&str);
1327                         var_append_value_name (ivar, val, &str);
1328
1329                         tab_text (t,
1330                                   1 + ivar_idx, heading_rows + n_cats * v + i,
1331                                   TAT_TITLE | TAB_LEFT,
1332                                   ds_cstr (&str)
1333                                   );
1334
1335                         ds_destroy (&str);
1336                       }
1337                   }
1338               }
1339
1340
1341             es = categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, i);
1342
1343
1344             total = es[v].missing + es[v].non_missing;
1345             tab_double (t,
1346                         heading_columns + 0,
1347                         heading_rows + n_cats * v + i,
1348                         0,
1349                         es[v].non_missing,
1350                         wfmt);
1351
1352
1353             tab_text_format (t,
1354                              heading_columns + 1,
1355                              heading_rows + n_cats * v + i,
1356                              0,
1357                              "%g%%",
1358                              100.0 * es[v].non_missing / total
1359                              );
1360
1361
1362             tab_double (t,
1363                         heading_columns + 2,
1364                         heading_rows + n_cats * v + i,
1365                         0,
1366                         es[v].missing,
1367                         wfmt);
1368
1369             tab_text_format (t,
1370                              heading_columns + 3,
1371                              heading_rows + n_cats * v + i,
1372                              0,
1373                              "%g%%",
1374                              100.0 * es[v].missing / total
1375                              );
1376             tab_double (t,
1377                         heading_columns + 4,
1378                         heading_rows + n_cats * v + i,
1379                         0,
1380                         total,
1381                         wfmt);
1382
1383             /* This can only be 100% can't it? */
1384             tab_text_format (t,
1385                              heading_columns + 5,
1386                              heading_rows + n_cats * v + i,
1387                              0,
1388                              "%g%%",
1389                              100.0 * (es[v].missing + es[v].non_missing)/ total
1390                              );
1391           }
1392         free (prev_values);
1393       }
1394
1395   tab_hline (t, TAL_1, heading_columns, nc - 1, 1);
1396   tab_hline (t, TAL_1, heading_columns, nc - 1, 2);
1397
1398   tab_submit (t);
1399 }
1400
1401
1402 /* Match a variable.
1403    If the match succeeds, the variable will be placed in VAR.
1404    Returns true if successful */
1405 static bool
1406 lex_match_variable (struct lexer *lexer,
1407                     const struct dictionary *dict, const struct variable **var)
1408 {
1409   if (lex_token (lexer) !=  T_ID)
1410
1411     return false;
1412
1413   *var = parse_variable_const  (lexer, dict);
1414
1415   if ( *var == NULL)
1416     return false;
1417   return true;
1418 }
1419
1420 /* Attempt to parse an interaction from LEXER */
1421 static struct interaction *
1422 parse_interaction (struct lexer *lexer, struct examine *ex)
1423 {
1424   const struct variable *v = NULL;
1425   struct interaction *iact = NULL;
1426
1427   if ( lex_match_variable (lexer, ex->dict, &v))
1428     {
1429       iact = interaction_create (v);
1430
1431       while (lex_match (lexer, T_BY))
1432         {
1433           if (!lex_match_variable (lexer, ex->dict, &v))
1434             {
1435               interaction_destroy (iact);
1436               return NULL;
1437             }
1438           interaction_add_variable (iact, v);
1439         }
1440       lex_match (lexer, T_COMMA);
1441     }
1442
1443   return iact;
1444 }
1445
1446
1447 static void *
1448 create_n (const void *aux1, void *aux2 UNUSED)
1449 {
1450   int v;
1451
1452   const struct examine *examine = aux1;
1453   struct exploratory_stats *es = pool_calloc (examine->pool, examine->n_dep_vars, sizeof (*es));
1454   struct subcase ordering;
1455   subcase_init (&ordering, 0, 0, SC_ASCEND);
1456
1457   for (v = 0; v < examine->n_dep_vars; v++)
1458     {
1459       es[v].sorted_writer = sort_create_writer (&ordering, examine->ex_proto);
1460       es[v].sorted_reader = NULL;
1461
1462       es[v].mom = moments_create (MOMENT_KURTOSIS);
1463       es[v].cmin = DBL_MAX;
1464
1465       es[v].maximum = -DBL_MAX;
1466       es[v].minimum =  DBL_MAX;
1467     }
1468
1469   subcase_destroy (&ordering);
1470   return es;
1471 }
1472
1473 static void
1474 update_n (const void *aux1, void *aux2 UNUSED, void *user_data,
1475           const struct ccase *c, double weight)
1476 {
1477   int v;
1478   const struct examine *examine = aux1;
1479   struct exploratory_stats *es = user_data;
1480
1481   for (v = 0; v < examine->n_dep_vars; v++)
1482     {
1483       struct ccase *outcase ;
1484       const struct variable *var = examine->dep_vars[v];
1485       const double x = case_data (c, var)->f;
1486
1487       if (var_is_value_missing (var, case_data (c, var), examine->exclude))
1488         {
1489           es[v].missing += weight;
1490           continue;
1491         }
1492
1493       outcase = case_create (examine->ex_proto);
1494
1495       if (x > es[v].maximum)
1496         es[v].maximum = x;
1497
1498       if (x < es[v].minimum)
1499         es[v].minimum =  x;
1500
1501       es[v].non_missing += weight;
1502
1503       moments_pass_one (es[v].mom, x, weight);
1504
1505       /* Save the value and the ID to the writer */
1506       assert (examine->id_idx != -1);
1507       case_data_rw_idx (outcase, EX_VAL)->f = x;
1508       value_copy (case_data_rw_idx (outcase, EX_ID),
1509                   case_data_idx (c, examine->id_idx), examine->id_width);
1510
1511       case_data_rw_idx (outcase, EX_WT)->f = weight;
1512
1513       es[v].cc += weight;
1514
1515       if (es[v].cmin > weight)
1516         es[v].cmin = weight;
1517
1518       casewriter_write (es[v].sorted_writer, outcase);
1519     }
1520 }
1521
1522 static void
1523 calculate_n (const void *aux1, void *aux2 UNUSED, void *user_data)
1524 {
1525   int v;
1526   const struct examine *examine = aux1;
1527   struct exploratory_stats *es = user_data;
1528
1529   for (v = 0; v < examine->n_dep_vars; v++)
1530     {
1531       int i;
1532       casenumber imin = 0;
1533       double imax = es[v].cc;
1534       struct casereader *reader;
1535       struct ccase *c;
1536       casenumber total_cases;
1537
1538       if (examine->histogram)
1539         {
1540           /* Sturges Rule */
1541           double bin_width = fabs (es[v].minimum - es[v].maximum)
1542             / (1 + log2 (es[v].cc))
1543             ;
1544
1545           bin_width = chart_rounded_tick (bin_width);
1546
1547           es[v].histogram =
1548             histogram_create (bin_width, es[v].minimum, es[v].maximum);
1549         }
1550
1551       es[v].sorted_reader = casewriter_make_reader (es[v].sorted_writer);
1552       total_cases = casereader_count_cases (es[v].sorted_reader);
1553       es[v].sorted_writer = NULL;
1554
1555       es[v].maxima = pool_calloc (examine->pool, examine->calc_extremes, sizeof (*es[v].maxima));
1556       es[v].minima = pool_calloc (examine->pool, examine->calc_extremes, sizeof (*es[v].minima));
1557       for (i = 0; i < examine->calc_extremes; ++i)
1558         {
1559           value_init_pool (examine->pool, &es[v].maxima[i].identity, examine->id_width) ;
1560           value_init_pool (examine->pool, &es[v].minima[i].identity, examine->id_width) ;
1561         }
1562
1563       for (reader = casereader_clone (es[v].sorted_reader);
1564            (c = casereader_read (reader)) != NULL; case_unref (c))
1565         {
1566           const double val = case_data_idx (c, EX_VAL)->f;
1567           const double wt = case_data_idx (c, EX_WT)->f;   /* FIXME: What about fractional weights ??? */
1568
1569           moments_pass_two (es[v].mom, val, wt);
1570
1571           if (es[v].histogram)
1572             histogram_add (es[v].histogram, val, wt);
1573
1574           if (imin < examine->calc_extremes)
1575             {
1576               int x;
1577               for (x = imin; x < examine->calc_extremes; ++x)
1578                 {
1579                   struct extremity *min = &es[v].minima[x];
1580                   min->val = val;
1581                   value_copy (&min->identity, case_data_idx (c, EX_ID), examine->id_width);
1582                 }
1583               imin += wt;
1584             }
1585
1586           imax -= wt;
1587           if (imax < examine->calc_extremes)
1588             {
1589               int x;
1590
1591               for (x = imax; x < imax + wt; ++x)
1592                 {
1593                   struct extremity *max;
1594
1595                   if (x >= examine->calc_extremes)
1596                     break;
1597
1598                   max = &es[v].maxima[x];
1599                   max->val = val;
1600                   value_copy (&max->identity, case_data_idx (c, EX_ID), examine->id_width);
1601                 }
1602             }
1603         }
1604       casereader_destroy (reader);
1605
1606       if (examine->calc_extremes > 0)
1607         {
1608           assert (es[v].minima[0].val == es[v].minimum);
1609           assert (es[v].maxima[0].val == es[v].maximum);
1610         }
1611
1612       {
1613         const int n_os = 5 + examine->n_percentiles;
1614         struct order_stats **os ;
1615         es[v].percentiles = pool_calloc (examine->pool, examine->n_percentiles, sizeof (*es[v].percentiles));
1616
1617         es[v].trimmed_mean = trimmed_mean_create (es[v].cc, 0.05);
1618
1619         os = xcalloc (n_os, sizeof *os);
1620         os[0] = &es[v].trimmed_mean->parent;
1621
1622         es[v].quartiles[0] = percentile_create (0.25, es[v].cc);
1623         es[v].quartiles[1] = percentile_create (0.5,  es[v].cc);
1624         es[v].quartiles[2] = percentile_create (0.75, es[v].cc);
1625
1626         os[1] = &es[v].quartiles[0]->parent;
1627         os[2] = &es[v].quartiles[1]->parent;
1628         os[3] = &es[v].quartiles[2]->parent;
1629
1630         es[v].hinges = tukey_hinges_create (es[v].cc, es[v].cmin);
1631         os[4] = &es[v].hinges->parent;
1632
1633         for (i = 0; i < examine->n_percentiles; ++i)
1634           {
1635             es[v].percentiles[i] = percentile_create (examine->ptiles[i] / 100.00, es[v].cc);
1636             os[5 + i] = &es[v].percentiles[i]->parent;
1637           }
1638
1639         order_stats_accumulate_idx (os, n_os,
1640                                     casereader_clone (es[v].sorted_reader),
1641                                     EX_WT, EX_VAL);
1642
1643         free (os);
1644       }
1645
1646       if (examine->boxplot)
1647         {
1648           struct order_stats *os;
1649
1650           es[v].box_whisker = box_whisker_create (es[v].hinges,
1651                                                   EX_ID, examine->id_var);
1652
1653           os = &es[v].box_whisker->parent;
1654           order_stats_accumulate_idx (&os, 1,
1655                                       casereader_clone (es[v].sorted_reader),
1656                                       EX_WT, EX_VAL);
1657         }
1658
1659       if (examine->npplot)
1660         {
1661           double n, mean, var;
1662           struct order_stats *os;
1663
1664           moments_calculate (es[v].mom, &n, &mean, &var, NULL, NULL);
1665
1666           es[v].np = np_create (n, mean, var);
1667
1668           os = &es[v].np->parent;
1669
1670           order_stats_accumulate_idx (&os, 1,
1671                                       casereader_clone (es[v].sorted_reader),
1672                                       EX_WT, EX_VAL);
1673         }
1674
1675     }
1676 }
1677
1678 static void
1679 cleanup_exploratory_stats (struct examine *cmd)
1680 {
1681   int i;
1682   for (i = 0; i < cmd->n_iacts; ++i)
1683     {
1684       int v;
1685       const size_t n_cats =  categoricals_n_count (cmd->cats, i);
1686
1687       for (v = 0; v < cmd->n_dep_vars; ++v)
1688         {
1689           int grp;
1690           for (grp = 0; grp < n_cats; ++grp)
1691             {
1692               int q;
1693               const struct exploratory_stats *es =
1694                 categoricals_get_user_data_by_category_real (cmd->cats, i, grp);
1695
1696               struct order_stats *os = &es[v].hinges->parent;
1697               struct statistic  *stat = &os->parent;
1698               stat->destroy (stat);
1699
1700               for (q = 0; q < 3 ; q++)
1701                 {
1702                   os = &es[v].quartiles[q]->parent;
1703                   stat = &os->parent;
1704                   stat->destroy (stat);
1705                 }
1706
1707               for (q = 0; q < cmd->n_percentiles ; q++)
1708                 {
1709                   os = &es[v].percentiles[q]->parent;
1710                   stat = &os->parent;
1711                   stat->destroy (stat);
1712                 }
1713
1714               os = &es[v].trimmed_mean->parent;
1715               stat = &os->parent;
1716               stat->destroy (stat);
1717
1718               os = &es[v].np->parent;
1719               if (os)
1720                 {
1721                   stat = &os->parent;
1722                   stat->destroy (stat);
1723                 }
1724
1725               statistic_destroy (&es[v].histogram->parent);
1726               moments_destroy (es[v].mom);
1727
1728               casereader_destroy (es[v].sorted_reader);
1729             }
1730         }
1731     }
1732 }
1733
1734
1735 static void
1736 run_examine (struct examine *cmd, struct casereader *input)
1737 {
1738   int i;
1739   struct ccase *c;
1740   struct casereader *reader;
1741
1742   struct payload payload;
1743   payload.create = create_n;
1744   payload.update = update_n;
1745   payload.destroy = calculate_n;
1746
1747   cmd->wv = dict_get_weight (cmd->dict);
1748
1749   cmd->cats
1750     = categoricals_create (cmd->iacts, cmd->n_iacts,
1751                            cmd->wv, cmd->exclude);
1752
1753   categoricals_set_payload (cmd->cats, &payload, cmd, NULL);
1754
1755   if (cmd->id_idx == -1)
1756     {
1757       struct ccase *c = casereader_peek (input,  0);
1758
1759       assert (cmd->id_var == NULL);
1760
1761       cmd->id_idx = case_get_value_cnt (c);
1762       input = casereader_create_arithmetic_sequence (input, 1.0, 1.0);
1763
1764       case_unref (c);
1765     }
1766
1767   /* FIXME: Filter out missing factor variables */
1768
1769   /* Remove cases on a listwise basis if requested */
1770   if ( cmd->missing_pw == false)
1771     input = casereader_create_filter_missing (input,
1772                                               cmd->dep_vars,
1773                                               cmd->n_dep_vars,
1774                                               cmd->exclude,
1775                                               NULL,
1776                                               NULL);
1777
1778   for (reader = input;
1779        (c = casereader_read (reader)) != NULL; case_unref (c))
1780     {
1781       categoricals_update (cmd->cats, c);
1782     }
1783   casereader_destroy (reader);
1784   categoricals_done (cmd->cats);
1785
1786   for (i = 0; i < cmd->n_iacts; ++i)
1787     {
1788       summary_report (cmd, i);
1789
1790       if (cmd->disp_extremes > 0)
1791         extremes_report (cmd, i);
1792
1793       if (cmd->n_percentiles > 0)
1794         percentiles_report (cmd, i);
1795
1796       if (cmd->boxplot)
1797         {
1798           switch (cmd->boxplot_mode)
1799             {
1800             case BP_GROUPS:
1801               show_boxplot_grouped (cmd, i);
1802               break;
1803             case BP_VARIABLES:
1804               show_boxplot_variabled (cmd, i);
1805               break;
1806             default:
1807               NOT_REACHED ();
1808               break;
1809             }
1810         }
1811
1812       if (cmd->histogram)
1813         show_histogram (cmd, i);
1814
1815       if (cmd->npplot)
1816         show_npplot (cmd, i);
1817
1818       if (cmd->descriptives)
1819         descriptives_report (cmd, i);
1820     }
1821
1822   cleanup_exploratory_stats (cmd);
1823   categoricals_destroy (cmd->cats);
1824 }
1825
1826
1827 int
1828 cmd_examine (struct lexer *lexer, struct dataset *ds)
1829 {
1830   int i;
1831   bool nototals_seen = false;
1832   bool totals_seen = false;
1833
1834   struct interaction **iacts_mem = NULL;
1835   struct examine examine;
1836   bool percentiles_seen = false;
1837
1838   examine.missing_pw = false;
1839   examine.disp_extremes = 0;
1840   examine.calc_extremes = 0;
1841   examine.descriptives = false;
1842   examine.conf = 0.95;
1843   examine.pc_alg = PC_HAVERAGE;
1844   examine.ptiles = NULL;
1845   examine.n_percentiles = 0;
1846   examine.id_idx = -1;
1847   examine.id_width = 0;
1848   examine.id_var = NULL;
1849   examine.boxplot_mode = BP_GROUPS;
1850
1851   examine.ex_proto = caseproto_create ();
1852
1853   examine.pool = pool_create ();
1854
1855   /* Allocate space for the first interaction.
1856      This is interaction is an empty one (for the totals).
1857      If no totals are requested, we will simply ignore this
1858      interaction.
1859   */
1860   examine.n_iacts = 1;
1861   examine.iacts = iacts_mem = pool_zalloc (examine.pool, sizeof (struct interaction *));
1862   examine.iacts[0] = interaction_create (NULL);
1863
1864   examine.exclude = MV_ANY;
1865   examine.histogram = false;
1866   examine.npplot = false;
1867   examine.boxplot = false;
1868
1869   examine.dict = dataset_dict (ds);
1870
1871   /* Accept an optional, completely pointless "/VARIABLES=" */
1872   lex_match (lexer, T_SLASH);
1873   if (lex_match_id  (lexer, "VARIABLES"))
1874     {
1875       if (! lex_force_match (lexer, T_EQUALS) )
1876         goto error;
1877     }
1878
1879   if (!parse_variables_const (lexer, examine.dict,
1880                               &examine.dep_vars, &examine.n_dep_vars,
1881                               PV_NO_DUPLICATE | PV_NUMERIC))
1882     goto error;
1883
1884   if (lex_match (lexer, T_BY))
1885     {
1886       struct interaction *iact = NULL;
1887       do
1888         {
1889           iact = parse_interaction (lexer, &examine);
1890           if (iact)
1891             {
1892               examine.n_iacts++;
1893               iacts_mem =
1894                 pool_nrealloc (examine.pool, iacts_mem,
1895                                examine.n_iacts,
1896                                sizeof (*iacts_mem));
1897
1898               iacts_mem[examine.n_iacts - 1] = iact;
1899             }
1900         }
1901       while (iact);
1902     }
1903
1904
1905   while (lex_token (lexer) != T_ENDCMD)
1906     {
1907       lex_match (lexer, T_SLASH);
1908
1909       if (lex_match_id (lexer, "STATISTICS"))
1910         {
1911           lex_match (lexer, T_EQUALS);
1912
1913           while (lex_token (lexer) != T_ENDCMD
1914                  && lex_token (lexer) != T_SLASH)
1915             {
1916               if (lex_match_id (lexer, "DESCRIPTIVES"))
1917                 {
1918                   examine.descriptives = true;
1919                 }
1920               else if (lex_match_id (lexer, "EXTREME"))
1921                 {
1922                   int extr = 5;
1923                   if (lex_match (lexer, T_LPAREN))
1924                     {
1925                       extr = lex_integer (lexer);
1926
1927                       if (extr < 0)
1928                         {
1929                           msg (MW, _("%s may not be negative. Using default value (%g)."), "EXTREME", 5.0);
1930                           extr = 5;
1931                         }
1932
1933                       lex_get (lexer);
1934                       if (! lex_force_match (lexer, T_RPAREN))
1935                         goto error;
1936                     }
1937                   examine.disp_extremes  = extr;
1938                 }
1939               else if (lex_match_id (lexer, "NONE"))
1940                 {
1941                 }
1942               else if (lex_match (lexer, T_ALL))
1943                 {
1944                   if (examine.disp_extremes == 0)
1945                     examine.disp_extremes = 5;
1946                 }
1947               else
1948                 {
1949                   lex_error (lexer, NULL);
1950                   goto error;
1951                 }
1952             }
1953         }
1954       else if (lex_match_id (lexer, "PERCENTILES"))
1955         {
1956           percentiles_seen = true;
1957           if (lex_match (lexer, T_LPAREN))
1958             {
1959               while (lex_is_number (lexer))
1960                 {
1961                   double p = lex_number (lexer);
1962
1963                   if ( p <= 0 || p >= 100.0)
1964                     {
1965                       lex_error (lexer,
1966                                  _("Percentiles must lie in the range (0, 100)"));
1967                       goto error;
1968                     }
1969
1970                   examine.n_percentiles++;
1971                   examine.ptiles =
1972                     xrealloc (examine.ptiles,
1973                               sizeof (*examine.ptiles) *
1974                               examine.n_percentiles);
1975
1976                   examine.ptiles[examine.n_percentiles - 1] = p;
1977
1978                   lex_get (lexer);
1979                   lex_match (lexer, T_COMMA);
1980                 }
1981               if (!lex_force_match (lexer, T_RPAREN))
1982                 goto error;
1983             }
1984
1985           lex_match (lexer, T_EQUALS);
1986
1987           while (lex_token (lexer) != T_ENDCMD
1988                  && lex_token (lexer) != T_SLASH)
1989             {
1990               if (lex_match_id (lexer, "HAVERAGE"))
1991                 {
1992                   examine.pc_alg = PC_HAVERAGE;
1993                 }
1994               else if (lex_match_id (lexer, "WAVERAGE"))
1995                 {
1996                   examine.pc_alg = PC_WAVERAGE;
1997                 }
1998               else if (lex_match_id (lexer, "ROUND"))
1999                 {
2000                   examine.pc_alg = PC_ROUND;
2001                 }
2002               else if (lex_match_id (lexer, "EMPIRICAL"))
2003                 {
2004                   examine.pc_alg = PC_EMPIRICAL;
2005                 }
2006               else if (lex_match_id (lexer, "AEMPIRICAL"))
2007                 {
2008                   examine.pc_alg = PC_AEMPIRICAL;
2009                 }
2010               else if (lex_match_id (lexer, "NONE"))
2011                 {
2012                   examine.pc_alg = PC_NONE;
2013                 }
2014               else
2015                 {
2016                   lex_error (lexer, NULL);
2017                   goto error;
2018                 }
2019             }
2020         }
2021       else if (lex_match_id (lexer, "TOTAL"))
2022         {
2023           totals_seen = true;
2024         }
2025       else if (lex_match_id (lexer, "NOTOTAL"))
2026         {
2027           nototals_seen = true;
2028         }
2029       else if (lex_match_id (lexer, "MISSING"))
2030         {
2031           lex_match (lexer, T_EQUALS);
2032
2033           while (lex_token (lexer) != T_ENDCMD
2034                  && lex_token (lexer) != T_SLASH)
2035             {
2036               if (lex_match_id (lexer, "LISTWISE"))
2037                 {
2038                   examine.missing_pw = false;
2039                 }
2040               else if (lex_match_id (lexer, "PAIRWISE"))
2041                 {
2042                   examine.missing_pw = true;
2043                 }
2044               else if (lex_match_id (lexer, "EXCLUDE"))
2045                 {
2046                   examine.exclude = MV_ANY;
2047                 }
2048               else if (lex_match_id (lexer, "INCLUDE"))
2049                 {
2050                   examine.exclude = MV_SYSTEM;
2051                 }
2052               else
2053                 {
2054                   lex_error (lexer, NULL);
2055                   goto error;
2056                 }
2057             }
2058         }
2059       else if (lex_match_id (lexer, "COMPARE"))
2060         {
2061           lex_match (lexer, T_EQUALS);
2062           if (lex_match_id (lexer, "VARIABLES"))
2063             {
2064               examine.boxplot_mode = BP_VARIABLES;
2065             }
2066           else if (lex_match_id (lexer, "GROUPS"))
2067             {
2068               examine.boxplot_mode = BP_GROUPS;
2069             }
2070           else
2071             {
2072               lex_error (lexer, NULL);
2073               goto error;
2074             }
2075         }
2076       else if (lex_match_id (lexer, "PLOT"))
2077         {
2078           lex_match (lexer, T_EQUALS);
2079
2080           while (lex_token (lexer) != T_ENDCMD
2081                  && lex_token (lexer) != T_SLASH)
2082             {
2083               if (lex_match_id (lexer, "BOXPLOT"))
2084                 {
2085                   examine.boxplot = true;
2086                 }
2087               else if (lex_match_id (lexer, "NPPLOT"))
2088                 {
2089                   examine.npplot = true;
2090                 }
2091               else if (lex_match_id (lexer, "HISTOGRAM"))
2092                 {
2093                   examine.histogram = true;
2094                 }
2095               else if (lex_match_id (lexer, "NONE"))
2096                 {
2097                   examine.histogram = false;
2098                   examine.npplot = false;
2099                   examine.boxplot = false;
2100                 }
2101               else if (lex_match (lexer, T_ALL))
2102                 {
2103                   examine.histogram = true;
2104                   examine.npplot = true;
2105                   examine.boxplot = true;
2106                 }
2107               else
2108                 {
2109                   lex_error (lexer, NULL);
2110                   goto error;
2111                 }
2112               lex_match (lexer, T_COMMA);
2113             }
2114         }
2115       else if (lex_match_id (lexer, "CINTERVAL"))
2116         {
2117           if ( !lex_force_num (lexer))
2118             goto error;
2119
2120           examine.conf = lex_number (lexer);
2121           lex_get (lexer);
2122         }
2123       else if (lex_match_id (lexer, "ID"))
2124         {
2125           lex_match (lexer, T_EQUALS);
2126
2127           examine.id_var = parse_variable_const (lexer, examine.dict);
2128         }
2129       else
2130         {
2131           lex_error (lexer, NULL);
2132           goto error;
2133         }
2134     }
2135
2136
2137   if ( totals_seen && nototals_seen)
2138     {
2139       msg (SE, _("%s and %s are mutually exclusive"),"TOTAL","NOTOTAL");
2140       goto error;
2141     }
2142
2143   /* If totals have been requested or if there are no factors
2144      in this analysis, then the totals need to be included. */
2145   if ( !nototals_seen || examine.n_iacts == 1)
2146     {
2147       examine.iacts = &iacts_mem[0];
2148     }
2149   else
2150     {
2151       examine.n_iacts--;
2152       examine.iacts = &iacts_mem[1];
2153     }
2154
2155
2156   if ( examine.id_var )
2157     {
2158       examine.id_idx = var_get_case_index (examine.id_var);
2159       examine.id_width = var_get_width (examine.id_var);
2160     }
2161
2162   examine.ex_proto = caseproto_add_width (examine.ex_proto, 0); /* value */
2163   examine.ex_proto = caseproto_add_width (examine.ex_proto, examine.id_width);   /* id */
2164   examine.ex_proto = caseproto_add_width (examine.ex_proto, 0); /* weight */
2165
2166
2167   if (examine.disp_extremes > 0)
2168     {
2169       examine.calc_extremes = examine.disp_extremes;
2170     }
2171
2172   if (examine.descriptives && examine.calc_extremes == 0)
2173     {
2174       /* Descriptives always displays the max and min */
2175       examine.calc_extremes = 1;
2176     }
2177
2178   if (percentiles_seen && examine.n_percentiles == 0)
2179     {
2180       examine.n_percentiles = 7;
2181       examine.ptiles = xcalloc (examine.n_percentiles,
2182                                 sizeof (*examine.ptiles));
2183
2184       examine.ptiles[0] = 5;
2185       examine.ptiles[1] = 10;
2186       examine.ptiles[2] = 25;
2187       examine.ptiles[3] = 50;
2188       examine.ptiles[4] = 75;
2189       examine.ptiles[5] = 90;
2190       examine.ptiles[6] = 95;
2191     }
2192
2193   assert (examine.calc_extremes >= examine.disp_extremes);
2194   {
2195     struct casegrouper *grouper;
2196     struct casereader *group;
2197     bool ok;
2198
2199     grouper = casegrouper_create_splits (proc_open (ds), examine.dict);
2200     while (casegrouper_get_next_group (grouper, &group))
2201       run_examine (&examine, group);
2202     ok = casegrouper_destroy (grouper);
2203     ok = proc_commit (ds) && ok;
2204   }
2205
2206   caseproto_unref (examine.ex_proto);
2207
2208   for (i = 0; i < examine.n_iacts; ++i)
2209     interaction_destroy (examine.iacts[i]);
2210
2211   free (examine.ptiles);
2212   free (examine.dep_vars);
2213   pool_destroy (examine.pool);
2214
2215   return CMD_SUCCESS;
2216
2217  error:
2218   caseproto_unref (examine.ex_proto);
2219   for (i = 0; i < examine.n_iacts; ++i)
2220     interaction_destroy (examine.iacts[i]);
2221   free (examine.dep_vars);
2222   free (examine.ptiles);
2223   pool_destroy (examine.pool);
2224
2225   return CMD_FAILURE;
2226 }