X-Git-Url: https://pintos-os.org/cgi-bin/gitweb.cgi?a=blobdiff_plain;f=src%2Fmath%2Fmoments.c;h=517219d54fcf5d3cc159b4b64c8f6454f64a5bd4;hb=4504e04f97804a96914ba12d9e450ddda49a7bab;hp=02208bca70d007d1bc9f61f07d518b9b985a621d;hpb=a9afcdd223cc1c1b7889601e405cff0bd47300d7;p=pspp diff --git a/src/math/moments.c b/src/math/moments.c index 02208bca70..517219d54f 100644 --- a/src/math/moments.c +++ b/src/math/moments.c @@ -1,5 +1,5 @@ /* PSPP - a program for statistical analysis. - Copyright (C) 1997-9, 2000 Free Software Foundation, Inc. + Copyright (C) 1997-9, 2000, 2010, 2011 Free Software Foundation, Inc. This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by @@ -15,19 +15,19 @@ along with this program. If not, see . */ #include -#include "moments.h" + +#include "math/moments.h" + #include -#include #include #include -#include -#include -#include -#include "xalloc.h" +#include "data/val-type.h" +#include "data/value.h" +#include "libpspp/misc.h" + +#include "gl/xalloc.h" -#include "gettext.h" -#define _(msgid) gettext (msgid) /* Calculates variance, skewness, and kurtosis into *VARIANCE, *SKEWNESS, and *KURTOSIS if they are non-null and not greater @@ -50,13 +50,13 @@ calc_moments (enum moment max_moment, /* From _SPSS Statistical Algorithms, 2nd ed., 0-918469-89-9, section "DESCRIPTIVES". */ - if (fabs (*variance) >= 1e-20) + if (fabs (s2) >= 1e-20) { if (max_moment >= MOMENT_SKEWNESS && skewness != NULL && w > 2.) { double s3 = s2 * sqrt (s2); double g1 = (w * d3) / ((w - 1.0) * (w - 2.0) * s3); - if (gsl_finite (g1)) + if (isfinite (g1)) *skewness = g1; } if (max_moment >= MOMENT_KURTOSIS && kurtosis != NULL && w > 3.) @@ -64,7 +64,7 @@ calc_moments (enum moment max_moment, double den = (w - 2.) * (w - 3.) * pow2 (s2); double g2 = (w * (w + 1) * d4 / (w - 1.) / den - 3. * pow2 (d2) / den); - if (gsl_finite (g2)) + if (isfinite (g2)) *kurtosis = g2; } } @@ -223,12 +223,8 @@ moments_calculate (const struct moments *m, } else { - /* After the second pass we can calculate any stat. We - don't support "online" computation during the second - pass, so As a simple self-check, the total weight for - the passes must agree. */ + /* After the second pass we can calculate any stat. */ assert (m->pass == 2); - assert (m->w1 == m->w2); if (m->w2 > 0.) { @@ -248,14 +244,14 @@ moments_destroy (struct moments *m) free (m); } -/* Calculates the requested moments on the CNT values in ARRAY. +/* Calculates the requested moments on the N values in ARRAY. Each value is given a weight of 1. The total weight is stored into *WEIGHT (trivially) and the mean, variance, skewness, and kurtosis are stored into *MEAN, *VARIANCE, *SKEWNESS, and *KURTOSIS, respectively. Any of the result pointers may be null, in which case no value is stored. */ void -moments_of_doubles (const double *array, size_t cnt, +moments_of_doubles (const double *array, size_t n, double *weight, double *mean, double *variance, double *skewness, double *kurtosis) @@ -274,21 +270,21 @@ moments_of_doubles (const double *array, size_t cnt, max_moment = MOMENT_MEAN; init_moments (&m, max_moment); - for (idx = 0; idx < cnt; idx++) + for (idx = 0; idx < n; idx++) moments_pass_one (&m, array[idx], 1.); - for (idx = 0; idx < cnt; idx++) + for (idx = 0; idx < n; idx++) moments_pass_two (&m, array[idx], 1.); moments_calculate (&m, weight, mean, variance, skewness, kurtosis); } -/* Calculates the requested moments on the CNT numeric values in +/* Calculates the requested moments on the N numeric values in ARRAY. Each value is given a weight of 1. The total weight is stored into *WEIGHT (trivially) and the mean, variance, skewness, and kurtosis are stored into *MEAN, *VARIANCE, *SKEWNESS, and *KURTOSIS, respectively. Any of the result pointers may be null, in which case no value is stored. */ void -moments_of_values (const union value *array, size_t cnt, +moments_of_values (const union value *array, size_t n, double *weight, double *mean, double *variance, double *skewness, double *kurtosis) @@ -307,9 +303,9 @@ moments_of_values (const union value *array, size_t cnt, max_moment = MOMENT_MEAN; init_moments (&m, max_moment); - for (idx = 0; idx < cnt; idx++) + for (idx = 0; idx < n; idx++) moments_pass_one (&m, array[idx].f, 1.); - for (idx = 0; idx < cnt; idx++) + for (idx = 0; idx < n; idx++) moments_pass_two (&m, array[idx].f, 1.); moments_calculate (&m, weight, mean, variance, skewness, kurtosis); } @@ -459,6 +455,13 @@ moments1_destroy (struct moments1 *m) free (m); } + +double +calc_semean (double variance, double W) +{ + return sqrt (variance / W); +} + /* Returns the standard error of the skewness for the given total weight W.