1 /* 2 * Copyright 2003-2004 The Apache Software Foundation. 3 * 4 * Licensed under the Apache License, Version 2.0 (the "License"); 5 * you may not use this file except in compliance with the License. 6 * You may obtain a copy of the License at 7 * 8 * http://www.apache.org/licenses/LICENSE-2.0 9 * 10 * Unless required by applicable law or agreed to in writing, software 11 * distributed under the License is distributed on an "AS IS" BASIS, 12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 13 * See the License for the specific language governing permissions and 14 * limitations under the License. 15 */ 16 package org.apache.commons.math.distribution; 17 18 import java.io.Serializable; 19 20 import org.apache.commons.math.MathException; 21 import org.apache.commons.math.special.Beta; 22 23 /** 24 * Default implementation of 25 * {@link org.apache.commons.math.distribution.FDistribution}. 26 * 27 * @version $Revision: 348519 $ $Date: 2005-11-23 12:12:18 -0700 (Wed, 23 Nov 2005) $ 28 */ 29 public class FDistributionImpl 30 extends AbstractContinuousDistribution 31 implements FDistribution, Serializable { 32 33 /** Serializable version identifier */ 34 private static final long serialVersionUID = -8516354193418641566L; 35 36 /** The numerator degrees of freedom*/ 37 private double numeratorDegreesOfFreedom; 38 39 /** The numerator degrees of freedom*/ 40 private double denominatorDegreesOfFreedom; 41 42 /** 43 * Create a F distribution using the given degrees of freedom. 44 * @param numeratorDegreesOfFreedom the numerator degrees of freedom. 45 * @param denominatorDegreesOfFreedom the denominator degrees of freedom. 46 */ 47 public FDistributionImpl(double numeratorDegreesOfFreedom, 48 double denominatorDegreesOfFreedom) { 49 super(); 50 setNumeratorDegreesOfFreedom(numeratorDegreesOfFreedom); 51 setDenominatorDegreesOfFreedom(denominatorDegreesOfFreedom); 52 } 53 54 /** 55 * For this disbution, X, this method returns P(X < x). 56 * 57 * The implementation of this method is based on: 58 * <ul> 59 * <li> 60 * <a href="http://mathworld.wolfram.com/F-Distribution.html"> 61 * F-Distribution</a>, equation (4).</li> 62 * </ul> 63 * 64 * @param x the value at which the CDF is evaluated. 65 * @return CDF for this distribution. 66 * @throws MathException if the cumulative probability can not be 67 * computed due to convergence or other numerical errors. 68 */ 69 public double cumulativeProbability(double x) throws MathException { 70 double ret; 71 if (x <= 0.0) { 72 ret = 0.0; 73 } else { 74 double n = getNumeratorDegreesOfFreedom(); 75 double m = getDenominatorDegreesOfFreedom(); 76 77 ret = Beta.regularizedBeta((n * x) / (m + n * x), 78 0.5 * n, 79 0.5 * m); 80 } 81 return ret; 82 } 83 84 /** 85 * For this distribution, X, this method returns the critical point x, such 86 * that P(X < x) = <code>p</code>. 87 * <p> 88 * Returns 0 for p=0 and <code>Double.POSITIVE_INFINITY</code> for p=1. 89 * 90 * @param p the desired probability 91 * @return x, such that P(X < x) = <code>p</code> 92 * @throws MathException if the inverse cumulative probability can not be 93 * computed due to convergence or other numerical errors. 94 * @throws IllegalArgumentException if <code>p</code> is not a valid 95 * probability. 96 */ 97 public double inverseCumulativeProbability(final double p) 98 throws MathException { 99 if (p == 0) { 100 return 0d; 101 } 102 if (p == 1) { 103 return Double.POSITIVE_INFINITY; 104 } 105 return super.inverseCumulativeProbability(p); 106 } 107 108 /** 109 * Access the domain value lower bound, based on <code>p</code>, used to 110 * bracket a CDF root. This method is used by 111 * {@link #inverseCumulativeProbability(double)} to find critical values. 112 * 113 * @param p the desired probability for the critical value 114 * @return domain value lower bound, i.e. 115 * P(X < <i>lower bound</i>) < <code>p</code> 116 */ 117 protected double getDomainLowerBound(double p) { 118 return 0.0; 119 } 120 121 /** 122 * Access the domain value upper bound, based on <code>p</code>, used to 123 * bracket a CDF root. This method is used by 124 * {@link #inverseCumulativeProbability(double)} to find critical values. 125 * 126 * @param p the desired probability for the critical value 127 * @return domain value upper bound, i.e. 128 * P(X < <i>upper bound</i>) > <code>p</code> 129 */ 130 protected double getDomainUpperBound(double p) { 131 return Double.MAX_VALUE; 132 } 133 134 /** 135 * Access the initial domain value, based on <code>p</code>, used to 136 * bracket a CDF root. This method is used by 137 * {@link #inverseCumulativeProbability(double)} to find critical values. 138 * 139 * @param p the desired probability for the critical value 140 * @return initial domain value 141 */ 142 protected double getInitialDomain(double p) { 143 return getDenominatorDegreesOfFreedom() / 144 (getDenominatorDegreesOfFreedom() - 2.0); 145 } 146 147 /** 148 * Modify the numerator degrees of freedom. 149 * @param degreesOfFreedom the new numerator degrees of freedom. 150 * @throws IllegalArgumentException if <code>degreesOfFreedom</code> is not 151 * positive. 152 */ 153 public void setNumeratorDegreesOfFreedom(double degreesOfFreedom) { 154 if (degreesOfFreedom <= 0.0) { 155 throw new IllegalArgumentException( 156 "degrees of freedom must be positive."); 157 } 158 this.numeratorDegreesOfFreedom = degreesOfFreedom; 159 } 160 161 /** 162 * Access the numerator degrees of freedom. 163 * @return the numerator degrees of freedom. 164 */ 165 public double getNumeratorDegreesOfFreedom() { 166 return numeratorDegreesOfFreedom; 167 } 168 169 /** 170 * Modify the denominator degrees of freedom. 171 * @param degreesOfFreedom the new denominator degrees of freedom. 172 * @throws IllegalArgumentException if <code>degreesOfFreedom</code> is not 173 * positive. 174 */ 175 public void setDenominatorDegreesOfFreedom(double degreesOfFreedom) { 176 if (degreesOfFreedom <= 0.0) { 177 throw new IllegalArgumentException( 178 "degrees of freedom must be positive."); 179 } 180 this.denominatorDegreesOfFreedom = degreesOfFreedom; 181 } 182 183 /** 184 * Access the denominator degrees of freedom. 185 * @return the denominator degrees of freedom. 186 */ 187 public double getDenominatorDegreesOfFreedom() { 188 return denominatorDegreesOfFreedom; 189 } 190 }