src/language/stats/descriptives.c

   1 /* PSPP - computes sample statistics.
   2    Copyright (C) 1997-9, 2000 Free Software Foundation, Inc.
   3    Written by Ben Pfaff <blp@gnu.org>.
   4
   5    This program is free software; you can redistribute it and/or
   6    modify it under the terms of the GNU General Public License as
   7    published by the Free Software Foundation; either version 2 of the
   8    License, or (at your option) any later version.
   9
  10    This program is distributed in the hope that it will be useful, but
  11    WITHOUT ANY WARRANTY; without even the implied warranty of
  12    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13    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, write to the Free Software
  17    Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
  18    02110-1301, USA. */
  19
  20 /* FIXME: Many possible optimizations. */
  21
  22 #include <config.h>
  23
  24 #include <limits.h>
  25 #include <math.h>
  26 #include <stdlib.h>
  27
  28 #include <data/case.h>
  29 #include <data/casefile.h>
  30 #include <data/dictionary.h>
  31 #include <data/procedure.h>
  32 #include <data/transformations.h>
  33 #include <data/variable.h>
  34 #include <language/command.h>
  35 #include <language/dictionary/split-file.h>
  36 #include <language/lexer/lexer.h>
  37 #include <language/lexer/variable-parser.h>
  38 #include <libpspp/alloc.h>
  39 #include <libpspp/array.h>
  40 #include <libpspp/compiler.h>
  41 #include <libpspp/magic.h>
  42 #include <libpspp/message.h>
  43 #include <libpspp/assertion.h>
  44 #include <math/moments.h>
  45 #include <output/manager.h>
  46 #include <output/table.h>
  47
  48 #include "gettext.h"
  49 #define _(msgid) gettext (msgid)
  50 #define N_(msgid) msgid
  51
  52 /* DESCRIPTIVES private data. */
  53
  54 struct dsc_proc;
  55
  56 /* Handling of missing values. */
  57 enum dsc_missing_type
  58   {
  59     DSC_VARIABLE,       /* Handle missing values on a per-variable basis. */
  60     DSC_LISTWISE        /* Discard entire case if any variable is missing. */
  61   };
  62
  63 /* Describes properties of a distribution for the purpose of
  64    calculating a Z-score. */
  65 struct dsc_z_score
  66   {
  67     int src_idx;                /* Source index into case data. */
  68     int dst_idx;                /* Destination index into case data. */
  69     double mean;                /* Distribution mean. */
  70     double std_dev;             /* Distribution standard deviation. */
  71     struct variable *v;         /* Variable on which z-score is based. */
  72   };
  73
  74 /* DESCRIPTIVES transformation (for calculating Z-scores). */
  75 struct dsc_trns
  76   {
  77     struct dsc_z_score *z_scores; /* Array of Z-scores. */
  78     int z_score_cnt;            /* Number of Z-scores. */
  79     struct variable **vars;     /* Variables for listwise missing checks. */
  80     size_t var_cnt;             /* Number of variables. */
  81     enum dsc_missing_type missing_type; /* Treatment of missing values. */
  82     int include_user_missing;   /* Nonzero to include user-missing values. */
  83   };
  84
  85 /* Statistics.  Used as bit indexes, so must be 32 or fewer. */
  86 enum dsc_statistic
  87   {
  88     DSC_MEAN = 0, DSC_SEMEAN, DSC_STDDEV, DSC_VARIANCE, DSC_KURTOSIS,
  89     DSC_SEKURT, DSC_SKEWNESS, DSC_SESKEW, DSC_RANGE, DSC_MIN,
  90     DSC_MAX, DSC_SUM, DSC_N_STATS,
  91
  92     /* Only valid as sort criteria. */
  93     DSC_NAME = -2,              /* Sort by name. */
  94     DSC_NONE = -1               /* Unsorted. */
  95   };
  96
  97 /* Describes one statistic. */
  98 struct dsc_statistic_info
  99   {
 100     const char *identifier;     /* Identifier. */
 101     const char *name;           /* Full name. */
 102     enum moment moment;         /* Highest moment needed to calculate. */
 103   };
 104
 105 /* Table of statistics, indexed by DSC_*. */
 106 static const struct dsc_statistic_info dsc_info[DSC_N_STATS] =
 107   {
 108     {"MEAN", N_("Mean"), MOMENT_MEAN},
 109     {"SEMEAN", N_("S E Mean"), MOMENT_VARIANCE},
 110     {"STDDEV", N_("Std Dev"), MOMENT_VARIANCE},
 111     {"VARIANCE", N_("Variance"), MOMENT_VARIANCE},
 112     {"KURTOSIS", N_("Kurtosis"), MOMENT_KURTOSIS},
 113     {"SEKURTOSIS", N_("S E Kurt"), MOMENT_NONE},
 114     {"SKEWNESS", N_("Skewness"), MOMENT_SKEWNESS},
 115     {"SESKEWNESS", N_("S E Skew"), MOMENT_NONE},
 116     {"RANGE", N_("Range"), MOMENT_NONE},
 117     {"MINIMUM", N_("Minimum"), MOMENT_NONE},
 118     {"MAXIMUM", N_("Maximum"), MOMENT_NONE},
 119     {"SUM", N_("Sum"), MOMENT_MEAN},
 120   };
 121
 122 /* Statistics calculated by default if none are explicitly
 123    requested. */
 124 #define DEFAULT_STATS                                                   \
 125         ((1ul << DSC_MEAN) | (1ul << DSC_STDDEV) | (1ul << DSC_MIN)     \
 126          | (1ul << DSC_MAX))
 127
 128 /* A variable specified on DESCRIPTIVES. */
 129 struct dsc_var
 130   {
 131     struct variable *v;         /* Variable to calculate on. */
 132     char z_name[LONG_NAME_LEN + 1]; /* Name for z-score variable. */
 133     double valid, missing;      /* Valid, missing counts. */
 134     struct moments *moments;    /* Moments. */
 135     double min, max;            /* Maximum and mimimum values. */
 136     double stats[DSC_N_STATS];  /* All the stats' values. */
 137   };
 138
 139 /* Output format. */
 140 enum dsc_format
 141   {
 142     DSC_LINE,           /* Abbreviated format. */
 143     DSC_SERIAL          /* Long format. */
 144   };
 145
 146 /* A DESCRIPTIVES procedure. */
 147 struct dsc_proc
 148   {
 149     /* Per-variable info. */
 150     struct dsc_var *vars;       /* Variables. */
 151     size_t var_cnt;             /* Number of variables. */
 152
 153     /* User options. */
 154     enum dsc_missing_type missing_type; /* Treatment of missing values. */
 155     int include_user_missing;   /* Nonzero to include user-missing values. */
 156     int show_var_labels;        /* Nonzero to show variable labels. */
 157     int show_index;             /* Nonzero to show variable index. */
 158     enum dsc_format format;     /* Output format. */
 159
 160     /* Accumulated results. */
 161     double missing_listwise;    /* Sum of weights of cases missing listwise. */
 162     double valid;               /* Sum of weights of valid cases. */
 163     bool bad_warn;               /* Warn if bad weight found. */
 164     enum dsc_statistic sort_by_stat; /* Statistic to sort by; -1: name. */
 165     int sort_ascending;         /* !0: ascending order; 0: descending. */
 166     unsigned long show_stats;   /* Statistics to display. */
 167     unsigned long calc_stats;   /* Statistics to calculate. */
 168     enum moment max_moment;     /* Highest moment needed for stats. */
 169   };
 170
 171 /* Parsing. */
 172 static enum dsc_statistic match_statistic (void);
 173 static void free_dsc_proc (struct dsc_proc *);
 174
 175 /* Z-score functions. */
 176 static bool try_name (struct dsc_proc *dsc, char *name);
 177 static bool generate_z_varname (struct dsc_proc *dsc, char *z_name,
 178                                const char *name, size_t *z_cnt);
 179 static void dump_z_table (struct dsc_proc *);
 180 static void setup_z_trns (struct dsc_proc *);
 181
 182 /* Procedure execution functions. */
 183 static bool calc_descriptives (const struct ccase *first,
 184                                const struct casefile *, void *dsc_);
 185 static void display (struct dsc_proc *dsc);
 186 \f
 187 /* Parser and outline. */
 188
 189 /* Handles DESCRIPTIVES. */
 190 int
 191 cmd_descriptives (void)
 192 {
 193   struct dsc_proc *dsc;
 194   struct variable **vars = NULL;
 195   size_t var_cnt = 0;
 196   int save_z_scores = 0;
 197   size_t z_cnt = 0;
 198   size_t i;
 199   bool ok;
 200
 201   /* Create and initialize dsc. */
 202   dsc = xmalloc (sizeof *dsc);
 203   dsc->vars = NULL;
 204   dsc->var_cnt = 0;
 205   dsc->missing_type = DSC_VARIABLE;
 206   dsc->include_user_missing = 0;
 207   dsc->show_var_labels = 1;
 208   dsc->show_index = 0;
 209   dsc->format = DSC_LINE;
 210   dsc->missing_listwise = 0.;
 211   dsc->valid = 0.;
 212   dsc->bad_warn = 1;
 213   dsc->sort_by_stat = DSC_NONE;
 214   dsc->sort_ascending = 1;
 215   dsc->show_stats = dsc->calc_stats = DEFAULT_STATS;
 216
 217   /* Parse DESCRIPTIVES. */
 218   while (token != '.')
 219     {
 220       if (lex_match_id ("MISSING"))
 221         {
 222           lex_match ('=');
 223           while (token != '.' && token != '/')
 224             {
 225               if (lex_match_id ("VARIABLE"))
 226                 dsc->missing_type = DSC_VARIABLE;
 227               else if (lex_match_id ("LISTWISE"))
 228                 dsc->missing_type = DSC_LISTWISE;
 229               else if (lex_match_id ("INCLUDE"))
 230                 dsc->include_user_missing = 1;
 231               else
 232                 {
 233                   lex_error (NULL);
 234                   goto error;
 235                 }
 236               lex_match (',');
 237             }
 238         }
 239       else if (lex_match_id ("SAVE"))
 240         save_z_scores = 1;
 241       else if (lex_match_id ("FORMAT"))
 242         {
 243           lex_match ('=');
 244           while (token != '.' && token != '/')
 245             {
 246               if (lex_match_id ("LABELS"))
 247                 dsc->show_var_labels = 1;
 248               else if (lex_match_id ("NOLABELS"))
 249                 dsc->show_var_labels = 0;
 250               else if (lex_match_id ("INDEX"))
 251                 dsc->show_index = 1;
 252               else if (lex_match_id ("NOINDEX"))
 253                 dsc->show_index = 0;
 254               else if (lex_match_id ("LINE"))
 255                 dsc->format = DSC_LINE;
 256               else if (lex_match_id ("SERIAL"))
 257                 dsc->format = DSC_SERIAL;
 258               else
 259                 {
 260                   lex_error (NULL);
 261                   goto error;
 262                 }
 263               lex_match (',');
 264             }
 265         }
 266       else if (lex_match_id ("STATISTICS"))
 267         {
 268           lex_match ('=');
 269           dsc->show_stats = 0;
 270           while (token != '.' && token != '/')
 271             {
 272               if (lex_match (T_ALL))
 273                 dsc->show_stats |= (1ul << DSC_N_STATS) - 1;
 274               else if (lex_match_id ("DEFAULT"))
 275                 dsc->show_stats |= DEFAULT_STATS;
 276               else
 277                 dsc->show_stats |= 1ul << (match_statistic ());
 278               lex_match (',');
 279             }
 280           if (dsc->show_stats == 0)
 281             dsc->show_stats = DEFAULT_STATS;
 282         }
 283       else if (lex_match_id ("SORT"))
 284         {
 285           lex_match ('=');
 286           if (lex_match_id ("NAME"))
 287             dsc->sort_by_stat = DSC_NAME;
 288           else
 289             {
 290               dsc->sort_by_stat = match_statistic ();
 291               if (dsc->sort_by_stat == DSC_NONE )
 292                 dsc->sort_by_stat = DSC_MEAN;
 293             }
 294           if (lex_match ('('))
 295             {
 296               if (lex_match_id ("A"))
 297                 dsc->sort_ascending = 1;
 298               else if (lex_match_id ("D"))
 299                 dsc->sort_ascending = 0;
 300               else
 301                 lex_error (NULL);
 302               lex_force_match (')');
 303             }
 304         }
 305       else if (var_cnt == 0)
 306         {
 307           if (lex_look_ahead () == '=')
 308             {
 309               lex_match_id ("VARIABLES");
 310               lex_match ('=');
 311             }
 312
 313           while (token != '.' && token != '/')
 314             {
 315               int i;
 316
 317               if (!parse_variables (dataset_dict (current_dataset), &vars, &var_cnt,
 318                                     PV_APPEND | PV_NO_DUPLICATE | PV_NUMERIC))
 319                 goto error;
 320
 321               dsc->vars = xnrealloc (dsc->vars, var_cnt, sizeof *dsc->vars);
 322               for (i = dsc->var_cnt; i < var_cnt; i++)
 323                 {
 324                   struct dsc_var *dv = &dsc->vars[i];
 325                   dv->v = vars[i];
 326                   dv->z_name[0] = '\0';
 327                   dv->moments = NULL;
 328                 }
 329               dsc->var_cnt = var_cnt;
 330
 331               if (lex_match ('('))
 332                 {
 333                   if (token != T_ID)
 334                     {
 335                       lex_error (NULL);
 336                       goto error;
 337                     }
 338                   if (try_name (dsc, tokid))
 339                     {
 340                       strcpy (dsc->vars[dsc->var_cnt - 1].z_name, tokid);
 341                       z_cnt++;
 342                     }
 343                   else
 344                     msg (SE, _("Z-score variable name %s would be"
 345                                " a duplicate variable name."), tokid);
 346                   lex_get ();
 347                   if (!lex_force_match (')'))
 348                     goto error;
 349                 }
 350             }
 351         }
 352       else
 353         {
 354           lex_error (NULL);
 355           goto error;
 356         }
 357
 358       lex_match ('/');
 359     }
 360   if (var_cnt == 0)
 361     {
 362       msg (SE, _("No variables specified."));
 363       goto error;
 364     }
 365
 366   /* Construct z-score varnames, show translation table. */
 367   if (z_cnt || save_z_scores)
 368     {
 369       if (save_z_scores)
 370         {
 371           size_t gen_cnt = 0;
 372
 373           for (i = 0; i < dsc->var_cnt; i++)
 374             if (dsc->vars[i].z_name[0] == 0)
 375               {
 376                 if (!generate_z_varname (dsc, dsc->vars[i].z_name,
 377                                          dsc->vars[i].v->name, &gen_cnt))
 378                   goto error;
 379                 z_cnt++;
 380               }
 381         }
 382       dump_z_table (dsc);
 383     }
 384
 385   /* Figure out statistics to display. */
 386   if (dsc->show_stats & (1ul << DSC_SKEWNESS))
 387     dsc->show_stats |= 1ul << DSC_SESKEW;
 388   if (dsc->show_stats & (1ul << DSC_KURTOSIS))
 389     dsc->show_stats |= 1ul << DSC_SEKURT;
 390
 391   /* Figure out which statistics to calculate. */
 392   dsc->calc_stats = dsc->show_stats;
 393   if (z_cnt > 0)
 394     dsc->calc_stats |= (1ul << DSC_MEAN) | (1ul << DSC_STDDEV);
 395   if (dsc->sort_by_stat >= 0)
 396     dsc->calc_stats |= 1ul << dsc->sort_by_stat;
 397   if (dsc->show_stats & (1ul << DSC_SESKEW))
 398     dsc->calc_stats |= 1ul << DSC_SKEWNESS;
 399   if (dsc->show_stats & (1ul << DSC_SEKURT))
 400     dsc->calc_stats |= 1ul << DSC_KURTOSIS;
 401
 402   /* Figure out maximum moment needed and allocate moments for
 403      the variables. */
 404   dsc->max_moment = MOMENT_NONE;
 405   for (i = 0; i < DSC_N_STATS; i++)
 406     if (dsc->calc_stats & (1ul << i) && dsc_info[i].moment > dsc->max_moment)
 407       dsc->max_moment = dsc_info[i].moment;
 408   if (dsc->max_moment != MOMENT_NONE)
 409     for (i = 0; i < dsc->var_cnt; i++)
 410       dsc->vars[i].moments = moments_create (dsc->max_moment);
 411
 412   /* Data pass. */
 413   ok = multipass_procedure_with_splits (current_dataset, calc_descriptives, dsc);
 414
 415   /* Z-scoring! */
 416   if (ok && z_cnt)
 417     setup_z_trns (dsc);
 418
 419   /* Done. */
 420   free (vars);
 421   free_dsc_proc (dsc);
 422   return ok ? CMD_SUCCESS : CMD_CASCADING_FAILURE;
 423
 424  error:
 425   free (vars);
 426   free_dsc_proc (dsc);
 427   return CMD_FAILURE;
 428 }
 429
 430 /* Returns the statistic named by the current token and skips past the token.
 431    Returns DSC_NONE if no statistic is given (e.g., subcommand with no
 432    specifiers). Emits an error if the current token ID does not name a
 433    statistic. */
 434 static enum dsc_statistic
 435 match_statistic (void)
 436 {
 437   if (token == T_ID)
 438     {
 439       enum dsc_statistic stat;
 440
 441       for (stat = 0; stat < DSC_N_STATS; stat++)
 442         if (lex_match_id (dsc_info[stat].identifier))
 443           return stat;
 444
 445       lex_get();
 446       lex_error (_("expecting statistic name: reverting to default"));
 447     }
 448
 449   return DSC_NONE;
 450 }
 451
 452 /* Frees DSC. */
 453 static void
 454 free_dsc_proc (struct dsc_proc *dsc)
 455 {
 456   size_t i;
 457
 458   if (dsc == NULL)
 459     return;
 460
 461   for (i = 0; i < dsc->var_cnt; i++)
 462     moments_destroy (dsc->vars[i].moments);
 463   free (dsc->vars);
 464   free (dsc);
 465 }
 466 \f
 467 /* Z scores. */
 468
 469 /* Returns false if NAME is a duplicate of any existing variable name or
 470    of any previously-declared z-var name; otherwise returns true. */
 471 static bool
 472 try_name (struct dsc_proc *dsc, char *name)
 473 {
 474   size_t i;
 475
 476   if (dict_lookup_var (dataset_dict (current_dataset), name) != NULL)
 477     return false;
 478   for (i = 0; i < dsc->var_cnt; i++)
 479     if (!strcasecmp (dsc->vars[i].z_name, name))
 480       return false;
 481   return true;
 482 }
 483
 484 /* Generates a name for a Z-score variable based on a variable
 485    named VAR_NAME, given that *Z_CNT generated variable names are
 486    known to already exist.  If successful, returns true and
 487    copies the new name into Z_NAME.  On failure, returns false. */
 488 static bool
 489 generate_z_varname (struct dsc_proc *dsc, char *z_name,
 490                     const char *var_name, size_t *z_cnt)
 491 {
 492   char name[LONG_NAME_LEN + 1];
 493
 494   /* Try a name based on the original variable name. */
 495   name[0] = 'Z';
 496   str_copy_trunc (name + 1, sizeof name - 1, var_name);
 497   if (try_name (dsc, name))
 498     {
 499       strcpy (z_name, name);
 500       return true;
 501     }
 502
 503   /* Generate a synthetic name. */
 504   for (;;)
 505     {
 506       (*z_cnt)++;
 507
 508       if (*z_cnt <= 99)
 509         sprintf (name, "ZSC%03d", *z_cnt);
 510       else if (*z_cnt <= 108)
 511         sprintf (name, "STDZ%02d", *z_cnt - 99);
 512       else if (*z_cnt <= 117)
 513         sprintf (name, "ZZZZ%02d", *z_cnt - 108);
 514       else if (*z_cnt <= 126)
 515         sprintf (name, "ZQZQ%02d", *z_cnt - 117);
 516       else
 517         {
 518           msg (SE, _("Ran out of generic names for Z-score variables.  "
 519                      "There are only 126 generic names: ZSC001-ZSC0999, "
 520                      "STDZ01-STDZ09, ZZZZ01-ZZZZ09, ZQZQ01-ZQZQ09."));
 521           return false;
 522         }
 523
 524       if (try_name (dsc, name))
 525         {
 526           strcpy (z_name, name);
 527           return true;
 528         }
 529     }
 530   NOT_REACHED();
 531 }
 532
 533 /* Outputs a table describing the mapping between source
 534    variables and Z-score variables. */
 535 static void
 536 dump_z_table (struct dsc_proc *dsc)
 537 {
 538   size_t cnt = 0;
 539   struct tab_table *t;
 540
 541   {
 542     size_t i;
 543
 544     for (i = 0; i < dsc->var_cnt; i++)
 545       if (dsc->vars[i].z_name[0] != '\0')
 546         cnt++;
 547   }
 548
 549   t = tab_create (2, cnt + 1, 0);
 550   tab_title (t, _("Mapping of variables to corresponding Z-scores."));
 551   tab_columns (t, SOM_COL_DOWN, 1);
 552   tab_headers (t, 0, 0, 1, 0);
 553   tab_box (t, TAL_1, TAL_1, TAL_0, TAL_1, 0, 0, 1, cnt);
 554   tab_hline (t, TAL_2, 0, 1, 1);
 555   tab_text (t, 0, 0, TAB_CENTER | TAT_TITLE, _("Source"));
 556   tab_text (t, 1, 0, TAB_CENTER | TAT_TITLE, _("Target"));
 557   tab_dim (t, tab_natural_dimensions);
 558
 559   {
 560     size_t i, y;
 561
 562     for (i = 0, y = 1; i < dsc->var_cnt; i++)
 563       if (dsc->vars[i].z_name[0] != '\0')
 564         {
 565           tab_text (t, 0, y, TAB_LEFT, dsc->vars[i].v->name);
 566           tab_text (t, 1, y++, TAB_LEFT, dsc->vars[i].z_name);
 567         }
 568   }
 569
 570   tab_submit (t);
 571 }
 572
 573 /* Transformation function to calculate Z-scores. Will return SYSMIS if any of
 574    the following are true: 1) mean or standard deviation is SYSMIS 2) score is
 575    SYSMIS 3) score is user missing and they were not included in the original
 576    analyis. 4) any of the variables in the original analysis were missing
 577    (either system or user-missing values that weren't included).
 578 */
 579 static int
 580 descriptives_trns_proc (void *trns_, struct ccase * c,
 581                         casenum_t case_idx UNUSED)
 582 {
 583   struct dsc_trns *t = trns_;
 584   struct dsc_z_score *z;
 585   struct variable **vars;
 586   int all_sysmis = 0;
 587
 588   if (t->missing_type == DSC_LISTWISE)
 589     {
 590       assert(t->vars);
 591       for (vars = t->vars; vars < t->vars + t->var_cnt; vars++)
 592         {
 593           double score = case_num (c, (*vars)->fv);
 594           if ( score == SYSMIS
 595                || (!t->include_user_missing
 596                    && mv_is_num_user_missing (&(*vars)->miss, score)))
 597             {
 598               all_sysmis = 1;
 599               break;
 600             }
 601         }
 602     }
 603
 604   for (z = t->z_scores; z < t->z_scores + t->z_score_cnt; z++)
 605     {
 606       double input = case_num (c, z->src_idx);
 607       double *output = &case_data_rw (c, z->dst_idx)->f;
 608
 609       if (z->mean == SYSMIS || z->std_dev == SYSMIS
 610           || all_sysmis || input == SYSMIS
 611           || (!t->include_user_missing
 612               && mv_is_num_user_missing (&z->v->miss, input)))
 613         *output = SYSMIS;
 614       else
 615         *output = (input - z->mean) / z->std_dev;
 616     }
 617   return TRNS_CONTINUE;
 618 }
 619
 620 /* Frees a descriptives_trns struct. */
 621 static bool
 622 descriptives_trns_free (void *trns_)
 623 {
 624   struct dsc_trns *t = trns_;
 625
 626   free (t->z_scores);
 627   assert((t->missing_type != DSC_LISTWISE) ^ (t->vars != NULL));
 628   free (t->vars);
 629   return true;
 630 }
 631
 632 /* Sets up a transformation to calculate Z scores. */
 633 static void
 634 setup_z_trns (struct dsc_proc *dsc)
 635 {
 636   struct dsc_trns *t;
 637   size_t cnt, i;
 638
 639   for (cnt = i = 0; i < dsc->var_cnt; i++)
 640     if (dsc->vars[i].z_name[0] != '\0')
 641       cnt++;
 642
 643   t = xmalloc (sizeof *t);
 644   t->z_scores = xnmalloc (cnt, sizeof *t->z_scores);
 645   t->z_score_cnt = cnt;
 646   t->missing_type = dsc->missing_type;
 647   t->include_user_missing = dsc->include_user_missing;
 648   if ( t->missing_type == DSC_LISTWISE )
 649     {
 650       t->var_cnt = dsc->var_cnt;
 651       t->vars = xnmalloc (t->var_cnt, sizeof *t->vars);
 652       for (i = 0; i < t->var_cnt; i++)
 653         t->vars[i] = dsc->vars[i].v;
 654     }
 655   else
 656     {
 657       t->var_cnt = 0;
 658       t->vars = NULL;
 659     }
 660
 661   for (cnt = i = 0; i < dsc->var_cnt; i++)
 662     {
 663       struct dsc_var *dv = &dsc->vars[i];
 664       if (dv->z_name[0] != '\0')
 665         {
 666           struct dsc_z_score *z;
 667           char *cp;
 668           struct variable *dst_var;
 669
 670           dst_var = dict_create_var_assert (dataset_dict (current_dataset), dv->z_name, 0);
 671           if (dv->v->label)
 672             {
 673               dst_var->label = xmalloc (strlen (dv->v->label) + 12);
 674               cp = stpcpy (dst_var->label, _("Z-score of "));
 675               strcpy (cp, dv->v->label);
 676             }
 677           else
 678             {
 679               dst_var->label = xmalloc (strlen (dv->v->name) + 12);
 680               cp = stpcpy (dst_var->label, _("Z-score of "));
 681               strcpy (cp, dv->v->name);
 682             }
 683
 684           z = &t->z_scores[cnt++];
 685           z->src_idx = dv->v->fv;
 686           z->dst_idx = dst_var->fv;
 687           z->mean = dv->stats[DSC_MEAN];
 688           z->std_dev = dv->stats[DSC_STDDEV];
 689           z->v = dv->v;
 690         }
 691     }
 692
 693   add_transformation (current_dataset,
 694                       descriptives_trns_proc, descriptives_trns_free, t);
 695 }
 696 \f
 697 /* Statistical calculation. */
 698
 699 static bool listwise_missing (struct dsc_proc *dsc, const struct ccase *c);
 700
 701 /* Calculates and displays descriptive statistics for the cases
 702    in CF. */
 703 static bool
 704 calc_descriptives (const struct ccase *first,
 705                    const struct casefile *cf, void *dsc_)
 706 {
 707   struct dsc_proc *dsc = dsc_;
 708   struct casereader *reader;
 709   struct ccase c;
 710   size_t i;
 711
 712   output_split_file_values (first);
 713
 714   for (i = 0; i < dsc->var_cnt; i++)
 715     {
 716       struct dsc_var *dv = &dsc->vars[i];
 717
 718       dv->valid = dv->missing = 0.0;
 719       if (dv->moments != NULL)
 720         moments_clear (dv->moments);
 721       dv->min = DBL_MAX;
 722       dv->max = -DBL_MAX;
 723     }
 724   dsc->missing_listwise = 0.;
 725   dsc->valid = 0.;
 726
 727   /* First pass to handle most of the work. */
 728   for (reader = casefile_get_reader (cf);
 729        casereader_read (reader, &c);
 730        case_destroy (&c))
 731     {
 732       double weight = dict_get_case_weight (dataset_dict (current_dataset), &c, &dsc->bad_warn);
 733       if (weight <= 0.0)
 734         continue;
 735
 736       /* Check for missing values. */
 737       if (listwise_missing (dsc, &c))
 738         {
 739           dsc->missing_listwise += weight;
 740           if (dsc->missing_type == DSC_LISTWISE)
 741             continue;
 742         }
 743       dsc->valid += weight;
 744
 745       for (i = 0; i < dsc->var_cnt; i++)
 746         {
 747           struct dsc_var *dv = &dsc->vars[i];
 748           double x = case_num (&c, dv->v->fv);
 749
 750           if (dsc->missing_type != DSC_LISTWISE
 751               && (x == SYSMIS
 752                   || (!dsc->include_user_missing
 753                       && mv_is_num_user_missing (&dv->v->miss, x))))
 754             {
 755               dv->missing += weight;
 756               continue;
 757             }
 758
 759           if (dv->moments != NULL)
 760             moments_pass_one (dv->moments, x, weight);
 761
 762           if (x < dv->min)
 763             dv->min = x;
 764           if (x > dv->max)
 765             dv->max = x;
 766         }
 767     }
 768   casereader_destroy (reader);
 769
 770   /* Second pass for higher-order moments. */
 771   if (dsc->max_moment > MOMENT_MEAN)
 772     {
 773       for (reader = casefile_get_reader (cf);
 774            casereader_read (reader, &c);
 775            case_destroy (&c))
 776         {
 777           double weight = dict_get_case_weight (dataset_dict (current_dataset), &c,
 778                                                 &dsc->bad_warn);
 779           if (weight <= 0.0)
 780             continue;
 781
 782           /* Check for missing values. */
 783           if (listwise_missing (dsc, &c)
 784               && dsc->missing_type == DSC_LISTWISE)
 785             continue;
 786
 787           for (i = 0; i < dsc->var_cnt; i++)
 788             {
 789               struct dsc_var *dv = &dsc->vars[i];
 790               double x = case_num (&c, dv->v->fv);
 791
 792               if (dsc->missing_type != DSC_LISTWISE
 793                   && (x == SYSMIS
 794                       || (!dsc->include_user_missing
 795                           && mv_is_num_user_missing (&dv->v->miss, x))))
 796                 continue;
 797
 798               if (dv->moments != NULL)
 799                 moments_pass_two (dv->moments, x, weight);
 800             }
 801         }
 802       casereader_destroy (reader);
 803     }
 804
 805   /* Calculate results. */
 806   for (i = 0; i < dsc->var_cnt; i++)
 807     {
 808       struct dsc_var *dv = &dsc->vars[i];
 809       double W;
 810       int j;
 811
 812       for (j = 0; j < DSC_N_STATS; j++)
 813         dv->stats[j] = SYSMIS;
 814
 815       dv->valid = W = dsc->valid - dv->missing;
 816
 817       if (dv->moments != NULL)
 818         moments_calculate (dv->moments, NULL,
 819                            &dv->stats[DSC_MEAN], &dv->stats[DSC_VARIANCE],
 820                            &dv->stats[DSC_SKEWNESS], &dv->stats[DSC_KURTOSIS]);
 821       if (dsc->calc_stats & (1ul << DSC_SEMEAN)
 822           && dv->stats[DSC_VARIANCE] != SYSMIS && W > 0.)
 823         dv->stats[DSC_SEMEAN] = sqrt (dv->stats[DSC_VARIANCE]) / sqrt (W);
 824       if (dsc->calc_stats & (1ul << DSC_STDDEV)
 825           && dv->stats[DSC_VARIANCE] != SYSMIS)
 826         dv->stats[DSC_STDDEV] = sqrt (dv->stats[DSC_VARIANCE]);
 827       if (dsc->calc_stats & (1ul << DSC_SEKURT))
 828         if (dv->stats[DSC_KURTOSIS] != SYSMIS)
 829             dv->stats[DSC_SEKURT] = calc_sekurt (W);
 830       if (dsc->calc_stats & (1ul << DSC_SESKEW)
 831           && dv->stats[DSC_SKEWNESS] != SYSMIS)
 832         dv->stats[DSC_SESKEW] = calc_seskew (W);
 833       dv->stats[DSC_RANGE] = ((dv->min == DBL_MAX || dv->max == -DBL_MAX)
 834                               ? SYSMIS : dv->max - dv->min);
 835       dv->stats[DSC_MIN] = dv->min == DBL_MAX ? SYSMIS : dv->min;
 836       dv->stats[DSC_MAX] = dv->max == -DBL_MAX ? SYSMIS : dv->max;
 837       if (dsc->calc_stats & (1ul << DSC_SUM))
 838         dv->stats[DSC_SUM] = W * dv->stats[DSC_MEAN];
 839     }
 840
 841   /* Output results. */
 842   display (dsc);
 843
 844   return true;
 845 }
 846
 847 /* Returns true if any of the descriptives variables in DSC's
 848    variable list have missing values in case C, false otherwise. */
 849 static bool
 850 listwise_missing (struct dsc_proc *dsc, const struct ccase *c)
 851 {
 852   size_t i;
 853
 854   for (i = 0; i < dsc->var_cnt; i++)
 855     {
 856       struct dsc_var *dv = &dsc->vars[i];
 857       double x = case_num (c, dv->v->fv);
 858
 859       if (x == SYSMIS
 860           || (!dsc->include_user_missing
 861               && mv_is_num_user_missing (&dv->v->miss, x)))
 862         return true;
 863     }
 864   return false;
 865 }
 866 \f
 867 /* Statistical display. */
 868
 869 static algo_compare_func descriptives_compare_dsc_vars;
 870
 871 /* Displays a table of descriptive statistics for DSC. */
 872 static void
 873 display (struct dsc_proc *dsc)
 874 {
 875   size_t i;
 876   int nc;
 877   struct tab_table *t;
 878
 879   nc = 1 + (dsc->format == DSC_SERIAL ? 2 : 1);
 880   for (i = 0; i < DSC_N_STATS; i++)
 881     if (dsc->show_stats & (1ul << i))
 882       nc++;
 883
 884   if (dsc->sort_by_stat != DSC_NONE)
 885     sort (dsc->vars, dsc->var_cnt, sizeof *dsc->vars,
 886           descriptives_compare_dsc_vars, dsc);
 887
 888   t = tab_create (nc, dsc->var_cnt + 1, 0);
 889   tab_headers (t, 1, 0, 1, 0);
 890   tab_box (t, TAL_1, TAL_1, -1, -1, 0, 0, nc - 1, dsc->var_cnt);
 891   tab_box (t, -1, -1, -1, TAL_1, 1, 0, nc - 1, dsc->var_cnt);
 892   tab_hline (t, TAL_2, 0, nc - 1, 1);
 893   tab_vline (t, TAL_2, 1, 0, dsc->var_cnt);
 894   tab_dim (t, tab_natural_dimensions);
 895
 896   nc = 0;
 897   tab_text (t, nc++, 0, TAB_LEFT | TAT_TITLE, _("Variable"));
 898   if (dsc->format == DSC_SERIAL)
 899     {
 900       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Valid N"));
 901       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Missing N"));
 902     }
 903   else
 904     tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, "N");
 905
 906   for (i = 0; i < DSC_N_STATS; i++)
 907     if (dsc->show_stats & (1ul << i))
 908       {
 909         const char *title = gettext (dsc_info[i].name);
 910         tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, title);
 911       }
 912
 913   for (i = 0; i < dsc->var_cnt; i++)
 914     {
 915       struct dsc_var *dv = &dsc->vars[i];
 916       size_t j;
 917
 918       nc = 0;
 919       tab_text (t, nc++, i + 1, TAB_LEFT, dv->v->name);
 920       tab_text (t, nc++, i + 1, TAT_PRINTF, "%g", dv->valid);
 921       if (dsc->format == DSC_SERIAL)
 922         tab_text (t, nc++, i + 1, TAT_PRINTF, "%g", dv->missing);
 923       for (j = 0; j < DSC_N_STATS; j++)
 924         if (dsc->show_stats & (1ul << j))
 925           tab_float (t, nc++, i + 1, TAB_NONE, dv->stats[j], 10, 3);
 926     }
 927
 928   tab_title (t, _("Valid cases = %g; cases with missing value(s) = %g."),
 929              dsc->valid, dsc->missing_listwise);
 930
 931   tab_submit (t);
 932 }
 933
 934 /* Compares `struct dsc_var's A and B according to the ordering
 935    specified by CMD. */
 936 static int
 937 descriptives_compare_dsc_vars (const void *a_, const void *b_, void *dsc_)
 938 {
 939   const struct dsc_var *a = a_;
 940   const struct dsc_var *b = b_;
 941   struct dsc_proc *dsc = dsc_;
 942
 943   int result;
 944
 945   if (dsc->sort_by_stat == DSC_NAME)
 946     result = strcasecmp (a->v->name, b->v->name);
 947   else
 948     {
 949       double as = a->stats[dsc->sort_by_stat];
 950       double bs = b->stats[dsc->sort_by_stat];
 951
 952       result = as < bs ? -1 : as > bs;
 953     }
 954
 955   if (!dsc->sort_ascending)
 956     result = -result;
 957
 958   return result;
 959 }