| Line 1... |
Line 1... |
| 1 |
/*
|
1 |
/*
|
| 2 |
* Mathlib : A C Library of Special Functions
|
2 |
* Mathlib : A C Library of Special Functions
|
| 3 |
* Copyright (C) 1998 Ross Ihaka
|
- |
|
| 4 |
* Copyright (C) 2000--2015 The R Core Team
|
3 |
* Copyright (C) 2000--2015 The R Core Team
|
| - |
|
4 |
* Copyright (C) 1998 Ross Ihaka
|
| 5 |
* Copyright (C) 2004--2015 The R Foundation
|
5 |
* Copyright (C) 2004--2015 The R Foundation
|
| 6 |
* based on AS 91 (C) 1979 Royal Statistical Society
|
6 |
* originally based on AS 91 (C) 1979 Royal Statistical Society
|
| 7 |
*
|
7 |
*
|
| 8 |
* This program is free software; you can redistribute it and/or modify
|
8 |
* This program is free software; you can redistribute it and/or modify
|
| 9 |
* it under the terms of the GNU General Public License as published by
|
9 |
* it under the terms of the GNU General Public License as published by
|
| 10 |
* the Free Software Foundation; either version 2 of the License, or
|
10 |
* the Free Software Foundation; either version 2 of the License, or
|
| 11 |
* (at your option) any later version.
|
11 |
* (at your option) any later version.
|
| Line 24... |
Line 24... |
| 24 |
* Compute the quantile function of the gamma distribution.
|
24 |
* Compute the quantile function of the gamma distribution.
|
| 25 |
*
|
25 |
*
|
| 26 |
* NOTES
|
26 |
* NOTES
|
| 27 |
*
|
27 |
*
|
| 28 |
* This function is based on the Applied Statistics
|
28 |
* This function is based on the Applied Statistics
|
| 29 |
* Algorithm AS 91 ("ppchi2") and via pgamma(.) AS 239.
|
29 |
* Algorithm AS 91 ("ppchi2") and originally (no longer since 2005) via pgamma(.) AS 239.
|
| 30 |
*
|
30 |
*
|
| 31 |
* R core improvements:
|
31 |
* R core improvements:
|
| 32 |
* o lower_tail, log_p
|
32 |
* o lower_tail, log_p
|
| 33 |
* o non-trivial result for p outside [0.000002, 0.999998]
|
33 |
* o non-trivial result for p outside [0.000002, 0.999998]
|
| 34 |
* o p ~ 1 no longer gives +Inf; final Newton step(s)
|
34 |
* o p ~ 1 no longer gives +Inf; final Newton step(s)
|
| Line 117... |
Line 117... |
| 117 |
/* shape = alpha */
|
117 |
/* shape = alpha */
|
| 118 |
{
|
118 |
{
|
| 119 |
#define EPS1 1e-2
|
119 |
#define EPS1 1e-2
|
| 120 |
#define EPS2 5e-7/* final precision of AS 91 */
|
120 |
#define EPS2 5e-7/* final precision of AS 91 */
|
| 121 |
#define EPS_N 1e-15/* precision of Newton step / iterations */
|
121 |
#define EPS_N 1e-15/* precision of Newton step / iterations */
|
| 122 |
#define LN_EPS -36.043653389117156 /* = log(.Machine$double.eps) iff IEEE_754 */
|
122 |
//define LN_EPS -36.043653389117156 /* = log(.Machine$double.eps) iff IEEE_754 */
|
| 123 |
|
123 |
|
| 124 |
#define MAXIT 1000/* was 20 */
|
124 |
#define MAXIT 1000/* was 20 */
|
| 125 |
|
125 |
|
| 126 |
#define pMIN 1e-100 /* was 0.000002 = 2e-6 */
|
126 |
#define pMIN 1e-100 /* was 0.000002 = 2e-6 */
|
| 127 |
#define pMAX (1-1e-14)/* was (1-1e-12) and 0.999998 = 1 - 2e-6 */
|
127 |
#define pMAX (1-1e-14)/* was (1-1e-12) and 0.999998 = 1 - 2e-6 */
|