3723 |
05 Feb 18 |
peter |
1 |
#ifndef _theplu_yat_statistics_chi2_ |
3723 |
05 Feb 18 |
peter |
2 |
#define _theplu_yat_statistics_chi2_ |
3723 |
05 Feb 18 |
peter |
3 |
|
3723 |
05 Feb 18 |
peter |
// $Id$ |
3723 |
05 Feb 18 |
peter |
5 |
|
3723 |
05 Feb 18 |
peter |
6 |
/* |
4207 |
26 Aug 22 |
peter |
Copyright (C) 2018, 2022 Peter Johansson |
3723 |
05 Feb 18 |
peter |
8 |
|
3723 |
05 Feb 18 |
peter |
This file is part of the yat library, http://dev.thep.lu.se/yat |
3723 |
05 Feb 18 |
peter |
10 |
|
3723 |
05 Feb 18 |
peter |
The yat library is free software; you can redistribute it and/or |
3723 |
05 Feb 18 |
peter |
modify it under the terms of the GNU General Public License as |
3723 |
05 Feb 18 |
peter |
published by the Free Software Foundation; either version 3 of the |
3723 |
05 Feb 18 |
peter |
License, or (at your option) any later version. |
3723 |
05 Feb 18 |
peter |
15 |
|
3723 |
05 Feb 18 |
peter |
The yat library is distributed in the hope that it will be useful, |
3723 |
05 Feb 18 |
peter |
but WITHOUT ANY WARRANTY; without even the implied warranty of |
3723 |
05 Feb 18 |
peter |
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
3723 |
05 Feb 18 |
peter |
General Public License for more details. |
3723 |
05 Feb 18 |
peter |
20 |
|
3723 |
05 Feb 18 |
peter |
You should have received a copy of the GNU General Public License |
3723 |
05 Feb 18 |
peter |
along with yat. If not, see <http://www.gnu.org/licenses/>. |
3723 |
05 Feb 18 |
peter |
23 |
*/ |
3723 |
05 Feb 18 |
peter |
24 |
|
3723 |
05 Feb 18 |
peter |
25 |
#include <yat/utility/Matrix.h> |
3723 |
05 Feb 18 |
peter |
26 |
|
3723 |
05 Feb 18 |
peter |
27 |
namespace theplu { |
3723 |
05 Feb 18 |
peter |
28 |
namespace yat { |
3723 |
05 Feb 18 |
peter |
29 |
namespace statistics { |
3723 |
05 Feb 18 |
peter |
30 |
|
3723 |
05 Feb 18 |
peter |
31 |
/** |
3723 |
05 Feb 18 |
peter |
\brief Chi squared test |
3723 |
05 Feb 18 |
peter |
33 |
|
3723 |
05 Feb 18 |
peter |
\since New in yat 0.16 |
3723 |
05 Feb 18 |
peter |
35 |
*/ |
3723 |
05 Feb 18 |
peter |
36 |
class Chi2 |
3723 |
05 Feb 18 |
peter |
37 |
{ |
3723 |
05 Feb 18 |
peter |
38 |
public: |
3723 |
05 Feb 18 |
peter |
39 |
/** |
3723 |
05 Feb 18 |
peter |
\return degrees of freedom |
3723 |
05 Feb 18 |
peter |
41 |
*/ |
3723 |
05 Feb 18 |
peter |
42 |
unsigned int dof(void) const; |
3723 |
05 Feb 18 |
peter |
43 |
|
3723 |
05 Feb 18 |
peter |
44 |
/** |
3723 |
05 Feb 18 |
peter |
The expected value in element ij is calculated as sum of values |
3723 |
05 Feb 18 |
peter |
in row i times sum of values in column j divided by total sum |
3723 |
05 Feb 18 |
peter |
of values. |
3723 |
05 Feb 18 |
peter |
48 |
|
3723 |
05 Feb 18 |
peter |
\return expected values if there is no association between row |
3723 |
05 Feb 18 |
peter |
and column |
3723 |
05 Feb 18 |
peter |
51 |
*/ |
3723 |
05 Feb 18 |
peter |
52 |
const utility::Matrix& expected(void) const; |
3723 |
05 Feb 18 |
peter |
53 |
|
3723 |
05 Feb 18 |
peter |
54 |
/** |
3723 |
05 Feb 18 |
peter |
\return two-sided p-value |
3723 |
05 Feb 18 |
peter |
56 |
*/ |
3723 |
05 Feb 18 |
peter |
57 |
double p_value(void) const; |
3723 |
05 Feb 18 |
peter |
58 |
|
3723 |
05 Feb 18 |
peter |
59 |
/** |
3725 |
08 Feb 18 |
peter |
Randomise the data such that row sums and column sums are |
3725 |
08 Feb 18 |
peter |
conserved and calculate chi-squared. Return fraction of time |
3725 |
08 Feb 18 |
peter |
the randomised chi-squared is greater (or equal) than the |
3725 |
08 Feb 18 |
peter |
chi-squared calculated on original data. |
3725 |
08 Feb 18 |
peter |
64 |
|
3725 |
08 Feb 18 |
peter |
\param N number of randomisations |
3725 |
08 Feb 18 |
peter |
66 |
|
3725 |
08 Feb 18 |
peter |
Complexity: N x M, where M is total number of samples |
3725 |
08 Feb 18 |
peter |
68 |
*/ |
3725 |
08 Feb 18 |
peter |
69 |
double p_value(size_t N) const; |
3725 |
08 Feb 18 |
peter |
70 |
|
3725 |
08 Feb 18 |
peter |
71 |
/** |
3723 |
05 Feb 18 |
peter |
\return Chi2 sum calculated previously calling operator(void) |
3723 |
05 Feb 18 |
peter |
73 |
*/ |
3723 |
05 Feb 18 |
peter |
74 |
double operator()(void); |
3723 |
05 Feb 18 |
peter |
75 |
|
3723 |
05 Feb 18 |
peter |
76 |
/** |
3723 |
05 Feb 18 |
peter |
Sum is calculated as |
3723 |
05 Feb 18 |
peter |
78 |
|
3723 |
05 Feb 18 |
peter |
\f$ sum_{ij} \frac{ (X_{ij}-E_{ij})^2}{E_{ij}} \f$ |
3723 |
05 Feb 18 |
peter |
80 |
|
3723 |
05 Feb 18 |
peter |
where is X is data and E is matrix with expected values (see |
3723 |
05 Feb 18 |
peter |
expected(void)); |
3723 |
05 Feb 18 |
peter |
83 |
|
3723 |
05 Feb 18 |
peter |
\return Chi2 sum |
3723 |
05 Feb 18 |
peter |
85 |
*/ |
4125 |
14 Jan 22 |
peter |
86 |
double operator()(const utility::MatrixBase& data); |
3723 |
05 Feb 18 |
peter |
87 |
|
3723 |
05 Feb 18 |
peter |
88 |
private: |
3723 |
05 Feb 18 |
peter |
89 |
double chi2_; |
3723 |
05 Feb 18 |
peter |
90 |
utility::Matrix expected_; |
3725 |
08 Feb 18 |
peter |
91 |
utility::Vector rowsum_; |
3725 |
08 Feb 18 |
peter |
92 |
utility::Vector colsum_; |
3725 |
08 Feb 18 |
peter |
93 |
unsigned long int total_; |
3725 |
08 Feb 18 |
peter |
94 |
|
4125 |
14 Jan 22 |
peter |
95 |
void calculate_expected(const utility::MatrixBase& data); |
4125 |
14 Jan 22 |
peter |
96 |
void calculate_rowsum(const utility::MatrixBase& data); |
4125 |
14 Jan 22 |
peter |
97 |
void calculate_colsum(const utility::MatrixBase& data); |
3725 |
08 Feb 18 |
peter |
98 |
|
3725 |
08 Feb 18 |
peter |
// Create a matrix with non-negative integers such that the sum of |
3725 |
08 Feb 18 |
peter |
// row r is rowsum_(r) and sum of column c is colsum_(c) |
3725 |
08 Feb 18 |
peter |
101 |
utility::Matrix generate_random_data(void) const; |
3725 |
08 Feb 18 |
peter |
102 |
|
4125 |
14 Jan 22 |
peter |
103 |
double chi2sum(const utility::MatrixBase& data, |
3725 |
08 Feb 18 |
peter |
104 |
const utility::Matrix& expected) const; |
3725 |
08 Feb 18 |
peter |
105 |
|
3725 |
08 Feb 18 |
peter |
106 |
size_t draw(const utility::Vector& count, unsigned long int N) const; |
3725 |
08 Feb 18 |
peter |
107 |
|
3725 |
08 Feb 18 |
peter |
108 |
struct Calculator |
3725 |
08 Feb 18 |
peter |
109 |
{ |
3725 |
08 Feb 18 |
peter |
110 |
double operator()(double x, double m) const; |
3725 |
08 Feb 18 |
peter |
111 |
}; |
3723 |
05 Feb 18 |
peter |
112 |
}; |
3723 |
05 Feb 18 |
peter |
113 |
|
3723 |
05 Feb 18 |
peter |
114 |
}}} // of namespace statistics, yat, and theplu |
3723 |
05 Feb 18 |
peter |
115 |
|
3723 |
05 Feb 18 |
peter |
116 |
#endif |