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

   1 /* PSPP - a program for statistical analysis.
   2    Copyright (C) 1997-2000, 2009-2012 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       proto = caseproto_create ();
 401       for (i = 0; i < 1 + 2 * z_cnt; i++)
 402         proto = caseproto_add_width (proto, 0);
 403       dsc->z_writer = autopaging_writer_create (proto);
 404       caseproto_unref (proto);
 405
 406       dump_z_table (dsc);
 407     }
 408
 409   /* Figure out statistics to display. */
 410   if (dsc->show_stats & (1ul << DSC_SKEWNESS))
 411     dsc->show_stats |= 1ul << DSC_SESKEW;
 412   if (dsc->show_stats & (1ul << DSC_KURTOSIS))
 413     dsc->show_stats |= 1ul << DSC_SEKURT;
 414
 415   /* Figure out which statistics to calculate. */
 416   dsc->calc_stats = dsc->show_stats;
 417   if (z_cnt > 0)
 418     dsc->calc_stats |= (1ul << DSC_MEAN) | (1ul << DSC_STDDEV);
 419   if (dsc->sort_by_stat >= 0)
 420     dsc->calc_stats |= 1ul << dsc->sort_by_stat;
 421   if (dsc->show_stats & (1ul << DSC_SESKEW))
 422     dsc->calc_stats |= 1ul << DSC_SKEWNESS;
 423   if (dsc->show_stats & (1ul << DSC_SEKURT))
 424     dsc->calc_stats |= 1ul << DSC_KURTOSIS;
 425
 426   /* Figure out maximum moment needed and allocate moments for
 427      the variables. */
 428   dsc->max_moment = MOMENT_NONE;
 429   for (i = 0; i < DSC_N_STATS; i++)
 430     if (dsc->calc_stats & (1ul << i) && dsc_info[i].moment > dsc->max_moment)
 431       dsc->max_moment = dsc_info[i].moment;
 432   if (dsc->max_moment != MOMENT_NONE)
 433     for (i = 0; i < dsc->var_cnt; i++)
 434       dsc->vars[i].moments = moments_create (dsc->max_moment);
 435
 436   /* Data pass. */
 437   grouper = casegrouper_create_splits (proc_open (ds), dict);
 438   while (casegrouper_get_next_group (grouper, &group))
 439     calc_descriptives (dsc, group, ds);
 440   ok = casegrouper_destroy (grouper);
 441   ok = proc_commit (ds) && ok;
 442
 443   /* Z-scoring! */
 444   if (ok && z_cnt)
 445     setup_z_trns (dsc, ds);
 446
 447   /* Done. */
 448   free (vars);
 449   free_dsc_proc (dsc);
 450   return ok ? CMD_SUCCESS : CMD_CASCADING_FAILURE;
 451
 452  error:
 453   free (vars);
 454   free_dsc_proc (dsc);
 455   return CMD_FAILURE;
 456 }
 457
 458 /* Returns the statistic named by the current token and skips past the token.
 459    Returns DSC_NONE if no statistic is given (e.g., subcommand with no
 460    specifiers). Emits an error if the current token ID does not name a
 461    statistic. */
 462 static enum dsc_statistic
 463 match_statistic (struct lexer *lexer)
 464 {
 465   if (lex_token (lexer) == T_ID)
 466     {
 467       enum dsc_statistic stat;
 468
 469       for (stat = 0; stat < DSC_N_STATS; stat++)
 470         if (lex_match_id (lexer, dsc_info[stat].identifier))
 471           return stat;
 472
 473       lex_get (lexer);
 474       lex_error (lexer, _("expecting statistic name: reverting to default"));
 475     }
 476
 477   return DSC_NONE;
 478 }
 479
 480 /* Frees DSC. */
 481 static void
 482 free_dsc_proc (struct dsc_proc *dsc)
 483 {
 484   size_t i;
 485
 486   if (dsc == NULL)
 487     return;
 488
 489   for (i = 0; i < dsc->var_cnt; i++)
 490     {
 491       struct dsc_var *dsc_var = &dsc->vars[i];
 492       free (dsc_var->z_name);
 493       moments_destroy (dsc_var->moments);
 494     }
 495   casewriter_destroy (dsc->z_writer);
 496   free (dsc->vars);
 497   free (dsc);
 498 }
 499 \f
 500 /* Z scores. */
 501
 502 /* Returns false if NAME is a duplicate of any existing variable name or
 503    of any previously-declared z-var name; otherwise returns true. */
 504 static bool
 505 try_name (const struct dictionary *dict, struct dsc_proc *dsc,
 506           const char *name)
 507 {
 508   size_t i;
 509
 510   if (dict_lookup_var (dict, name) != NULL)
 511     return false;
 512   for (i = 0; i < dsc->var_cnt; i++)
 513     {
 514       struct dsc_var *dsc_var = &dsc->vars[i];
 515       if (dsc_var->z_name != NULL && !utf8_strcasecmp (dsc_var->z_name, name))
 516         return false;
 517     }
 518   return true;
 519 }
 520
 521 /* Generates a name for a Z-score variable based on a variable
 522    named VAR_NAME, given that *Z_CNT generated variable names are
 523    known to already exist.  If successful, returns the new name
 524    as a dynamically allocated string.  On failure, returns NULL. */
 525 static char *
 526 generate_z_varname (const struct dictionary *dict, struct dsc_proc *dsc,
 527                     const char *var_name, int *z_cnt)
 528 {
 529   char *z_name, *trunc_name;
 530
 531   /* Try a name based on the original variable name. */
 532   z_name = xasprintf ("Z%s", var_name);
 533   trunc_name = utf8_encoding_trunc (z_name, dict_get_encoding (dict),
 534                                     ID_MAX_LEN);
 535   free (z_name);
 536   if (try_name (dict, dsc, trunc_name))
 537     return trunc_name;
 538   free (trunc_name);
 539
 540   /* Generate a synthetic name. */
 541   for (;;)
 542     {
 543       char name[8];
 544
 545       (*z_cnt)++;
 546
 547       if (*z_cnt <= 99)
 548         sprintf (name, "ZSC%03d", *z_cnt);
 549       else if (*z_cnt <= 108)
 550         sprintf (name, "STDZ%02d", *z_cnt - 99);
 551       else if (*z_cnt <= 117)
 552         sprintf (name, "ZZZZ%02d", *z_cnt - 108);
 553       else if (*z_cnt <= 126)
 554         sprintf (name, "ZQZQ%02d", *z_cnt - 117);
 555       else
 556         {
 557           msg (SE, _("Ran out of generic names for Z-score variables.  "
 558                      "There are only 126 generic names: ZSC001-ZSC0999, "
 559                      "STDZ01-STDZ09, ZZZZ01-ZZZZ09, ZQZQ01-ZQZQ09."));
 560           return NULL;
 561         }
 562
 563       if (try_name (dict, dsc, name))
 564         return xstrdup (name);
 565     }
 566   NOT_REACHED();
 567 }
 568
 569 /* Outputs a table describing the mapping between source
 570    variables and Z-score variables. */
 571 static void
 572 dump_z_table (struct dsc_proc *dsc)
 573 {
 574   size_t cnt = 0;
 575   struct tab_table *t;
 576
 577   {
 578     size_t i;
 579
 580     for (i = 0; i < dsc->var_cnt; i++)
 581       if (dsc->vars[i].z_name != NULL)
 582         cnt++;
 583   }
 584
 585   t = tab_create (2, cnt + 1);
 586   tab_title (t, _("Mapping of variables to corresponding Z-scores."));
 587   tab_headers (t, 0, 0, 1, 0);
 588   tab_box (t, TAL_1, TAL_1, TAL_0, TAL_1, 0, 0, 1, cnt);
 589   tab_hline (t, TAL_2, 0, 1, 1);
 590   tab_text (t, 0, 0, TAB_CENTER | TAT_TITLE, _("Source"));
 591   tab_text (t, 1, 0, TAB_CENTER | TAT_TITLE, _("Target"));
 592
 593   {
 594     size_t i, y;
 595
 596     for (i = 0, y = 1; i < dsc->var_cnt; i++)
 597       if (dsc->vars[i].z_name != NULL)
 598         {
 599           tab_text (t, 0, y, TAB_LEFT, var_to_string (dsc->vars[i].v));
 600           tab_text (t, 1, y++, TAB_LEFT, dsc->vars[i].z_name);
 601         }
 602   }
 603
 604   tab_submit (t);
 605 }
 606
 607 /* Transformation function to calculate Z-scores. Will return SYSMIS if any of
 608    the following are true: 1) mean or standard deviation is SYSMIS 2) score is
 609    SYSMIS 3) score is user missing and they were not included in the original
 610    analyis. 4) any of the variables in the original analysis were missing
 611    (either system or user-missing values that weren't included).
 612 */
 613 static int
 614 descriptives_trns_proc (void *trns_, struct ccase **c,
 615                         casenumber case_idx UNUSED)
 616 {
 617   struct dsc_trns *t = trns_;
 618   struct dsc_z_score *z;
 619   const struct variable **vars;
 620   int all_sysmis = 0;
 621
 622   if (t->count <= 0)
 623     {
 624       struct ccase *z_case;
 625
 626       z_case = casereader_read (t->z_reader);
 627       if (z_case)
 628         {
 629           size_t z_idx = 0;
 630
 631           t->count = case_num_idx (z_case, z_idx++);
 632           for (z = t->z_scores; z < t->z_scores + t->z_score_cnt; z++)
 633             {
 634               z->mean = case_num_idx (z_case, z_idx++);
 635               z->std_dev = case_num_idx (z_case, z_idx++);
 636             }
 637           case_unref (z_case);
 638         }
 639       else
 640         {
 641           if (t->ok)
 642             {
 643               msg (SE, _("Internal error processing Z scores"));
 644               t->ok = false;
 645             }
 646           for (z = t->z_scores; z < t->z_scores + t->z_score_cnt; z++)
 647             z->mean = z->std_dev = SYSMIS;
 648         }
 649     }
 650   t->count--;
 651
 652   if (t->missing_type == DSC_LISTWISE)
 653     {
 654       assert(t->vars);
 655       for (vars = t->vars; vars < t->vars + t->var_cnt; vars++)
 656         {
 657           double score = case_num (*c, *vars);
 658           if (var_is_num_missing (*vars, score, t->exclude))
 659             {
 660               all_sysmis = 1;
 661               break;
 662             }
 663         }
 664     }
 665
 666   *c = case_unshare (*c);
 667   for (z = t->z_scores; z < t->z_scores + t->z_score_cnt; z++)
 668     {
 669       double input = case_num (*c, z->src_var);
 670       double *output = &case_data_rw (*c, z->z_var)->f;
 671
 672       if (z->mean == SYSMIS || z->std_dev == SYSMIS || all_sysmis
 673           || var_is_num_missing (z->src_var, input, t->exclude))
 674         *output = SYSMIS;
 675       else
 676         *output = (input - z->mean) / z->std_dev;
 677     }
 678   return TRNS_CONTINUE;
 679 }
 680
 681 /* Frees a descriptives_trns struct. */
 682 static bool
 683 descriptives_trns_free (void *trns_)
 684 {
 685   struct dsc_trns *t = trns_;
 686   bool ok = t->ok && !casereader_error (t->z_reader);
 687
 688   free (t->z_scores);
 689   casereader_destroy (t->z_reader);
 690   assert((t->missing_type != DSC_LISTWISE) ^ (t->vars != NULL));
 691   free (t->vars);
 692   free (t);
 693
 694   return ok;
 695 }
 696
 697 /* Sets up a transformation to calculate Z scores. */
 698 static void
 699 setup_z_trns (struct dsc_proc *dsc, struct dataset *ds)
 700 {
 701   struct dsc_trns *t;
 702   size_t cnt, i;
 703
 704   for (cnt = i = 0; i < dsc->var_cnt; i++)
 705     if (dsc->vars[i].z_name != NULL)
 706       cnt++;
 707
 708   t = xmalloc (sizeof *t);
 709   t->z_scores = xnmalloc (cnt, sizeof *t->z_scores);
 710   t->z_score_cnt = cnt;
 711   t->missing_type = dsc->missing_type;
 712   t->exclude = dsc->exclude;
 713   if ( t->missing_type == DSC_LISTWISE )
 714     {
 715       t->var_cnt = dsc->var_cnt;
 716       t->vars = xnmalloc (t->var_cnt, sizeof *t->vars);
 717       for (i = 0; i < t->var_cnt; i++)
 718         t->vars[i] = dsc->vars[i].v;
 719     }
 720   else
 721     {
 722       t->var_cnt = 0;
 723       t->vars = NULL;
 724     }
 725   t->z_reader = casewriter_make_reader (dsc->z_writer);
 726   t->count = 0;
 727   t->ok = true;
 728   dsc->z_writer = NULL;
 729
 730   for (cnt = i = 0; i < dsc->var_cnt; i++)
 731     {
 732       struct dsc_var *dv = &dsc->vars[i];
 733       if (dv->z_name != NULL)
 734         {
 735           struct dsc_z_score *z;
 736           struct variable *dst_var;
 737           char *label;
 738
 739           dst_var = dict_create_var_assert (dataset_dict (ds), dv->z_name, 0);
 740
 741           label = xasprintf (_("Z-score of %s"),var_to_string (dv->v));
 742           var_set_label (dst_var, label, false);
 743           free (label);
 744
 745           z = &t->z_scores[cnt++];
 746           z->src_var = dv->v;
 747           z->z_var = dst_var;
 748         }
 749     }
 750
 751   add_transformation (ds,
 752                       descriptives_trns_proc, descriptives_trns_free, t);
 753 }
 754 \f
 755 /* Statistical calculation. */
 756
 757 static bool listwise_missing (struct dsc_proc *dsc, const struct ccase *c);
 758
 759 /* Calculates and displays descriptive statistics for the cases
 760    in CF. */
 761 static void
 762 calc_descriptives (struct dsc_proc *dsc, struct casereader *group,
 763                    struct dataset *ds)
 764 {
 765   struct casereader *pass1, *pass2;
 766   casenumber count;
 767   struct ccase *c;
 768   size_t z_idx;
 769   size_t i;
 770
 771   c = casereader_peek (group, 0);
 772   if (c == NULL)
 773     {
 774       casereader_destroy (group);
 775       return;
 776     }
 777   output_split_file_values (ds, c);
 778   case_unref (c);
 779
 780   group = casereader_create_filter_weight (group, dataset_dict (ds),
 781                                            NULL, NULL);
 782
 783   pass1 = group;
 784   pass2 = dsc->max_moment <= MOMENT_MEAN ? NULL : casereader_clone (pass1);
 785
 786   for (i = 0; i < dsc->var_cnt; i++)
 787     {
 788       struct dsc_var *dv = &dsc->vars[i];
 789
 790       dv->valid = dv->missing = 0.0;
 791       if (dv->moments != NULL)
 792         moments_clear (dv->moments);
 793       dv->min = DBL_MAX;
 794       dv->max = -DBL_MAX;
 795     }
 796   dsc->missing_listwise = 0.;
 797   dsc->valid = 0.;
 798
 799   /* First pass to handle most of the work. */
 800   count = 0;
 801   for (; (c = casereader_read (pass1)) != NULL; case_unref (c))
 802     {
 803       double weight = dict_get_case_weight (dataset_dict (ds), c, NULL);
 804
 805       /* Check for missing values. */
 806       if (listwise_missing (dsc, c))
 807         {
 808           dsc->missing_listwise += weight;
 809           if (dsc->missing_type == DSC_LISTWISE)
 810             continue;
 811         }
 812       dsc->valid += weight;
 813
 814       for (i = 0; i < dsc->var_cnt; i++)
 815         {
 816           struct dsc_var *dv = &dsc->vars[i];
 817           double x = case_num (c, dv->v);
 818
 819           if (var_is_num_missing (dv->v, x, dsc->exclude))
 820             {
 821               dv->missing += weight;
 822               continue;
 823             }
 824
 825           if (dv->moments != NULL)
 826             moments_pass_one (dv->moments, x, weight);
 827
 828           if (x < dv->min)
 829             dv->min = x;
 830           if (x > dv->max)
 831             dv->max = x;
 832         }
 833
 834       count++;
 835     }
 836   if (!casereader_destroy (pass1))
 837     {
 838       casereader_destroy (pass2);
 839       return;
 840     }
 841
 842   /* Second pass for higher-order moments. */
 843   if (dsc->max_moment > MOMENT_MEAN)
 844     {
 845       for (; (c = casereader_read (pass2)) != NULL; case_unref (c))
 846         {
 847           double weight = dict_get_case_weight (dataset_dict (ds), c, NULL);
 848
 849           /* Check for missing values. */
 850           if (dsc->missing_type == DSC_LISTWISE && listwise_missing (dsc, c))
 851             continue;
 852
 853           for (i = 0; i < dsc->var_cnt; i++)
 854             {
 855               struct dsc_var *dv = &dsc->vars[i];
 856               double x = case_num (c, dv->v);
 857
 858               if (var_is_num_missing (dv->v, x, dsc->exclude))
 859                 continue;
 860
 861               if (dv->moments != NULL)
 862                 moments_pass_two (dv->moments, x, weight);
 863             }
 864         }
 865       if (!casereader_destroy (pass2))
 866         return;
 867     }
 868
 869   /* Calculate results. */
 870   if (dsc->z_writer)
 871     {
 872       c = case_create (casewriter_get_proto (dsc->z_writer));
 873       z_idx = 0;
 874       case_data_rw_idx (c, z_idx++)->f = count;
 875     }
 876   else
 877     c = NULL;
 878
 879   for (i = 0; i < dsc->var_cnt; i++)
 880     {
 881       struct dsc_var *dv = &dsc->vars[i];
 882       double W;
 883       int j;
 884
 885       for (j = 0; j < DSC_N_STATS; j++)
 886         dv->stats[j] = SYSMIS;
 887
 888       dv->valid = W = dsc->valid - dv->missing;
 889
 890       if (dv->moments != NULL)
 891         moments_calculate (dv->moments, NULL,
 892                            &dv->stats[DSC_MEAN], &dv->stats[DSC_VARIANCE],
 893                            &dv->stats[DSC_SKEWNESS], &dv->stats[DSC_KURTOSIS]);
 894       if (dsc->calc_stats & (1ul << DSC_SEMEAN)
 895           && dv->stats[DSC_VARIANCE] != SYSMIS && W > 0.)
 896         dv->stats[DSC_SEMEAN] = sqrt (dv->stats[DSC_VARIANCE]) / sqrt (W);
 897       if (dsc->calc_stats & (1ul << DSC_STDDEV)
 898           && dv->stats[DSC_VARIANCE] != SYSMIS)
 899         dv->stats[DSC_STDDEV] = sqrt (dv->stats[DSC_VARIANCE]);
 900       if (dsc->calc_stats & (1ul << DSC_SEKURT))
 901         if (dv->stats[DSC_KURTOSIS] != SYSMIS)
 902             dv->stats[DSC_SEKURT] = calc_sekurt (W);
 903       if (dsc->calc_stats & (1ul << DSC_SESKEW)
 904           && dv->stats[DSC_SKEWNESS] != SYSMIS)
 905         dv->stats[DSC_SESKEW] = calc_seskew (W);
 906       dv->stats[DSC_RANGE] = ((dv->min == DBL_MAX || dv->max == -DBL_MAX)
 907                               ? SYSMIS : dv->max - dv->min);
 908       dv->stats[DSC_MIN] = dv->min == DBL_MAX ? SYSMIS : dv->min;
 909       dv->stats[DSC_MAX] = dv->max == -DBL_MAX ? SYSMIS : dv->max;
 910       if (dsc->calc_stats & (1ul << DSC_SUM))
 911         dv->stats[DSC_SUM] = W * dv->stats[DSC_MEAN];
 912
 913       if (dv->z_name)
 914         {
 915           case_data_rw_idx (c, z_idx++)->f = dv->stats[DSC_MEAN];
 916           case_data_rw_idx (c, z_idx++)->f = dv->stats[DSC_STDDEV];
 917         }
 918     }
 919
 920   if (c != NULL)
 921     casewriter_write (dsc->z_writer, c);
 922
 923   /* Output results. */
 924   display (dsc);
 925 }
 926
 927 /* Returns true if any of the descriptives variables in DSC's
 928    variable list have missing values in case C, false otherwise. */
 929 static bool
 930 listwise_missing (struct dsc_proc *dsc, const struct ccase *c)
 931 {
 932   size_t i;
 933
 934   for (i = 0; i < dsc->var_cnt; i++)
 935     {
 936       struct dsc_var *dv = &dsc->vars[i];
 937       double x = case_num (c, dv->v);
 938
 939       if (var_is_num_missing (dv->v, x, dsc->exclude))
 940         return true;
 941     }
 942   return false;
 943 }
 944 \f
 945 /* Statistical display. */
 946
 947 static algo_compare_func descriptives_compare_dsc_vars;
 948
 949 /* Displays a table of descriptive statistics for DSC. */
 950 static void
 951 display (struct dsc_proc *dsc)
 952 {
 953   size_t i;
 954   int nc;
 955   struct tab_table *t;
 956
 957   nc = 1 + (dsc->format == DSC_SERIAL ? 2 : 1);
 958   for (i = 0; i < DSC_N_STATS; i++)
 959     if (dsc->show_stats & (1ul << i))
 960       nc++;
 961
 962   if (dsc->sort_by_stat != DSC_NONE)
 963     sort (dsc->vars, dsc->var_cnt, sizeof *dsc->vars,
 964           descriptives_compare_dsc_vars, dsc);
 965
 966   t = tab_create (nc, dsc->var_cnt + 1);
 967   tab_headers (t, 1, 0, 1, 0);
 968   tab_box (t, TAL_1, TAL_1, -1, -1, 0, 0, nc - 1, dsc->var_cnt);
 969   tab_box (t, -1, -1, -1, TAL_1, 1, 0, nc - 1, dsc->var_cnt);
 970   tab_hline (t, TAL_2, 0, nc - 1, 1);
 971   tab_vline (t, TAL_2, 1, 0, dsc->var_cnt);
 972
 973   nc = 0;
 974   tab_text (t, nc++, 0, TAB_LEFT | TAT_TITLE, _("Variable"));
 975   if (dsc->format == DSC_SERIAL)
 976     {
 977       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Valid N"));
 978       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Missing N"));
 979     }
 980   else
 981     tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, "N");
 982
 983   for (i = 0; i < DSC_N_STATS; i++)
 984     if (dsc->show_stats & (1ul << i))
 985       {
 986         const char *title = gettext (dsc_info[i].name);
 987         tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, title);
 988       }
 989
 990   for (i = 0; i < dsc->var_cnt; i++)
 991     {
 992       struct dsc_var *dv = &dsc->vars[i];
 993       size_t j;
 994
 995       nc = 0;
 996       tab_text (t, nc++, i + 1, TAB_LEFT, var_to_string (dv->v));
 997       tab_text_format (t, nc++, i + 1, 0, "%g", dv->valid);
 998       if (dsc->format == DSC_SERIAL)
 999         tab_text_format (t, nc++, i + 1, 0, "%g", dv->missing);
1000
1001       for (j = 0; j < DSC_N_STATS; j++)
1002         if (dsc->show_stats & (1ul << j))
1003           tab_double (t, nc++, i + 1, TAB_NONE, dv->stats[j], NULL);
1004     }
1005
1006   tab_title (t, _("Valid cases = %g; cases with missing value(s) = %g."),
1007              dsc->valid, dsc->missing_listwise);
1008
1009   tab_submit (t);
1010 }
1011
1012 /* Compares `struct dsc_var's A and B according to the ordering
1013    specified by CMD. */
1014 static int
1015 descriptives_compare_dsc_vars (const void *a_, const void *b_, const void *dsc_)
1016 {
1017   const struct dsc_var *a = a_;
1018   const struct dsc_var *b = b_;
1019   const struct dsc_proc *dsc = dsc_;
1020
1021   int result;
1022
1023   if (dsc->sort_by_stat == DSC_NAME)
1024     result = utf8_strcasecmp (var_get_name (a->v), var_get_name (b->v));
1025   else
1026     {
1027       double as = a->stats[dsc->sort_by_stat];
1028       double bs = b->stats[dsc->sort_by_stat];
1029
1030       result = as < bs ? -1 : as > bs;
1031     }
1032
1033   if (!dsc->sort_ascending)
1034     result = -result;
1035
1036   return result;
1037 }