/
eministat_analysis.erl
191 lines (172 loc) · 7.96 KB
/
eministat_analysis.erl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
-module(eministat_analysis).
-include("eministat.hrl").
-export([outlier_variance/3]).
-export([relative/3]).
%% -- OUTLIER VARIANCE ------
%%
%% @doc outlier_variance/3 computes the severity of the outlier variance
%%
%% http://www.ellipticgroup.com/misc/article_supplement.pdf
%%
%% @end
outlier_variance(_, Sigma, _) when Sigma < 0.000000000000001 -> {0.0, unaffected};
outlier_variance(Mu, Sigma, A) ->
MinBy = fun(F, Q, R) -> min(F(Q), F(R)) end,
MuA = Mu / A,
MugMin = MuA / 2,
SigmaG = min(MugMin / 4, Sigma / math:sqrt(A)),
SigmaG2 = SigmaG * SigmaG,
Sigma2 = Sigma * Sigma,
VarOut = fun(C) ->
AC = A - C,
(AC / A) * (Sigma2 - AC * SigmaG2)
end,
CMax = fun(X) ->
K = MuA - X,
D = K * K,
AD = A * D,
K0 = -A * AD,
K1 = Sigma2 - A * SigmaG2 + AD,
Det = K1 * K1 - 4 * SigmaG2 * K0,
trunc(-2 * K0 / (K1 + math:sqrt(Det)))
end,
VarOutMin = MinBy(VarOut, 1, (MinBy(CMax, 0, MugMin))) / Sigma2,
case VarOutMin of
K when K < 0.01 -> {K, unaffected};
K when K < 0.1 -> {K, slight};
K when K < 0.5 -> {K, moderate};
K when K > 0.5 -> {K, severe}
end.
%% -- RELATIVE #dataset{} COMPARISONS -----
relative(#dataset { n = DsN } = Ds, #dataset { n = RsN } = Rs, ConfIdx) ->
I = DsN + RsN - 2,
T = element(ConfIdx, student_lookup(I)),
Spool1 = (DsN - 1) * eministat_ds:variance(Ds) + (RsN - 1) * eministat_ds:variance(Rs),
Spool = math:sqrt(Spool1 / I),
S = Spool * math:sqrt(1.0 / DsN + 1.0 / RsN),
D = eministat_ds:mean(Ds) - eministat_ds:mean(Rs),
E = T * S,
case abs(D) > E of
false ->
{no_difference, element(ConfIdx, student_pct())};
true ->
{difference,
#{ confidence_level => element(ConfIdx, student_pct()),
difference => {D, E},
difference_pct => {
D * 100 / eministat_ds:mean(Rs),
E * 100 / eministat_ds:mean(Rs)},
pooled_s => Spool
}}
end.
%% -- STUDENT's T TABLES -----
%% Constant tables, represented as tuples for O(1) lookup speeds
student_pct() -> {80.0, 90.0, 95.0, 98.0, 99.0, 99.5}.
student_inf() ->
{ 1.282, 1.645, 1.960, 2.326, 2.576, 3.090 }. %% inf
-define(NSTUDENT, 100). %% Number of elements in the students distribution lookup table
student() ->
{
{ 3.078, 6.314, 12.706, 31.821, 63.657, 318.313 }, %% 1.
{ 1.886, 2.920, 4.303, 6.965, 9.925, 22.327 }, %% 2.
{ 1.638, 2.353, 3.182, 4.541, 5.841, 10.215 }, %% 3.
{ 1.533, 2.132, 2.776, 3.747, 4.604, 7.173 }, %% 4.
{ 1.476, 2.015, 2.571, 3.365, 4.032, 5.893 }, %% 5.
{ 1.440, 1.943, 2.447, 3.143, 3.707, 5.208 }, %% 6.
{ 1.415, 1.895, 2.365, 2.998, 3.499, 4.782 }, %% 7.
{ 1.397, 1.860, 2.306, 2.896, 3.355, 4.499 }, %% 8.
{ 1.383, 1.833, 2.262, 2.821, 3.250, 4.296 }, %% 9.
{ 1.372, 1.812, 2.228, 2.764, 3.169, 4.143 }, %% 10.
{ 1.363, 1.796, 2.201, 2.718, 3.106, 4.024 }, %% 11.
{ 1.356, 1.782, 2.179, 2.681, 3.055, 3.929 }, %% 12.
{ 1.350, 1.771, 2.160, 2.650, 3.012, 3.852 }, %% 13.
{ 1.345, 1.761, 2.145, 2.624, 2.977, 3.787 }, %% 14.
{ 1.341, 1.753, 2.131, 2.602, 2.947, 3.733 }, %% 15.
{ 1.337, 1.746, 2.120, 2.583, 2.921, 3.686 }, %% 16.
{ 1.333, 1.740, 2.110, 2.567, 2.898, 3.646 }, %% 17.
{ 1.330, 1.734, 2.101, 2.552, 2.878, 3.610 }, %% 18.
{ 1.328, 1.729, 2.093, 2.539, 2.861, 3.579 }, %% 19.
{ 1.325, 1.725, 2.086, 2.528, 2.845, 3.552 }, %% 20.
{ 1.323, 1.721, 2.080, 2.518, 2.831, 3.527 }, %% 21.
{ 1.321, 1.717, 2.074, 2.508, 2.819, 3.505 }, %% 22.
{ 1.319, 1.714, 2.069, 2.500, 2.807, 3.485 }, %% 23.
{ 1.318, 1.711, 2.064, 2.492, 2.797, 3.467 }, %% 24.
{ 1.316, 1.708, 2.060, 2.485, 2.787, 3.450 }, %% 25.
{ 1.315, 1.706, 2.056, 2.479, 2.779, 3.435 }, %% 26.
{ 1.314, 1.703, 2.052, 2.473, 2.771, 3.421 }, %% 27.
{ 1.313, 1.701, 2.048, 2.467, 2.763, 3.408 }, %% 28.
{ 1.311, 1.699, 2.045, 2.462, 2.756, 3.396 }, %% 29.
{ 1.310, 1.697, 2.042, 2.457, 2.750, 3.385 }, %% 30.
{ 1.309, 1.696, 2.040, 2.453, 2.744, 3.375 }, %% 31.
{ 1.309, 1.694, 2.037, 2.449, 2.738, 3.365 }, %% 32.
{ 1.308, 1.692, 2.035, 2.445, 2.733, 3.356 }, %% 33.
{ 1.307, 1.691, 2.032, 2.441, 2.728, 3.348 }, %% 34.
{ 1.306, 1.690, 2.030, 2.438, 2.724, 3.340 }, %% 35.
{ 1.306, 1.688, 2.028, 2.434, 2.719, 3.333 }, %% 36.
{ 1.305, 1.687, 2.026, 2.431, 2.715, 3.326 }, %% 37.
{ 1.304, 1.686, 2.024, 2.429, 2.712, 3.319 }, %% 38.
{ 1.304, 1.685, 2.023, 2.426, 2.708, 3.313 }, %% 39.
{ 1.303, 1.684, 2.021, 2.423, 2.704, 3.307 }, %% 40.
{ 1.303, 1.683, 2.020, 2.421, 2.701, 3.301 }, %% 41.
{ 1.302, 1.682, 2.018, 2.418, 2.698, 3.296 }, %% 42.
{ 1.302, 1.681, 2.017, 2.416, 2.695, 3.291 }, %% 43.
{ 1.301, 1.680, 2.015, 2.414, 2.692, 3.286 }, %% 44.
{ 1.301, 1.679, 2.014, 2.412, 2.690, 3.281 }, %% 45.
{ 1.300, 1.679, 2.013, 2.410, 2.687, 3.277 }, %% 46.
{ 1.300, 1.678, 2.012, 2.408, 2.685, 3.273 }, %% 47.
{ 1.299, 1.677, 2.011, 2.407, 2.682, 3.269 }, %% 48.
{ 1.299, 1.677, 2.010, 2.405, 2.680, 3.265 }, %% 49.
{ 1.299, 1.676, 2.009, 2.403, 2.678, 3.261 }, %% 50.
{ 1.298, 1.675, 2.008, 2.402, 2.676, 3.258 }, %% 51.
{ 1.298, 1.675, 2.007, 2.400, 2.674, 3.255 }, %% 52.
{ 1.298, 1.674, 2.006, 2.399, 2.672, 3.251 }, %% 53.
{ 1.297, 1.674, 2.005, 2.397, 2.670, 3.248 }, %% 54.
{ 1.297, 1.673, 2.004, 2.396, 2.668, 3.245 }, %% 55.
{ 1.297, 1.673, 2.003, 2.395, 2.667, 3.242 }, %% 56.
{ 1.297, 1.672, 2.002, 2.394, 2.665, 3.239 }, %% 57.
{ 1.296, 1.672, 2.002, 2.392, 2.663, 3.237 }, %% 58.
{ 1.296, 1.671, 2.001, 2.391, 2.662, 3.234 }, %% 59.
{ 1.296, 1.671, 2.000, 2.390, 2.660, 3.232 }, %% 60.
{ 1.296, 1.670, 2.000, 2.389, 2.659, 3.229 }, %% 61.
{ 1.295, 1.670, 1.999, 2.388, 2.657, 3.227 }, %% 62.
{ 1.295, 1.669, 1.998, 2.387, 2.656, 3.225 }, %% 63.
{ 1.295, 1.669, 1.998, 2.386, 2.655, 3.223 }, %% 64.
{ 1.295, 1.669, 1.997, 2.385, 2.654, 3.220 }, %% 65.
{ 1.295, 1.668, 1.997, 2.384, 2.652, 3.218 }, %% 66.
{ 1.294, 1.668, 1.996, 2.383, 2.651, 3.216 }, %% 67.
{ 1.294, 1.668, 1.995, 2.382, 2.650, 3.214 }, %% 68.
{ 1.294, 1.667, 1.995, 2.382, 2.649, 3.213 }, %% 69.
{ 1.294, 1.667, 1.994, 2.381, 2.648, 3.211 }, %% 70.
{ 1.294, 1.667, 1.994, 2.380, 2.647, 3.209 }, %% 71.
{ 1.293, 1.666, 1.993, 2.379, 2.646, 3.207 }, %% 72.
{ 1.293, 1.666, 1.993, 2.379, 2.645, 3.206 }, %% 73.
{ 1.293, 1.666, 1.993, 2.378, 2.644, 3.204 }, %% 74.
{ 1.293, 1.665, 1.992, 2.377, 2.643, 3.202 }, %% 75.
{ 1.293, 1.665, 1.992, 2.376, 2.642, 3.201 }, %% 76.
{ 1.293, 1.665, 1.991, 2.376, 2.641, 3.199 }, %% 77.
{ 1.292, 1.665, 1.991, 2.375, 2.640, 3.198 }, %% 78.
{ 1.292, 1.664, 1.990, 2.374, 2.640, 3.197 }, %% 79.
{ 1.292, 1.664, 1.990, 2.374, 2.639, 3.195 }, %% 80.
{ 1.292, 1.664, 1.990, 2.373, 2.638, 3.194 }, %% 81.
{ 1.292, 1.664, 1.989, 2.373, 2.637, 3.193 }, %% 82.
{ 1.292, 1.663, 1.989, 2.372, 2.636, 3.191 }, %% 83.
{ 1.292, 1.663, 1.989, 2.372, 2.636, 3.190 }, %% 84.
{ 1.292, 1.663, 1.988, 2.371, 2.635, 3.189 }, %% 85.
{ 1.291, 1.663, 1.988, 2.370, 2.634, 3.188 }, %% 86.
{ 1.291, 1.663, 1.988, 2.370, 2.634, 3.187 }, %% 87.
{ 1.291, 1.662, 1.987, 2.369, 2.633, 3.185 }, %% 88.
{ 1.291, 1.662, 1.987, 2.369, 2.632, 3.184 }, %% 89.
{ 1.291, 1.662, 1.987, 2.368, 2.632, 3.183 }, %% 90.
{ 1.291, 1.662, 1.986, 2.368, 2.631, 3.182 }, %% 91.
{ 1.291, 1.662, 1.986, 2.368, 2.630, 3.181 }, %% 92.
{ 1.291, 1.661, 1.986, 2.367, 2.630, 3.180 }, %% 93.
{ 1.291, 1.661, 1.986, 2.367, 2.629, 3.179 }, %% 94.
{ 1.291, 1.661, 1.985, 2.366, 2.629, 3.178 }, %% 95.
{ 1.290, 1.661, 1.985, 2.366, 2.628, 3.177 }, %% 96.
{ 1.290, 1.661, 1.985, 2.365, 2.627, 3.176 }, %% 97.
{ 1.290, 1.661, 1.984, 2.365, 2.627, 3.175 }, %% 98.
{ 1.290, 1.660, 1.984, 2.365, 2.626, 3.175 }, %% 99.
{ 1.290, 1.660, 1.984, 2.364, 2.626, 3.174 } %% 100.
}.
student_lookup(I) when I > ?NSTUDENT -> student_inf();
student_lookup(I) -> element(I, student()).