aGrUM 3.2.0
a C++ library for (probabilistic) graphical models
chi2.cpp
Go to the documentation of this file.
1/****************************************************************************
2 * This file is part of the aGrUM/pyAgrum library. *
3 * *
4 * Copyright (c) 2005-2026 by *
5 * - Pierre-Henri WUILLEMIN(_at_LIP6) *
6 * - Christophe GONZALES(_at_AMU) *
7 * *
8 * The aGrUM/pyAgrum library is free software; you can redistribute it *
9 * and/or modify it under the terms of either : *
10 * *
11 * - the GNU Lesser General Public License as published by *
12 * the Free Software Foundation, either version 3 of the License, *
13 * or (at your option) any later version, *
14 * - the MIT license (MIT), *
15 * - or both in dual license, as here. *
16 * *
17 * (see https://agrum.gitlab.io/articles/dual-licenses-lgplv3mit.html) *
18 * *
19 * This aGrUM/pyAgrum library is distributed in the hope that it will be *
20 * useful, but WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, *
21 * INCLUDING BUT NOT LIMITED TO THE WARRANTIES MERCHANTABILITY or FITNESS *
22 * FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE *
23 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER *
24 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, *
25 * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR *
26 * OTHER DEALINGS IN THE SOFTWARE. *
27 * *
28 * See LICENCES for more details. *
29 * *
30 * SPDX-FileCopyrightText: Copyright 2005-2026 *
31 * - Pierre-Henri WUILLEMIN(_at_LIP6) *
32 * - Christophe GONZALES(_at_AMU) *
33 * SPDX-License-Identifier: LGPL-3.0-or-later OR MIT *
34 * *
35 * Contact : info_at_agrum_dot_org *
36 * homepage : http://agrum.gitlab.io *
37 * gitlab : https://gitlab.com/agrumery/agrum *
38 * *
39 ****************************************************************************/
40
41
51
52#include <agrum/base/core/math/chi2.h>
53
54#ifndef DOXYGEN_SHOULD_SKIP_THIS
55
56// constants used by Gary Perlman for his code for computing chi2 critical
57// values
58# define GUM_Z_MAX 6.0 // maximum meaningful z value
59# define GUM_CHI_EPSILON 0.000001 // accuracy of critchi approximation
60# define GUM_CHI_MAX 99999.0 // maximum chi square value
61# define GUM_LOG_SQRT_PI 0.5723649429247000870717135 // std::log (std::sqrt (pi))
62# define GUM_I_SQRT_PI 0.5641895835477562869480795 // 1 / std::sqrt (pi)
63# define GUM_BIGX 20.0 // max value to represent exp (x)
64# define _gum_ex(x) (((x) < -GUM_BIGX) ? 0.0 : std::exp(x))
65
66#endif /* DOXYGEN_SHOULD_SKIP_THIS */
67
68
69namespace gum {
70
71 // computes the probability of normal z value (used by the cache)
72 double Chi2::_probaZValue_(double z) {
73 double y, x, w;
74
75 if (z == 0.0) x = 0.0;
76 else {
77 y = 0.5 * std::fabs(z);
78
79 if (y >= (GUM_Z_MAX * 0.5)) x = 1.0;
80 else if (y < 1.0) {
81 w = y * y;
82 x = ((((((((0.000124818987 * w - 0.001075204047) * w + 0.005198775019) * w - 0.019198292004)
83 * w
84 + 0.059054035642)
85 * w
86 - 0.151968751364)
87 * w
88 + 0.319152932694)
89 * w
90 - 0.531923007300)
91 * w
92 + 0.797884560593)
93 * y * 2.0;
94 } else {
95 y -= 2.0;
96 x = (((((((((((((-0.000045255659 * y + 0.000152529290) * y - 0.000019538132) * y
97 - 0.000676904986)
98 * y
99 + 0.001390604284)
100 * y
101 - 0.000794620820)
102 * y
103 - 0.002034254874)
104 * y
105 + 0.006549791214)
106 * y
107 - 0.010557625006)
108 * y
109 + 0.011630447319)
110 * y
111 - 0.009279453341)
112 * y
113 + 0.005353579108)
114 * y
115 - 0.002141268741)
116 * y
117 + 0.000535310849)
118 * y
119 + 0.999936657524;
120 }
121 }
122
123 return (z > 0.0 ? ((x + 1.0) * 0.5) : ((1.0 - x) * 0.5));
124 }
125
126 // computes the probability of chi2 value (used by the cache)
127 double Chi2::probaChi2(double x, Size df) {
128 double a, y = 0, s;
129 double e, c, z;
130 int even; /* true if df is an even number */
131
132 if ((x <= 0.0) || (df < 1)) return (1.0);
133
134 a = 0.5 * x;
135
136 even = (2 * (df / 2)) == df;
137
138 if (df > 1) y = _gum_ex(-a);
139
140 s = (even ? y : (2.0 * _probaZValue_(-std::sqrt(x))));
141
142 if (df > 2) {
143 x = 0.5 * (df - 1.0);
144 z = (even ? 1.0 : 0.5);
145
146 if (a > GUM_BIGX) {
147 e = (even ? 0.0 : GUM_LOG_SQRT_PI);
148 c = std::log(a);
149
150 while (z <= x) {
151 e = std::log(z) + e;
152 s += _gum_ex(c * z - a - e);
153 z += 1.0;
154 }
155
156 return (s);
157 } else {
158 e = (even ? 1.0 : (GUM_I_SQRT_PI / std::sqrt(a)));
159 c = 0.0;
160
161 while (z <= x) {
162 e = e * (a / z);
163 c = c + e;
164 z += 1.0;
165 }
166
167 return (c * y + s);
168 }
169 } else return (s);
170 }
171
172 // computes the critical value of a given chi2 test
173 double Chi2::criticalValue(double proba, Size df) {
174 double minchisq = 0.0;
175 double maxchisq = GUM_CHI_MAX;
176 double chisqval;
177
178 if (df == 0) return (0.0);
179
180 if (proba <= 0.0) return (maxchisq);
181 else if (proba >= 1.0) return (0.0);
182
183 chisqval = df / std::sqrt(proba); /* fair first value */
184
185 while (maxchisq - minchisq > GUM_CHI_EPSILON) {
186 if (probaChi2(chisqval, df) < proba) maxchisq = chisqval;
187 else minchisq = chisqval;
188
189 chisqval = (maxchisq + minchisq) * 0.5;
190 }
191
192 return (chisqval);
193 }
194
195} /* namespace gum */
chi2.h
The class that represents the chi2 distribution.
gum::Chi2::criticalValue
static double criticalValue(double proba, Size df)
Computes the critical chi2 value for a given confidence probability and number of degrees of freedom.
Definition chi2.cpp:173
gum::Chi2::_probaZValue_
static double _probaZValue_(double z)
Computes the probability of normal z value.
Definition chi2.cpp:72
gum::Chi2::probaChi2
static double probaChi2(double x, Size df)
Computes the probability of chi2 value.
Definition chi2.cpp:127
gum::Size
std::size_t Size
In aGrUM, hashed values are unsigned long int.
Definition types.h:74
gum
gum is the global namespace for all aGrUM entities
Definition agrum.h:46