mirror of
				https://github.com/saitohirga/WSJT-X.git
				synced 2025-10-31 04:50:34 -04:00 
			
		
		
		
	
		
			
	
	
		
			137 lines
		
	
	
		
			4.8 KiB
		
	
	
	
		
			Plaintext
		
	
	
	
	
	
		
		
			
		
	
	
			137 lines
		
	
	
		
			4.8 KiB
		
	
	
	
		
			Plaintext
		
	
	
	
	
	
|  | [section:owens_t Owen's T function] | ||
|  | 
 | ||
|  | [h4 Synopsis] | ||
|  | 
 | ||
|  | `` | ||
|  | #include <boost/math/special_functions/owens_t.hpp> | ||
|  | `` | ||
|  | 
 | ||
|  |   namespace boost{ namespace math{ | ||
|  |    | ||
|  |   template <class T> | ||
|  |   ``__sf_result`` owens_t(T h, T a); | ||
|  |    | ||
|  |   template <class T, class ``__Policy``> | ||
|  |   ``__sf_result`` owens_t(T h, T a, const ``__Policy``&); | ||
|  |    | ||
|  |   }} // namespaces | ||
|  |    | ||
|  | [h4 Description] | ||
|  | 
 | ||
|  | Returns the | ||
|  | [@http://en.wikipedia.org/wiki/Owen%27s_T_function Owens_t function] | ||
|  | of ['h] and ['a]. | ||
|  | 
 | ||
|  | [optional_policy] | ||
|  | 
 | ||
|  | [sixemspace][sixemspace][equation owens_t] | ||
|  | 
 | ||
|  | [$../graphs/plot_owens_t.png] | ||
|  | 
 | ||
|  | The function `owens_t(h, a)` gives the probability | ||
|  | of the event ['(X > h and 0 < Y < a * X)], | ||
|  | where ['X] and ['Y] are independent standard normal random variables. | ||
|  | 
 | ||
|  | For h and a > 0, T(h,a), | ||
|  | gives the volume of an uncorrelated bivariate normal distribution | ||
|  | with zero means and unit variances over the area between | ||
|  | ['y = ax] and ['y = 0] and to the right of ['x = h]. | ||
|  | 
 | ||
|  | That is the area shaded in the figure below (Owens 1956). | ||
|  | 
 | ||
|  | [graph owens_integration_area] | ||
|  | 
 | ||
|  | and is also illustrated by a 3D plot. | ||
|  | 
 | ||
|  | [$../graphs/plot_owens_3d_xyp.png] | ||
|  | 
 | ||
|  | This function is used in the computation of the __skew_normal_distrib. | ||
|  | It is also used in the computation of bivariate and | ||
|  | multivariate normal distribution probabilities. | ||
|  | The return type of this function is computed using the __arg_promotion_rules: | ||
|  | the result is of type `double` when T is an integer type, and type T otherwise. | ||
|  | 
 | ||
|  | Owen's original paper (page 1077) provides some additional corner cases. | ||
|  | 
 | ||
|  | [: ['T(h, 0) = 0]] | ||
|  | 
 | ||
|  | [:['T(0, a) = [frac12][pi] arctan(a)]] | ||
|  | 
 | ||
|  | [:['T(h, 1) = [frac12] G(h) \[1 - G(h)\]]] | ||
|  | 
 | ||
|  | [:['T(h, [infin]) = G(|h|)]] | ||
|  | 
 | ||
|  | where G(h) is the univariate normal with zero mean and unit variance integral from -[infin] to h.   | ||
|  | 
 | ||
|  | [h4 Accuracy] | ||
|  | 
 | ||
|  | Over the built-in types and range tested, | ||
|  | errors are less than 10 * std::numeric_limits<RealType>::epsilon(). | ||
|  | 
 | ||
|  | [table_owens_t] | ||
|  | 
 | ||
|  | [h4 Testing] | ||
|  | 
 | ||
|  | Test data was generated by Patefield and Tandy algorithms T1 and T4, | ||
|  | and also the suggested reference routine T7. | ||
|  | 
 | ||
|  | * T1 was rejected if the result was too small compared to `atan(a)` (ie cancellation), | ||
|  | * T4 was rejected if there was no convergence, | ||
|  | * Both were rejected if they didn't agree. | ||
|  | 
 | ||
|  | Over the built-in types and range tested, | ||
|  | errors are less than 10 std::numeric_limits<RealType>::epsilon(). | ||
|  | 
 | ||
|  | However, that there was a whole domain (large ['h], small ['a]) | ||
|  | where it was not possible to generate any reliable test values | ||
|  | (all the methods got rejected for one reason or another). | ||
|  | 
 | ||
|  | There are also two sets of sanity tests: spot values are computed using __Mathematica and __R. | ||
|  | 
 | ||
|  | 
 | ||
|  | [h4 Implementation] | ||
|  | 
 | ||
|  | The function was proposed and evaluated by | ||
|  | [@http://projecteuclid.org/DPubS?service=UI&version=1.0&verb=Display&handle=euclid.aoms/1177728074 | ||
|  | Donald. B. Owen, Tables for computing bivariate normal probabilities,  | ||
|  | Ann. Math. Statist., 27, 1075-1090 (1956)]. | ||
|  | 
 | ||
|  | The algorithms of Patefield, M. and Tandy, D. | ||
|  | "Fast and accurate Calculation of Owen's T-Function", Journal of Statistical Software, 5 (5), 1 - 25 (2000) | ||
|  | are adapted for C++ with arbitrary RealType. | ||
|  | 
 | ||
|  | The Patefield-Tandy algorithm provides six methods of evalualution (T1 to T6); | ||
|  | the best method is selected according to the values of ['a] and ['h]. | ||
|  | See the original paper and the source in | ||
|  | [@../../../../boost/math/special_functions/owens_t.hpp owens_t.hpp] for details. | ||
|  | 
 | ||
|  | The Patefield-Tandy algorithm is accurate to approximately 20 decimal places, so for | ||
|  | types with greater precision we use: | ||
|  | 
 | ||
|  | * A modified version of T1 which folds the calculation of ['atan(h)] into the T1 series | ||
|  | (to avoid subtracting two values similar in magnitude), and then accelerates the | ||
|  | resulting alternating series using method 1 from H. Cohen, F. Rodriguez Villegas, D. Zagier,  | ||
|  | "Convergence acceleration of alternating series", Bonn, (1991).  The result is valid everywhere, | ||
|  | but doesn't always converge, or may become too divergent in the first few terms to sum accurately. | ||
|  | This is used for ['ah < 1]. | ||
|  | * A modified version of T2 which is accelerated in the same manner as T1.  This is used for ['h > 1]. | ||
|  | * A version of T4 only when both T1 and T2 have failed to produce an accurate answer. | ||
|  | * Fallback to the Patefiled Tandy algorithm when all the above methods fail: this happens not at all | ||
|  | for our test data at 100 decimal digits precision.  However, there is a difficult area when  | ||
|  | ['a] is very close to 1 and the precision increases which may cause this to happen in very exceptional  | ||
|  | circumstances. | ||
|  | 
 | ||
|  | Using the above algorithm and a 100-decimal digit type, results accurate to 80 decimal places were obtained | ||
|  | in the difficult area where ['a] is close to 1, and greater than 95 decimal places elsewhere. | ||
|  | 
 | ||
|  | [endsect] [/section:owens_t The owens_t Function] | ||
|  | 
 | ||
|  | [/  | ||
|  |   Copyright 2012 Bejamin Sobotta, 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). | ||
|  | ] | ||
|  | 
 | ||
|  | 
 |