001/*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements.  See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License.  You may obtain a copy of the License at
008 *
009 *      http://www.apache.org/licenses/LICENSE-2.0
010 *
011 * Unless required by applicable law or agreed to in writing, software
012 * distributed under the License is distributed on an "AS IS" BASIS,
013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014 * See the License for the specific language governing permissions and
015 * limitations under the License.
016 */
017package org.apache.commons.statistics.distribution;
018
019import org.apache.commons.numbers.gamma.RegularizedBeta;
020import org.apache.commons.numbers.gamma.LogGamma;
021
022/**
023 * Implementation of <a href='http://en.wikipedia.org/wiki/Student&apos;s_t-distribution'>Student's t-distribution</a>.
024 */
025public class TDistribution extends AbstractContinuousDistribution {
026    /** The degrees of freedom. */
027    private final double degreesOfFreedom;
028    /** degreesOfFreedom / 2. */
029    private final double dofOver2;
030    /** Cached value. */
031    private final double factor;
032    /** Cached value. */
033    private final double mean;
034    /** Cached value. */
035    private final double variance;
036
037    /**
038     * Creates a distribution.
039     *
040     * @param degreesOfFreedom Degrees of freedom.
041     * @throws IllegalArgumentException if {@code degreesOfFreedom <= 0}
042     */
043    public TDistribution(double degreesOfFreedom) {
044        if (degreesOfFreedom <= 0) {
045            throw new DistributionException(DistributionException.NEGATIVE,
046                                            degreesOfFreedom);
047        }
048        this.degreesOfFreedom = degreesOfFreedom;
049
050        dofOver2 = 0.5 * degreesOfFreedom;
051        factor = LogGamma.value(dofOver2 + 0.5) -
052                 0.5 * (Math.log(Math.PI) + Math.log(degreesOfFreedom)) -
053                 LogGamma.value(dofOver2);
054        mean = degreesOfFreedom > 1 ? 0 :
055            Double.NaN;
056        variance = degreesOfFreedom > 2 ? degreesOfFreedom / (degreesOfFreedom - 2) :
057            degreesOfFreedom > 1 && degreesOfFreedom <= 2 ? Double.POSITIVE_INFINITY :
058            Double.NaN;
059    }
060
061    /**
062     * Access the degrees of freedom.
063     *
064     * @return the degrees of freedom.
065     */
066    public double getDegreesOfFreedom() {
067        return degreesOfFreedom;
068    }
069
070    /** {@inheritDoc} */
071    @Override
072    public double density(double x) {
073        return Math.exp(logDensity(x));
074    }
075
076    /** {@inheritDoc} */
077    @Override
078    public double logDensity(double x) {
079        final double nPlus1Over2 = dofOver2 + 0.5;
080        return factor - nPlus1Over2 * Math.log1p(x * x / degreesOfFreedom);
081    }
082
083    /** {@inheritDoc} */
084    @Override
085    public double cumulativeProbability(double x) {
086        if (x == 0) {
087            return 0.5;
088        } else {
089            final double t =
090                RegularizedBeta.value(degreesOfFreedom / (degreesOfFreedom + (x * x)),
091                                      dofOver2,
092                                      0.5);
093            return x < 0 ?
094                0.5 * t :
095                1 - 0.5 * t;
096        }
097    }
098
099    /**
100     * {@inheritDoc}
101     *
102     * For degrees of freedom parameter {@code df}, the mean is
103     * <ul>
104     *  <li>zero if {@code df > 1}, and</li>
105     *  <li>undefined ({@code Double.NaN}) otherwise.</li>
106     * </ul>
107     */
108    @Override
109    public double getMean() {
110        return mean;
111    }
112
113    /**
114     * {@inheritDoc}
115     *
116     * For degrees of freedom parameter {@code df}, the variance is
117     * <ul>
118     *  <li>{@code df / (df - 2)} if {@code df > 2},</li>
119     *  <li>infinite ({@code Double.POSITIVE_INFINITY}) if {@code 1 < df <= 2}, and</li>
120     *  <li>undefined ({@code Double.NaN}) otherwise.</li>
121     * </ul>
122     */
123    @Override
124    public double getVariance() {
125        return variance;
126    }
127
128    /**
129     * {@inheritDoc}
130     *
131     * The lower bound of the support is always negative infinity..
132     *
133     * @return lower bound of the support (always
134     * {@code Double.NEGATIVE_INFINITY})
135     */
136    @Override
137    public double getSupportLowerBound() {
138        return Double.NEGATIVE_INFINITY;
139    }
140
141    /**
142     * {@inheritDoc}
143     *
144     * The upper bound of the support is always positive infinity.
145     *
146     * @return upper bound of the support (always
147     * {@code Double.POSITIVE_INFINITY})
148     */
149    @Override
150    public double getSupportUpperBound() {
151        return Double.POSITIVE_INFINITY;
152    }
153
154    /**
155     * {@inheritDoc}
156     *
157     * The support of this distribution is connected.
158     *
159     * @return {@code true}
160     */
161    @Override
162    public boolean isSupportConnected() {
163        return true;
164    }
165}