pintos-os.org Git - pspp/blob - src/language/stats/descriptives.c

   1 /* PSPP - a program for statistical analysis.
   2    Copyright (C) 1997-2000, 2009-2014 Free Software Foundation, Inc.
   3
   4    This program is free software: you can redistribute it and/or modify
   5    it under the terms of the GNU General Public License as published by
   6    the Free Software Foundation, either version 3 of the License, or
   7    (at your option) any later version.
   8
   9    This program is distributed in the hope that it will be useful,
  10    but WITHOUT ANY WARRANTY; without even the implied warranty of
  11    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  12    GNU General Public License for more details.
  13
  14    You should have received a copy of the GNU General Public License
  15    along with this program.  If not, see <http://www.gnu.org/licenses/>. */
  16
  17 #include <config.h>
  18
  19 #include <float.h>
  20 #include <limits.h>
  21 #include <math.h>
  22 #include <stdlib.h>
  23
  24 #include "data/casegrouper.h"
  25 #include "data/casereader.h"
  26 #include "data/casewriter.h"
  27 #include "data/dataset.h"
  28 #include "data/dictionary.h"
  29 #include "data/transformations.h"
  30 #include "data/variable.h"
  31 #include "language/command.h"
  32 #include "language/dictionary/split-file.h"
  33 #include "language/lexer/lexer.h"
  34 #include "language/lexer/variable-parser.h"
  35 #include "libpspp/array.h"
  36 #include "libpspp/assertion.h"
  37 #include "libpspp/compiler.h"
  38 #include "libpspp/i18n.h"
  39 #include "libpspp/message.h"
  40 #include "math/moments.h"
  41 #include "output/tab.h"
  42
  43 #include "gl/xalloc.h"
  44
  45 #include "gettext.h"
  46 #define _(msgid) gettext (msgid)
  47 #define N_(msgid) msgid
  48
  49 /* DESCRIPTIVES private data. */
  50
  51 struct dsc_proc;
  52
  53 /* Handling of missing values. */
  54 enum dsc_missing_type
  55   {
  56     DSC_VARIABLE,       /* Handle missing values on a per-variable basis. */
  57     DSC_LISTWISE        /* Discard entire case if any variable is missing. */
  58   };
  59
  60 /* Describes properties of a distribution for the purpose of
  61    calculating a Z-score. */
  62 struct dsc_z_score
  63   {
  64     const struct variable *src_var;   /* Variable on which z-score is based. */
  65     struct variable *z_var;     /* New z-score variable. */
  66     double mean;                /* Distribution mean. */
  67     double std_dev;             /* Distribution standard deviation. */
  68   };
  69
  70 /* DESCRIPTIVES transformation (for calculating Z-scores). */
  71 struct dsc_trns
  72   {
  73     struct dsc_z_score *z_scores; /* Array of Z-scores. */
  74     int z_score_cnt;            /* Number of Z-scores. */
  75     const struct variable **vars;     /* Variables for listwise missing checks. */
  76     size_t var_cnt;             /* Number of variables. */
  77     enum dsc_missing_type missing_type; /* Treatment of missing values. */
  78     enum mv_class exclude;      /* Classes of missing values to exclude. */
  79     struct casereader *z_reader; /* Reader for count, mean, stddev. */
  80     casenumber count;            /* Number left in this SPLIT FILE group.*/
  81     bool ok;
  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     const struct variable *v;         /* Variable to calculate on. */
 131     char *z_name;                     /* 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     enum mv_class exclude;      /* Classes of missing values to exclude. */
 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     /* Z scores. */
 170     struct casewriter *z_writer; /* Mean and stddev per SPLIT FILE group. */
 171   };
 172
 173 /* Parsing. */
 174 static enum dsc_statistic match_statistic (struct lexer *);
 175 static void free_dsc_proc (struct dsc_proc *);
 176
 177 /* Z-score functions. */
 178 static bool try_name (const struct dictionary *dict,
 179                       struct dsc_proc *dsc, const char *name);
 180 static char *generate_z_varname (const struct dictionary *dict,
 181                                  struct dsc_proc *dsc,
 182                                  const char *name, int *z_cnt);
 183 static void dump_z_table (struct dsc_proc *);
 184 static void setup_z_trns (struct dsc_proc *, struct dataset *);
 185
 186 /* Procedure execution functions. */
 187 static void calc_descriptives (struct dsc_proc *, struct casereader *,
 188                                struct dataset *);
 189 static void display (struct dsc_proc *dsc);
 190 \f
 191 /* Parser and outline. */
 192
 193 /* Handles DESCRIPTIVES. */
 194 int
 195 cmd_descriptives (struct lexer *lexer, struct dataset *ds)
 196 {
 197   struct dictionary *dict = dataset_dict (ds);
 198   struct dsc_proc *dsc;
 199   const struct variable **vars = NULL;
 200   size_t var_cnt = 0;
 201   int save_z_scores = 0;
 202   int z_cnt = 0;
 203   size_t i;
 204   bool ok;
 205
 206   struct casegrouper *grouper;
 207   struct casereader *group;
 208
 209   /* Create and initialize dsc. */
 210   dsc = xmalloc (sizeof *dsc);
 211   dsc->vars = NULL;
 212   dsc->var_cnt = 0;
 213   dsc->missing_type = DSC_VARIABLE;
 214   dsc->exclude = MV_ANY;
 215   dsc->show_var_labels = 1;
 216   dsc->show_index = 0;
 217   dsc->format = DSC_LINE;
 218   dsc->missing_listwise = 0.;
 219   dsc->valid = 0.;
 220   dsc->bad_warn = 1;
 221   dsc->sort_by_stat = DSC_NONE;
 222   dsc->sort_ascending = 1;
 223   dsc->show_stats = dsc->calc_stats = DEFAULT_STATS;
 224   dsc->z_writer = NULL;
 225
 226   /* Parse DESCRIPTIVES. */
 227   while (lex_token (lexer) != T_ENDCMD)
 228     {
 229       if (lex_match_id (lexer, "MISSING"))
 230         {
 231           lex_match (lexer, T_EQUALS);
 232           while (lex_token (lexer) != T_ENDCMD && lex_token (lexer) != T_SLASH)
 233             {
 234               if (lex_match_id (lexer, "VARIABLE"))
 235                 dsc->missing_type = DSC_VARIABLE;
 236               else if (lex_match_id (lexer, "LISTWISE"))
 237                 dsc->missing_type = DSC_LISTWISE;
 238               else if (lex_match_id (lexer, "INCLUDE"))
 239                 dsc->exclude = MV_SYSTEM;
 240               else
 241                 {
 242                   lex_error (lexer, NULL);
 243                   goto error;
 244                 }
 245               lex_match (lexer, T_COMMA);
 246             }
 247         }
 248       else if (lex_match_id (lexer, "SAVE"))
 249         save_z_scores = 1;
 250       else if (lex_match_id (lexer, "FORMAT"))
 251         {
 252           lex_match (lexer, T_EQUALS);
 253           while (lex_token (lexer) != T_ENDCMD && lex_token (lexer) != T_SLASH)
 254             {
 255               if (lex_match_id (lexer, "LABELS"))
 256                 dsc->show_var_labels = 1;
 257               else if (lex_match_id (lexer, "NOLABELS"))
 258                 dsc->show_var_labels = 0;
 259               else if (lex_match_id (lexer, "INDEX"))
 260                 dsc->show_index = 1;
 261               else if (lex_match_id (lexer, "NOINDEX"))
 262                 dsc->show_index = 0;
 263               else if (lex_match_id (lexer, "LINE"))
 264                 dsc->format = DSC_LINE;
 265               else if (lex_match_id (lexer, "SERIAL"))
 266                 dsc->format = DSC_SERIAL;
 267               else
 268                 {
 269                   lex_error (lexer, NULL);
 270                   goto error;
 271                 }
 272               lex_match (lexer, T_COMMA);
 273             }
 274         }
 275       else if (lex_match_id (lexer, "STATISTICS"))
 276         {
 277           lex_match (lexer, T_EQUALS);
 278           dsc->show_stats = 0;
 279           while (lex_token (lexer) != T_ENDCMD && lex_token (lexer) != T_SLASH)
 280             {
 281               if (lex_match (lexer, T_ALL))
 282                 dsc->show_stats |= (1ul << DSC_N_STATS) - 1;
 283               else if (lex_match_id (lexer, "DEFAULT"))
 284                 dsc->show_stats |= DEFAULT_STATS;
 285               else
 286                 dsc->show_stats |= 1ul << (match_statistic (lexer));
 287               lex_match (lexer, T_COMMA);
 288             }
 289           if (dsc->show_stats == 0)
 290             dsc->show_stats = DEFAULT_STATS;
 291         }
 292       else if (lex_match_id (lexer, "SORT"))
 293         {
 294           lex_match (lexer, T_EQUALS);
 295           if (lex_match_id (lexer, "NAME"))
 296             dsc->sort_by_stat = DSC_NAME;
 297           else
 298             {
 299               dsc->sort_by_stat = match_statistic (lexer);
 300               if (dsc->sort_by_stat == DSC_NONE )
 301                 dsc->sort_by_stat = DSC_MEAN;
 302             }
 303           if (lex_match (lexer, T_LPAREN))
 304             {
 305               if (lex_match_id (lexer, "A"))
 306                 dsc->sort_ascending = 1;
 307               else if (lex_match_id (lexer, "D"))
 308                 dsc->sort_ascending = 0;
 309               else
 310                 lex_error (lexer, NULL);
 311               lex_force_match (lexer, T_RPAREN);
 312             }
 313         }
 314       else if (var_cnt == 0)
 315         {
 316           if (lex_next_token (lexer, 1) == T_EQUALS)
 317             {
 318               lex_match_id (lexer, "VARIABLES");
 319               lex_match (lexer, T_EQUALS);
 320             }
 321
 322           while (lex_token (lexer) != T_ENDCMD && lex_token (lexer) != T_SLASH)
 323             {
 324               int i;
 325
 326               if (!parse_variables_const (lexer, dict, &vars, &var_cnt,
 327                                     PV_APPEND | PV_NO_DUPLICATE | PV_NUMERIC))
 328                 goto error;
 329
 330               dsc->vars = xnrealloc ((void *)dsc->vars, var_cnt, sizeof *dsc->vars);
 331               for (i = dsc->var_cnt; i < var_cnt; i++)
 332                 {
 333                   struct dsc_var *dv = &dsc->vars[i];
 334                   dv->v = vars[i];
 335                   dv->z_name = NULL;
 336                   dv->moments = NULL;
 337                 }
 338               dsc->var_cnt = var_cnt;
 339
 340               if (lex_match (lexer, T_LPAREN))
 341                 {
 342                   if (lex_token (lexer) != T_ID)
 343                     {
 344                       lex_error (lexer, NULL);
 345                       goto error;
 346                     }
 347                   if (try_name (dict, dsc, lex_tokcstr (lexer)))
 348                     {
 349                       struct dsc_var *dsc_var = &dsc->vars[dsc->var_cnt - 1];
 350                       dsc_var->z_name = xstrdup (lex_tokcstr (lexer));
 351                       z_cnt++;
 352                     }
 353                   else
 354                     msg (SE, _("Z-score variable name %s would be"
 355                                " a duplicate variable name."), lex_tokcstr (lexer));
 356                   lex_get (lexer);
 357                   if (!lex_force_match (lexer, T_RPAREN))
 358                     goto error;
 359                 }
 360             }
 361         }
 362       else
 363         {
 364           lex_error (lexer, NULL);
 365           goto error;
 366         }
 367
 368       lex_match (lexer, T_SLASH);
 369     }
 370   if (var_cnt == 0)
 371     {
 372       msg (SE, _("No variables specified."));
 373       goto error;
 374     }
 375
 376   /* Construct z-score varnames, show translation table. */
 377   if (z_cnt || save_z_scores)
 378     {
 379       struct caseproto *proto;
 380
 381       if (save_z_scores)
 382         {
 383           int gen_cnt = 0;
 384
 385           for (i = 0; i < dsc->var_cnt; i++)
 386             {
 387               struct dsc_var *dsc_var = &dsc->vars[i];
 388               if (dsc_var->z_name == NULL)
 389                 {
 390                   const char *name = var_get_name (dsc_var->v);
 391                   dsc_var->z_name = generate_z_varname (dict, dsc, name,
 392                                                         &gen_cnt);
 393                   if (dsc_var->z_name == NULL)
 394                     goto error;
 395
 396                   z_cnt++;
 397                 }
 398             }
 399         }
 400
 401       /* It would be better to handle Z scores correctly (however we define
 402          that) when TEMPORARY is in effect, but in the meantime this at least
 403          prevents a use-after-free error.  See bug #38786.  */
 404       if (proc_make_temporary_transformations_permanent (ds))
 405         msg (SW, _("DESCRIPTIVES with Z scores ignores TEMPORARY.  "
 406                    "Temporary transformations will be made permanent."));
 407
 408       proto = caseproto_create ();
 409       for (i = 0; i < 1 + 2 * z_cnt; i++)
 410         proto = caseproto_add_width (proto, 0);
 411       dsc->z_writer = autopaging_writer_create (proto);
 412       caseproto_unref (proto);
 413
 414       dump_z_table (dsc);
 415     }
 416
 417   /* Figure out statistics to display. */
 418   if (dsc->show_stats & (1ul << DSC_SKEWNESS))
 419     dsc->show_stats |= 1ul << DSC_SESKEW;
 420   if (dsc->show_stats & (1ul << DSC_KURTOSIS))
 421     dsc->show_stats |= 1ul << DSC_SEKURT;
 422
 423   /* Figure out which statistics to calculate. */
 424   dsc->calc_stats = dsc->show_stats;
 425   if (z_cnt > 0)
 426     dsc->calc_stats |= (1ul << DSC_MEAN) | (1ul << DSC_STDDEV);
 427   if (dsc->sort_by_stat >= 0)
 428     dsc->calc_stats |= 1ul << dsc->sort_by_stat;
 429   if (dsc->show_stats & (1ul << DSC_SESKEW))
 430     dsc->calc_stats |= 1ul << DSC_SKEWNESS;
 431   if (dsc->show_stats & (1ul << DSC_SEKURT))
 432     dsc->calc_stats |= 1ul << DSC_KURTOSIS;
 433
 434   /* Figure out maximum moment needed and allocate moments for
 435      the variables. */
 436   dsc->max_moment = MOMENT_NONE;
 437   for (i = 0; i < DSC_N_STATS; i++)
 438     if (dsc->calc_stats & (1ul << i) && dsc_info[i].moment > dsc->max_moment)
 439       dsc->max_moment = dsc_info[i].moment;
 440   if (dsc->max_moment != MOMENT_NONE)
 441     for (i = 0; i < dsc->var_cnt; i++)
 442       dsc->vars[i].moments = moments_create (dsc->max_moment);
 443
 444   /* Data pass. */
 445   grouper = casegrouper_create_splits (proc_open_filtering (ds, z_cnt == 0),
 446                                        dict);
 447   while (casegrouper_get_next_group (grouper, &group))
 448     calc_descriptives (dsc, group, ds);
 449   ok = casegrouper_destroy (grouper);
 450   ok = proc_commit (ds) && ok;
 451
 452   /* Z-scoring! */
 453   if (ok && z_cnt)
 454     setup_z_trns (dsc, ds);
 455
 456   /* Done. */
 457   free (vars);
 458   free_dsc_proc (dsc);
 459   return ok ? CMD_SUCCESS : CMD_CASCADING_FAILURE;
 460
 461  error:
 462   free (vars);
 463   free_dsc_proc (dsc);
 464   return CMD_FAILURE;
 465 }
 466
 467 /* Returns the statistic named by the current token and skips past the token.
 468    Returns DSC_NONE if no statistic is given (e.g., subcommand with no
 469    specifiers). Emits an error if the current token ID does not name a
 470    statistic. */
 471 static enum dsc_statistic
 472 match_statistic (struct lexer *lexer)
 473 {
 474   if (lex_token (lexer) == T_ID)
 475     {
 476       enum dsc_statistic stat;
 477
 478       for (stat = 0; stat < DSC_N_STATS; stat++)
 479         if (lex_match_id (lexer, dsc_info[stat].identifier))
 480           return stat;
 481
 482       lex_get (lexer);
 483       lex_error (lexer, _("expecting statistic name: reverting to default"));
 484     }
 485
 486   return DSC_NONE;
 487 }
 488
 489 /* Frees DSC. */
 490 static void
 491 free_dsc_proc (struct dsc_proc *dsc)
 492 {
 493   size_t i;
 494
 495   if (dsc == NULL)
 496     return;
 497
 498   for (i = 0; i < dsc->var_cnt; i++)
 499     {
 500       struct dsc_var *dsc_var = &dsc->vars[i];
 501       free (dsc_var->z_name);
 502       moments_destroy (dsc_var->moments);
 503     }
 504   casewriter_destroy (dsc->z_writer);
 505   free (dsc->vars);
 506   free (dsc);
 507 }
 508 \f
 509 /* Z scores. */
 510
 511 /* Returns false if NAME is a duplicate of any existing variable name or
 512    of any previously-declared z-var name; otherwise returns true. */
 513 static bool
 514 try_name (const struct dictionary *dict, struct dsc_proc *dsc,
 515           const char *name)
 516 {
 517   size_t i;
 518
 519   if (dict_lookup_var (dict, name) != NULL)
 520     return false;
 521   for (i = 0; i < dsc->var_cnt; i++)
 522     {
 523       struct dsc_var *dsc_var = &dsc->vars[i];
 524       if (dsc_var->z_name != NULL && !utf8_strcasecmp (dsc_var->z_name, name))
 525         return false;
 526     }
 527   return true;
 528 }
 529
 530 /* Generates a name for a Z-score variable based on a variable
 531    named VAR_NAME, given that *Z_CNT generated variable names are
 532    known to already exist.  If successful, returns the new name
 533    as a dynamically allocated string.  On failure, returns NULL. */
 534 static char *
 535 generate_z_varname (const struct dictionary *dict, struct dsc_proc *dsc,
 536                     const char *var_name, int *z_cnt)
 537 {
 538   char *z_name, *trunc_name;
 539
 540   /* Try a name based on the original variable name. */
 541   z_name = xasprintf ("Z%s", var_name);
 542   trunc_name = utf8_encoding_trunc (z_name, dict_get_encoding (dict),
 543                                     ID_MAX_LEN);
 544   free (z_name);
 545   if (try_name (dict, dsc, trunc_name))
 546     return trunc_name;
 547   free (trunc_name);
 548
 549   /* Generate a synthetic name. */
 550   for (;;)
 551     {
 552       char name[8];
 553
 554       (*z_cnt)++;
 555
 556       if (*z_cnt <= 99)
 557         sprintf (name, "ZSC%03d", *z_cnt);
 558       else if (*z_cnt <= 108)
 559         sprintf (name, "STDZ%02d", *z_cnt - 99);
 560       else if (*z_cnt <= 117)
 561         sprintf (name, "ZZZZ%02d", *z_cnt - 108);
 562       else if (*z_cnt <= 126)
 563         sprintf (name, "ZQZQ%02d", *z_cnt - 117);
 564       else
 565         {
 566           msg (SE, _("Ran out of generic names for Z-score variables.  "
 567                      "There are only 126 generic names: ZSC001-ZSC0999, "
 568                      "STDZ01-STDZ09, ZZZZ01-ZZZZ09, ZQZQ01-ZQZQ09."));
 569           return NULL;
 570         }
 571
 572       if (try_name (dict, dsc, name))
 573         return xstrdup (name);
 574     }
 575   NOT_REACHED();
 576 }
 577
 578 /* Outputs a table describing the mapping between source
 579    variables and Z-score variables. */
 580 static void
 581 dump_z_table (struct dsc_proc *dsc)
 582 {
 583   size_t cnt = 0;
 584   struct tab_table *t;
 585
 586   {
 587     size_t i;
 588
 589     for (i = 0; i < dsc->var_cnt; i++)
 590       if (dsc->vars[i].z_name != NULL)
 591         cnt++;
 592   }
 593
 594   t = tab_create (2, cnt + 1);
 595   tab_title (t, _("Mapping of variables to corresponding Z-scores."));
 596   tab_headers (t, 0, 0, 1, 0);
 597   tab_box (t, TAL_1, TAL_1, TAL_0, TAL_1, 0, 0, 1, cnt);
 598   tab_hline (t, TAL_2, 0, 1, 1);
 599   tab_text (t, 0, 0, TAB_CENTER | TAT_TITLE, _("Source"));
 600   tab_text (t, 1, 0, TAB_CENTER | TAT_TITLE, _("Target"));
 601
 602   {
 603     size_t i, y;
 604
 605     for (i = 0, y = 1; i < dsc->var_cnt; i++)
 606       if (dsc->vars[i].z_name != NULL)
 607         {
 608           tab_text (t, 0, y, TAB_LEFT, var_to_string (dsc->vars[i].v));
 609           tab_text (t, 1, y++, TAB_LEFT, dsc->vars[i].z_name);
 610         }
 611   }
 612
 613   tab_submit (t);
 614 }
 615
 616 /* Transformation function to calculate Z-scores. Will return SYSMIS if any of
 617    the following are true: 1) mean or standard deviation is SYSMIS 2) score is
 618    SYSMIS 3) score is user missing and they were not included in the original
 619    analyis. 4) any of the variables in the original analysis were missing
 620    (either system or user-missing values that weren't included).
 621 */
 622 static int
 623 descriptives_trns_proc (void *trns_, struct ccase **c,
 624                         casenumber case_idx UNUSED)
 625 {
 626   struct dsc_trns *t = trns_;
 627   struct dsc_z_score *z;
 628   const struct variable **vars;
 629   int all_sysmis = 0;
 630
 631   if (t->count <= 0)
 632     {
 633       struct ccase *z_case;
 634
 635       z_case = casereader_read (t->z_reader);
 636       if (z_case)
 637         {
 638           size_t z_idx = 0;
 639
 640           t->count = case_num_idx (z_case, z_idx++);
 641           for (z = t->z_scores; z < t->z_scores + t->z_score_cnt; z++)
 642             {
 643               z->mean = case_num_idx (z_case, z_idx++);
 644               z->std_dev = case_num_idx (z_case, z_idx++);
 645             }
 646           case_unref (z_case);
 647         }
 648       else
 649         {
 650           if (t->ok)
 651             {
 652               msg (SE, _("Internal error processing Z scores"));
 653               t->ok = false;
 654             }
 655           for (z = t->z_scores; z < t->z_scores + t->z_score_cnt; z++)
 656             z->mean = z->std_dev = SYSMIS;
 657         }
 658     }
 659   t->count--;
 660
 661   if (t->missing_type == DSC_LISTWISE)
 662     {
 663       assert(t->vars);
 664       for (vars = t->vars; vars < t->vars + t->var_cnt; vars++)
 665         {
 666           double score = case_num (*c, *vars);
 667           if (var_is_num_missing (*vars, score, t->exclude))
 668             {
 669               all_sysmis = 1;
 670               break;
 671             }
 672         }
 673     }
 674
 675   *c = case_unshare (*c);
 676   for (z = t->z_scores; z < t->z_scores + t->z_score_cnt; z++)
 677     {
 678       double input = case_num (*c, z->src_var);
 679       double *output = &case_data_rw (*c, z->z_var)->f;
 680
 681       if (z->mean == SYSMIS || z->std_dev == SYSMIS || all_sysmis
 682           || var_is_num_missing (z->src_var, input, t->exclude))
 683         *output = SYSMIS;
 684       else
 685         *output = (input - z->mean) / z->std_dev;
 686     }
 687   return TRNS_CONTINUE;
 688 }
 689
 690 /* Frees a descriptives_trns struct. */
 691 static bool
 692 descriptives_trns_free (void *trns_)
 693 {
 694   struct dsc_trns *t = trns_;
 695   bool ok = t->ok && !casereader_error (t->z_reader);
 696
 697   free (t->z_scores);
 698   casereader_destroy (t->z_reader);
 699   assert((t->missing_type != DSC_LISTWISE) ^ (t->vars != NULL));
 700   free (t->vars);
 701   free (t);
 702
 703   return ok;
 704 }
 705
 706 /* Sets up a transformation to calculate Z scores. */
 707 static void
 708 setup_z_trns (struct dsc_proc *dsc, struct dataset *ds)
 709 {
 710   struct dsc_trns *t;
 711   size_t cnt, i;
 712
 713   for (cnt = i = 0; i < dsc->var_cnt; i++)
 714     if (dsc->vars[i].z_name != NULL)
 715       cnt++;
 716
 717   t = xmalloc (sizeof *t);
 718   t->z_scores = xnmalloc (cnt, sizeof *t->z_scores);
 719   t->z_score_cnt = cnt;
 720   t->missing_type = dsc->missing_type;
 721   t->exclude = dsc->exclude;
 722   if ( t->missing_type == DSC_LISTWISE )
 723     {
 724       t->var_cnt = dsc->var_cnt;
 725       t->vars = xnmalloc (t->var_cnt, sizeof *t->vars);
 726       for (i = 0; i < t->var_cnt; i++)
 727         t->vars[i] = dsc->vars[i].v;
 728     }
 729   else
 730     {
 731       t->var_cnt = 0;
 732       t->vars = NULL;
 733     }
 734   t->z_reader = casewriter_make_reader (dsc->z_writer);
 735   t->count = 0;
 736   t->ok = true;
 737   dsc->z_writer = NULL;
 738
 739   for (cnt = i = 0; i < dsc->var_cnt; i++)
 740     {
 741       struct dsc_var *dv = &dsc->vars[i];
 742       if (dv->z_name != NULL)
 743         {
 744           struct dsc_z_score *z;
 745           struct variable *dst_var;
 746           char *label;
 747
 748           dst_var = dict_create_var_assert (dataset_dict (ds), dv->z_name, 0);
 749
 750           label = xasprintf (_("Z-score of %s"),var_to_string (dv->v));
 751           var_set_label (dst_var, label, false);
 752           free (label);
 753
 754           z = &t->z_scores[cnt++];
 755           z->src_var = dv->v;
 756           z->z_var = dst_var;
 757         }
 758     }
 759
 760   add_transformation (ds,
 761                       descriptives_trns_proc, descriptives_trns_free, t);
 762 }
 763 \f
 764 /* Statistical calculation. */
 765
 766 static bool listwise_missing (struct dsc_proc *dsc, const struct ccase *c);
 767
 768 /* Calculates and displays descriptive statistics for the cases
 769    in CF. */
 770 static void
 771 calc_descriptives (struct dsc_proc *dsc, struct casereader *group,
 772                    struct dataset *ds)
 773 {
 774   struct casereader *pass1, *pass2;
 775   casenumber count;
 776   struct ccase *c;
 777   size_t z_idx;
 778   size_t i;
 779
 780   c = casereader_peek (group, 0);
 781   if (c == NULL)
 782     {
 783       casereader_destroy (group);
 784       return;
 785     }
 786   output_split_file_values (ds, c);
 787   case_unref (c);
 788
 789   group = casereader_create_filter_weight (group, dataset_dict (ds),
 790                                            NULL, NULL);
 791
 792   pass1 = group;
 793   pass2 = dsc->max_moment <= MOMENT_MEAN ? NULL : casereader_clone (pass1);
 794
 795   for (i = 0; i < dsc->var_cnt; i++)
 796     {
 797       struct dsc_var *dv = &dsc->vars[i];
 798
 799       dv->valid = dv->missing = 0.0;
 800       if (dv->moments != NULL)
 801         moments_clear (dv->moments);
 802       dv->min = DBL_MAX;
 803       dv->max = -DBL_MAX;
 804     }
 805   dsc->missing_listwise = 0.;
 806   dsc->valid = 0.;
 807
 808   /* First pass to handle most of the work. */
 809   count = 0;
 810   for (; (c = casereader_read (pass1)) != NULL; case_unref (c))
 811     {
 812       double weight = dict_get_case_weight (dataset_dict (ds), c, NULL);
 813
 814       /* Check for missing values. */
 815       if (listwise_missing (dsc, c))
 816         {
 817           dsc->missing_listwise += weight;
 818           if (dsc->missing_type == DSC_LISTWISE)
 819             continue;
 820         }
 821       dsc->valid += weight;
 822
 823       for (i = 0; i < dsc->var_cnt; i++)
 824         {
 825           struct dsc_var *dv = &dsc->vars[i];
 826           double x = case_num (c, dv->v);
 827
 828           if (var_is_num_missing (dv->v, x, dsc->exclude))
 829             {
 830               dv->missing += weight;
 831               continue;
 832             }
 833
 834           if (dv->moments != NULL)
 835             moments_pass_one (dv->moments, x, weight);
 836
 837           if (x < dv->min)
 838             dv->min = x;
 839           if (x > dv->max)
 840             dv->max = x;
 841         }
 842
 843       count++;
 844     }
 845   if (!casereader_destroy (pass1))
 846     {
 847       casereader_destroy (pass2);
 848       return;
 849     }
 850
 851   /* Second pass for higher-order moments. */
 852   if (dsc->max_moment > MOMENT_MEAN)
 853     {
 854       for (; (c = casereader_read (pass2)) != NULL; case_unref (c))
 855         {
 856           double weight = dict_get_case_weight (dataset_dict (ds), c, NULL);
 857
 858           /* Check for missing values. */
 859           if (dsc->missing_type == DSC_LISTWISE && listwise_missing (dsc, c))
 860             continue;
 861
 862           for (i = 0; i < dsc->var_cnt; i++)
 863             {
 864               struct dsc_var *dv = &dsc->vars[i];
 865               double x = case_num (c, dv->v);
 866
 867               if (var_is_num_missing (dv->v, x, dsc->exclude))
 868                 continue;
 869
 870               if (dv->moments != NULL)
 871                 moments_pass_two (dv->moments, x, weight);
 872             }
 873         }
 874       if (!casereader_destroy (pass2))
 875         return;
 876     }
 877
 878   /* Calculate results. */
 879   if (dsc->z_writer)
 880     {
 881       c = case_create (casewriter_get_proto (dsc->z_writer));
 882       z_idx = 0;
 883       case_data_rw_idx (c, z_idx++)->f = count;
 884     }
 885   else
 886     c = NULL;
 887
 888   for (i = 0; i < dsc->var_cnt; i++)
 889     {
 890       struct dsc_var *dv = &dsc->vars[i];
 891       double W;
 892       int j;
 893
 894       for (j = 0; j < DSC_N_STATS; j++)
 895         dv->stats[j] = SYSMIS;
 896
 897       dv->valid = W = dsc->valid - dv->missing;
 898
 899       if (dv->moments != NULL)
 900         moments_calculate (dv->moments, NULL,
 901                            &dv->stats[DSC_MEAN], &dv->stats[DSC_VARIANCE],
 902                            &dv->stats[DSC_SKEWNESS], &dv->stats[DSC_KURTOSIS]);
 903       if (dsc->calc_stats & (1ul << DSC_SEMEAN)
 904           && dv->stats[DSC_VARIANCE] != SYSMIS && W > 0.)
 905         dv->stats[DSC_SEMEAN] = sqrt (dv->stats[DSC_VARIANCE]) / sqrt (W);
 906       if (dsc->calc_stats & (1ul << DSC_STDDEV)
 907           && dv->stats[DSC_VARIANCE] != SYSMIS)
 908         dv->stats[DSC_STDDEV] = sqrt (dv->stats[DSC_VARIANCE]);
 909       if (dsc->calc_stats & (1ul << DSC_SEKURT))
 910         if (dv->stats[DSC_KURTOSIS] != SYSMIS)
 911             dv->stats[DSC_SEKURT] = calc_sekurt (W);
 912       if (dsc->calc_stats & (1ul << DSC_SESKEW)
 913           && dv->stats[DSC_SKEWNESS] != SYSMIS)
 914         dv->stats[DSC_SESKEW] = calc_seskew (W);
 915       dv->stats[DSC_RANGE] = ((dv->min == DBL_MAX || dv->max == -DBL_MAX)
 916                               ? SYSMIS : dv->max - dv->min);
 917       dv->stats[DSC_MIN] = dv->min == DBL_MAX ? SYSMIS : dv->min;
 918       dv->stats[DSC_MAX] = dv->max == -DBL_MAX ? SYSMIS : dv->max;
 919       if (dsc->calc_stats & (1ul << DSC_SUM))
 920         dv->stats[DSC_SUM] = W * dv->stats[DSC_MEAN];
 921
 922       if (dv->z_name)
 923         {
 924           case_data_rw_idx (c, z_idx++)->f = dv->stats[DSC_MEAN];
 925           case_data_rw_idx (c, z_idx++)->f = dv->stats[DSC_STDDEV];
 926         }
 927     }
 928
 929   if (c != NULL)
 930     casewriter_write (dsc->z_writer, c);
 931
 932   /* Output results. */
 933   display (dsc);
 934 }
 935
 936 /* Returns true if any of the descriptives variables in DSC's
 937    variable list have missing values in case C, false otherwise. */
 938 static bool
 939 listwise_missing (struct dsc_proc *dsc, const struct ccase *c)
 940 {
 941   size_t i;
 942
 943   for (i = 0; i < dsc->var_cnt; i++)
 944     {
 945       struct dsc_var *dv = &dsc->vars[i];
 946       double x = case_num (c, dv->v);
 947
 948       if (var_is_num_missing (dv->v, x, dsc->exclude))
 949         return true;
 950     }
 951   return false;
 952 }
 953 \f
 954 /* Statistical display. */
 955
 956 static algo_compare_func descriptives_compare_dsc_vars;
 957
 958 /* Displays a table of descriptive statistics for DSC. */
 959 static void
 960 display (struct dsc_proc *dsc)
 961 {
 962   size_t i;
 963   int nc;
 964   struct tab_table *t;
 965
 966   nc = 1 + (dsc->format == DSC_SERIAL ? 2 : 1);
 967   for (i = 0; i < DSC_N_STATS; i++)
 968     if (dsc->show_stats & (1ul << i))
 969       nc++;
 970
 971   if (dsc->sort_by_stat != DSC_NONE)
 972     sort (dsc->vars, dsc->var_cnt, sizeof *dsc->vars,
 973           descriptives_compare_dsc_vars, dsc);
 974
 975   t = tab_create (nc, dsc->var_cnt + 1);
 976   tab_headers (t, 1, 0, 1, 0);
 977   tab_box (t, TAL_1, TAL_1, -1, -1, 0, 0, nc - 1, dsc->var_cnt);
 978   tab_box (t, -1, -1, -1, TAL_1, 1, 0, nc - 1, dsc->var_cnt);
 979   tab_hline (t, TAL_2, 0, nc - 1, 1);
 980   tab_vline (t, TAL_2, 1, 0, dsc->var_cnt);
 981
 982   nc = 0;
 983   tab_text (t, nc++, 0, TAB_LEFT | TAT_TITLE, _("Variable"));
 984   if (dsc->format == DSC_SERIAL)
 985     {
 986       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Valid N"));
 987       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Missing N"));
 988     }
 989   else
 990     tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, "N");
 991
 992   for (i = 0; i < DSC_N_STATS; i++)
 993     if (dsc->show_stats & (1ul << i))
 994       {
 995         const char *title = gettext (dsc_info[i].name);
 996         tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, title);
 997       }
 998
 999   for (i = 0; i < dsc->var_cnt; i++)
1000     {
1001       struct dsc_var *dv = &dsc->vars[i];
1002       size_t j;
1003
1004       nc = 0;
1005       tab_text (t, nc++, i + 1, TAB_LEFT, var_to_string (dv->v));
1006       tab_text_format (t, nc++, i + 1, 0, "%.*g", DBL_DIG + 1, dv->valid);
1007       if (dsc->format == DSC_SERIAL)
1008         tab_text_format (t, nc++, i + 1, 0, "%.*g", DBL_DIG + 1, dv->missing);
1009
1010       for (j = 0; j < DSC_N_STATS; j++)
1011         if (dsc->show_stats & (1ul << j))
1012           tab_double (t, nc++, i + 1, TAB_NONE, dv->stats[j], NULL);
1013     }
1014
1015   tab_title (t, _("Valid cases = %.*g; cases with missing value(s) = %.*g."),
1016              DBL_DIG + 1, dsc->valid,
1017              DBL_DIG + 1, dsc->missing_listwise);
1018
1019   tab_submit (t);
1020 }
1021
1022 /* Compares `struct dsc_var's A and B according to the ordering
1023    specified by CMD. */
1024 static int
1025 descriptives_compare_dsc_vars (const void *a_, const void *b_, const void *dsc_)
1026 {
1027   const struct dsc_var *a = a_;
1028   const struct dsc_var *b = b_;
1029   const struct dsc_proc *dsc = dsc_;
1030
1031   int result;
1032
1033   if (dsc->sort_by_stat == DSC_NAME)
1034     result = utf8_strcasecmp (var_get_name (a->v), var_get_name (b->v));
1035   else
1036     {
1037       double as = a->stats[dsc->sort_by_stat];
1038       double bs = b->stats[dsc->sort_by_stat];
1039
1040       result = as < bs ? -1 : as > bs;
1041     }
1042
1043   if (!dsc->sort_ascending)
1044     result = -result;
1045
1046   return result;
1047 }