    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     */
    017    package org.apache.commons.math.distribution;
    018
    019    import java.io.Serializable;
    020
    021    import org.apache.commons.math.MathException;
    022    import org.apache.commons.math.MathRuntimeException;
    023    import org.apache.commons.math.exception.util.LocalizedFormats;
    024    import org.apache.commons.math.util.FastMath;
    025
    026    /**
    027     * The default implementation of {@link ExponentialDistribution}.
    028     *
    029     * @version $Revision: 1055914 $ $Date: 2011-01-06 16:34:34 +0100 (jeu. 06 janv. 2011) $
    030     */
    031    public class ExponentialDistributionImpl extends AbstractContinuousDistribution
    032        implements ExponentialDistribution, Serializable {
    033
    034        /**
    035         * Default inverse cumulative probability accuracy
    036         * @since 2.1
    037         */
    038        public static final double DEFAULT_INVERSE_ABSOLUTE_ACCURACY = 1e-9;
    039
    040        /** Serializable version identifier */
    041        private static final long serialVersionUID = 2401296428283614780L;
    042
    043        /** The mean of this distribution. */
    044        private double mean;
    045
    046        /** Inverse cumulative probability accuracy */
    047        private final double solverAbsoluteAccuracy;
    048
    049        /**
    050         * Create a exponential distribution with the given mean.
    051         * @param mean mean of this distribution.
    052         */
    053        public ExponentialDistributionImpl(double mean) {
    054            this(mean, DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
    055        }
    056
    057        /**
    058         * Create a exponential distribution with the given mean.
    059         * @param mean mean of this distribution.
    060         * @param inverseCumAccuracy the maximum absolute error in inverse cumulative probability estimates
    061         * (defaults to {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY})
    062         * @since 2.1
    063         */
    064        public ExponentialDistributionImpl(double mean, double inverseCumAccuracy) {
    065            super();
    066            setMeanInternal(mean);
    067            solverAbsoluteAccuracy = inverseCumAccuracy;
    068        }
    069
    070        /**
    071         * Modify the mean.
    072         * @param mean the new mean.
    073         * @throws IllegalArgumentException if <code>mean</code> is not positive.
    074         * @deprecated as of 2.1 (class will become immutable in 3.0)
    075         */
    076        @Deprecated
    077        public void setMean(double mean) {
    078            setMeanInternal(mean);
    079        }
    080        /**
    081         * Modify the mean.
    082         * @param newMean the new mean.
    083         * @throws IllegalArgumentException if <code>newMean</code> is not positive.
    084         */
    085        private void setMeanInternal(double newMean) {
    086            if (newMean <= 0.0) {
    087                throw MathRuntimeException.createIllegalArgumentException(
    088                      LocalizedFormats.NOT_POSITIVE_MEAN, newMean);
    089            }
    090            this.mean = newMean;
    091        }
    092
    093        /**
    094         * Access the mean.
    095         * @return the mean.
    096         */
    097        public double getMean() {
    098            return mean;
    099        }
    100
    101        /**
    102         * Return the probability density for a particular point.
    103         *
    104         * @param x The point at which the density should be computed.
    105         * @return The pdf at point x.
    106         * @deprecated - use density(double)
    107         */
    108        @Deprecated
    109        public double density(Double x) {
    110            return density(x.doubleValue());
    111        }
    112
    113        /**
    114         * Return the probability density for a particular point.
    115         *
    116         * @param x The point at which the density should be computed.
    117         * @return The pdf at point x.
    118         * @since 2.1
    119         */
    120        @Override
    121        public double density(double x) {
    122            if (x < 0) {
    123                return 0;
    124            }
    125            return FastMath.exp(-x / mean) / mean;
    126        }
    127
    128        /**
    129         * For this distribution, X, this method returns P(X &lt; x).
    130         *
    131         * The implementation of this method is based on:
    132         * <ul>
    133         * <li>
    134         * <a href="http://mathworld.wolfram.com/ExponentialDistribution.html">
    135         * Exponential Distribution</a>, equation (1).</li>
    136         * </ul>
    137         *
    138         * @param x the value at which the CDF is evaluated.
    139         * @return CDF for this distribution.
    140         * @throws MathException if the cumulative probability can not be
    141         *            computed due to convergence or other numerical errors.
    142         */
    143        public double cumulativeProbability(double x) throws MathException{
    144            double ret;
    145            if (x <= 0.0) {
    146                ret = 0.0;
    147            } else {
    148                ret = 1.0 - FastMath.exp(-x / mean);
    149            }
    150            return ret;
    151        }
    152
    153        /**
    154         * For this distribution, X, this method returns the critical point x, such
    155         * that P(X &lt; x) = <code>p</code>.
    156         * <p>
    157         * Returns 0 for p=0 and <code>Double.POSITIVE_INFINITY</code> for p=1.</p>
    158         *
    159         * @param p the desired probability
    160         * @return x, such that P(X &lt; x) = <code>p</code>
    161         * @throws MathException if the inverse cumulative probability can not be
    162         *            computed due to convergence or other numerical errors.
    163         * @throws IllegalArgumentException if p < 0 or p > 1.
    164         */
    165        @Override
    166        public double inverseCumulativeProbability(double p) throws MathException {
    167            double ret;
    168
    169            if (p < 0.0 || p > 1.0) {
    170                throw MathRuntimeException.createIllegalArgumentException(
    171                      LocalizedFormats.OUT_OF_RANGE_SIMPLE, p, 0.0, 1.0);
    172            } else if (p == 1.0) {
    173                ret = Double.POSITIVE_INFINITY;
    174            } else {
    175                ret = -mean * FastMath.log(1.0 - p);
    176            }
    177
    178            return ret;
    179        }
    180
    181        /**
    182         * Generates a random value sampled from this distribution.
    183         *
    184         * <p><strong>Algorithm Description</strong>: Uses the <a
    185         * href="http://www.jesus.ox.ac.uk/~clifford/a5/chap1/node5.html"> Inversion
    186         * Method</a> to generate exponentially distributed random values from
    187         * uniform deviates. </p>
    188         *
    189         * @return random value
    190         * @since 2.2
    191         * @throws MathException if an error occurs generating the random value
    192         */
    193        @Override
    194        public double sample() throws MathException {
    195            return randomData.nextExponential(mean);
    196        }
    197
    198        /**
    199         * Access the domain value lower bound, based on <code>p</code>, used to
    200         * bracket a CDF root.
    201         *
    202         * @param p the desired probability for the critical value
    203         * @return domain value lower bound, i.e.
    204         *         P(X &lt; <i>lower bound</i>) &lt; <code>p</code>
    205         */
    206        @Override
    207        protected double getDomainLowerBound(double p) {
    208            return 0;
    209        }
    210
    211        /**
    212         * Access the domain value upper bound, based on <code>p</code>, used to
    213         * bracket a CDF root.
    214         *
    215         * @param p the desired probability for the critical value
    216         * @return domain value upper bound, i.e.
    217         *         P(X &lt; <i>upper bound</i>) &gt; <code>p</code>
    218         */
    219        @Override
    220        protected double getDomainUpperBound(double p) {
    221            // NOTE: exponential is skewed to the left
    222            // NOTE: therefore, P(X < &mu;) > .5
    223
    224            if (p < .5) {
    225                // use mean
    226                return mean;
    227            } else {
    228                // use max
    229                return Double.MAX_VALUE;
    230            }
    231        }
    232
    233        /**
    234         * Access the initial domain value, based on <code>p</code>, used to
    235         * bracket a CDF root.
    236         *
    237         * @param p the desired probability for the critical value
    238         * @return initial domain value
    239         */
    240        @Override
    241        protected double getInitialDomain(double p) {
    242            // TODO: try to improve on this estimate
    243            // TODO: what should really happen here is not derive from AbstractContinuousDistribution
    244            // TODO: because the inverse cumulative distribution is simple.
    245            // Exponential is skewed to the left, therefore, P(X < &mu;) > .5
    246            if (p < .5) {
    247                // use 1/2 mean
    248                return mean * .5;
    249            } else {
    250                // use mean
    251                return mean;
    252            }
    253        }
    254
    255        /**
    256         * Return the absolute accuracy setting of the solver used to estimate
    257         * inverse cumulative probabilities.
    258         *
    259         * @return the solver absolute accuracy
    260         * @since 2.1
    261         */
    262        @Override
    263        protected double getSolverAbsoluteAccuracy() {
    264            return solverAbsoluteAccuracy;
    265        }
    266
    267        /**
    268         * Returns the lower bound of the support for the distribution.
    269         *
    270         * The lower bound of the support is always 0, regardless of the mean.
    271         *
    272         * @return lower bound of the support (always 0)
    273         * @since 2.2
    274         */
    275        public double getSupportLowerBound() {
    276            return 0;
    277        }
    278
    279        /**
    280         * Returns the upper bound of the support for the distribution.
    281         *
    282         * The upper bound of the support is always positive infinity,
    283         * regardless of the mean.
    284         *
    285         * @return upper bound of the support (always Double.POSITIVE_INFINITY)
    286         * @since 2.2
    287         */
    288        public double getSupportUpperBound() {
    289            return Double.POSITIVE_INFINITY;
    290        }
    291
    292        /**
    293         * Returns the mean of the distribution.
    294         *
    295         * For mean parameter <code>k</code>, the mean is
    296         * <code>k</code>
    297         *
    298         * @return the mean
    299         * @since 2.2
    300         */
    301        public double getNumericalMean() {
    302            return getMean();
    303        }
    304
    305        /**
    306         * Returns the variance of the distribution.
    307         *
    308         * For mean parameter <code>k</code>, the variance is
    309         * <code>k^2</code>
    310         *
    311         * @return the variance
    312         * @since 2.2
    313         */
    314        public double getNumericalVariance() {
    315            final double m = getMean();
    316            return m * m;
    317        }
    318
    319    }