1 /*
2 * Mathlib : A C Library of Special Functions
3 * Copyright (C) 1998 Ross Ihaka
4 * Copyright (C) 2000-8 The R Core Team
5 * Copyright (C) 2005 The R Foundation
6 *
7 * This program is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2 of the License, or
10 * (at your option) any later version.
11 *
12 * This program is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with this program; if not, a copy is available at
19 * http://www.r-project.org/Licenses/
20 *
21 * DESCRIPTION
22 *
23 * The quantile function of the F distribution.
24 */
25
26 #include "nmath.h"
27 #include "dpq.h"
28
qf(double p,double df1,double df2,int lower_tail,int log_p)29 double qf(double p, double df1, double df2, int lower_tail, int log_p)
30 {
31 #ifdef IEEE_754
32 if (ISNAN(p) || ISNAN(df1) || ISNAN(df2))
33 return p + df1 + df2;
34 #endif
35 if (df1 <= 0. || df2 <= 0.) ML_ERR_return_NAN;
36
37 R_Q_P01_boundaries(p, 0, ML_POSINF);
38
39 /* fudge the extreme DF cases -- qbeta doesn't do this well.
40 But we still need to fudge the infinite ones.
41 */
42
43 if (df1 <= df2 && df2 > 4e5) {
44 if(!R_FINITE(df1)) /* df1 == df2 == Inf : */
45 return 1.;
46 /* else */
47 return qchisq(p, df1, lower_tail, log_p) / df1;
48 }
49 if (df1 > 4e5) { /* and so df2 < df1 */
50 return df2 / qchisq(p, df2, !lower_tail, log_p);
51 }
52
53 p = (1. / qbeta(p, df2/2, df1/2, !lower_tail, log_p) - 1.) * (df2 / df1);
54 return ML_VALID(p) ? p : ML_NAN;
55 }
56