src/descript.c

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