+++ /dev/null
-/* PSPP - a program for statistical analysis.
- Copyright (C) 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
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- This program is distributed in the hope that it will be useful,
- but WITHOUT ANY WARRANTY; without even the implied warranty of
- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- GNU General Public License for more details.
-
- You should have received a copy of the GNU General Public License
- along with this program. If not, see <http://www.gnu.org/licenses/>.
-*/
-
-#include <config.h>
-#include "median.h"
-
-#include <gsl/gsl_cdf.h>
-
-#include "data/format.h"
-
-
-#include "data/variable.h"
-#include "data/case.h"
-#include "data/dictionary.h"
-#include "data/dataset.h"
-#include "data/casereader.h"
-#include "data/casewriter.h"
-#include "data/subcase.h"
-#include "data/value.h"
-
-#include "math/percentiles.h"
-#include "math/sort.h"
-
-#include "libpspp/cast.h"
-#include "libpspp/hmap.h"
-#include "libpspp/array.h"
-#include "libpspp/str.h"
-#include "libpspp/misc.h"
-
-#include "output/pivot-table.h"
-
-#include "gettext.h"
-#define N_(msgid) msgid
-#define _(msgid) gettext (msgid)
-
-
-struct val_node
-{
- struct hmap_node node;
- union value val;
- casenumber le;
- casenumber gt;
-};
-
-struct results
-{
- const struct variable *var;
- struct val_node **sorted_array;
- double n;
- double median;
- double chisq;
-};
-
-
-
-static int
-val_node_cmp_3way (const void *a_, const void *b_, const void *aux)
-{
- const struct variable *indep_var = aux;
- const struct val_node *const *a = a_;
- const struct val_node *const *b = b_;
-
- return value_compare_3way (&(*a)->val, &(*b)->val, var_get_width (indep_var));
-}
-
-static void
-show_frequencies (const struct n_sample_test *nst, const struct results *results, int n_vals, const struct dictionary *);
-
-static void
-show_test_statistics (const struct n_sample_test *nst, const struct results *results, int, const struct dictionary *);
-
-
-static struct val_node *
-find_value (const struct hmap *map, const union value *val,
- const struct variable *var)
-{
- struct val_node *foo = NULL;
- size_t hash = value_hash (val, var_get_width (var), 0);
- HMAP_FOR_EACH_WITH_HASH (foo, struct val_node, node, hash, map)
- if (value_equal (val, &foo->val, var_get_width (var)))
- break;
-
- return foo;
-}
-
-void
-median_execute (const struct dataset *ds,
- struct casereader *input,
- enum mv_class exclude,
- const struct npar_test *test,
- bool exact UNUSED,
- double timer UNUSED)
-{
- const struct dictionary *dict = dataset_dict (ds);
- const struct variable *wvar = dict_get_weight (dict);
- bool warn = true;
- int v;
- const struct median_test *mt = UP_CAST (test, const struct median_test,
- parent.parent);
-
- const struct n_sample_test *nst = UP_CAST (test, const struct n_sample_test,
- parent);
-
- const bool n_sample_test = (value_compare_3way (&nst->val2, &nst->val1,
- var_get_width (nst->indep_var)) > 0);
-
- struct results *results = XCALLOC (nst->n_vars, struct results);
- int n_vals = 0;
- for (v = 0; v < nst->n_vars; ++v)
- {
- double count = 0;
- double cc = 0;
- double median = mt->median;
- const struct variable *var = nst->vars[v];
- struct ccase *c;
- struct hmap map = HMAP_INITIALIZER (map);
- struct casereader *r = casereader_clone (input);
-
-
-
- if (n_sample_test == false)
- {
- struct val_node *vn = XZALLOC (struct val_node);
- value_clone (&vn->val, &nst->val1, var_get_width (nst->indep_var));
- hmap_insert (&map, &vn->node, value_hash (&nst->val1,
- var_get_width (nst->indep_var), 0));
-
- vn = xzalloc (sizeof *vn);
- value_clone (&vn->val, &nst->val2, var_get_width (nst->indep_var));
- hmap_insert (&map, &vn->node, value_hash (&nst->val2,
- var_get_width (nst->indep_var), 0));
- }
-
- if (median == SYSMIS)
- {
- struct percentile *ptl;
- struct order_stats *os;
-
- struct casereader *rr;
- struct subcase sc;
- struct casewriter *writer;
- subcase_init_var (&sc, var, SC_ASCEND);
- rr = casereader_clone (r);
- writer = sort_create_writer (&sc, casereader_get_proto (rr));
-
- for (; (c = casereader_read (rr)) != NULL;)
- {
- if (var_is_value_missing (var, case_data (c, var)) & exclude)
- {
- case_unref (c);
- continue;
- }
-
- cc += dict_get_case_weight (dict, c, &warn);
- casewriter_write (writer, c);
- }
- subcase_uninit (&sc);
- casereader_destroy (rr);
-
- rr = casewriter_make_reader (writer);
-
- ptl = percentile_create (0.5, cc);
- os = &ptl->parent;
-
- order_stats_accumulate (&os, 1,
- rr,
- wvar,
- var,
- exclude);
-
- median = percentile_calculate (ptl, PC_HAVERAGE);
- statistic_destroy (&ptl->parent.parent);
- }
-
- results[v].median = median;
-
-
- for (; (c = casereader_read (r)) != NULL; case_unref (c))
- {
- struct val_node *vn ;
- const double weight = dict_get_case_weight (dict, c, &warn);
- const union value *val = case_data (c, var);
- const union value *indep_val = case_data (c, nst->indep_var);
-
- if (var_is_value_missing (var, case_data (c, var)) & exclude)
- {
- continue;
- }
-
- if (n_sample_test)
- {
- int width = var_get_width (nst->indep_var);
- /* Ignore out of range values */
- if (
- value_compare_3way (indep_val, &nst->val1, width) < 0
- ||
- value_compare_3way (indep_val, &nst->val2, width) > 0
- )
- {
- continue;
- }
- }
-
- vn = find_value (&map, indep_val, nst->indep_var);
- if (vn == NULL)
- {
- if (n_sample_test == true)
- {
- int width = var_get_width (nst->indep_var);
- vn = xzalloc (sizeof *vn);
- value_clone (&vn->val, indep_val, width);
-
- hmap_insert (&map, &vn->node, value_hash (indep_val, width, 0));
- }
- else
- {
- continue;
- }
- }
-
- if (val->f <= median)
- vn->le += weight;
- else
- vn->gt += weight;
-
- count += weight;
- }
- casereader_destroy (r);
-
- {
- int x = 0;
- struct val_node *vn = NULL;
- double r_0 = 0;
- double r_1 = 0;
- HMAP_FOR_EACH (vn, struct val_node, node, &map)
- {
- r_0 += vn->le;
- r_1 += vn->gt;
- }
-
- results[v].n = count;
- results[v].sorted_array = XCALLOC (hmap_count (&map), struct val_node *);
- results[v].var = var;
-
- HMAP_FOR_EACH (vn, struct val_node, node, &map)
- {
- double e_0j = r_0 * (vn->le + vn->gt) / count;
- double e_1j = r_1 * (vn->le + vn->gt) / count;
-
- results[v].chisq += pow2 (vn->le - e_0j) / e_0j;
- results[v].chisq += pow2 (vn->gt - e_1j) / e_1j;
-
- results[v].sorted_array[x++] = vn;
- }
-
- n_vals = x;
- hmap_destroy (&map);
-
- sort (results[v].sorted_array, x, sizeof (*results[v].sorted_array),
- val_node_cmp_3way, nst->indep_var);
-
- }
- }
-
- casereader_destroy (input);
-
- show_frequencies (nst, results, n_vals, dict);
- show_test_statistics (nst, results, n_vals, dict);
-
- for (v = 0; v < nst->n_vars; ++v)
- {
- int i;
- const struct results *rs = results + v;
-
- for (i = 0; i < n_vals; ++i)
- {
- struct val_node *vn = rs->sorted_array[i];
- value_destroy (&vn->val, var_get_width (nst->indep_var));
- free (vn);
- }
- free (rs->sorted_array);
- }
- free (results);
-}
-
-
-
-static void
-show_frequencies (const struct n_sample_test *nst, const struct results *results, int n_vals, const struct dictionary *dict)
-{
- struct pivot_table *table = pivot_table_create (N_("Frequencies"));
- pivot_table_set_weight_var (table, dict_get_weight (dict));
-
- struct pivot_dimension *indep = pivot_dimension_create__ (
- table, PIVOT_AXIS_COLUMN, pivot_value_new_variable (nst->indep_var));
- indep->root->show_label = true;
- for (int i = 0; i < n_vals; ++i)
- pivot_category_create_leaf_rc (
- indep->root, pivot_value_new_var_value (
- nst->indep_var, &results->sorted_array[i]->val), PIVOT_RC_COUNT);
- pivot_dimension_create (table, PIVOT_AXIS_ROW, N_("Statistics"),
- N_("> Median"), N_("≤ Median"));
-
- struct pivot_dimension *dep = pivot_dimension_create (
- table, PIVOT_AXIS_ROW, N_("Dependent Variables"));
-
- for (int v = 0; v < nst->n_vars; ++v)
- {
- const struct results *rs = &results[v];
-
- int dep_idx = pivot_category_create_leaf (
- dep->root, pivot_value_new_variable (rs->var));
-
- for (int indep_idx = 0; indep_idx < n_vals; indep_idx++)
- {
- const struct val_node *vn = rs->sorted_array[indep_idx];
- pivot_table_put3 (table, indep_idx, 0, dep_idx,
- pivot_value_new_number (vn->gt));
- pivot_table_put3 (table, indep_idx, 1, dep_idx,
- pivot_value_new_number (vn->le));
- }
- }
-
- pivot_table_submit (table);
-}
-
-static void
-show_test_statistics (const struct n_sample_test *nst,
- const struct results *results,
- int n_vals,
- const struct dictionary *dict)
-{
- struct pivot_table *table = pivot_table_create (N_("Test Statistics"));
- pivot_table_set_weight_var (table, dict_get_weight (dict));
-
- pivot_dimension_create (table, PIVOT_AXIS_COLUMN, N_("Statistics"),
- N_("N"), PIVOT_RC_COUNT,
- N_("Median"),
- N_("Chi-Square"), PIVOT_RC_OTHER,
- N_("df"), PIVOT_RC_COUNT,
- N_("Asymp. Sig."), PIVOT_RC_SIGNIFICANCE);
-
- struct pivot_dimension *variables = pivot_dimension_create (
- table, PIVOT_AXIS_ROW, N_("Variables"));
-
- for (int v = 0; v < nst->n_vars; ++v)
- {
- double df = n_vals - 1;
- const struct results *rs = &results[v];
-
- int var_idx = pivot_category_create_leaf (
- variables->root, pivot_value_new_variable (rs->var));
-
- double entries[] = {
- rs->n,
- rs->median,
- rs->chisq,
- df,
- gsl_cdf_chisq_Q (rs->chisq, df),
- };
- for (size_t i = 0; i < sizeof entries / sizeof *entries; i++)
- {
- struct pivot_value *value
- = pivot_value_new_number (entries[i]);
- if (i == 1)
- value->numeric.format = *var_get_print_format (rs->var);
- pivot_table_put2 (table, i, var_idx, value);
- }
- }
-
- pivot_table_submit (table);
-}