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

   1 /* PSPP - a program for statistical analysis.
   2    Copyright (C) 2011, 2012, 2013, 2019 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 "data/case.h"
  20 #include "data/format.h"
  21 #include "data/variable.h"
  22
  23 #include "libpspp/bt.h"
  24 #include "libpspp/hmap.h"
  25 #include "libpspp/misc.h"
  26 #include "libpspp/pool.h"
  27
  28 #include "math/moments.h"
  29
  30 #include <math.h>
  31
  32 #include "means.h"
  33
  34 #include "gettext.h"
  35 #define _(msgid) gettext (msgid)
  36 #define N_(msgid) (msgid)
  37
  38 struct per_var_data
  39 {
  40 };
  41
  42 struct per_var_data_simple
  43 {
  44   struct per_var_data parent;
  45   double acc;
  46 };
  47
  48 struct per_var_data_moment
  49 {
  50   struct per_var_data parent;
  51   struct moments1 *mom;
  52 };
  53
  54 static struct per_var_data *
  55 default_create (struct pool *pool)
  56 {
  57   struct per_var_data_moment *pvd = pool_alloc (pool, sizeof *pvd);
  58
  59   pvd->mom = moments1_create (MOMENT_KURTOSIS);
  60
  61   return (struct per_var_data *) pvd;
  62 }
  63
  64 static void
  65 default_update (struct per_var_data *stat, double w, double x)
  66 {
  67   struct per_var_data_moment *pvd = (struct per_var_data_moment *)stat;
  68
  69   moments1_add (pvd->mom, x, w);
  70 }
  71
  72 \f
  73
  74 /* HARMONIC MEAN: The reciprocal of the sum of the reciprocals:
  75    1 / ( 1/(x_0) + 1/(x_1) + ... + 1/(x_{n-1}) ) */
  76
  77 struct harmonic_mean
  78 {
  79   struct per_var_data parent;
  80   double rsum;
  81   double n;
  82 };
  83
  84 static struct per_var_data *
  85 harmonic_create (struct pool *pool)
  86 {
  87   struct harmonic_mean *hm = pool_alloc (pool, sizeof *hm);
  88
  89   hm->rsum = 0;
  90   hm->n = 0;
  91
  92   return (struct per_var_data *) hm;
  93 }
  94
  95
  96 static void
  97 harmonic_update (struct per_var_data *stat, double w, double x)
  98 {
  99   struct harmonic_mean *hm = (struct harmonic_mean *)stat;
 100   hm->rsum  += w / x;
 101   hm->n += w;
 102 }
 103
 104
 105 static double
 106 harmonic_get (const struct per_var_data *pvd)
 107 {
 108   const struct harmonic_mean *hm = (const struct harmonic_mean *) pvd;
 109
 110   return hm->n / hm->rsum;
 111 }
 112
 113 \f
 114
 115 /* GEOMETRIC MEAN:  The nth root of the product of all n observations
 116    pow ((x_0 * x_1 * ... x_{n - 1}), 1/n)  */
 117 struct geometric_mean
 118 {
 119   struct per_var_data parent;
 120   double prod;
 121   double n;
 122 };
 123
 124 static struct per_var_data *
 125 geometric_create (struct pool *pool)
 126 {
 127   struct geometric_mean *gm = pool_alloc (pool, sizeof *gm);
 128
 129   gm->prod = 1.0;
 130   gm->n = 0;
 131
 132   return (struct per_var_data *) gm;
 133 }
 134
 135 static void
 136 geometric_update (struct per_var_data  *pvd, double w, double x)
 137 {
 138   struct geometric_mean *gm = (struct geometric_mean *)pvd;
 139   gm->prod  *=  pow (x, w);
 140   gm->n += w;
 141 }
 142
 143
 144 static double
 145 geometric_get (const struct per_var_data *pvd)
 146 {
 147   const struct geometric_mean *gm = (const struct geometric_mean *)pvd;
 148   return pow (gm->prod, 1.0 / gm->n);
 149 }
 150
 151 \f
 152
 153 static double
 154 sum_get (const struct per_var_data *pvd)
 155 {
 156   double n, mean;
 157
 158   moments1_calculate (((struct per_var_data_moment *)pvd)->mom, &n, &mean, 0, 0, 0);
 159
 160   return mean * n;
 161 }
 162
 163
 164 static double
 165 n_get (const struct per_var_data *pvd)
 166 {
 167   double n;
 168
 169   moments1_calculate (((struct per_var_data_moment *)pvd)->mom, &n, 0, 0, 0, 0);
 170
 171   return n;
 172 }
 173
 174 static double
 175 arithmean_get (const struct per_var_data *pvd)
 176 {
 177   double n, mean;
 178
 179   moments1_calculate (((struct per_var_data_moment *)pvd)->mom, &n, &mean, 0, 0, 0);
 180
 181   return mean;
 182 }
 183
 184 static double
 185 variance_get (const struct per_var_data *pvd)
 186 {
 187   double n, mean, variance;
 188
 189   moments1_calculate (((struct per_var_data_moment *)pvd)->mom, &n, &mean, &variance, 0, 0);
 190
 191   return variance;
 192 }
 193
 194
 195 static double
 196 stddev_get (const struct per_var_data *pvd)
 197 {
 198   return sqrt (variance_get (pvd));
 199 }
 200
 201
 202 \f
 203
 204 static double
 205 skew_get (const struct per_var_data *pvd)
 206 {
 207   double skew;
 208
 209   moments1_calculate (((struct per_var_data_moment *)pvd)->mom, NULL, NULL, NULL, &skew, 0);
 210
 211   return skew;
 212 }
 213
 214 static double
 215 sekurt_get (const struct per_var_data *pvd)
 216 {
 217   double n;
 218
 219   moments1_calculate (((struct per_var_data_moment *)pvd)->mom, &n, NULL, NULL, NULL, NULL);
 220
 221   return calc_sekurt (n);
 222 }
 223
 224 static double
 225 seskew_get (const struct per_var_data *pvd)
 226 {
 227   double n;
 228
 229   moments1_calculate (((struct per_var_data_moment *)pvd)->mom, &n, NULL, NULL, NULL, NULL);
 230
 231   return calc_seskew (n);
 232 }
 233
 234 static double
 235 kurt_get (const struct per_var_data *pvd)
 236 {
 237   double kurt;
 238
 239   moments1_calculate (((struct per_var_data_moment *)pvd)->mom, NULL, NULL, NULL, NULL, &kurt);
 240
 241   return kurt;
 242 }
 243
 244 static double
 245 semean_get (const struct per_var_data *pvd)
 246 {
 247   double n, var;
 248
 249   moments1_calculate (((struct per_var_data_moment *)pvd)->mom, &n, NULL, &var, NULL, NULL);
 250
 251   return sqrt (var / n);
 252 }
 253
 254 \f
 255
 256 /* MIN: The smallest (closest to minus infinity) value. */
 257
 258 static struct per_var_data *
 259 min_create (struct pool *pool)
 260 {
 261   struct per_var_data_simple *pvd = pool_alloc (pool, sizeof *pvd);
 262
 263   pvd->acc = DBL_MAX;
 264
 265   return (struct per_var_data *) pvd;
 266 }
 267
 268 static void
 269 min_update (struct per_var_data *pvd, double w UNUSED, double x)
 270 {
 271   double *r = &((struct per_var_data_simple *)pvd)->acc;
 272
 273   if (x < *r)
 274     *r = x;
 275 }
 276
 277 static double
 278 min_get (const struct per_var_data *pvd)
 279 {
 280   double *r = &((struct per_var_data_simple *)pvd)->acc;
 281
 282   return *r;
 283 }
 284
 285 /* MAX: The largest (closest to plus infinity) value. */
 286
 287 static struct per_var_data *
 288 max_create (struct pool *pool)
 289 {
 290   struct per_var_data_simple *pvd = pool_alloc (pool, sizeof *pvd);
 291
 292   pvd->acc = -DBL_MAX;
 293
 294   return (struct per_var_data *) pvd;
 295 }
 296
 297 static void
 298 max_update (struct per_var_data *pvd, double w UNUSED, double x)
 299 {
 300   double *r = &((struct per_var_data_simple *)pvd)->acc;
 301
 302   if (x > *r)
 303     *r = x;
 304 }
 305
 306 static double
 307 max_get (const struct per_var_data *pvd)
 308 {
 309   double *r = &((struct per_var_data_simple *)pvd)->acc;
 310
 311   return *r;
 312 }
 313
 314 \f
 315
 316 struct range
 317 {
 318   struct per_var_data parent;
 319   double min;
 320   double max;
 321 };
 322
 323 static struct per_var_data *
 324 range_create (struct pool *pool)
 325 {
 326   struct range *r = pool_alloc (pool, sizeof *r);
 327
 328   r->min = DBL_MAX;
 329   r->max = -DBL_MAX;
 330
 331   return (struct per_var_data *) r;
 332 }
 333
 334 static void
 335 range_update (struct per_var_data *pvd, double w UNUSED, double x)
 336 {
 337   struct range *r = (struct range *) pvd;
 338
 339   if (x > r->max)
 340     r->max = x;
 341
 342   if (x < r->min)
 343     r->min = x;
 344 }
 345
 346 static double
 347 range_get (const struct per_var_data *pvd)
 348 {
 349   const struct range *r = (struct range *) pvd;
 350
 351   return r->max - r->min;
 352 }
 353
 354 \f
 355
 356 /* LAST: The last value (the one closest to the end of the file).  */
 357
 358 static struct per_var_data *
 359 last_create (struct pool *pool)
 360 {
 361   struct per_var_data_simple *pvd = pool_alloc (pool, sizeof *pvd);
 362
 363   return (struct per_var_data *) pvd;
 364 }
 365
 366 static void
 367 last_update (struct per_var_data *pvd, double w UNUSED, double x)
 368 {
 369   struct per_var_data_simple *stat = (struct per_var_data_simple *) pvd;
 370
 371   stat->acc = x;
 372 }
 373
 374 static double
 375 last_get (const struct per_var_data *pvd)
 376 {
 377   const struct per_var_data_simple *stat = (struct per_var_data_simple *) pvd;
 378
 379   return stat->acc;
 380 }
 381
 382 /* FIRST: The first value (the one closest to the start of the file).  */
 383
 384 static struct per_var_data *
 385 first_create (struct pool *pool)
 386 {
 387   struct per_var_data_simple *pvd = pool_alloc (pool, sizeof *pvd);
 388
 389   pvd->acc = SYSMIS;
 390
 391   return (struct per_var_data *) pvd;
 392 }
 393
 394 static void
 395 first_update (struct per_var_data *pvd, double w UNUSED, double x)
 396 {
 397   struct per_var_data_simple *stat = (struct per_var_data_simple *) pvd;
 398
 399   if (stat->acc == SYSMIS)
 400     stat->acc = x;
 401 }
 402
 403 static double
 404 first_get (const struct per_var_data *pvd)
 405 {
 406   const struct per_var_data_simple *stat = (struct per_var_data_simple *) pvd;
 407
 408   return stat->acc;
 409 }
 410
 411 /* Table of cell_specs */
 412 const struct cell_spec cell_spec[n_MEANS_STATISTICS] = {
 413   {N_("Mean"), "MEAN", default_create, default_update, arithmean_get},
 414   {N_("N"), "COUNT", default_create, default_update, n_get},
 415   {N_("Std. Deviation"), "STDDEV", default_create, default_update, stddev_get},
 416 #if 0
 417   {N_("Median"), "MEDIAN", default_create, default_update, NULL},
 418   {N_("Group Median"), "GMEDIAN", default_create, default_update, NULL},
 419 #endif
 420   {N_("S.E. Mean"), "SEMEAN", default_create, default_update, semean_get},
 421   {N_("Sum"), "SUM", default_create, default_update, sum_get},
 422   {N_("Minimum"), "MIN", min_create, min_update, min_get},
 423   {N_("Maximum"), "MAX", max_create, max_update, max_get},
 424   {N_("Range"), "RANGE", range_create, range_update, range_get},
 425   {N_("Variance"), "VARIANCE", default_create, default_update, variance_get},
 426   {N_("Kurtosis"), "KURT", default_create, default_update, kurt_get},
 427   {N_("S.E. Kurt"), "SEKURT", default_create, default_update, sekurt_get},
 428   {N_("Skewness"), "SKEW", default_create, default_update, skew_get},
 429   {N_("S.E. Skew"), "SESKEW", default_create, default_update, seskew_get},
 430   {N_("First"), "FIRST", first_create, first_update, first_get},
 431   {N_("Last"), "LAST", last_create, last_update, last_get},
 432 #if 0
 433   {N_("Percent N"), "NPCT", default_create, default_update, NULL},
 434   {N_("Percent Sum"), "SPCT", default_create, default_update, NULL},
 435 #endif
 436   {N_("Harmonic Mean"), "HARMONIC", harmonic_create, harmonic_update, harmonic_get},
 437   {N_("Geom. Mean"), "GEOMETRIC", geometric_create, geometric_update, geometric_get}
 438 };