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   return ok;
 693 }
 694
 695 /* Sets up a transformation to calculate Z scores. */
 696 static void
 697 setup_z_trns (struct dsc_proc *dsc, struct dataset *ds)
 698 {
 699   struct dsc_trns *t;
 700   size_t cnt, i;
 701
 702   for (cnt = i = 0; i < dsc->var_cnt; i++)
 703     if (dsc->vars[i].z_name != NULL)
 704       cnt++;
 705
 706   t = xmalloc (sizeof *t);
 707   t->z_scores = xnmalloc (cnt, sizeof *t->z_scores);
 708   t->z_score_cnt = cnt;
 709   t->missing_type = dsc->missing_type;
 710   t->exclude = dsc->exclude;
 711   if ( t->missing_type == DSC_LISTWISE )
 712     {
 713       t->var_cnt = dsc->var_cnt;
 714       t->vars = xnmalloc (t->var_cnt, sizeof *t->vars);
 715       for (i = 0; i < t->var_cnt; i++)
 716         t->vars[i] = dsc->vars[i].v;
 717     }
 718   else
 719     {
 720       t->var_cnt = 0;
 721       t->vars = NULL;
 722     }
 723   t->z_reader = casewriter_make_reader (dsc->z_writer);
 724   t->count = 0;
 725   t->ok = true;
 726   dsc->z_writer = NULL;
 727
 728   for (cnt = i = 0; i < dsc->var_cnt; i++)
 729     {
 730       struct dsc_var *dv = &dsc->vars[i];
 731       if (dv->z_name != NULL)
 732         {
 733           struct dsc_z_score *z;
 734           struct variable *dst_var;
 735
 736           dst_var = dict_create_var_assert (dataset_dict (ds), dv->z_name, 0);
 737           var_set_label (dst_var,
 738                          xasprintf (_("Z-score of %s"),var_to_string (dv->v)),
 739                          false);
 740
 741           z = &t->z_scores[cnt++];
 742           z->src_var = dv->v;
 743           z->z_var = dst_var;
 744         }
 745     }
 746
 747   add_transformation (ds,
 748                       descriptives_trns_proc, descriptives_trns_free, t);
 749 }
 750 \f
 751 /* Statistical calculation. */
 752
 753 static bool listwise_missing (struct dsc_proc *dsc, const struct ccase *c);
 754
 755 /* Calculates and displays descriptive statistics for the cases
 756    in CF. */
 757 static void
 758 calc_descriptives (struct dsc_proc *dsc, struct casereader *group,
 759                    struct dataset *ds)
 760 {
 761   struct casereader *pass1, *pass2;
 762   casenumber count;
 763   struct ccase *c;
 764   size_t z_idx;
 765   size_t i;
 766
 767   c = casereader_peek (group, 0);
 768   if (c == NULL)
 769     {
 770       casereader_destroy (group);
 771       return;
 772     }
 773   output_split_file_values (ds, c);
 774   case_unref (c);
 775
 776   group = casereader_create_filter_weight (group, dataset_dict (ds),
 777                                            NULL, NULL);
 778
 779   pass1 = group;
 780   pass2 = dsc->max_moment <= MOMENT_MEAN ? NULL : casereader_clone (pass1);
 781
 782   for (i = 0; i < dsc->var_cnt; i++)
 783     {
 784       struct dsc_var *dv = &dsc->vars[i];
 785
 786       dv->valid = dv->missing = 0.0;
 787       if (dv->moments != NULL)
 788         moments_clear (dv->moments);
 789       dv->min = DBL_MAX;
 790       dv->max = -DBL_MAX;
 791     }
 792   dsc->missing_listwise = 0.;
 793   dsc->valid = 0.;
 794
 795   /* First pass to handle most of the work. */
 796   count = 0;
 797   for (; (c = casereader_read (pass1)) != NULL; case_unref (c))
 798     {
 799       double weight = dict_get_case_weight (dataset_dict (ds), c, NULL);
 800
 801       /* Check for missing values. */
 802       if (listwise_missing (dsc, c))
 803         {
 804           dsc->missing_listwise += weight;
 805           if (dsc->missing_type == DSC_LISTWISE)
 806             continue;
 807         }
 808       dsc->valid += weight;
 809
 810       for (i = 0; i < dsc->var_cnt; i++)
 811         {
 812           struct dsc_var *dv = &dsc->vars[i];
 813           double x = case_num (c, dv->v);
 814
 815           if (var_is_num_missing (dv->v, x, dsc->exclude))
 816             {
 817               dv->missing += weight;
 818               continue;
 819             }
 820
 821           if (dv->moments != NULL)
 822             moments_pass_one (dv->moments, x, weight);
 823
 824           if (x < dv->min)
 825             dv->min = x;
 826           if (x > dv->max)
 827             dv->max = x;
 828         }
 829
 830       count++;
 831     }
 832   if (!casereader_destroy (pass1))
 833     {
 834       casereader_destroy (pass2);
 835       return;
 836     }
 837
 838   /* Second pass for higher-order moments. */
 839   if (dsc->max_moment > MOMENT_MEAN)
 840     {
 841       for (; (c = casereader_read (pass2)) != NULL; case_unref (c))
 842         {
 843           double weight = dict_get_case_weight (dataset_dict (ds), c, NULL);
 844
 845           /* Check for missing values. */
 846           if (dsc->missing_type == DSC_LISTWISE && listwise_missing (dsc, c))
 847             continue;
 848
 849           for (i = 0; i < dsc->var_cnt; i++)
 850             {
 851               struct dsc_var *dv = &dsc->vars[i];
 852               double x = case_num (c, dv->v);
 853
 854               if (var_is_num_missing (dv->v, x, dsc->exclude))
 855                 continue;
 856
 857               if (dv->moments != NULL)
 858                 moments_pass_two (dv->moments, x, weight);
 859             }
 860         }
 861       if (!casereader_destroy (pass2))
 862         return;
 863     }
 864
 865   /* Calculate results. */
 866   if (dsc->z_writer)
 867     {
 868       c = case_create (casewriter_get_proto (dsc->z_writer));
 869       z_idx = 0;
 870       case_data_rw_idx (c, z_idx++)->f = count;
 871     }
 872   else
 873     c = NULL;
 874
 875   for (i = 0; i < dsc->var_cnt; i++)
 876     {
 877       struct dsc_var *dv = &dsc->vars[i];
 878       double W;
 879       int j;
 880
 881       for (j = 0; j < DSC_N_STATS; j++)
 882         dv->stats[j] = SYSMIS;
 883
 884       dv->valid = W = dsc->valid - dv->missing;
 885
 886       if (dv->moments != NULL)
 887         moments_calculate (dv->moments, NULL,
 888                            &dv->stats[DSC_MEAN], &dv->stats[DSC_VARIANCE],
 889                            &dv->stats[DSC_SKEWNESS], &dv->stats[DSC_KURTOSIS]);
 890       if (dsc->calc_stats & (1ul << DSC_SEMEAN)
 891           && dv->stats[DSC_VARIANCE] != SYSMIS && W > 0.)
 892         dv->stats[DSC_SEMEAN] = sqrt (dv->stats[DSC_VARIANCE]) / sqrt (W);
 893       if (dsc->calc_stats & (1ul << DSC_STDDEV)
 894           && dv->stats[DSC_VARIANCE] != SYSMIS)
 895         dv->stats[DSC_STDDEV] = sqrt (dv->stats[DSC_VARIANCE]);
 896       if (dsc->calc_stats & (1ul << DSC_SEKURT))
 897         if (dv->stats[DSC_KURTOSIS] != SYSMIS)
 898             dv->stats[DSC_SEKURT] = calc_sekurt (W);
 899       if (dsc->calc_stats & (1ul << DSC_SESKEW)
 900           && dv->stats[DSC_SKEWNESS] != SYSMIS)
 901         dv->stats[DSC_SESKEW] = calc_seskew (W);
 902       dv->stats[DSC_RANGE] = ((dv->min == DBL_MAX || dv->max == -DBL_MAX)
 903                               ? SYSMIS : dv->max - dv->min);
 904       dv->stats[DSC_MIN] = dv->min == DBL_MAX ? SYSMIS : dv->min;
 905       dv->stats[DSC_MAX] = dv->max == -DBL_MAX ? SYSMIS : dv->max;
 906       if (dsc->calc_stats & (1ul << DSC_SUM))
 907         dv->stats[DSC_SUM] = W * dv->stats[DSC_MEAN];
 908
 909       if (dv->z_name)
 910         {
 911           case_data_rw_idx (c, z_idx++)->f = dv->stats[DSC_MEAN];
 912           case_data_rw_idx (c, z_idx++)->f = dv->stats[DSC_STDDEV];
 913         }
 914     }
 915
 916   if (c != NULL)
 917     casewriter_write (dsc->z_writer, c);
 918
 919   /* Output results. */
 920   display (dsc);
 921 }
 922
 923 /* Returns true if any of the descriptives variables in DSC's
 924    variable list have missing values in case C, false otherwise. */
 925 static bool
 926 listwise_missing (struct dsc_proc *dsc, const struct ccase *c)
 927 {
 928   size_t i;
 929
 930   for (i = 0; i < dsc->var_cnt; i++)
 931     {
 932       struct dsc_var *dv = &dsc->vars[i];
 933       double x = case_num (c, dv->v);
 934
 935       if (var_is_num_missing (dv->v, x, dsc->exclude))
 936         return true;
 937     }
 938   return false;
 939 }
 940 \f
 941 /* Statistical display. */
 942
 943 static algo_compare_func descriptives_compare_dsc_vars;
 944
 945 /* Displays a table of descriptive statistics for DSC. */
 946 static void
 947 display (struct dsc_proc *dsc)
 948 {
 949   size_t i;
 950   int nc;
 951   struct tab_table *t;
 952
 953   nc = 1 + (dsc->format == DSC_SERIAL ? 2 : 1);
 954   for (i = 0; i < DSC_N_STATS; i++)
 955     if (dsc->show_stats & (1ul << i))
 956       nc++;
 957
 958   if (dsc->sort_by_stat != DSC_NONE)
 959     sort (dsc->vars, dsc->var_cnt, sizeof *dsc->vars,
 960           descriptives_compare_dsc_vars, dsc);
 961
 962   t = tab_create (nc, dsc->var_cnt + 1);
 963   tab_headers (t, 1, 0, 1, 0);
 964   tab_box (t, TAL_1, TAL_1, -1, -1, 0, 0, nc - 1, dsc->var_cnt);
 965   tab_box (t, -1, -1, -1, TAL_1, 1, 0, nc - 1, dsc->var_cnt);
 966   tab_hline (t, TAL_2, 0, nc - 1, 1);
 967   tab_vline (t, TAL_2, 1, 0, dsc->var_cnt);
 968
 969   nc = 0;
 970   tab_text (t, nc++, 0, TAB_LEFT | TAT_TITLE, _("Variable"));
 971   if (dsc->format == DSC_SERIAL)
 972     {
 973       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Valid N"));
 974       tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, _("Missing N"));
 975     }
 976   else
 977     tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, "N");
 978
 979   for (i = 0; i < DSC_N_STATS; i++)
 980     if (dsc->show_stats & (1ul << i))
 981       {
 982         const char *title = gettext (dsc_info[i].name);
 983         tab_text (t, nc++, 0, TAB_CENTER | TAT_TITLE, title);
 984       }
 985
 986   for (i = 0; i < dsc->var_cnt; i++)
 987     {
 988       struct dsc_var *dv = &dsc->vars[i];
 989       size_t j;
 990
 991       nc = 0;
 992       tab_text (t, nc++, i + 1, TAB_LEFT, var_to_string (dv->v));
 993       tab_text_format (t, nc++, i + 1, 0, "%g", dv->valid);
 994       if (dsc->format == DSC_SERIAL)
 995         tab_text_format (t, nc++, i + 1, 0, "%g", dv->missing);
 996
 997       for (j = 0; j < DSC_N_STATS; j++)
 998         if (dsc->show_stats & (1ul << j))
 999           tab_double (t, nc++, i + 1, TAB_NONE, dv->stats[j], NULL);
1000     }
1001
1002   tab_title (t, _("Valid cases = %g; cases with missing value(s) = %g."),
1003              dsc->valid, dsc->missing_listwise);
1004
1005   tab_submit (t);
1006 }
1007
1008 /* Compares `struct dsc_var's A and B according to the ordering
1009    specified by CMD. */
1010 static int
1011 descriptives_compare_dsc_vars (const void *a_, const void *b_, const void *dsc_)
1012 {
1013   const struct dsc_var *a = a_;
1014   const struct dsc_var *b = b_;
1015   const struct dsc_proc *dsc = dsc_;
1016
1017   int result;
1018
1019   if (dsc->sort_by_stat == DSC_NAME)
1020     result = utf8_strcasecmp (var_get_name (a->v), var_get_name (b->v));
1021   else
1022     {
1023       double as = a->stats[dsc->sort_by_stat];
1024       double bs = b->stats[dsc->sort_by_stat];
1025
1026       result = as < bs ? -1 : as > bs;
1027     }
1028
1029   if (!dsc->sort_ascending)
1030     result = -result;
1031
1032   return result;
1033 }