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

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