View Javadoc

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