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., 59 Temple Place - Suite 330, Boston, MA
  18    02111-1307, 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 "lexer.h"
  33 #include "error.h"
  34 #include "magic.h"
  35 #include "moments.h"
  36 #include "som.h"
  37 #include "tab.h"
  38 #include "var.h"
  39 #include "vfm.h"
  40
  41 /* DESCRIPTIVES private data. */
  42
  43 struct dsc_proc;
  44
  45 /* Handling of missing values. */
  46 enum dsc_missing_type
  47   {
  48     DSC_VARIABLE,       /* Handle missing values on a per-variable basis. */
  49     DSC_LISTWISE        /* Discard entire case if any variable is missing. */
  50   };
  51
  52 /* Describes properties of a distribution for the purpose of
  53    calculating a Z-score. */
  54 struct dsc_z_score
  55   {
  56     int src_idx;                /* Source index into case data. */
  57     int dst_idx;                /* Destination index into case data. */
  58     double mean;                /* Distribution mean. */
  59     double std_dev;             /* Distribution standard deviation. */
  60     struct variable *v;         /* Variable on which z-score is based. */
  61   };
  62
  63 /* DESCRIPTIVES transformation (for calculating Z-scores). */
  64 struct dsc_trns
  65   {
  66     struct trns_header h;
  67     struct dsc_z_score *z_scores; /* Array of Z-scores. */
  68     int z_score_cnt;            /* Number of Z-scores. */
  69     struct variable **vars;     /* Variables for listwise missing checks. */
  70     size_t var_cnt;             /* Number of variables. */
  71     enum dsc_missing_type missing_type; /* Treatment of missing values. */
  72     int include_user_missing;   /* Nonzero to include user-missing values. */
  73   };
  74
  75 /* Statistics.  Used as bit indexes, so must be 32 or fewer. */
  76 enum dsc_statistic
  77   {
  78     DSC_MEAN = 0, DSC_SEMEAN, DSC_STDDEV, DSC_VARIANCE, DSC_KURTOSIS,
  79     DSC_SEKURT, DSC_SKEWNESS, DSC_SESKEW, DSC_RANGE, DSC_MIN,
  80     DSC_MAX, DSC_SUM, DSC_N_STATS,
  81
  82     /* Only valid as sort criteria. */
  83     DSC_NAME = -2,              /* Sort by name. */
  84     DSC_NONE = -1               /* Unsorted. */
  85   };
  86
  87 /* Describes one statistic. */
  88 struct dsc_statistic_info
  89   {
  90     const char *identifier;     /* Identifier. */
  91     const char *name;           /* Full name. */
  92     enum moment moment;         /* Highest moment needed to calculate. */
  93   };
  94
  95 /* Table of statistics, indexed by DSC_*. */
  96 static const struct dsc_statistic_info dsc_info[DSC_N_STATS] =
  97   {
  98     {"MEAN", N_("Mean"), MOMENT_MEAN},
  99     {"SEMEAN", N_("S E Mean"), MOMENT_VARIANCE},
 100     {"STDDEV", N_("Std Dev"), MOMENT_VARIANCE},
 101     {"VARIANCE", N_("Variance"), MOMENT_VARIANCE},
 102     {"KURTOSIS", N_("Kurtosis"), MOMENT_KURTOSIS},
 103     {"SEKURTOSIS", N_("S E Kurt"), MOMENT_NONE},
 104     {"SKEWNESS", N_("Skewness"), MOMENT_SKEWNESS},
 105     {"SESKEWNESS", N_("S E Skew"), MOMENT_NONE},
 106     {"RANGE", N_("Range"), MOMENT_NONE},
 107     {"MINIMUM", N_("Minimum"), MOMENT_NONE},
 108     {"MAXIMUM", N_("Maximum"), MOMENT_NONE},
 109     {"SUM", N_("Sum"), MOMENT_MEAN},
 110   };
 111
 112 /* Statistics calculated by default if none are explicitly
 113    requested. */
 114 #define DEFAULT_STATS                                                   \
 115         ((1ul << DSC_MEAN) | (1ul << DSC_STDDEV) | (1ul << DSC_MIN)     \
 116          | (1ul << DSC_MAX))
 117
 118 /* A variable specified on DESCRIPTIVES. */
 119 struct dsc_var
 120   {
 121     struct variable *v;         /* Variable to calculate on. */
 122     char z_name[9];             /* Name for z-score variable. */
 123     double valid, missing;      /* Valid, missing counts. */
 124     struct moments *moments;    /* Moments. */
 125     double min, max;            /* Maximum and mimimum values. */
 126     double stats[DSC_N_STATS];  /* All the stats' values. */
 127   };
 128
 129 /* Output format. */
 130 enum dsc_format
 131   {
 132     DSC_LINE,           /* Abbreviated format. */
 133     DSC_SERIAL          /* Long format. */
 134   };
 135
 136 /* A DESCRIPTIVES procedure. */
 137 struct dsc_proc
 138   {
 139     /* Per-variable info. */
 140     struct dsc_var *vars;       /* Variables. */
 141     size_t var_cnt;             /* Number of variables. */
 142
 143     /* User options. */
 144     enum dsc_missing_type missing_type; /* Treatment of missing values. */
 145     int include_user_missing;   /* Nonzero to include user-missing values. */
 146     int show_var_labels;        /* Nonzero to show variable labels. */
 147     int show_index;             /* Nonzero to show variable index. */
 148     enum dsc_format format;     /* Output format. */
 149
 150     /* Accumulated results. */
 151     double missing_listwise;    /* Sum of weights of cases missing listwise. */
 152     double valid;               /* Sum of weights of valid cases. */
 153     int bad_warn;               /* Warn if bad weight found. */
 154     enum dsc_statistic sort_by_stat; /* Statistic to sort by; -1: name. */
 155     int sort_ascending;         /* !0: ascending order; 0: descending. */
 156     unsigned long show_stats;   /* Statistics to display. */
 157     unsigned long calc_stats;   /* Statistics to calculate. */
 158     enum moment max_moment;     /* Highest moment needed for stats. */
 159   };
 160
 161 /* Parsing. */
 162 static enum dsc_statistic match_statistic (void);
 163 static void free_dsc_proc (struct dsc_proc *);
 164
 165 /* Z-score functions. */
 166 static int try_name (struct dsc_proc *dsc, char *name);
 167 static int generate_z_varname (struct dsc_proc *dsc, char *z_name,
 168                                const char *name, int *z_cnt);
 169 static void dump_z_table (struct dsc_proc *);
 170 static void setup_z_trns (struct dsc_proc *);
 171
 172 /* Procedure execution functions. */
 173 static void calc_descriptives (const struct casefile *, void *dsc_);
 174 static void display (struct dsc_proc *dsc);
 175 \f
 176 /* Parser and outline. */
 177
 178 /* Handles DESCRIPTIVES. */
 179 int
 180 cmd_descriptives (void)
 181 {
 182   struct dsc_proc *dsc;
 183   struct variable **vars = NULL;
 184   int var_cnt = 0;
 185   int save_z_scores = 0;
 186   int z_cnt = 0;
 187   int i;
 188
 189   /* Create and initialize dsc. */
 190   dsc = xmalloc (sizeof *dsc);
 191   dsc->vars = NULL;
 192   dsc->var_cnt = 0;
 193   dsc->missing_type = DSC_VARIABLE;
 194   dsc->include_user_missing = 0;
 195   dsc->show_var_labels = 1;
 196   dsc->show_index = 0;
 197   dsc->format = DSC_LINE;
 198   dsc->missing_listwise = 0.;
 199   dsc->valid = 0.;
 200   dsc->bad_warn = 1;
 201   dsc->sort_by_stat = DSC_NONE;
 202   dsc->sort_ascending = 1;
 203   dsc->show_stats = dsc->calc_stats = DEFAULT_STATS;
 204
 205   /* Parse DESCRIPTIVES. */
 206   while (token != '.')
 207     {
 208       if (lex_match_id ("MISSING"))
 209         {
 210           lex_match ('=');
 211           while (token != '.' && token != '/')
 212             {
 213               if (lex_match_id ("VARIABLE"))
 214                 dsc->missing_type = DSC_VARIABLE;
 215               else if (lex_match_id ("LISTWISE"))
 216                 dsc->missing_type = DSC_LISTWISE;
 217               else if (lex_match_id ("INCLUDE"))
 218                 dsc->include_user_missing = 1;
 219               else
 220                 {
 221                   lex_error (NULL);
 222                   goto error;
 223                 }
 224               lex_match (',');
 225             }
 226         }
 227       else if (lex_match_id ("SAVE"))
 228         save_z_scores = 1;
 229       else if (lex_match_id ("FORMAT"))
 230         {
 231           lex_match ('=');
 232           while (token != '.' && token != '/')
 233             {
 234               if (lex_match_id ("LABELS"))
 235                 dsc->show_var_labels = 1;
 236               else if (lex_match_id ("NOLABELS"))
 237                 dsc->show_var_labels = 0;
 238               else if (lex_match_id ("INDEX"))
 239                 dsc->show_index = 1;
 240               else if (lex_match_id ("NOINDEX"))
 241                 dsc->show_index = 0;
 242               else if (lex_match_id ("LINE"))
 243                 dsc->format = DSC_LINE;
 244               else if (lex_match_id ("SERIAL"))
 245                 dsc->format = DSC_SERIAL;
 246               else
 247                 {
 248                   lex_error (NULL);
 249                   goto error;
 250                 }
 251               lex_match (',');
 252             }
 253         }
 254       else if (lex_match_id ("STATISTICS"))
 255         {
 256           lex_match ('=');
 257           dsc->show_stats = 0;
 258           while (token != '.' && token != '/')
 259             {
 260               if (lex_match (T_ALL))
 261                 dsc->show_stats |= (1ul << DSC_N_STATS) - 1;
 262               else if (lex_match_id ("DEFAULT"))
 263                 dsc->show_stats |= DEFAULT_STATS;
 264               else
 265                 {
 266                   dsc->show_stats |= 1ul << (match_statistic ());
 267                   if (dsc->show_stats == DSC_NONE)
 268                     dsc->show_stats = DEFAULT_STATS;
 269                 }
 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 = xrealloc (dsc->vars, sizeof *dsc->vars * var_cnt);
 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           int 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   int 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 (!strcmp (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, int *z_cnt)
 483 {
 484   char name[10];
 485
 486   /* Try a name based on the original variable name. */
 487   name[0] = 'Z';
 488   strcpy (name + 1, var_name);
 489   name[8] = '\0';
 490   if (try_name (dsc, name))
 491     {
 492       strcpy (z_name, name);
 493       return 1;
 494     }
 495
 496   /* Generate a synthetic name. */
 497   for (;;)
 498     {
 499       (*z_cnt)++;
 500
 501       if (*z_cnt <= 99)
 502         sprintf (name, "ZSC%03d", *z_cnt);
 503       else if (*z_cnt <= 108)
 504         sprintf (name, "STDZ%02d", *z_cnt - 99);
 505       else if (*z_cnt <= 117)
 506         sprintf (name, "ZZZZ%02d", *z_cnt - 108);
 507       else if (*z_cnt <= 126)
 508         sprintf (name, "ZQZQ%02d", *z_cnt - 117);
 509       else
 510         {
 511           msg (SE, _("Ran out of generic names for Z-score variables.  "
 512                      "There are only 126 generic names: ZSC001-ZSC0999, "
 513                      "STDZ01-STDZ09, ZZZZ01-ZZZZ09, ZQZQ01-ZQZQ09."));
 514           return 0;
 515         }
 516
 517       if (try_name (dsc, name))
 518         {
 519           strcpy (z_name, name);
 520           return 1;
 521         }
 522     }
 523 }
 524
 525 /* Outputs a table describing the mapping between source
 526    variables and Z-score variables. */
 527 static void
 528 dump_z_table (struct dsc_proc *dsc)
 529 {
 530   int cnt = 0;
 531   struct tab_table *t;
 532
 533   {
 534     int i;
 535
 536     for (i = 0; i < dsc->var_cnt; i++)
 537       if (dsc->vars[i].z_name[0] != '\0')
 538         cnt++;
 539   }
 540
 541   t = tab_create (2, cnt + 1, 0);
 542   tab_title (t, 0, _("Mapping of variables to corresponding Z-scores."));
 543   tab_columns (t, SOM_COL_DOWN, 1);
 544   tab_headers (t, 0, 0, 1, 0);
 545   tab_box (t, TAL_1, TAL_1, TAL_0, TAL_1, 0, 0, 1, cnt);
 546   tab_hline (t, TAL_2, 0, 1, 1);
 547   tab_text (t, 0, 0, TAB_CENTER | TAT_TITLE, _("Source"));
 548   tab_text (t, 1, 0, TAB_CENTER | TAT_TITLE, _("Target"));
 549   tab_dim (t, tab_natural_dimensions);
 550
 551   {
 552     int i, y;
 553
 554     for (i = 0, y = 1; i < dsc->var_cnt; i++)
 555       if (dsc->vars[i].z_name[0] != '\0')
 556         {
 557           tab_text (t, 0, y, TAB_LEFT, dsc->vars[i].v->name);
 558           tab_text (t, 1, y++, TAB_LEFT, dsc->vars[i].z_name);
 559         }
 560   }
 561
 562   tab_submit (t);
 563 }
 564
 565 /* Transformation function to calculate Z-scores. Will return SYSMIS if any of
 566    the following are true: 1) mean or standard deviation is SYSMIS 2) score is
 567    SYSMIS 3) score is user missing and they were not included in the original
 568    analyis. 4) any of the variables in the original analysis were missing
 569    (either system or user-missing values that weren't included).
 570 */
 571 static int
 572 descriptives_trns_proc (struct trns_header *trns, struct ccase * c,
 573                         int case_idx UNUSED)
 574 {
 575   struct dsc_trns *t = (struct dsc_trns *) trns;
 576   struct dsc_z_score *z;
 577   struct variable **vars;
 578   int all_sysmis = 0;
 579
 580   if (t->missing_type == DSC_LISTWISE)
 581     {
 582       assert(t->vars);
 583       for (vars = t->vars; vars < t->vars + t->var_cnt; vars++)
 584         {
 585           double score = case_num (c, (*vars)->fv);
 586           if ( score == SYSMIS || (!t->include_user_missing
 587                                    && is_num_user_missing(score, *vars)) )
 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 && is_num_user_missing(input, z->v)))
 603         *output = SYSMIS;
 604       else
 605         *output = (input - z->mean) / z->std_dev;
 606     }
 607   return -1;
 608 }
 609
 610 /* Frees a descriptives_trns struct. */
 611 static void
 612 descriptives_trns_free (struct trns_header * trns)
 613 {
 614   struct dsc_trns *t = (struct dsc_trns *) trns;
 615
 616   free (t->z_scores);
 617   assert((t->missing_type != DSC_LISTWISE) ^ (t->vars != NULL));
 618   free (t->vars);
 619 }
 620
 621 /* Sets up a transformation to calculate Z scores. */
 622 static void
 623 setup_z_trns (struct dsc_proc *dsc)
 624 {
 625   struct dsc_trns *t;
 626   int cnt, i;
 627
 628   for (cnt = i = 0; i < dsc->var_cnt; i++)
 629     if (dsc->vars[i].z_name[0] != '\0')
 630       cnt++;
 631
 632   t = xmalloc (sizeof *t);
 633   t->h.proc = descriptives_trns_proc;
 634   t->h.free = descriptives_trns_free;
 635   t->z_scores = xmalloc (cnt * sizeof *t->z_scores);
 636   t->z_score_cnt = cnt;
 637   t->missing_type = dsc->missing_type;
 638   t->include_user_missing = dsc->include_user_missing;
 639   if ( t->missing_type == DSC_LISTWISE )
 640     {
 641       t->var_cnt = dsc->var_cnt;
 642       t->vars = xmalloc(t->var_cnt * sizeof *t->vars);
 643       for (i = 0; i < t->var_cnt; i++)
 644         t->vars[i] = dsc->vars[i].v;
 645     }
 646   else
 647     {
 648       t->var_cnt = 0;
 649       t->vars = NULL;
 650     }
 651
 652
 653   for (cnt = i = 0; i < dsc->var_cnt; i++)
 654     {
 655       struct dsc_var *dv = &dsc->vars[i];
 656       if (dv->z_name[0] != '\0')
 657         {
 658           struct dsc_z_score *z;
 659           char *cp;
 660           struct variable *dst_var;
 661
 662           dst_var = dict_create_var_assert (default_dict, dv->z_name, 0);
 663           dst_var->init = 0;
 664           if (dv->v->label)
 665             {
 666               dst_var->label = xmalloc (strlen (dv->v->label) + 12);
 667               cp = stpcpy (dst_var->label, _("Z-score of "));
 668               strcpy (cp, dv->v->label);
 669             }
 670           else
 671             {
 672               dst_var->label = xmalloc (strlen (dv->v->name) + 12);
 673               cp = stpcpy (dst_var->label, _("Z-score of "));
 674               strcpy (cp, dv->v->name);
 675             }
 676
 677           z = &t->z_scores[cnt++];
 678           z->src_idx = dv->v->fv;
 679           z->dst_idx = dst_var->fv;
 680           z->mean = dv->stats[DSC_MEAN];
 681           z->std_dev = dv->stats[DSC_STDDEV];
 682           z->v = dv->v;
 683         }
 684     }
 685
 686   add_transformation ((struct trns_header *) t);
 687 }
 688 \f
 689 /* Statistical calculation. */
 690
 691 static int listwise_missing (struct dsc_proc *dsc, const struct ccase *c);
 692
 693 /* Calculates and displays descriptive statistics for the cases
 694    in CF. */
 695 static void
 696 calc_descriptives (const struct casefile *cf, void *dsc_)
 697 {
 698   struct dsc_proc *dsc = dsc_;
 699   struct casereader *reader;
 700   struct ccase c;
 701   int i;
 702
 703   for (i = 0; i < dsc->var_cnt; i++)
 704     {
 705       struct dsc_var *dv = &dsc->vars[i];
 706
 707       dv->valid = dv->missing = 0.0;
 708       if (dv->moments != NULL)
 709         moments_clear (dv->moments);
 710       dv->min = DBL_MAX;
 711       dv->max = -DBL_MAX;
 712     }
 713   dsc->missing_listwise = 0.;
 714   dsc->valid = 0.;
 715
 716   /* First pass to handle most of the work. */
 717   for (reader = casefile_get_reader (cf);
 718        casereader_read (reader, &c);
 719        case_destroy (&c))
 720     {
 721       double weight = dict_get_case_weight (default_dict, &c, &dsc->bad_warn);
 722       if (weight <= 0.0)
 723         continue;
 724
 725       /* Check for missing values. */
 726       if (listwise_missing (dsc, &c))
 727         {
 728           dsc->missing_listwise += weight;
 729           if (dsc->missing_type == DSC_LISTWISE)
 730             continue;
 731         }
 732       dsc->valid += weight;
 733
 734       for (i = 0; i < dsc->var_cnt; i++)
 735         {
 736           struct dsc_var *dv = &dsc->vars[i];
 737           double x = case_num (&c, dv->v->fv);
 738
 739           if (dsc->missing_type != DSC_LISTWISE
 740               && (x == SYSMIS
 741                   || (!dsc->include_user_missing
 742                       && is_num_user_missing (x, dv->v))))
 743             {
 744               dv->missing += weight;
 745               continue;
 746             }
 747
 748           if (dv->moments != NULL)
 749             moments_pass_one (dv->moments, x, weight);
 750
 751           if (x < dv->min)
 752             dv->min = x;
 753           if (x > dv->max)
 754             dv->max = x;
 755         }
 756     }
 757   casereader_destroy (reader);
 758
 759   /* Second pass for higher-order moments. */
 760   if (dsc->max_moment > MOMENT_MEAN)
 761     {
 762       for (reader = casefile_get_reader (cf);
 763            casereader_read (reader, &c);
 764            case_destroy (&c))
 765         {
 766           double weight = dict_get_case_weight (default_dict, &c,
 767                                                 &dsc->bad_warn);
 768           if (weight <= 0.0)
 769             continue;
 770
 771           /* Check for missing values. */
 772           if (listwise_missing (dsc, &c)
 773               && dsc->missing_type == DSC_LISTWISE)
 774             continue;
 775
 776           for (i = 0; i < dsc->var_cnt; i++)
 777             {
 778               struct dsc_var *dv = &dsc->vars[i];
 779               double x = case_num (&c, dv->v->fv);
 780
 781               if (dsc->missing_type != DSC_LISTWISE
 782                   && (x == SYSMIS
 783                       || (!dsc->include_user_missing
 784                           && is_num_user_missing (x, dv->v))))
 785                 continue;
 786
 787               if (dv->moments != NULL)
 788                 moments_pass_two (dv->moments, x, weight);
 789             }
 790         }
 791       casereader_destroy (reader);
 792     }
 793
 794   /* Calculate results. */
 795   for (i = 0; i < dsc->var_cnt; i++)
 796     {
 797       struct dsc_var *dv = &dsc->vars[i];
 798       double W;
 799       int j;
 800
 801       for (j = 0; j < DSC_N_STATS; j++)
 802         dv->stats[j] = SYSMIS;
 803
 804       dv->valid = W = dsc->valid - dv->missing;
 805
 806       if (dv->moments != NULL)
 807         moments_calculate (dv->moments, NULL,
 808                            &dv->stats[DSC_MEAN], &dv->stats[DSC_VARIANCE],
 809                            &dv->stats[DSC_SKEWNESS], &dv->stats[DSC_KURTOSIS]);
 810       if (dsc->calc_stats & (1ul << DSC_SEMEAN)
 811           && dv->stats[DSC_VARIANCE] != SYSMIS && W > 0.)
 812         dv->stats[DSC_SEMEAN] = sqrt (dv->stats[DSC_VARIANCE]) / sqrt (W);
 813       if (dsc->calc_stats & (1ul << DSC_STDDEV)
 814           && dv->stats[DSC_VARIANCE] != SYSMIS)
 815         dv->stats[DSC_STDDEV] = sqrt (dv->stats[DSC_VARIANCE]);
 816       if (dsc->calc_stats & (1ul << DSC_SEKURT))
 817         if (dv->stats[DSC_KURTOSIS] != SYSMIS)
 818             dv->stats[DSC_SEKURT] = calc_sekurt (W);
 819       if (dsc->calc_stats & (1ul << DSC_SESKEW)
 820           && dv->stats[DSC_SKEWNESS] != SYSMIS)
 821         dv->stats[DSC_SESKEW] = calc_seskew (W);
 822       dv->stats[DSC_RANGE] = ((dv->min == DBL_MAX || dv->max == -DBL_MAX)
 823                               ? SYSMIS : dv->max - dv->min);
 824       dv->stats[DSC_MIN] = dv->min == DBL_MAX ? SYSMIS : dv->min;
 825       dv->stats[DSC_MAX] = dv->max == -DBL_MAX ? SYSMIS : dv->max;
 826       if (dsc->calc_stats & (1ul << DSC_SUM))
 827         dv->stats[DSC_SUM] = W * dv->stats[DSC_MEAN];
 828     }
 829
 830   /* Output results. */
 831   display (dsc);
 832 }
 833
 834 /* Returns nonzero if any of the descriptives variables in DSC's
 835    variable list have missing values in case C, zero otherwise. */
 836 static int
 837 listwise_missing (struct dsc_proc *dsc, const struct ccase *c)
 838 {
 839   int i;
 840
 841   for (i = 0; i < dsc->var_cnt; i++)
 842     {
 843       struct dsc_var *dv = &dsc->vars[i];
 844       double x = case_num (c, dv->v->fv);
 845
 846       if (x == SYSMIS
 847           || (!dsc->include_user_missing && is_num_user_missing (x, dv->v)))
 848         return 1;
 849     }
 850   return 0;
 851 }
 852 \f
 853 /* Statistical display. */
 854
 855 static algo_compare_func descriptives_compare_dsc_vars;
 856
 857 /* Displays a table of descriptive statistics for DSC. */
 858 static void
 859 display (struct dsc_proc *dsc)
 860 {
 861   int i, j;
 862   int nc;
 863   struct tab_table *t;
 864
 865   nc = 1 + (dsc->format == DSC_SERIAL ? 2 : 1);
 866   for (i = 0; i < DSC_N_STATS; i++)
 867     if (dsc->show_stats & (1ul << i))
 868       nc++;
 869
 870   if (dsc->sort_by_stat != DSC_NONE)
 871     sort (dsc->vars, dsc->var_cnt, sizeof *dsc->vars,
 872           descriptives_compare_dsc_vars, dsc);
 873
 874   t = tab_create (nc, dsc->var_cnt + 1, 0);
 875   tab_headers (t, 1, 0, 1, 0);
 876   tab_box (t, TAL_1, TAL_1, -1, -1, 0, 0, nc - 1, dsc->var_cnt);
 877   tab_box (t, -1, -1, -1, TAL_1, 1, 0, nc - 1, dsc->var_cnt);
 878   tab_hline (t, TAL_2, 0, nc - 1, 1);
 879   tab_vline (t, TAL_2, 1, 0, dsc->var_cnt);
 880   tab_dim (t, tab_natural_dimensions);
 881
 882   nc = 0;
 883   tab_text (t, nc++, 0, TAB_LEFT | TAT_TITLE, _("Variable"));
 884   if (dsc->format == DSC_SERIAL)
 885     {
 886       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Valid N"));
 887       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Missing N"));
 888     }
 889   else
 890     tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, "N");
 891
 892   for (i = 0; i < DSC_N_STATS; i++)
 893     if (dsc->show_stats & (1ul << i))
 894       {
 895         const char *title = gettext (dsc_info[i].name);
 896         tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, title);
 897       }
 898
 899   for (i = 0; i < dsc->var_cnt; i++)
 900     {
 901       struct dsc_var *dv = &dsc->vars[i];
 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 = strcmp (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 }