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