Squashed 'third_party/boostorg/math/' content from commit 0e9549f
Change-Id: I7c2a13cb6a5beea4a471341510d8364cedd71613
git-subtree-dir: third_party/boostorg/math
git-subtree-split: 0e9549ff2f854e6edafaf4627d65026f2f533a18
diff --git a/doc/distributions/nc_chi_squared.qbk b/doc/distributions/nc_chi_squared.qbk
new file mode 100644
index 0000000..a1c2a32
--- /dev/null
+++ b/doc/distributions/nc_chi_squared.qbk
@@ -0,0 +1,273 @@
+[section:nc_chi_squared_dist Noncentral Chi-Squared Distribution]
+
+``#include <boost/math/distributions/non_central_chi_squared.hpp>``
+
+ namespace boost{ namespace math{
+
+ template <class RealType = double,
+ class ``__Policy`` = ``__policy_class`` >
+ class non_central_chi_squared_distribution;
+
+ typedef non_central_chi_squared_distribution<> non_central_chi_squared;
+
+ template <class RealType, class ``__Policy``>
+ class non_central_chi_squared_distribution
+ {
+ public:
+ typedef RealType value_type;
+ typedef Policy policy_type;
+
+ // Constructor:
+ non_central_chi_squared_distribution(RealType v, RealType lambda);
+
+ // Accessor to degrees of freedom parameter v:
+ RealType degrees_of_freedom()const;
+
+ // Accessor to non centrality parameter lambda:
+ RealType non_centrality()const;
+
+ // Parameter finders:
+ static RealType find_degrees_of_freedom(RealType lambda, RealType x, RealType p);
+ template <class A, class B, class C>
+ static RealType find_degrees_of_freedom(const complemented3_type<A,B,C>& c);
+
+ static RealType find_non_centrality(RealType v, RealType x, RealType p);
+ template <class A, class B, class C>
+ static RealType find_non_centrality(const complemented3_type<A,B,C>& c);
+ };
+
+ }} // namespaces
+
+The noncentral chi-squared distribution is a generalization of the
+__chi_squared_distrib. If X[sub i] are [nu] independent, normally
+distributed random variables with means [mu][sub i] and variances
+[sigma][sub i][super 2], then the random variable
+
+[equation nc_chi_squ_ref1]
+
+is distributed according to the noncentral chi-squared distribution.
+
+The noncentral chi-squared distribution has two parameters:
+[nu] which specifies the number of degrees of freedom
+(i.e. the number of X[sub i]), and [lambda] which is related to the
+mean of the random variables X[sub i] by:
+
+[equation nc_chi_squ_ref2]
+
+(Note that some references define [lambda] as one half of the above sum).
+
+This leads to a PDF of:
+
+[equation nc_chi_squ_ref3]
+
+where ['f(x;k)] is the central chi-squared distribution PDF, and
+['I[sub v](x)] is a modified Bessel function of the first kind.
+
+The following graph illustrates how the distribution changes
+for different values of [lambda]:
+
+[graph nccs_pdf]
+
+[h4 Member Functions]
+
+ non_central_chi_squared_distribution(RealType v, RealType lambda);
+
+Constructs a Chi-Squared distribution with /v/ degrees of freedom
+and non-centrality parameter /lambda/.
+
+Requires v > 0 and lambda >= 0, 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 parameter /lambda/ from which this object was constructed.
+
+
+ static RealType find_degrees_of_freedom(RealType lambda, RealType x, RealType p);
+
+This function returns the number of degrees of freedom /v/ such that:
+`cdf(non_central_chi_squared<RealType, Policy>(v, lambda), x) == p`
+
+ template <class A, class B, class C>
+ static RealType find_degrees_of_freedom(const complemented3_type<A,B,C>& c);
+
+When called with argument `boost::math::complement(lambda, x, q)`
+this function returns the number of degrees of freedom /v/ such that:
+
+`cdf(complement(non_central_chi_squared<RealType, Policy>(v, lambda), x)) == q`.
+
+ static RealType find_non_centrality(RealType v, RealType x, RealType p);
+
+This function returns the non centrality parameter /lambda/ such that:
+
+`cdf(non_central_chi_squared<RealType, Policy>(v, lambda), x) == p`
+
+ template <class A, class B, class C>
+ static RealType find_non_centrality(const complemented3_type<A,B,C>& c);
+
+When called with argument `boost::math::complement(v, x, q)`
+this function returns the non centrality parameter /lambda/ such that:
+
+`cdf(complement(non_central_chi_squared<RealType, Policy>(v, lambda), x)) == q`.
+
+[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 \[0, +[infin]\].
+
+[h4 Examples]
+
+There is a
+[link math_toolkit.stat_tut.weg.nccs_eg worked example]
+for the noncentral chi-squared distribution.
+
+[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.
+The failures in the comparison to the [@http://www.r-project.org/ R Math library],
+seem to be mostly in the corner cases when the probablity would be very small.
+Unless otherwise specified any floating-point type that is narrower
+than the one shown will have __zero_error.
+
+[table_non_central_chi_squared_CDF]
+
+[table_non_central_chi_squared_CDF_complement]
+
+Error rates for the quantile
+functions are broadly similar. Special mention should go to
+the `mode` function: there is no closed form for this function,
+so it is evaluated numerically by finding the maxima of the PDF:
+in principal this can not produce an accuracy greater than the
+square root of the machine epsilon.
+
+[h4 Tests]
+
+There are two sets of test data used to verify this implementation:
+firstly we can compare with published data, for example with
+Table 6 of "Self-Validating Computations of Probabilities for
+Selected Central and Noncentral Univariate Probability Functions",
+Morgan C. Wang and William J. Kennedy,
+Journal of the American Statistical Association,
+Vol. 89, No. 427. (Sep., 1994), pp. 878-887.
+Secondly, we have tables of test data, computed with this
+implementation and using interval arithmetic - this data should
+be accurate to at least 50 decimal digits - and is the used for
+our accuracy tests.
+
+[h4 Implementation]
+
+The CDF and its complement are evaluated as follows:
+
+First we determine which of the two values (the CDF or its
+complement) is likely to be the smaller: for this we can use the
+relation due to Temme (see "Asymptotic and Numerical Aspects of the
+Noncentral Chi-Square Distribution", N. M. Temme, Computers Math. Applic.
+Vol 25, No. 5, 55-63, 1993) that:
+
+F([nu],[lambda];[nu]+[lambda]) [asymp] 0.5
+
+and so compute the CDF when the random variable is less than
+[nu]+[lambda], and its complement when the random variable is
+greater than [nu]+[lambda]. If necessary the computed result
+is then subtracted from 1 to give the desired result (the CDF or its
+complement).
+
+For small values of the non centrality parameter, the CDF is computed
+using the method of Ding (see "Algorithm AS 275: Computing the Non-Central
+#2 Distribution Function", Cherng G. Ding, Applied Statistics, Vol. 41,
+No. 2. (1992), pp. 478-482). This uses the following series representation:
+
+[equation nc_chi_squ_ref4]
+
+which requires just one call to __gamma_p_derivative with the subsequent
+terms being computed by recursion as shown above.
+
+For larger values of the non-centrality parameter, Ding's method can take
+an unreasonable number of terms before convergence is achieved. Furthermore,
+the largest term is not the first term, so in extreme cases the first term may
+be zero, leading to a zero result, even though the true value may be non-zero.
+
+Therefore, when the non-centrality parameter is greater than 200, the method due
+to Krishnamoorthy (see "Computing discrete mixtures of continuous distributions:
+noncentral chisquare, noncentral t and the distribution of the
+square of the sample multiple correlation coefficient",
+Denise Benton and K. Krishnamoorthy, Computational Statistics &
+Data Analysis, 43, (2003), 249-267) is used.
+
+This method uses the well known sum:
+
+[equation nc_chi_squ_ref5]
+
+Where P[sub a](x) is the incomplete gamma function.
+
+The method starts at the [lambda]th term, which is where the Poisson weighting
+function achieves its maximum value, although this is not necessarily
+the largest overall term. Subsequent terms are calculated via the normal
+recurrence relations for the incomplete gamma function, and iteration proceeds
+both forwards and backwards until sufficient precision has been achieved. It
+should be noted that recurrence in the forwards direction of P[sub a](x) is
+numerically unstable. However, since we always start /after/ the largest
+term in the series, numeric instability is introduced more slowly than the
+series converges.
+
+Computation of the complement of the CDF uses an extension of Krishnamoorthy's
+method, given that:
+
+[equation nc_chi_squ_ref6]
+
+we can again start at the [lambda]'th term and proceed in both directions from
+there until the required precision is achieved. This time it is backwards
+recursion on the incomplete gamma function Q[sub a](x) which is unstable.
+However, as long as we start well /before/ the largest term, this is not an
+issue in practice.
+
+The PDF is computed directly using the relation:
+
+[equation nc_chi_squ_ref3]
+
+Where ['f(x; v)] is the PDF of the central __chi_squared_distrib and
+['I[sub v](x)] is a modified Bessel function, see __cyl_bessel_i.
+For small values of the
+non-centrality parameter the relation in terms of __cyl_bessel_i
+is used. However, this method fails for large values of the
+non-centrality parameter, so in that case the infinite sum is
+evaluated using the method of Benton and Krishnamoorthy, and
+the usual recurrence relations for successive terms.
+
+The quantile functions are computed by numeric inversion of the CDF.
+An improve starting quess is from
+Thomas Luu,
+[@http://discovery.ucl.ac.uk/1482128/, Fast and accurate parallel computation of quantile functions for random number generation, Doctorial Thesis, 2016].
+
+There is no [@http://en.wikipedia.org/wiki/Closed_form closed form]
+for the mode of the noncentral chi-squared
+distribution: it is computed numerically by finding the maximum
+of the PDF. Likewise, the median is computed numerically via
+the quantile.
+
+The remaining non-member functions use the following formulas:
+
+[equation nc_chi_squ_ref7]
+
+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_chi_squared_dist]
+
+[/ nc_chi_squared.qbk
+ Copyright 2008 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).
+]
+