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