File: nc_t.qbk

package info (click to toggle)
boost1.88 1.88.0-1
  • links: PTS, VCS
  • area: main
  • in suites: trixie
  • size: 576,932 kB
  • sloc: cpp: 4,149,234; xml: 136,789; ansic: 35,092; python: 33,910; asm: 5,698; sh: 4,604; ada: 1,681; makefile: 1,633; pascal: 1,139; perl: 1,124; sql: 640; yacc: 478; ruby: 271; java: 77; lisp: 24; csh: 6
file content (223 lines) | stat: -rw-r--r-- 7,821 bytes parent folder | download | duplicates (9)
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
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
[section:nc_t_dist Noncentral T Distribution]

``#include <boost/math/distributions/non_central_t.hpp>``

   namespace boost{ namespace math{

   template <class RealType = double,
             class ``__Policy``   = ``__policy_class`` >
   class non_central_t_distribution;

   typedef non_central_t_distribution<> non_central_t;

   template <class RealType, class ``__Policy``>
   class non_central_t_distribution
   {
   public:
      typedef RealType  value_type;
      typedef Policy    policy_type;

      // Constructor:
      non_central_t_distribution(RealType v, RealType delta);

      // Accessor to degrees_of_freedom parameter v:
      RealType degrees_of_freedom()const;

      // Accessor to non-centrality parameter delta:
      RealType non_centrality()const;
   };

   }} // namespaces

The noncentral T distribution is a generalization of the __students_t_distrib.
Let X have a normal distribution with mean [delta] and variance 1, and let
['[nu] S[super 2]] have
a chi-squared distribution with degrees of freedom [nu]. Assume that
X and S[super 2] are independent.
The distribution of [role serif_italic t[sub [nu]]([delta])=X/S] is called a
noncentral t distribution with degrees of freedom [nu] and noncentrality parameter [delta].

This gives the following PDF:

[equation nc_t_ref1]

where [role serif_italic [sub 1]F[sub 1](a;b;x)] is a confluent hypergeometric function.

The following graph illustrates how the distribution changes
for different values of [nu] and [delta]:

[graph nc_t_pdf]
[graph nc_t_cdf]

[h4 Member Functions]

      non_central_t_distribution(RealType v, RealType delta);

Constructs a non-central t distribution with degrees of freedom
parameter /v/ and non-centrality parameter /delta/.

Requires /v/ > 0 (including positive infinity) and finite /delta/, otherwise calls __domain_error.

      RealType degrees_of_freedom()const;

Returns the parameter /v/ from which this object was constructed.

      RealType non_centrality()const;

Returns the non-centrality parameter /delta/ from which this object was constructed.

[h4 Non-member Accessors]

All the [link math_toolkit.dist_ref.nmp usual non-member accessor functions]
that are generic to all distributions are supported: __usual_accessors.

The domain of the random variable is \[-[infin], +[infin]\].

[h4 Accuracy]

The following table shows the peak errors
(in units of [@http://en.wikipedia.org/wiki/Machine_epsilon epsilon])
found on various platforms with various floating-point types.
Unless otherwise specified, any floating-point type that is narrower
than the one shown will have __zero_error.

[table_non_central_t_CDF]

[table_non_central_t_CDF_complement]

[caution The complexity of the current algorithm is dependent upon
[delta][super 2]: consequently the time taken to evaluate the CDF
increases rapidly for [delta] > 500, likewise the accuracy decreases
rapidly for very large [delta].]

Accuracy for the quantile and PDF functions should be broadly similar.
The /mode/ is determined numerically and cannot
in principal be more accurate than the square root of
floating-point type FPT epsilon, accessed using `boost::math::tools::epsilon<FPT>()`.
For 64-bit `double`, epsilon is about 1e-16, so the fractional accuracy is limited to 1e-8.

[h4 Tests]

There are two sets of tests of this distribution:

Basic sanity checks compare this implementation to the test values given in
"Computing discrete mixtures of continuous
distributions: noncentral chisquare, noncentral t
and the distribution of the square of the sample
multiple correlation coefficient."
Denise Benton, K. Krishnamoorthy,
Computational Statistics & Data Analysis 43 (2003) 249-267.

Accuracy checks use test data computed with this
implementation and arbitrary precision interval arithmetic:
this test data is believed to be accurate to at least 50
decimal places.

The cases of large (or infinite) [nu] and/or large [delta] has received special
treatment to avoid catastrophic loss of accuracy.
New tests have been added to confirm the improvement achieved.

From Boost 1.52, degrees of freedom [nu] can be +[infin]
when the normal distribution located at [delta]
(equivalent to the central Student's t distribution)
is used in place for accuracy and speed.

[h4 Implementation]

The CDF is computed using a modification of the method
described in
"Computing discrete mixtures of continuous
distributions: noncentral chisquare, noncentral t
and the distribution of the square of the sample
multiple correlation coefficient."
Denise Benton, K. Krishnamoorthy,
Computational Statistics & Data Analysis 43 (2003) 249-267.

This uses the following formula for the CDF:

[equation nc_t_ref2]

Where I[sub x](a,b) is the incomplete beta function, and
[Phi](x) is the normal CDF at x.

Iteration starts at the largest of the Poisson weighting terms
(at i = [delta][super 2] / 2) and then proceeds in both directions
as per Benton and Krishnamoorthy's paper.

Alternatively, by considering what happens when t = [infin], we have
x = 1, and therefore I[sub x](a,b) = 1 and:

[equation nc_t_ref3]

From this we can easily show that:

[equation nc_t_ref4]

and therefore we have a means to compute either the probability or its
complement directly without the risk of cancellation error.  The
crossover criterion for choosing whether to calculate the CDF or
its complement is the same as for the
__non_central_beta_distrib.

The PDF can be computed by a very similar method using:

[equation nc_t_ref5]

Where I[sub x][super '](a,b) is the derivative of the incomplete beta function.

For both the PDF and CDF we switch to approximating the distribution by a
Student's t distribution centred on [delta] when [nu] is very large.
The crossover location appears to be when [delta]/(4[nu]) < [epsilon],
this location was estimated by inspection of equation 2.6 in
"A Comparison of Approximations To Percentiles of the
Noncentral t-Distribution".  H. Sahai and M. M. Ojeda,
Revista Investigacion Operacional Vol 21, No 2, 2000, page 123.

Equation 2.6 is a Fisher-Cornish expansion by Eeden and Johnson.
The second term includes the ratio [delta]/(4[nu]),
so when this term become negligible, this and following terms can be ignored,
leaving just Student's t distribution centred on [delta].

This was also confirmed by experimental testing.

See also

* "Some Approximations to the Percentage Points of the Noncentral
t-Distribution". C. van Eeden. International Statistical Review, 29, 4-31.

* "Continuous Univariate Distributions".  N.L. Johnson, S. Kotz and
N. Balkrishnan. 1995. John Wiley and Sons New York.

The quantile is calculated via the usual
__root_finding_without_derivatives method
with the initial guess taken as the quantile of a normal approximation
to the noncentral T.

There is no closed form for the mode, so this is computed via
functional maximisation of the PDF.

The remaining functions (mean, variance etc) are implemented
using the formulas given in
Weisstein, Eric W. "Noncentral Student's t-Distribution."
From MathWorld--A Wolfram Web Resource.
[@http://mathworld.wolfram.com/NoncentralStudentst-Distribution.html
http://mathworld.wolfram.com/NoncentralStudentst-Distribution.html]
and in the
[@http://reference.wolfram.com/mathematica/ref/NoncentralStudentTDistribution.html
Mathematica documentation].

Some analytic properties of noncentral distributions
(particularly unimodality, and monotonicity of their modes)
are surveyed and summarized by:

Andrea van Aubel & Wolfgang Gawronski, Applied Mathematics and Computation, 141 (2003) 3-12.

[endsect] [/section:nc_t_dist]

[/ nc_t.qbk
  Copyright 2008, 2012 John Maddock and Paul A. Bristow.
  Distributed under the Boost Software License, Version 1.0.
  (See accompanying file LICENSE_1_0.txt or copy at
  http://www.boost.org/LICENSE_1_0.txt).
]