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.exception.NotStrictlyPositiveException;
022 import org.apache.commons.math.exception.OutOfRangeException;
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 $Id: ExponentialDistributionImpl.java 1178295 2011-10-03 04:36:27Z psteitz $
030 */
031 public class ExponentialDistributionImpl extends AbstractContinuousDistribution
032 implements ExponentialDistribution, Serializable {
033 /**
034 * Default inverse cumulative probability accuracy.
035 * @since 2.1
036 */
037 public static final double DEFAULT_INVERSE_ABSOLUTE_ACCURACY = 1e-9;
038 /** Serializable version identifier */
039 private static final long serialVersionUID = 2401296428283614780L;
040 /** The mean of this distribution. */
041 private final double mean;
042 /** Inverse cumulative probability accuracy. */
043 private final double solverAbsoluteAccuracy;
044
045 /**
046 * Create a exponential distribution with the given mean.
047 * @param mean mean of this distribution.
048 */
049 public ExponentialDistributionImpl(double mean) {
050 this(mean, DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
051 }
052
053 /**
054 * Create a exponential distribution with the given mean.
055 *
056 * @param mean Mean of this distribution.
057 * @param inverseCumAccuracy Maximum absolute error in inverse
058 * cumulative probability estimates (defaults to
059 * {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY}).
060 * @throws NotStrictlyPositiveException if {@code mean <= 0}.
061 * @since 2.1
062 */
063 public ExponentialDistributionImpl(double mean, double inverseCumAccuracy) {
064 if (mean <= 0) {
065 throw new NotStrictlyPositiveException(LocalizedFormats.MEAN, mean);
066 }
067 this.mean = mean;
068 solverAbsoluteAccuracy = inverseCumAccuracy;
069 }
070
071 /**
072 * {@inheritDoc}
073 */
074 public double getMean() {
075 return mean;
076 }
077
078 /**
079 * {@inheritDoc}
080 */
081 @Override
082 public double density(double x) {
083 if (x < 0) {
084 return 0;
085 }
086 return FastMath.exp(-x / mean) / mean;
087 }
088
089 /**
090 * For this distribution, X, this method returns P(X < x).
091 *
092 * The implementation of this method is based on:
093 * <ul>
094 * <li>
095 * <a href="http://mathworld.wolfram.com/ExponentialDistribution.html">
096 * Exponential Distribution</a>, equation (1).</li>
097 * </ul>
098 *
099 * @param x Value at which the CDF is evaluated.
100 * @return the CDF for this distribution.
101 */
102 public double cumulativeProbability(double x) {
103 double ret;
104 if (x <= 0.0) {
105 ret = 0.0;
106 } else {
107 ret = 1.0 - FastMath.exp(-x / mean);
108 }
109 return ret;
110 }
111
112 /**
113 * For this distribution, X, this method returns the critical point x, such
114 * that {@code P(X < x) = p}.
115 * It will return 0 when p = 0 and {@code Double.POSITIVE_INFINITY}
116 * when p = 1.
117 *
118 * @param p Desired probability.
119 * @return {@code x}, such that {@code P(X < x) = p}.
120 * @throws OutOfRangeException if {@code p < 0} or {@code p > 1}.
121 */
122 @Override
123 public double inverseCumulativeProbability(double p) {
124 double ret;
125
126 if (p < 0.0 || p > 1.0) {
127 throw new OutOfRangeException(p, 0.0, 1.0);
128 } else if (p == 1.0) {
129 ret = Double.POSITIVE_INFINITY;
130 } else {
131 ret = -mean * FastMath.log(1.0 - p);
132 }
133
134 return ret;
135 }
136
137 /**
138 * Generates a random value sampled from this distribution.
139 *
140 * <p><strong>Algorithm Description</strong>: Uses the <a
141 * href="http://www.jesus.ox.ac.uk/~clifford/a5/chap1/node5.html"> Inversion
142 * Method</a> to generate exponentially distributed random values from
143 * uniform deviates.</p>
144 *
145 * @return a random value.
146 * @since 2.2
147 */
148 @Override
149 public double sample() {
150 return randomData.nextExponential(mean);
151 }
152
153 /**
154 * Access the domain value lower bound, based on {@code p}, used to
155 * bracket a CDF root.
156 *
157 * @param p Desired probability for the critical value.
158 * @return the domain value lower bound, i.e. {@code P(X < 'lower bound') < p}.
159 */
160 @Override
161 protected double getDomainLowerBound(double p) {
162 return 0;
163 }
164
165 /**
166 * Access the domain value upper bound, based on {@code p}, used to
167 * bracket a CDF root.
168 *
169 * @param p Desired probability for the critical value.
170 * @return the domain value upper bound, i.e. {@code P(X < 'upper bound') > p}.
171 */
172 @Override
173 protected double getDomainUpperBound(double p) {
174 // NOTE: exponential is skewed to the left
175 // NOTE: therefore, P(X < μ) > .5
176
177 if (p < 0.5) {
178 // use mean
179 return mean;
180 } else {
181 // use max
182 return Double.MAX_VALUE;
183 }
184 }
185
186 /**
187 * Access the initial domain value, based on {@code p}, used to
188 * bracket a CDF root.
189 *
190 * @param p Desired probability for the critical value.
191 * @return the initial domain value.
192 */
193 @Override
194 protected double getInitialDomain(double p) {
195 // TODO: try to improve on this estimate
196 // TODO: what should really happen here is not derive from AbstractContinuousDistribution
197 // TODO: because the inverse cumulative distribution is simple.
198 // Exponential is skewed to the left, therefore, P(X < μ) > .5
199 if (p < 0.5) {
200 // use 1/2 mean
201 return mean * 0.5;
202 } else {
203 // use mean
204 return mean;
205 }
206 }
207
208 /**
209 * Return the absolute accuracy setting of the solver used to estimate
210 * inverse cumulative probabilities.
211 *
212 * @return the solver absolute accuracy.
213 * @since 2.1
214 */
215 @Override
216 protected double getSolverAbsoluteAccuracy() {
217 return solverAbsoluteAccuracy;
218 }
219
220 /**
221 * {@inheritDoc}
222 *
223 * The lower bound of the support is always 0 no matter the mean parameter.
224 *
225 * @return lower bound of the support (always 0)
226 */
227 @Override
228 public double getSupportLowerBound() {
229 return 0;
230 }
231
232 /**
233 * {@inheritDoc}
234 *
235 * The upper bound of the support is always positive infinity
236 * no matter the mean parameter.
237 *
238 * @return upper bound of the support (always Double.POSITIVE_INFINITY)
239 */
240 @Override
241 public double getSupportUpperBound() {
242 return Double.POSITIVE_INFINITY;
243 }
244
245 /**
246 * {@inheritDoc}
247 *
248 * For mean parameter <code>k</code>, the mean is
249 * <code>k</code>
250 *
251 * @return {@inheritDoc}
252 */
253 @Override
254 protected double calculateNumericalMean() {
255 return getMean();
256 }
257
258 /**
259 * {@inheritDoc}
260 *
261 * For mean parameter <code>k</code>, the variance is
262 * <code>k^2</code>
263 *
264 * @return {@inheritDoc}
265 */
266 @Override
267 protected double calculateNumericalVariance() {
268 final double m = getMean();
269 return m * m;
270 }
271
272 /**
273 * {@inheritDoc}
274 */
275 @Override
276 public boolean isSupportLowerBoundInclusive() {
277 return true;
278 }
279
280 /**
281 * {@inheritDoc}
282 */
283 @Override
284 public boolean isSupportUpperBoundInclusive() {
285 return false;
286 }
287 }