/* Pspp - a program for statistical analysis.
- Copyright (C) 2008 Free Software Foundation, Inc.
+ Copyright (C) 2008, 2009 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
#include <config.h>
+
#include "wilcoxon.h"
-#include <data/variable.h>
+
+#include <gsl/gsl_cdf.h>
+#include <math.h>
+#include <signal.h>
+#include <unistd.h>
+
#include <data/casereader.h>
#include <data/casewriter.h>
-#include <data/case-ordering.h>
-#include <math/sort.h>
-#include <libpspp/message.h>
-#include <xalloc.h>
-#include <output/table.h>
-#include <data/procedure.h>
#include <data/dictionary.h>
-#include <misc/wx-mp-sr.h>
+#include <data/format.h>
+#include <data/procedure.h>
+#include <data/subcase.h>
+#include <data/variable.h>
#include <gsl/gsl_cdf.h>
-#include <unistd.h>
-#include <signal.h>
#include <libpspp/assertion.h>
+#include <libpspp/message.h>
+#include <libpspp/misc.h>
+#include <math/sort.h>
+#include <math/wilcoxon-sig.h>
+#include <output/table.h>
+#include <signal.h>
+#include <unistd.h>
-static double timed_wilcoxon_significance (double w, long int n, double timer);
-
+#include "minmax.h"
+#include "xalloc.h"
static double
append_difference (const struct ccase *c, casenumber n UNUSED, void *aux)
}
static void show_ranks_box (const struct wilcoxon_state *,
- const struct two_sample_test *);
+ const struct two_sample_test *,
+ const struct dictionary *);
static void show_tests_box (const struct wilcoxon_state *,
const struct two_sample_test *,
struct wilcoxon_state *ws = xcalloc (sizeof (*ws), t2s->n_pairs);
const struct variable *weight = dict_get_weight (dict);
- struct variable *weightx = var_create_internal (WEIGHT_IDX);
+ struct variable *weightx = var_create_internal (WEIGHT_IDX, 0);
+ struct caseproto *proto;
input =
casereader_create_filter_weight (input, dict, &warn, NULL);
+ proto = caseproto_create ();
+ proto = caseproto_add_width (proto, 0);
+ proto = caseproto_add_width (proto, 0);
+ if (weight != NULL)
+ proto = caseproto_add_width (proto, 0);
+
for (i = 0 ; i < t2s->n_pairs; ++i )
{
struct casereader *r = casereader_clone (input);
struct casewriter *writer;
- struct ccase c;
- struct case_ordering *ordering = case_ordering_create ();
+ struct ccase *c;
+ struct subcase ordering;
variable_pair *vp = &t2s->pairs[i];
- const int reader_width = weight ? 3 : 2;
-
- ws[i].sign = var_create_internal (0);
- ws[i].absdiff = var_create_internal (1);
-
- case_ordering_add_var (ordering, ws[i].absdiff, SRT_ASCEND);
-
+ ws[i].sign = var_create_internal (0, 0);
+ ws[i].absdiff = var_create_internal (1, 0);
r = casereader_create_filter_missing (r, *vp, 2,
exclude,
NULL, NULL);
- writer = sort_create_writer (ordering, reader_width);
- while (casereader_read (r, &c))
- {
- struct ccase output;
- double d = append_difference (&c, 0, vp);
+ subcase_init_var (&ordering, ws[i].absdiff, SC_ASCEND);
+ writer = sort_create_writer (&ordering, proto);
+ subcase_destroy (&ordering);
- case_create (&output, reader_width);
+ for (; (c = casereader_read (r)) != NULL; case_unref (c))
+ {
+ struct ccase *output = case_create (proto);
+ double d = append_difference (c, 0, vp);
if (d > 0)
{
- case_data_rw (&output, ws[i].sign)->f = 1.0;
+ case_data_rw (output, ws[i].sign)->f = 1.0;
}
else if (d < 0)
{
- case_data_rw (&output, ws[i].sign)->f = -1.0;
+ case_data_rw (output, ws[i].sign)->f = -1.0;
}
else
{
double w = 1.0;
if (weight)
- w = case_data (&c, weight)->f;
+ w = case_data (c, weight)->f;
/* Central point values should be dropped */
ws[i].n_zeros += w;
- case_destroy (&c);
- continue;
+ case_unref (output);
+ continue;
}
- case_data_rw (&output, ws[i].absdiff)->f = fabs (d);
+ case_data_rw (output, ws[i].absdiff)->f = fabs (d);
if (weight)
- case_data_rw (&output, weightx)->f = case_data (&c, weight)->f;
+ case_data_rw (output, weightx)->f = case_data (c, weight)->f;
- casewriter_write (writer, &output);
- case_destroy (&c);
+ casewriter_write (writer, output);
}
casereader_destroy (r);
ws[i].reader = casewriter_make_reader (writer);
}
+ caseproto_unref (proto);
for (i = 0 ; i < t2s->n_pairs; ++i )
{
struct casereader *rr ;
- struct ccase c;
+ struct ccase *c;
enum rank_error err = 0;
rr = casereader_create_append_rank (ws[i].reader, ws[i].absdiff,
distinct_callback, &ws[i]
);
- while (casereader_read (rr, &c))
+ for (; (c = casereader_read (rr)) != NULL; case_unref (c))
{
- double sign = case_data (&c, ws[i].sign)->f;
- double rank = case_data_idx (&c, weight ? 3 : 2)->f;
+ double sign = case_data (c, ws[i].sign)->f;
+ double rank = case_data_idx (c, weight ? 3 : 2)->f;
double w = 1.0;
if (weight)
- w = case_data (&c, weightx)->f;
+ w = case_data (c, weightx)->f;
if ( sign > 0 )
{
}
else
NOT_REACHED ();
-
- case_destroy (&c);
}
casereader_destroy (rr);
var_destroy (weightx);
- show_ranks_box (ws, t2s);
+ show_ranks_box (ws, t2s, dict);
show_tests_box (ws, t2s, exact, timer);
for (i = 0 ; i < t2s->n_pairs; ++i )
#define _(msgid) gettext (msgid)
static void
-show_ranks_box (const struct wilcoxon_state *ws, const struct two_sample_test *t2s)
+show_ranks_box (const struct wilcoxon_state *ws,
+ const struct two_sample_test *t2s,
+ const struct dictionary *dict)
{
size_t i;
+
+ const struct variable *wv = dict_get_weight (dict);
+ const struct fmt_spec *wfmt = wv ? var_get_print_format (wv) : & F_8_0;
+
struct tab_table *table = tab_create (5, 1 + 4 * t2s->n_pairs, 0);
- tab_dim (table, tab_natural_dimensions);
+ tab_dim (table, tab_natural_dimensions, NULL);
tab_title (table, _("Ranks"));
/* N */
- tab_float (table, 2, 1 + i * 4, TAB_RIGHT, ws[i].negatives.n, 8, 0);
- tab_float (table, 2, 2 + i * 4, TAB_RIGHT, ws[i].positives.n, 8, 0);
- tab_float (table, 2, 3 + i * 4, TAB_RIGHT, ws[i].n_zeros, 8, 0);
+ tab_double (table, 2, 1 + i * 4, TAB_RIGHT, ws[i].negatives.n, wfmt);
+ tab_double (table, 2, 2 + i * 4, TAB_RIGHT, ws[i].positives.n, wfmt);
+ tab_double (table, 2, 3 + i * 4, TAB_RIGHT, ws[i].n_zeros, wfmt);
- tab_float (table, 2, 4 + i * 4, TAB_RIGHT,
- ws[i].n_zeros + ws[i].positives.n + ws[i].negatives.n, 8, 0);
+ tab_double (table, 2, 4 + i * 4, TAB_RIGHT,
+ ws[i].n_zeros + ws[i].positives.n + ws[i].negatives.n, wfmt);
/* Sums */
- tab_float (table, 4, 1 + i * 4, TAB_RIGHT, ws[i].negatives.sum, 8, 2);
- tab_float (table, 4, 2 + i * 4, TAB_RIGHT, ws[i].positives.sum, 8, 2);
+ tab_double (table, 4, 1 + i * 4, TAB_RIGHT, ws[i].negatives.sum, NULL);
+ tab_double (table, 4, 2 + i * 4, TAB_RIGHT, ws[i].positives.sum, NULL);
/* Means */
- tab_float (table, 3, 1 + i * 4, TAB_RIGHT,
- ws[i].negatives.sum / (double) ws[i].negatives.n, 8, 2);
+ tab_double (table, 3, 1 + i * 4, TAB_RIGHT,
+ ws[i].negatives.sum / (double) ws[i].negatives.n, NULL);
- tab_float (table, 3, 2 + i * 4, TAB_RIGHT,
- ws[i].positives.sum / (double) ws[i].positives.n, 8, 2);
+ tab_double (table, 3, 2 + i * 4, TAB_RIGHT,
+ ws[i].positives.sum / (double) ws[i].positives.n, NULL);
}
show_tests_box (const struct wilcoxon_state *ws,
const struct two_sample_test *t2s,
bool exact,
- double timer
+ double timer UNUSED
)
{
size_t i;
struct tab_table *table = tab_create (1 + t2s->n_pairs, exact ? 5 : 3, 0);
- tab_dim (table, tab_natural_dimensions);
+ tab_dim (table, tab_natural_dimensions, NULL);
tab_title (table, _("Test Statistics"));
tab_text (table, 0, 1, TAB_LEFT, _("Z"));
- tab_text (table, 0, 2, TAB_LEFT, _("Asymp. Sig (2-tailed)"));
+ tab_text (table, 0, 2, TAB_LEFT, _("Asymp. Sig. (2-tailed)"));
if ( exact )
{
- tab_text (table, 0, 3, TAB_LEFT, _("Exact Sig (2-tailed)"));
- tab_text (table, 0, 4, TAB_LEFT, _("Exact Sig (1-tailed)"));
+ tab_text (table, 0, 3, TAB_LEFT, _("Exact Sig. (2-tailed)"));
+ tab_text (table, 0, 4, TAB_LEFT, _("Exact Sig. (1-tailed)"));
#if 0
tab_text (table, 0, 5, TAB_LEFT, _("Point Probability"));
z /= sqrt (n * (n + 1) * (2*n + 1)/24.0 - ws[i].tiebreaker / 48.0);
- tab_float (table, 1 + i, 1, TAB_RIGHT, z, 8, 3);
+ tab_double (table, 1 + i, 1, TAB_RIGHT, z, NULL);
- tab_float (table, 1 + i, 2, TAB_RIGHT,
+ tab_double (table, 1 + i, 2, TAB_RIGHT,
2.0 * gsl_cdf_ugaussian_P (z),
- 8, 3);
+ NULL);
if (exact)
{
- double p =
- timed_wilcoxon_significance (ws[i].positives.sum,
- n,
- timer );
-
- if ( p == SYSMIS)
+ double p = LevelOfSignificanceWXMPSR (ws[i].positives.sum, n);
+ if (p < 0)
{
- msg (MW, _("Exact significance was not calculated after %.2f minutes. Skipping test."), timer);
+ msg (MW, ("Too many pairs to calculate exact significance."));
}
else
{
- tab_float (table, 1 + i, 3, TAB_RIGHT, p, 8, 3);
- tab_float (table, 1 + i, 4, TAB_RIGHT, p / 2.0, 8, 3);
+ tab_double (table, 1 + i, 3, TAB_RIGHT, p, NULL);
+ tab_double (table, 1 + i, 4, TAB_RIGHT, p / 2.0, NULL);
}
}
}
tab_submit (table);
}
-
-\f
-
-#include <setjmp.h>
-
-static sigjmp_buf env;
-
-static void
-give_up_callback (int signal UNUSED)
-{
- siglongjmp (env, 1);
-}
-
-static double
-timed_wilcoxon_significance (double w, long int n, double timer)
-{
- double p = SYSMIS;
-
- sigset_t set;
-
- struct sigaction timeout_action;
- struct sigaction old_action;
-
- if (timer <= 0 )
- return LevelOfSignificanceWXMPSR (w, n);
-
- sigemptyset (&set);
-
- timeout_action.sa_mask = set;
- timeout_action.sa_flags = 0;
-
- timeout_action.sa_handler = give_up_callback;
-
- if ( 0 == sigsetjmp (env, 1))
- {
- sigaction (SIGALRM, &timeout_action, &old_action);
- alarm (timer * 60.0);
-
- p = LevelOfSignificanceWXMPSR (w, n);
- }
-
- sigaction (SIGALRM, &old_action, NULL);
-
- return p;
-}
projects / pspp-builds.git / blobdiff