FDistribution.java

  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.  *      https://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.  */

  18. /*
  19.  * This is not the original file distributed by the Apache Software Foundation
  20.  * It has been modified by the Hipparchus project
  21.  */

  23. package org.hipparchus.distribution.continuous;

  25. import org.hipparchus.exception.LocalizedCoreFormats;
  26. import org.hipparchus.exception.MathIllegalArgumentException;
  27. import org.hipparchus.special.Beta;
  28. import org.hipparchus.util.FastMath;

  30. /**
  31.  * Implementation of the F-distribution.
  32.  *
  33.  * @see <a href="http://en.wikipedia.org/wiki/F-distribution">F-distribution (Wikipedia)</a>
  34.  * @see <a href="http://mathworld.wolfram.com/F-Distribution.html">F-distribution (MathWorld)</a>
  35.  */
  36. public class FDistribution extends AbstractRealDistribution {
  37.     /** Serializable version identifier. */
  38.     private static final long serialVersionUID = 20160320L;
  39.     /** The numerator degrees of freedom. */
  40.     private final double numeratorDegreesOfFreedom;
  41.     /** The numerator degrees of freedom. */
  42.     private final double denominatorDegreesOfFreedom;
  43.     /** Cached numerical variance */
  44.     private final double numericalVariance;

  46.     /**
  47.      * Creates an F distribution using the given degrees of freedom.
  48.      *
  49.      * @param numeratorDegreesOfFreedom Numerator degrees of freedom.
  50.      * @param denominatorDegreesOfFreedom Denominator degrees of freedom.
  51.      * @throws MathIllegalArgumentException if
  52.      * {@code numeratorDegreesOfFreedom <= 0} or
  53.      * {@code denominatorDegreesOfFreedom <= 0}.
  54.      */
  55.     public FDistribution(double numeratorDegreesOfFreedom,
  56.                          double denominatorDegreesOfFreedom)
  57.         throws MathIllegalArgumentException {
  58.         this(numeratorDegreesOfFreedom, denominatorDegreesOfFreedom,
  59.              DEFAULT_SOLVER_ABSOLUTE_ACCURACY);
  60.     }


  63.     /**
  64.      * Creates an F distribution.
  65.      *
  66.      * @param numeratorDegreesOfFreedom Numerator degrees of freedom.
  67.      * @param denominatorDegreesOfFreedom Denominator degrees of freedom.
  68.      * @param inverseCumAccuracy the maximum absolute error in inverse
  69.      * cumulative probability estimates.
  70.      * @throws MathIllegalArgumentException if {@code numeratorDegreesOfFreedom <= 0} or
  71.      * {@code denominatorDegreesOfFreedom <= 0}.
  72.      */
  73.     public FDistribution(double numeratorDegreesOfFreedom,
  74.                          double denominatorDegreesOfFreedom,
  75.                          double inverseCumAccuracy)
  76.         throws MathIllegalArgumentException {
  77.         super(inverseCumAccuracy);

  79.         if (numeratorDegreesOfFreedom <= 0) {
  80.             throw new MathIllegalArgumentException(LocalizedCoreFormats.DEGREES_OF_FREEDOM,
  81.                                                    numeratorDegreesOfFreedom);
  82.         }
  83.         if (denominatorDegreesOfFreedom <= 0) {
  84.             throw new MathIllegalArgumentException(LocalizedCoreFormats.DEGREES_OF_FREEDOM,
  85.                                                    denominatorDegreesOfFreedom);
  86.         }
  87.         this.numeratorDegreesOfFreedom   = numeratorDegreesOfFreedom;
  88.         this.denominatorDegreesOfFreedom = denominatorDegreesOfFreedom;
  89.         this.numericalVariance           = calculateNumericalVariance();
  90.     }

  92.     /**
  93.      * {@inheritDoc}
  94.      */
  95.     @Override
  96.     public double density(double x) {
  97.         return FastMath.exp(logDensity(x));
  98.     }

  100.     /** {@inheritDoc} **/
  101.     @Override
  102.     public double logDensity(double x) {
  103.         final double nhalf = numeratorDegreesOfFreedom / 2;
  104.         final double mhalf = denominatorDegreesOfFreedom / 2;
  105.         final double logx = FastMath.log(x);
  106.         final double logn = FastMath.log(numeratorDegreesOfFreedom);
  107.         final double logm = FastMath.log(denominatorDegreesOfFreedom);
  108.         final double lognxm = FastMath.log(numeratorDegreesOfFreedom * x +
  109.                                            denominatorDegreesOfFreedom);
  110.         return nhalf * logn + nhalf * logx - logx +
  111.                mhalf * logm - nhalf * lognxm - mhalf * lognxm -
  112.                Beta.logBeta(nhalf, mhalf);
  113.     }

  115.     /**
  116.      * {@inheritDoc}
  117.      *
  118.      * The implementation of this method is based on
  119.      * <ul>
  120.      *  <li>
  121.      *   <a href="http://mathworld.wolfram.com/F-Distribution.html">
  122.      *   F-Distribution</a>, equation (4).
  123.      *  </li>
  124.      * </ul>
  125.      */
  126.     @Override
  127.     public double cumulativeProbability(double x)  {
  128.         double ret;
  129.         if (x <= 0) {
  130.             ret = 0;
  131.         } else {
  132.             double n = numeratorDegreesOfFreedom;
  133.             double m = denominatorDegreesOfFreedom;

  135.             ret = Beta.regularizedBeta((n * x) / (m + n * x),
  136.                 0.5 * n,
  137.                 0.5 * m);
  138.         }
  139.         return ret;
  140.     }

  142.     /**
  143.      * Access the numerator degrees of freedom.
  144.      *
  145.      * @return the numerator degrees of freedom.
  146.      */
  147.     public double getNumeratorDegreesOfFreedom() {
  148.         return numeratorDegreesOfFreedom;
  149.     }

  151.     /**
  152.      * Access the denominator degrees of freedom.
  153.      *
  154.      * @return the denominator degrees of freedom.
  155.      */
  156.     public double getDenominatorDegreesOfFreedom() {
  157.         return denominatorDegreesOfFreedom;
  158.     }

  160.     /**
  161.      * {@inheritDoc}
  162.      *
  163.      * For denominator degrees of freedom parameter {@code b}, the mean is
  164.      * <ul>
  165.      *  <li>if {@code b > 2} then {@code b / (b - 2)},</li>
  166.      *  <li>else undefined ({@code Double.NaN}).
  167.      * </ul>
  168.      */
  169.     @Override
  170.     public double getNumericalMean() {
  171.         final double denominatorDF = getDenominatorDegreesOfFreedom();

  173.         if (denominatorDF > 2) {
  174.             return denominatorDF / (denominatorDF - 2);
  175.         }

  177.         return Double.NaN;
  178.     }

  180.     /**
  181.      * {@inheritDoc}
  182.      *
  183.      * For numerator degrees of freedom parameter {@code a} and denominator
  184.      * degrees of freedom parameter {@code b}, the variance is
  185.      * <ul>
  186.      *  <li>
  187.      *    if {@code b > 4} then
  188.      *    {@code [2 * b^2 * (a + b - 2)] / [a * (b - 2)^2 * (b - 4)]},
  189.      *  </li>
  190.      *  <li>else undefined ({@code Double.NaN}).
  191.      * </ul>
  192.      */
  193.     @Override
  194.     public double getNumericalVariance() {
  195.         return numericalVariance;
  196.     }

  198.     /**
  199.      * Calculates the numerical variance.
  200.      *
  201.      * @return the variance of this distribution
  202.      */
  203.     private double calculateNumericalVariance() {
  204.         final double denominatorDF = getDenominatorDegreesOfFreedom();

  206.         if (denominatorDF > 4) {
  207.             final double numeratorDF = getNumeratorDegreesOfFreedom();
  208.             final double denomDFMinusTwo = denominatorDF - 2;

  210.             return ( 2 * (denominatorDF * denominatorDF) * (numeratorDF + denominatorDF - 2) ) /
  211.                    ( (numeratorDF * (denomDFMinusTwo * denomDFMinusTwo) * (denominatorDF - 4)) );
  212.         }

  214.         return Double.NaN;
  215.     }

  217.     /**
  218.      * {@inheritDoc}
  219.      *
  220.      * The lower bound of the support is always 0 no matter the parameters.
  221.      *
  222.      * @return lower bound of the support (always 0)
  223.      */
  224.     @Override
  225.     public double getSupportLowerBound() {
  226.         return 0;
  227.     }

  229.     /**
  230.      * {@inheritDoc}
  231.      *
  232.      * The upper bound of the support is always positive infinity
  233.      * no matter the parameters.
  234.      *
  235.      * @return upper bound of the support (always Double.POSITIVE_INFINITY)
  236.      */
  237.     @Override
  238.     public double getSupportUpperBound() {
  239.         return Double.POSITIVE_INFINITY;
  240.     }

  242.     /**
  243.      * {@inheritDoc}
  244.      *
  245.      * The support of this distribution is connected.
  246.      *
  247.      * @return {@code true}
  248.      */
  249.     @Override
  250.     public boolean isSupportConnected() {
  251.         return true;
  252.     }
  253. }
Created with JaCoCo 0.8.12.202403310830