1 /* PSPP - a program for statistical analysis. -*-c-*-
2 Copyright (C) 2010 Free Software Foundation, Inc.
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.
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.
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/>.
22 #include <gsl/gsl_cdf.h>
25 #include <data/format.h>
27 #include <libpspp/misc.h>
28 #include <libpspp/message.h>
29 #include <data/procedure.h>
30 #include <data/casereader.h>
31 #include <data/dictionary.h>
32 #include <data/variable.h>
36 #define _(msgid) gettext (msgid)
45 const struct dictionary *dict;
48 static void show_ranks_box (const struct one_sample_test *ost,
49 const struct friedman *fr);
51 static void show_sig_box (const struct one_sample_test *ost,
52 const struct friedman *fr);
61 cmp_x (const void *a_, const void *b_)
63 const struct datum *a = a_;
64 const struct datum *b = b_;
73 cmp_posn (const void *a_, const void *b_)
75 const struct datum *a = a_;
76 const struct datum *b = b_;
78 if (a->posn < b->posn)
81 return (a->posn > b->posn);
85 friedman_execute (const struct dataset *ds,
86 struct casereader *input,
87 enum mv_class exclude,
88 const struct npar_test *test,
92 double numerator = 0.0;
93 double denominator = 0.0;
96 const struct dictionary *dict = dataset_dict (ds);
97 const struct variable *weight = dict_get_weight (dict);
99 struct one_sample_test *ost = UP_CAST (test, struct one_sample_test, parent);
100 struct friedman_test *ft = UP_CAST (ost, struct friedman_test, parent);
103 double sigma_t = 0.0;
104 struct datum *row = xcalloc (ost->n_vars, sizeof *row);
107 fr.rank_sum = xcalloc (ost->n_vars, sizeof *fr.rank_sum);
110 for (v = 0; v < ost->n_vars; ++v)
113 fr.rank_sum[v] = 0.0;
116 input = casereader_create_filter_weight (input, dict, &warn, NULL);
117 input = casereader_create_filter_missing (input,
118 ost->vars, ost->n_vars,
121 for (; (c = casereader_read (input)); case_unref (c))
123 double prev_x = SYSMIS;
126 const double w = weight ? case_data (c, weight)->f: 1.0;
130 for (v = 0; v < ost->n_vars; ++v)
132 const struct variable *var = ost->vars[v];
133 const union value *val = case_data (c, var);
137 qsort (row, ost->n_vars, sizeof *row, cmp_x);
138 for (v = 0; v < ost->n_vars; ++v)
141 /* Replace value by the Rank */
147 for (i = v - run_length; i < v; ++i)
149 row[i].x *= run_length ;
151 row[i].x /= run_length + 1;
153 row[v].x = row[v-1].x;
160 double t = run_length + 1;
161 sigma_t += w * (pow3 (t) - t);
169 double t = run_length + 1;
170 sigma_t += w * (pow3 (t) - t );
173 qsort (row, ost->n_vars, sizeof *row, cmp_posn);
175 for (v = 0; v < ost->n_vars; ++v)
176 fr.rank_sum[v] += row[v].x * w;
178 casereader_destroy (input);
182 for (v = 0; v < ost->n_vars; ++v)
184 numerator += pow2 (fr.rank_sum[v]);
189 numerator *= 12.0 / (fr.cc * ost->n_vars * ( ost->n_vars + 1));
190 numerator -= 3 * fr.cc * ( ost->n_vars + 1);
192 denominator = 1 - sigma_t / ( fr.cc * ost->n_vars * ( pow2 (ost->n_vars) - 1));
194 fr.chi_sq = numerator / denominator;
199 fr.w -= 3 * pow2 (fr.cc) *
200 ost->n_vars * pow2 (ost->n_vars + 1);
202 fr.w /= pow2 (fr.cc) * (pow3 (ost->n_vars) - ost->n_vars)
208 show_ranks_box (ost, &fr);
209 show_sig_box (ost, &fr);
216 #include <output/tab.h>
219 show_ranks_box (const struct one_sample_test *ost, const struct friedman *fr)
222 const int row_headers = 1;
223 const int column_headers = 1;
224 struct tab_table *table =
225 tab_create (row_headers + 1, column_headers + ost->n_vars);
227 tab_headers (table, row_headers, 0, column_headers, 0);
229 tab_title (table, _("Ranks"));
231 /* Vertical lines inside the box */
232 tab_box (table, 1, 0, -1, TAL_1,
233 row_headers, 0, tab_nc (table) - 1, tab_nr (table) - 1 );
235 /* Box around the table */
236 tab_box (table, TAL_2, TAL_2, -1, -1,
237 0, 0, tab_nc (table) - 1, tab_nr (table) - 1 );
240 tab_text (table, 1, 0, 0, _("Mean Rank"));
242 tab_hline (table, TAL_2, 0, tab_nc (table) - 1, column_headers);
243 tab_vline (table, TAL_2, row_headers, 0, tab_nr (table) - 1);
245 for (i = 0 ; i < ost->n_vars ; ++i)
247 tab_text (table, 0, row_headers + i,
248 TAB_LEFT, var_to_string (ost->vars[i]));
250 tab_double (table, 1, row_headers + i,
251 0, fr->rank_sum[i] / fr->cc, 0);
259 show_sig_box (const struct one_sample_test *ost, const struct friedman *fr)
261 const struct friedman_test *ft = UP_CAST (ost, const struct friedman_test, parent);
264 const struct variable *weight = dict_get_weight (fr->dict);
265 const struct fmt_spec *wfmt = weight ? var_get_print_format (weight) : &F_8_0;
267 const int row_headers = 1;
268 const int column_headers = 0;
269 struct tab_table *table =
270 tab_create (row_headers + 1, column_headers + (ft->kendalls_w ? 5 : 4));
272 tab_headers (table, row_headers, 0, column_headers, 0);
274 tab_title (table, _("Test Statistics"));
276 tab_text (table, 0, column_headers + row++,
277 TAT_TITLE | TAB_LEFT , _("N"));
280 tab_text (table, 0, column_headers + row++,
281 TAT_TITLE | TAB_LEFT , _("Kendall's W"));
283 tab_text (table, 0, column_headers + row++,
284 TAT_TITLE | TAB_LEFT , _("Chi-Square"));
286 tab_text (table, 0, column_headers + row++,
287 TAT_TITLE | TAB_LEFT, _("df"));
289 tab_text (table, 0, column_headers + row++,
290 TAT_TITLE | TAB_LEFT, _("Asymp. Sig."));
292 /* Box around the table */
293 tab_box (table, TAL_2, TAL_2, -1, -1,
294 0, 0, tab_nc (table) - 1, tab_nr (table) - 1 );
297 tab_hline (table, TAL_2, 0, tab_nc (table) -1, column_headers);
298 tab_vline (table, TAL_2, row_headers, 0, tab_nr (table) - 1);
301 tab_double (table, 1, column_headers + row++,
305 tab_double (table, 1, column_headers + row++,
308 tab_double (table, 1, column_headers + row++,
311 tab_double (table, 1, column_headers + row++,
312 0, ost->n_vars - 1, &F_8_0);
314 tab_double (table, 1, column_headers + row++,
315 0, gsl_cdf_chisq_Q (fr->chi_sq, ost->n_vars - 1),