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