source: trunk/test/kolmogorov_smirnov_test.cc @ 1658

Last change on this file since 1658 was 1626, checked in by Peter, 13 years ago

closes #457 - kolmogorov_smirnov_test fails on Mac OS 10.4

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 3.9 KB
Line 
1// $Id: kolmogorov_smirnov_test.cc 1626 2008-11-15 16:01:46Z peter $
2
3/*
4  Copyright (C) 2008 Peter Johansson
5
6  This file is part of the yat library, http://dev.thep.lu.se/yat
7
8  The yat library is free software; you can redistribute it and/or
9  modify it under the terms of the GNU General Public License as
10  published by the Free Software Foundation; either version 3 of the
11  License, or (at your option) any later version.
12
13  The yat library is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16  General Public License for more details.
17
18  You should have received a copy of the GNU General Public License
19  along with yat. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22#include "Suite.h"
23
24#include "yat/statistics/Averager.h"
25#include "yat/statistics/KolmogorovSmirnov.h"
26#include "yat/random/random.h"
27
28#include <cmath>
29#include <iostream>
30
31using namespace theplu::yat;
32
33void test_one_sample(test::Suite&);
34void test_two_sample(test::Suite&);
35void test_p_value(test::Suite&);
36void test_reset(test::Suite&);
37
38int main(int argc, char* argv[])
39{ 
40  test::Suite suite(argc, argv);
41
42  test_one_sample(suite);
43  test_two_sample(suite);
44  test_p_value(suite);
45  test_reset(suite);
46
47  return suite.return_value();
48}
49
50void test_one_sample(test::Suite& suite)
51{
52  std::vector<double> correct(11);
53  for (size_t i=0; i<correct.size(); ++i) {
54    double s1 = 1.0 - i/10.0;
55    double s2 = 0.0-i/10.0;
56    if (std::abs(s1)>std::abs(s2))
57      correct[i] = s1;
58    else
59      correct[i] = s2;
60  }
61
62  for (size_t i=0; i<11; ++i) {
63    statistics::KolmogorovSmirnov ks;
64    for (size_t j=0; j<11; ++j) {
65      ks.add(j, i==j);
66    }     
67    double score = ks.signed_score();
68    if (!suite.add(suite.equal(score, correct[i]))) {
69      suite.err() << "signed_score(void) failed\n";
70    }
71  }
72
73  statistics::KolmogorovSmirnov ks;
74  for (size_t i=0; i<11; ++i) {
75    ks.add(i, i==0);
76  }     
77  size_t n=110000;
78  double p = ks.p_value(n);
79  double p_correct = 2.0/11.0;
80  double margin = 10*std::sqrt(p_correct*(1-p_correct)/n);
81  if (p>p_correct+margin || p<p_correct-margin) {
82    suite.err() << "Error: p-value: " << p << "\n"
83                << "expected approximately: " << p_correct << "\n"
84                << "and at most " << margin << "deviation\n";
85    suite.add(false);
86  } 
87}
88
89
90void test_two_sample(test::Suite& suite)
91{
92  suite.err() << "testing two sample\n";
93  statistics::KolmogorovSmirnov ks;
94  for (size_t i=0; i<5; ++i)
95    ks.add(i+0.5, i<2);
96  suite.add(suite.equal(ks.score(), 1.0));
97  size_t n=100000;
98  double p = ks.p_value(n);
99  double p_correct = 0.2;
100  double margin=10*std::sqrt(p_correct*(1-p_correct)/n);
101  if (std::abs(p-p_correct)>margin) {
102    suite.add(false);
103    suite.err() << "Error: p = " << p << "\n"
104                << "correct p would be: " << p_correct << "\n"
105                << "expected a difference less than " << margin << "margin\n";
106    suite.err() << p << std::endl;
107  }
108}
109
110
111void test_p_value(test::Suite& suite)
112{
113  statistics::KolmogorovSmirnov ks;
114  for (size_t i=0; i<100; ++i) {
115    ks.add(i, true);
116    ks.add(i+14.5, false);
117  } 
118  suite.add(suite.equal(ks.score(), 0.15, 10));
119
120  statistics::Averager a;
121  for (size_t n=0; n<100; ++n) {
122    a.add(ks.p_value(100));
123  }
124  double margin = 5 * a.standard_error(); 
125  double p_approx = ks.p_value();
126  if (std::abs(a.mean()-p_approx)>margin) {
127    suite.add(false);
128    suite.err() << "Error: unexpected large deviation between p_values\n"
129                << "permutation p-value: " << a.mean() << "\n"
130                << "analytical approximation: " << p_approx << "\n"
131                << "expected deviation to be smaller than " << margin << "\n";
132  }
133 
134}
135
136void test_reset(test::Suite& suite)
137{
138  suite.err() << "testing reset\n";
139  statistics::KolmogorovSmirnov ks;
140  ks.add(1.0, true);
141  ks.add(2.0, false);
142  ks.add(3.0, true);
143  double score = ks.score();
144  double p = ks.p_value();
145  ks.reset();
146  ks.add(1.0, true);
147  ks.add(2.0, false);
148  ks.add(3.0, true);
149  suite.add(suite.equal(ks.score(), score));
150  suite.add(suite.equal(ks.p_value(), p));
151}
Note: See TracBrowser for help on using the repository browser.