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