src/language/stats/descriptives.c

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