pintos-os.org Git - pspp/blob - src/descript.c

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