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    *      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   */
17  
18  /*
19   * This is not the original file distributed by the Apache Software Foundation
20   * It has been modified by the Hipparchus project
21   */
22  package org.hipparchus.distribution.continuous;
23  
24  import org.hipparchus.exception.LocalizedCoreFormats;
25  import org.hipparchus.exception.MathIllegalArgumentException;
26  import org.hipparchus.special.Beta;
27  import org.hipparchus.special.Gamma;
28  import org.hipparchus.util.FastMath;
29  
30  /**
31   * Implementation of Student's t-distribution.
32   */
33  public class TDistribution extends AbstractRealDistribution {
34      /** Serializable version identifier */
35      private static final long serialVersionUID = 20160320L;
36      /** The degrees of freedom. */
37      private final double degreesOfFreedom;
38      /** Static computation factor based on degreesOfFreedom. */
39      private final double factor;
40  
41      /**
42       * Create a t distribution using the given degrees of freedom.
43       *
44       * @param degreesOfFreedom Degrees of freedom.
45       * @throws MathIllegalArgumentException if {@code degreesOfFreedom <= 0}
46       */
47      public TDistribution(double degreesOfFreedom)
48          throws MathIllegalArgumentException {
49          this(degreesOfFreedom, DEFAULT_SOLVER_ABSOLUTE_ACCURACY);
50      }
51  
52      /**
53       * Create a t distribution using the given degrees of freedom and the
54       * specified inverse cumulative probability absolute accuracy.
55       *
56       * @param degreesOfFreedom Degrees of freedom.
57       * @param inverseCumAccuracy the maximum absolute error in inverse
58       * cumulative probability estimates
59       * (defaults to {@link #DEFAULT_SOLVER_ABSOLUTE_ACCURACY}).
60       * @throws MathIllegalArgumentException if {@code degreesOfFreedom <= 0}
61       */
62      public TDistribution(double degreesOfFreedom, double inverseCumAccuracy)
63          throws MathIllegalArgumentException {
64          super(inverseCumAccuracy);
65  
66          if (degreesOfFreedom <= 0) {
67              throw new MathIllegalArgumentException(LocalizedCoreFormats.DEGREES_OF_FREEDOM,
68                                                     degreesOfFreedom);
69          }
70          this.degreesOfFreedom = degreesOfFreedom;
71  
72          final double n = degreesOfFreedom;
73          final double nPlus1Over2 = (n + 1) / 2;
74          factor = Gamma.logGamma(nPlus1Over2) -
75                   0.5 * (FastMath.log(FastMath.PI) + FastMath.log(n)) -
76                   Gamma.logGamma(n / 2);
77      }
78  
79      /**
80       * Access the degrees of freedom.
81       *
82       * @return the degrees of freedom.
83       */
84      public double getDegreesOfFreedom() {
85          return degreesOfFreedom;
86      }
87  
88      /** {@inheritDoc} */
89      @Override
90      public double density(double x) {
91          return FastMath.exp(logDensity(x));
92      }
93  
94      /** {@inheritDoc} */
95      @Override
96      public double logDensity(double x) {
97          final double n = degreesOfFreedom;
98          final double nPlus1Over2 = (n + 1) / 2;
99          return factor - nPlus1Over2 * FastMath.log(1 + x * x / n);
100     }
101 
102     /** {@inheritDoc} */
103     @Override
104     public double cumulativeProbability(double x) {
105         double ret;
106         if (x == 0) {
107             ret = 0.5;
108         } else {
109             double t =
110                 Beta.regularizedBeta(
111                     degreesOfFreedom / (degreesOfFreedom + (x * x)),
112                     0.5 * degreesOfFreedom,
113                     0.5);
114             if (x < 0.0) {
115                 ret = 0.5 * t;
116             } else {
117                 ret = 1.0 - 0.5 * t;
118             }
119         }
120 
121         return ret;
122     }
123 
124     /**
125      * {@inheritDoc}
126      *
127      * For degrees of freedom parameter {@code df}, the mean is
128      * <ul>
129      *  <li>if {@code df > 1} then {@code 0},</li>
130      * <li>else undefined ({@code Double.NaN}).</li>
131      * </ul>
132      */
133     @Override
134     public double getNumericalMean() {
135         final double df = getDegreesOfFreedom();
136 
137         if (df > 1) {
138             return 0;
139         }
140 
141         return Double.NaN;
142     }
143 
144     /**
145      * {@inheritDoc}
146      *
147      * For degrees of freedom parameter {@code df}, the variance is
148      * <ul>
149      *  <li>if {@code df > 2} then {@code df / (df - 2)},</li>
150      *  <li>if {@code 1 < df <= 2} then positive infinity
151      *  ({@code Double.POSITIVE_INFINITY}),</li>
152      *  <li>else undefined ({@code Double.NaN}).</li>
153      * </ul>
154      */
155     @Override
156     public double getNumericalVariance() {
157         final double df = getDegreesOfFreedom();
158 
159         if (df > 2) {
160             return df / (df - 2);
161         }
162 
163         if (df > 1 && df <= 2) {
164             return Double.POSITIVE_INFINITY;
165         }
166 
167         return Double.NaN;
168     }
169 
170     /**
171      * {@inheritDoc}
172      *
173      * The lower bound of the support is always negative infinity no matter the
174      * parameters.
175      *
176      * @return lower bound of the support (always
177      * {@code Double.NEGATIVE_INFINITY})
178      */
179     @Override
180     public double getSupportLowerBound() {
181         return Double.NEGATIVE_INFINITY;
182     }
183 
184     /**
185      * {@inheritDoc}
186      *
187      * The upper bound of the support is always positive infinity no matter the
188      * parameters.
189      *
190      * @return upper bound of the support (always
191      * {@code Double.POSITIVE_INFINITY})
192      */
193     @Override
194     public double getSupportUpperBound() {
195         return Double.POSITIVE_INFINITY;
196     }
197 
198     /**
199      * {@inheritDoc}
200      *
201      * The support of this distribution is connected.
202      *
203      * @return {@code true}
204      */
205     @Override
206     public boolean isSupportConnected() {
207         return true;
208     }
209 }