http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/distribution/PascalDistribution.java ---------------------------------------------------------------------- diff --git a/src/main/java/org/apache/commons/math3/distribution/PascalDistribution.java b/src/main/java/org/apache/commons/math3/distribution/PascalDistribution.java deleted file mode 100644 index bf39933..0000000 --- a/src/main/java/org/apache/commons/math3/distribution/PascalDistribution.java +++ /dev/null @@ -1,248 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math3.distribution; - -import org.apache.commons.math3.exception.NotStrictlyPositiveException; -import org.apache.commons.math3.exception.OutOfRangeException; -import org.apache.commons.math3.exception.util.LocalizedFormats; -import org.apache.commons.math3.random.RandomGenerator; -import org.apache.commons.math3.random.Well19937c; -import org.apache.commons.math3.special.Beta; -import org.apache.commons.math3.util.CombinatoricsUtils; -import org.apache.commons.math3.util.FastMath; - -/** - * <p> - * Implementation of the Pascal distribution. The Pascal distribution is a - * special case of the Negative Binomial distribution where the number of - * successes parameter is an integer. - * </p> - * <p> - * There are various ways to express the probability mass and distribution - * functions for the Pascal distribution. The present implementation represents - * the distribution of the number of failures before {@code r} successes occur. - * This is the convention adopted in e.g. - * <a href="http://mathworld.wolfram.com/NegativeBinomialDistribution.html";>MathWorld</a>, - * but <em>not</em> in - * <a href="http://en.wikipedia.org/wiki/Negative_binomial_distribution";>Wikipedia</a>. - * </p> - * <p> - * For a random variable {@code X} whose values are distributed according to this - * distribution, the probability mass function is given by<br/> - * {@code P(X = k) = C(k + r - 1, r - 1) * p^r * (1 - p)^k,}<br/> - * where {@code r} is the number of successes, {@code p} is the probability of - * success, and {@code X} is the total number of failures. {@code C(n, k)} is - * the binomial coefficient ({@code n} choose {@code k}). The mean and variance - * of {@code X} are<br/> - * {@code E(X) = (1 - p) * r / p, var(X) = (1 - p) * r / p^2.}<br/> - * Finally, the cumulative distribution function is given by<br/> - * {@code P(X <= k) = I(p, r, k + 1)}, - * where I is the regularized incomplete Beta function. - * </p> - * - * @see <a href="http://en.wikipedia.org/wiki/Negative_binomial_distribution";> - * Negative binomial distribution (Wikipedia)</a> - * @see <a href="http://mathworld.wolfram.com/NegativeBinomialDistribution.html";> - * Negative binomial distribution (MathWorld)</a> - * @since 1.2 (changed to concrete class in 3.0) - */ -public class PascalDistribution extends AbstractIntegerDistribution { - /** Serializable version identifier. */ - private static final long serialVersionUID = 6751309484392813623L; - /** The number of successes. */ - private final int numberOfSuccesses; - /** The probability of success. */ - private final double probabilityOfSuccess; - /** The value of {@code log(p)}, where {@code p} is the probability of success, - * stored for faster computation. */ - private final double logProbabilityOfSuccess; - /** The value of {@code log(1-p)}, where {@code p} is the probability of success, - * stored for faster computation. */ - private final double log1mProbabilityOfSuccess; - - /** - * Create a Pascal distribution with the given number of successes and - * probability of success. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param r Number of successes. - * @param p Probability of success. - * @throws NotStrictlyPositiveException if the number of successes is not positive - * @throws OutOfRangeException if the probability of success is not in the - * range {@code [0, 1]}. - */ - public PascalDistribution(int r, double p) - throws NotStrictlyPositiveException, OutOfRangeException { - this(new Well19937c(), r, p); - } - - /** - * Create a Pascal distribution with the given number of successes and - * probability of success. - * - * @param rng Random number generator. - * @param r Number of successes. - * @param p Probability of success. - * @throws NotStrictlyPositiveException if the number of successes is not positive - * @throws OutOfRangeException if the probability of success is not in the - * range {@code [0, 1]}. - * @since 3.1 - */ - public PascalDistribution(RandomGenerator rng, - int r, - double p) - throws NotStrictlyPositiveException, OutOfRangeException { - super(rng); - - if (r <= 0) { - throw new NotStrictlyPositiveException(LocalizedFormats.NUMBER_OF_SUCCESSES, - r); - } - if (p < 0 || p > 1) { - throw new OutOfRangeException(p, 0, 1); - } - - numberOfSuccesses = r; - probabilityOfSuccess = p; - logProbabilityOfSuccess = FastMath.log(p); - log1mProbabilityOfSuccess = FastMath.log1p(-p); - } - - /** - * Access the number of successes for this distribution. - * - * @return the number of successes. - */ - public int getNumberOfSuccesses() { - return numberOfSuccesses; - } - - /** - * Access the probability of success for this distribution. - * - * @return the probability of success. - */ - public double getProbabilityOfSuccess() { - return probabilityOfSuccess; - } - - /** {@inheritDoc} */ - public double probability(int x) { - double ret; - if (x < 0) { - ret = 0.0; - } else { - ret = CombinatoricsUtils.binomialCoefficientDouble(x + - numberOfSuccesses - 1, numberOfSuccesses - 1) * - FastMath.pow(probabilityOfSuccess, numberOfSuccesses) * - FastMath.pow(1.0 - probabilityOfSuccess, x); - } - return ret; - } - - /** {@inheritDoc} */ - @Override - public double logProbability(int x) { - double ret; - if (x < 0) { - ret = Double.NEGATIVE_INFINITY; - } else { - ret = CombinatoricsUtils.binomialCoefficientLog(x + - numberOfSuccesses - 1, numberOfSuccesses - 1) + - logProbabilityOfSuccess * numberOfSuccesses + - log1mProbabilityOfSuccess * x; - } - return ret; - } - - /** {@inheritDoc} */ - public double cumulativeProbability(int x) { - double ret; - if (x < 0) { - ret = 0.0; - } else { - ret = Beta.regularizedBeta(probabilityOfSuccess, - numberOfSuccesses, x + 1.0); - } - return ret; - } - - /** - * {@inheritDoc} - * - * For number of successes {@code r} and probability of success {@code p}, - * the mean is {@code r * (1 - p) / p}. - */ - public double getNumericalMean() { - final double p = getProbabilityOfSuccess(); - final double r = getNumberOfSuccesses(); - return (r * (1 - p)) / p; - } - - /** - * {@inheritDoc} - * - * For number of successes {@code r} and probability of success {@code p}, - * the variance is {@code r * (1 - p) / p^2}. - */ - public double getNumericalVariance() { - final double p = getProbabilityOfSuccess(); - final double r = getNumberOfSuccesses(); - return r * (1 - p) / (p * p); - } - - /** - * {@inheritDoc} - * - * The lower bound of the support is always 0 no matter the parameters. - * - * @return lower bound of the support (always 0) - */ - public int getSupportLowerBound() { - return 0; - } - - /** - * {@inheritDoc} - * - * The upper bound of the support is always positive infinity no matter the - * parameters. Positive infinity is symbolized by {@code Integer.MAX_VALUE}. - * - * @return upper bound of the support (always {@code Integer.MAX_VALUE} - * for positive infinity) - */ - public int getSupportUpperBound() { - return Integer.MAX_VALUE; - } - - /** - * {@inheritDoc} - * - * The support of this distribution is connected. - * - * @return {@code true} - */ - public boolean isSupportConnected() { - return true; - } -}
http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/distribution/PoissonDistribution.java ---------------------------------------------------------------------- diff --git a/src/main/java/org/apache/commons/math3/distribution/PoissonDistribution.java b/src/main/java/org/apache/commons/math3/distribution/PoissonDistribution.java deleted file mode 100644 index 0d885e9..0000000 --- a/src/main/java/org/apache/commons/math3/distribution/PoissonDistribution.java +++ /dev/null @@ -1,395 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math3.distribution; - -import org.apache.commons.math3.exception.NotStrictlyPositiveException; -import org.apache.commons.math3.exception.util.LocalizedFormats; -import org.apache.commons.math3.random.RandomGenerator; -import org.apache.commons.math3.random.Well19937c; -import org.apache.commons.math3.special.Gamma; -import org.apache.commons.math3.util.CombinatoricsUtils; -import org.apache.commons.math3.util.FastMath; -import org.apache.commons.math3.util.MathUtils; - -/** - * Implementation of the Poisson distribution. - * - * @see <a href="http://en.wikipedia.org/wiki/Poisson_distribution";>Poisson distribution (Wikipedia)</a> - * @see <a href="http://mathworld.wolfram.com/PoissonDistribution.html";>Poisson distribution (MathWorld)</a> - */ -public class PoissonDistribution extends AbstractIntegerDistribution { - /** - * Default maximum number of iterations for cumulative probability calculations. - * @since 2.1 - */ - public static final int DEFAULT_MAX_ITERATIONS = 10000000; - /** - * Default convergence criterion. - * @since 2.1 - */ - public static final double DEFAULT_EPSILON = 1e-12; - /** Serializable version identifier. */ - private static final long serialVersionUID = -3349935121172596109L; - /** Distribution used to compute normal approximation. */ - private final NormalDistribution normal; - /** Distribution needed for the {@link #sample()} method. */ - private final ExponentialDistribution exponential; - /** Mean of the distribution. */ - private final double mean; - - /** - * Maximum number of iterations for cumulative probability. Cumulative - * probabilities are estimated using either Lanczos series approximation - * of {@link Gamma#regularizedGammaP(double, double, double, int)} - * or continued fraction approximation of - * {@link Gamma#regularizedGammaQ(double, double, double, int)}. - */ - private final int maxIterations; - - /** Convergence criterion for cumulative probability. */ - private final double epsilon; - - /** - * Creates a new Poisson distribution with specified mean. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param p the Poisson mean - * @throws NotStrictlyPositiveException if {@code p <= 0}. - */ - public PoissonDistribution(double p) throws NotStrictlyPositiveException { - this(p, DEFAULT_EPSILON, DEFAULT_MAX_ITERATIONS); - } - - /** - * Creates a new Poisson distribution with specified mean, convergence - * criterion and maximum number of iterations. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param p Poisson mean. - * @param epsilon Convergence criterion for cumulative probabilities. - * @param maxIterations the maximum number of iterations for cumulative - * probabilities. - * @throws NotStrictlyPositiveException if {@code p <= 0}. - * @since 2.1 - */ - public PoissonDistribution(double p, double epsilon, int maxIterations) - throws NotStrictlyPositiveException { - this(new Well19937c(), p, epsilon, maxIterations); - } - - /** - * Creates a new Poisson distribution with specified mean, convergence - * criterion and maximum number of iterations. - * - * @param rng Random number generator. - * @param p Poisson mean. - * @param epsilon Convergence criterion for cumulative probabilities. - * @param maxIterations the maximum number of iterations for cumulative - * probabilities. - * @throws NotStrictlyPositiveException if {@code p <= 0}. - * @since 3.1 - */ - public PoissonDistribution(RandomGenerator rng, - double p, - double epsilon, - int maxIterations) - throws NotStrictlyPositiveException { - super(rng); - - if (p <= 0) { - throw new NotStrictlyPositiveException(LocalizedFormats.MEAN, p); - } - mean = p; - this.epsilon = epsilon; - this.maxIterations = maxIterations; - - // Use the same RNG instance as the parent class. - normal = new NormalDistribution(rng, p, FastMath.sqrt(p), - NormalDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY); - exponential = new ExponentialDistribution(rng, 1, - ExponentialDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY); - } - - /** - * Creates a new Poisson distribution with the specified mean and - * convergence criterion. - * - * @param p Poisson mean. - * @param epsilon Convergence criterion for cumulative probabilities. - * @throws NotStrictlyPositiveException if {@code p <= 0}. - * @since 2.1 - */ - public PoissonDistribution(double p, double epsilon) - throws NotStrictlyPositiveException { - this(p, epsilon, DEFAULT_MAX_ITERATIONS); - } - - /** - * Creates a new Poisson distribution with the specified mean and maximum - * number of iterations. - * - * @param p Poisson mean. - * @param maxIterations Maximum number of iterations for cumulative - * probabilities. - * @since 2.1 - */ - public PoissonDistribution(double p, int maxIterations) { - this(p, DEFAULT_EPSILON, maxIterations); - } - - /** - * Get the mean for the distribution. - * - * @return the mean for the distribution. - */ - public double getMean() { - return mean; - } - - /** {@inheritDoc} */ - public double probability(int x) { - final double logProbability = logProbability(x); - return logProbability == Double.NEGATIVE_INFINITY ? 0 : FastMath.exp(logProbability); - } - - /** {@inheritDoc} */ - @Override - public double logProbability(int x) { - double ret; - if (x < 0 || x == Integer.MAX_VALUE) { - ret = Double.NEGATIVE_INFINITY; - } else if (x == 0) { - ret = -mean; - } else { - ret = -SaddlePointExpansion.getStirlingError(x) - - SaddlePointExpansion.getDeviancePart(x, mean) - - 0.5 * FastMath.log(MathUtils.TWO_PI) - 0.5 * FastMath.log(x); - } - return ret; - } - - /** {@inheritDoc} */ - public double cumulativeProbability(int x) { - if (x < 0) { - return 0; - } - if (x == Integer.MAX_VALUE) { - return 1; - } - return Gamma.regularizedGammaQ((double) x + 1, mean, epsilon, - maxIterations); - } - - /** - * Calculates the Poisson distribution function using a normal - * approximation. The {@code N(mean, sqrt(mean))} distribution is used - * to approximate the Poisson distribution. The computation uses - * "half-correction" (evaluating the normal distribution function at - * {@code x + 0.5}). - * - * @param x Upper bound, inclusive. - * @return the distribution function value calculated using a normal - * approximation. - */ - public double normalApproximateProbability(int x) { - // calculate the probability using half-correction - return normal.cumulativeProbability(x + 0.5); - } - - /** - * {@inheritDoc} - * - * For mean parameter {@code p}, the mean is {@code p}. - */ - public double getNumericalMean() { - return getMean(); - } - - /** - * {@inheritDoc} - * - * For mean parameter {@code p}, the variance is {@code p}. - */ - public double getNumericalVariance() { - return getMean(); - } - - /** - * {@inheritDoc} - * - * The lower bound of the support is always 0 no matter the mean parameter. - * - * @return lower bound of the support (always 0) - */ - public int getSupportLowerBound() { - return 0; - } - - /** - * {@inheritDoc} - * - * The upper bound of the support is positive infinity, - * regardless of the parameter values. There is no integer infinity, - * so this method returns {@code Integer.MAX_VALUE}. - * - * @return upper bound of the support (always {@code Integer.MAX_VALUE} for - * positive infinity) - */ - public int getSupportUpperBound() { - return Integer.MAX_VALUE; - } - - /** - * {@inheritDoc} - * - * The support of this distribution is connected. - * - * @return {@code true} - */ - public boolean isSupportConnected() { - return true; - } - - /** - * {@inheritDoc} - * <p> - * <strong>Algorithm Description</strong>: - * <ul> - * <li>For small means, uses simulation of a Poisson process - * using Uniform deviates, as described - * <a href="http://irmi.epfl.ch/cmos/Pmmi/interactive/rng7.htm";> here</a>. - * The Poisson process (and hence value returned) is bounded by 1000 * mean. - * </li> - * <li>For large means, uses the rejection algorithm described in - * <blockquote> - * Devroye, Luc. (1981).<i>The Computer Generation of Poisson Random Variables</i><br> - * <strong>Computing</strong> vol. 26 pp. 197-207.<br> - * </blockquote> - * </li> - * </ul> - * </p> - * - * @return a random value. - * @since 2.2 - */ - @Override - public int sample() { - return (int) FastMath.min(nextPoisson(mean), Integer.MAX_VALUE); - } - - /** - * @param meanPoisson Mean of the Poisson distribution. - * @return the next sample. - */ - private long nextPoisson(double meanPoisson) { - final double pivot = 40.0d; - if (meanPoisson < pivot) { - double p = FastMath.exp(-meanPoisson); - long n = 0; - double r = 1.0d; - double rnd = 1.0d; - - while (n < 1000 * meanPoisson) { - rnd = random.nextDouble(); - r *= rnd; - if (r >= p) { - n++; - } else { - return n; - } - } - return n; - } else { - final double lambda = FastMath.floor(meanPoisson); - final double lambdaFractional = meanPoisson - lambda; - final double logLambda = FastMath.log(lambda); - final double logLambdaFactorial = CombinatoricsUtils.factorialLog((int) lambda); - final long y2 = lambdaFractional < Double.MIN_VALUE ? 0 : nextPoisson(lambdaFractional); - final double delta = FastMath.sqrt(lambda * FastMath.log(32 * lambda / FastMath.PI + 1)); - final double halfDelta = delta / 2; - final double twolpd = 2 * lambda + delta; - final double a1 = FastMath.sqrt(FastMath.PI * twolpd) * FastMath.exp(1 / (8 * lambda)); - final double a2 = (twolpd / delta) * FastMath.exp(-delta * (1 + delta) / twolpd); - final double aSum = a1 + a2 + 1; - final double p1 = a1 / aSum; - final double p2 = a2 / aSum; - final double c1 = 1 / (8 * lambda); - - double x = 0; - double y = 0; - double v = 0; - int a = 0; - double t = 0; - double qr = 0; - double qa = 0; - for (;;) { - final double u = random.nextDouble(); - if (u <= p1) { - final double n = random.nextGaussian(); - x = n * FastMath.sqrt(lambda + halfDelta) - 0.5d; - if (x > delta || x < -lambda) { - continue; - } - y = x < 0 ? FastMath.floor(x) : FastMath.ceil(x); - final double e = exponential.sample(); - v = -e - (n * n / 2) + c1; - } else { - if (u > p1 + p2) { - y = lambda; - break; - } else { - x = delta + (twolpd / delta) * exponential.sample(); - y = FastMath.ceil(x); - v = -exponential.sample() - delta * (x + 1) / twolpd; - } - } - a = x < 0 ? 1 : 0; - t = y * (y + 1) / (2 * lambda); - if (v < -t && a == 0) { - y = lambda + y; - break; - } - qr = t * ((2 * y + 1) / (6 * lambda) - 1); - qa = qr - (t * t) / (3 * (lambda + a * (y + 1))); - if (v < qa) { - y = lambda + y; - break; - } - if (v > qr) { - continue; - } - if (v < y * logLambda - CombinatoricsUtils.factorialLog((int) (y + lambda)) + logLambdaFactorial) { - y = lambda + y; - break; - } - } - return y2 + (long) y; - } - } -} http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/distribution/RealDistribution.java ---------------------------------------------------------------------- diff --git a/src/main/java/org/apache/commons/math3/distribution/RealDistribution.java b/src/main/java/org/apache/commons/math3/distribution/RealDistribution.java deleted file mode 100644 index cb2045d..0000000 --- a/src/main/java/org/apache/commons/math3/distribution/RealDistribution.java +++ /dev/null @@ -1,196 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math3.distribution; - -import org.apache.commons.math3.exception.NumberIsTooLargeException; -import org.apache.commons.math3.exception.OutOfRangeException; - -/** - * Base interface for distributions on the reals. - * - * @since 3.0 - */ -public interface RealDistribution { - /** - * For a random variable {@code X} whose values are distributed according - * to this distribution, this method returns {@code P(X = x)}. In other - * words, this method represents the probability mass function (PMF) - * for the distribution. - * - * @param x the point at which the PMF is evaluated - * @return the value of the probability mass function at point {@code x} - */ - double probability(double x); - - /** - * Returns the probability density function (PDF) of this distribution - * evaluated at the specified point {@code x}. In general, the PDF is - * the derivative of the {@link #cumulativeProbability(double) CDF}. - * If the derivative does not exist at {@code x}, then an appropriate - * replacement should be returned, e.g. {@code Double.POSITIVE_INFINITY}, - * {@code Double.NaN}, or the limit inferior or limit superior of the - * difference quotient. - * - * @param x the point at which the PDF is evaluated - * @return the value of the probability density function at point {@code x} - */ - double density(double x); - - /** - * For a random variable {@code X} whose values are distributed according - * to this distribution, this method returns {@code P(X <= x)}. In other - * words, this method represents the (cumulative) distribution function - * (CDF) for this distribution. - * - * @param x the point at which the CDF is evaluated - * @return the probability that a random variable with this - * distribution takes a value less than or equal to {@code x} - */ - double cumulativeProbability(double x); - - /** - * For a random variable {@code X} whose values are distributed according - * to this distribution, this method returns {@code P(x0 < X <= x1)}. - * - * @param x0 the exclusive lower bound - * @param x1 the inclusive upper bound - * @return the probability that a random variable with this distribution - * takes a value between {@code x0} and {@code x1}, - * excluding the lower and including the upper endpoint - * @throws NumberIsTooLargeException if {@code x0 > x1} - * - * @deprecated As of 3.1. In 4.0, this method will be renamed - * {@code probability(double x0, double x1)}. - */ - @Deprecated - double cumulativeProbability(double x0, double x1) throws NumberIsTooLargeException; - - /** - * Computes the quantile function of this distribution. For a random - * variable {@code X} distributed according to this distribution, the - * returned value is - * <ul> - * <li><code>inf{x in R | P(X<=x) >= p}</code> for {@code 0 < p <= 1},</li> - * <li><code>inf{x in R | P(X<=x) > 0}</code> for {@code p = 0}.</li> - * </ul> - * - * @param p the cumulative probability - * @return the smallest {@code p}-quantile of this distribution - * (largest 0-quantile for {@code p = 0}) - * @throws OutOfRangeException if {@code p < 0} or {@code p > 1} - */ - double inverseCumulativeProbability(double p) throws OutOfRangeException; - - /** - * Use this method to get the numerical value of the mean of this - * distribution. - * - * @return the mean or {@code Double.NaN} if it is not defined - */ - double getNumericalMean(); - - /** - * Use this method to get the numerical value of the variance of this - * distribution. - * - * @return the variance (possibly {@code Double.POSITIVE_INFINITY} as - * for certain cases in {@link TDistribution}) or {@code Double.NaN} if it - * is not defined - */ - double getNumericalVariance(); - - /** - * Access the lower bound of the support. This method must return the same - * value as {@code inverseCumulativeProbability(0)}. In other words, this - * method must return - * <p><code>inf {x in R | P(X <= x) > 0}</code>.</p> - * - * @return lower bound of the support (might be - * {@code Double.NEGATIVE_INFINITY}) - */ - double getSupportLowerBound(); - - /** - * Access the upper bound of the support. This method must return the same - * value as {@code inverseCumulativeProbability(1)}. In other words, this - * method must return - * <p><code>inf {x in R | P(X <= x) = 1}</code>.</p> - * - * @return upper bound of the support (might be - * {@code Double.POSITIVE_INFINITY}) - */ - double getSupportUpperBound(); - - /** - * Whether or not the lower bound of support is in the domain of the density - * function. Returns true iff {@code getSupporLowerBound()} is finite and - * {@code density(getSupportLowerBound())} returns a non-NaN, non-infinite - * value. - * - * @return true if the lower bound of support is finite and the density - * function returns a non-NaN, non-infinite value there - * @deprecated to be removed in 4.0 - */ - @Deprecated - boolean isSupportLowerBoundInclusive(); - - /** - * Whether or not the upper bound of support is in the domain of the density - * function. Returns true iff {@code getSupportUpperBound()} is finite and - * {@code density(getSupportUpperBound())} returns a non-NaN, non-infinite - * value. - * - * @return true if the upper bound of support is finite and the density - * function returns a non-NaN, non-infinite value there - * @deprecated to be removed in 4.0 - */ - @Deprecated - boolean isSupportUpperBoundInclusive(); - - /** - * Use this method to get information about whether the support is connected, - * i.e. whether all values between the lower and upper bound of the support - * are included in the support. - * - * @return whether the support is connected or not - */ - boolean isSupportConnected(); - - /** - * Reseed the random generator used to generate samples. - * - * @param seed the new seed - */ - void reseedRandomGenerator(long seed); - - /** - * Generate a random value sampled from this distribution. - * - * @return a random value. - */ - double sample(); - - /** - * Generate a random sample from the distribution. - * - * @param sampleSize the number of random values to generate - * @return an array representing the random sample - * @throws org.apache.commons.math3.exception.NotStrictlyPositiveException - * if {@code sampleSize} is not positive - */ - double[] sample(int sampleSize); -} http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/distribution/SaddlePointExpansion.java ---------------------------------------------------------------------- diff --git a/src/main/java/org/apache/commons/math3/distribution/SaddlePointExpansion.java b/src/main/java/org/apache/commons/math3/distribution/SaddlePointExpansion.java deleted file mode 100644 index 63b3293..0000000 --- a/src/main/java/org/apache/commons/math3/distribution/SaddlePointExpansion.java +++ /dev/null @@ -1,200 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math3.distribution; - -import org.apache.commons.math3.special.Gamma; -import org.apache.commons.math3.util.FastMath; -import org.apache.commons.math3.util.MathUtils; - -/** - * <p> - * Utility class used by various distributions to accurately compute their - * respective probability mass functions. The implementation for this class is - * based on the Catherine Loader's <a target="_blank" - * href="http://www.herine.net/stat/software/dbinom.html";>dbinom</a> routines. - * </p> - * <p> - * This class is not intended to be called directly. - * </p> - * <p> - * References: - * <ol> - * <li>Catherine Loader (2000). "Fast and Accurate Computation of Binomial - * Probabilities.". <a target="_blank" - * href="http://www.herine.net/stat/papers/dbinom.pdf";> - * http://www.herine.net/stat/papers/dbinom.pdf</a></li> - * </ol> - * </p> - * - * @since 2.1 - */ -final class SaddlePointExpansion { - - /** 1/2 * log(2 π). */ - private static final double HALF_LOG_2_PI = 0.5 * FastMath.log(MathUtils.TWO_PI); - - /** exact Stirling expansion error for certain values. */ - private static final double[] EXACT_STIRLING_ERRORS = { 0.0, /* 0.0 */ - 0.1534264097200273452913848, /* 0.5 */ - 0.0810614667953272582196702, /* 1.0 */ - 0.0548141210519176538961390, /* 1.5 */ - 0.0413406959554092940938221, /* 2.0 */ - 0.03316287351993628748511048, /* 2.5 */ - 0.02767792568499833914878929, /* 3.0 */ - 0.02374616365629749597132920, /* 3.5 */ - 0.02079067210376509311152277, /* 4.0 */ - 0.01848845053267318523077934, /* 4.5 */ - 0.01664469118982119216319487, /* 5.0 */ - 0.01513497322191737887351255, /* 5.5 */ - 0.01387612882307074799874573, /* 6.0 */ - 0.01281046524292022692424986, /* 6.5 */ - 0.01189670994589177009505572, /* 7.0 */ - 0.01110455975820691732662991, /* 7.5 */ - 0.010411265261972096497478567, /* 8.0 */ - 0.009799416126158803298389475, /* 8.5 */ - 0.009255462182712732917728637, /* 9.0 */ - 0.008768700134139385462952823, /* 9.5 */ - 0.008330563433362871256469318, /* 10.0 */ - 0.007934114564314020547248100, /* 10.5 */ - 0.007573675487951840794972024, /* 11.0 */ - 0.007244554301320383179543912, /* 11.5 */ - 0.006942840107209529865664152, /* 12.0 */ - 0.006665247032707682442354394, /* 12.5 */ - 0.006408994188004207068439631, /* 13.0 */ - 0.006171712263039457647532867, /* 13.5 */ - 0.005951370112758847735624416, /* 14.0 */ - 0.005746216513010115682023589, /* 14.5 */ - 0.005554733551962801371038690 /* 15.0 */ - }; - - /** - * Default constructor. - */ - private SaddlePointExpansion() { - super(); - } - - /** - * Compute the error of Stirling's series at the given value. - * <p> - * References: - * <ol> - * <li>Eric W. Weisstein. "Stirling's Series." From MathWorld--A Wolfram Web - * Resource. <a target="_blank" - * href="http://mathworld.wolfram.com/StirlingsSeries.html";> - * http://mathworld.wolfram.com/StirlingsSeries.html</a></li> - * </ol> - * </p> - * - * @param z the value. - * @return the Striling's series error. - */ - static double getStirlingError(double z) { - double ret; - if (z < 15.0) { - double z2 = 2.0 * z; - if (FastMath.floor(z2) == z2) { - ret = EXACT_STIRLING_ERRORS[(int) z2]; - } else { - ret = Gamma.logGamma(z + 1.0) - (z + 0.5) * FastMath.log(z) + - z - HALF_LOG_2_PI; - } - } else { - double z2 = z * z; - ret = (0.083333333333333333333 - - (0.00277777777777777777778 - - (0.00079365079365079365079365 - - (0.000595238095238095238095238 - - 0.0008417508417508417508417508 / - z2) / z2) / z2) / z2) / z; - } - return ret; - } - - /** - * A part of the deviance portion of the saddle point approximation. - * <p> - * References: - * <ol> - * <li>Catherine Loader (2000). "Fast and Accurate Computation of Binomial - * Probabilities.". <a target="_blank" - * href="http://www.herine.net/stat/papers/dbinom.pdf";> - * http://www.herine.net/stat/papers/dbinom.pdf</a></li> - * </ol> - * </p> - * - * @param x the x value. - * @param mu the average. - * @return a part of the deviance. - */ - static double getDeviancePart(double x, double mu) { - double ret; - if (FastMath.abs(x - mu) < 0.1 * (x + mu)) { - double d = x - mu; - double v = d / (x + mu); - double s1 = v * d; - double s = Double.NaN; - double ej = 2.0 * x * v; - v *= v; - int j = 1; - while (s1 != s) { - s = s1; - ej *= v; - s1 = s + ej / ((j * 2) + 1); - ++j; - } - ret = s1; - } else { - ret = x * FastMath.log(x / mu) + mu - x; - } - return ret; - } - - /** - * Compute the logarithm of the PMF for a binomial distribution - * using the saddle point expansion. - * - * @param x the value at which the probability is evaluated. - * @param n the number of trials. - * @param p the probability of success. - * @param q the probability of failure (1 - p). - * @return log(p(x)). - */ - static double logBinomialProbability(int x, int n, double p, double q) { - double ret; - if (x == 0) { - if (p < 0.1) { - ret = -getDeviancePart(n, n * q) - n * p; - } else { - ret = n * FastMath.log(q); - } - } else if (x == n) { - if (q < 0.1) { - ret = -getDeviancePart(n, n * p) - n * q; - } else { - ret = n * FastMath.log(p); - } - } else { - ret = getStirlingError(n) - getStirlingError(x) - - getStirlingError(n - x) - getDeviancePart(x, n * p) - - getDeviancePart(n - x, n * q); - double f = (MathUtils.TWO_PI * x * (n - x)) / n; - ret = -0.5 * FastMath.log(f) + ret; - } - return ret; - } -} http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/distribution/TDistribution.java ---------------------------------------------------------------------- diff --git a/src/main/java/org/apache/commons/math3/distribution/TDistribution.java b/src/main/java/org/apache/commons/math3/distribution/TDistribution.java deleted file mode 100644 index 0c0e1cb..0000000 --- a/src/main/java/org/apache/commons/math3/distribution/TDistribution.java +++ /dev/null @@ -1,272 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math3.distribution; - -import org.apache.commons.math3.exception.NotStrictlyPositiveException; -import org.apache.commons.math3.exception.util.LocalizedFormats; -import org.apache.commons.math3.random.RandomGenerator; -import org.apache.commons.math3.random.Well19937c; -import org.apache.commons.math3.special.Beta; -import org.apache.commons.math3.special.Gamma; -import org.apache.commons.math3.util.FastMath; - -/** - * Implementation of Student's t-distribution. - * - * @see "<a href='http://en.wikipedia.org/wiki/Student's_t-distribution'>Student's t-distribution (Wikipedia)</a>" - * @see "<a href='http://mathworld.wolfram.com/Studentst-Distribution.html'>Student's t-distribution (MathWorld)</a>" - */ -public class TDistribution extends AbstractRealDistribution { - /** - * Default inverse cumulative probability accuracy. - * @since 2.1 - */ - public static final double DEFAULT_INVERSE_ABSOLUTE_ACCURACY = 1e-9; - /** Serializable version identifier */ - private static final long serialVersionUID = -5852615386664158222L; - /** The degrees of freedom. */ - private final double degreesOfFreedom; - /** Inverse cumulative probability accuracy. */ - private final double solverAbsoluteAccuracy; - /** Static computation factor based on degreesOfFreedom. */ - private final double factor; - - /** - * Create a t distribution using the given degrees of freedom. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param degreesOfFreedom Degrees of freedom. - * @throws NotStrictlyPositiveException if {@code degreesOfFreedom <= 0} - */ - public TDistribution(double degreesOfFreedom) - throws NotStrictlyPositiveException { - this(degreesOfFreedom, DEFAULT_INVERSE_ABSOLUTE_ACCURACY); - } - - /** - * Create a t distribution using the given degrees of freedom and the - * specified inverse cumulative probability absolute accuracy. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param degreesOfFreedom Degrees of freedom. - * @param inverseCumAccuracy the maximum absolute error in inverse - * cumulative probability estimates - * (defaults to {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY}). - * @throws NotStrictlyPositiveException if {@code degreesOfFreedom <= 0} - * @since 2.1 - */ - public TDistribution(double degreesOfFreedom, double inverseCumAccuracy) - throws NotStrictlyPositiveException { - this(new Well19937c(), degreesOfFreedom, inverseCumAccuracy); - } - - /** - * Creates a t distribution. - * - * @param rng Random number generator. - * @param degreesOfFreedom Degrees of freedom. - * @throws NotStrictlyPositiveException if {@code degreesOfFreedom <= 0} - * @since 3.3 - */ - public TDistribution(RandomGenerator rng, double degreesOfFreedom) - throws NotStrictlyPositiveException { - this(rng, degreesOfFreedom, DEFAULT_INVERSE_ABSOLUTE_ACCURACY); - } - - /** - * Creates a t distribution. - * - * @param rng Random number generator. - * @param degreesOfFreedom Degrees of freedom. - * @param inverseCumAccuracy the maximum absolute error in inverse - * cumulative probability estimates - * (defaults to {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY}). - * @throws NotStrictlyPositiveException if {@code degreesOfFreedom <= 0} - * @since 3.1 - */ - public TDistribution(RandomGenerator rng, - double degreesOfFreedom, - double inverseCumAccuracy) - throws NotStrictlyPositiveException { - super(rng); - - if (degreesOfFreedom <= 0) { - throw new NotStrictlyPositiveException(LocalizedFormats.DEGREES_OF_FREEDOM, - degreesOfFreedom); - } - this.degreesOfFreedom = degreesOfFreedom; - solverAbsoluteAccuracy = inverseCumAccuracy; - - final double n = degreesOfFreedom; - final double nPlus1Over2 = (n + 1) / 2; - factor = Gamma.logGamma(nPlus1Over2) - - 0.5 * (FastMath.log(FastMath.PI) + FastMath.log(n)) - - Gamma.logGamma(n / 2); - } - - /** - * Access the degrees of freedom. - * - * @return the degrees of freedom. - */ - public double getDegreesOfFreedom() { - return degreesOfFreedom; - } - - /** {@inheritDoc} */ - public double density(double x) { - return FastMath.exp(logDensity(x)); - } - - /** {@inheritDoc} */ - @Override - public double logDensity(double x) { - final double n = degreesOfFreedom; - final double nPlus1Over2 = (n + 1) / 2; - return factor - nPlus1Over2 * FastMath.log(1 + x * x / n); - } - - /** {@inheritDoc} */ - public double cumulativeProbability(double x) { - double ret; - if (x == 0) { - ret = 0.5; - } else { - double t = - Beta.regularizedBeta( - degreesOfFreedom / (degreesOfFreedom + (x * x)), - 0.5 * degreesOfFreedom, - 0.5); - if (x < 0.0) { - ret = 0.5 * t; - } else { - ret = 1.0 - 0.5 * t; - } - } - - return ret; - } - - /** {@inheritDoc} */ - @Override - protected double getSolverAbsoluteAccuracy() { - return solverAbsoluteAccuracy; - } - - /** - * {@inheritDoc} - * - * For degrees of freedom parameter {@code df}, the mean is - * <ul> - * <li>if {@code df > 1} then {@code 0},</li> - * <li>else undefined ({@code Double.NaN}).</li> - * </ul> - */ - public double getNumericalMean() { - final double df = getDegreesOfFreedom(); - - if (df > 1) { - return 0; - } - - return Double.NaN; - } - - /** - * {@inheritDoc} - * - * For degrees of freedom parameter {@code df}, the variance is - * <ul> - * <li>if {@code df > 2} then {@code df / (df - 2)},</li> - * <li>if {@code 1 < df <= 2} then positive infinity - * ({@code Double.POSITIVE_INFINITY}),</li> - * <li>else undefined ({@code Double.NaN}).</li> - * </ul> - */ - public double getNumericalVariance() { - final double df = getDegreesOfFreedom(); - - if (df > 2) { - return df / (df - 2); - } - - if (df > 1 && df <= 2) { - return Double.POSITIVE_INFINITY; - } - - return Double.NaN; - } - - /** - * {@inheritDoc} - * - * The lower bound of the support is always negative infinity no matter the - * parameters. - * - * @return lower bound of the support (always - * {@code Double.NEGATIVE_INFINITY}) - */ - public double getSupportLowerBound() { - return Double.NEGATIVE_INFINITY; - } - - /** - * {@inheritDoc} - * - * The upper bound of the support is always positive infinity no matter the - * parameters. - * - * @return upper bound of the support (always - * {@code Double.POSITIVE_INFINITY}) - */ - public double getSupportUpperBound() { - return Double.POSITIVE_INFINITY; - } - - /** {@inheritDoc} */ - public boolean isSupportLowerBoundInclusive() { - return false; - } - - /** {@inheritDoc} */ - public boolean isSupportUpperBoundInclusive() { - return false; - } - - /** - * {@inheritDoc} - * - * The support of this distribution is connected. - * - * @return {@code true} - */ - public boolean isSupportConnected() { - return true; - } -} http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/distribution/TriangularDistribution.java ---------------------------------------------------------------------- diff --git a/src/main/java/org/apache/commons/math3/distribution/TriangularDistribution.java b/src/main/java/org/apache/commons/math3/distribution/TriangularDistribution.java deleted file mode 100644 index 1a81ecb..0000000 --- a/src/main/java/org/apache/commons/math3/distribution/TriangularDistribution.java +++ /dev/null @@ -1,282 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -package org.apache.commons.math3.distribution; - -import org.apache.commons.math3.exception.NumberIsTooLargeException; -import org.apache.commons.math3.exception.NumberIsTooSmallException; -import org.apache.commons.math3.exception.OutOfRangeException; -import org.apache.commons.math3.exception.util.LocalizedFormats; -import org.apache.commons.math3.random.RandomGenerator; -import org.apache.commons.math3.random.Well19937c; -import org.apache.commons.math3.util.FastMath; - -/** - * Implementation of the triangular real distribution. - * - * @see <a href="http://en.wikipedia.org/wiki/Triangular_distribution";> - * Triangular distribution (Wikipedia)</a> - * - * @since 3.0 - */ -public class TriangularDistribution extends AbstractRealDistribution { - /** Serializable version identifier. */ - private static final long serialVersionUID = 20120112L; - /** Lower limit of this distribution (inclusive). */ - private final double a; - /** Upper limit of this distribution (inclusive). */ - private final double b; - /** Mode of this distribution. */ - private final double c; - /** Inverse cumulative probability accuracy. */ - private final double solverAbsoluteAccuracy; - - /** - * Creates a triangular real distribution using the given lower limit, - * upper limit, and mode. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param a Lower limit of this distribution (inclusive). - * @param b Upper limit of this distribution (inclusive). - * @param c Mode of this distribution. - * @throws NumberIsTooLargeException if {@code a >= b} or if {@code c > b}. - * @throws NumberIsTooSmallException if {@code c < a}. - */ - public TriangularDistribution(double a, double c, double b) - throws NumberIsTooLargeException, NumberIsTooSmallException { - this(new Well19937c(), a, c, b); - } - - /** - * Creates a triangular distribution. - * - * @param rng Random number generator. - * @param a Lower limit of this distribution (inclusive). - * @param b Upper limit of this distribution (inclusive). - * @param c Mode of this distribution. - * @throws NumberIsTooLargeException if {@code a >= b} or if {@code c > b}. - * @throws NumberIsTooSmallException if {@code c < a}. - * @since 3.1 - */ - public TriangularDistribution(RandomGenerator rng, - double a, - double c, - double b) - throws NumberIsTooLargeException, NumberIsTooSmallException { - super(rng); - - if (a >= b) { - throw new NumberIsTooLargeException( - LocalizedFormats.LOWER_BOUND_NOT_BELOW_UPPER_BOUND, - a, b, false); - } - if (c < a) { - throw new NumberIsTooSmallException( - LocalizedFormats.NUMBER_TOO_SMALL, c, a, true); - } - if (c > b) { - throw new NumberIsTooLargeException( - LocalizedFormats.NUMBER_TOO_LARGE, c, b, true); - } - - this.a = a; - this.c = c; - this.b = b; - solverAbsoluteAccuracy = FastMath.max(FastMath.ulp(a), FastMath.ulp(b)); - } - - /** - * Returns the mode {@code c} of this distribution. - * - * @return the mode {@code c} of this distribution - */ - public double getMode() { - return c; - } - - /** - * {@inheritDoc} - * - * <p> - * For this distribution, the returned value is not really meaningful, - * since exact formulas are implemented for the computation of the - * {@link #inverseCumulativeProbability(double)} (no solver is invoked). - * </p> - * <p> - * For lower limit {@code a} and upper limit {@code b}, the current - * implementation returns {@code max(ulp(a), ulp(b)}. - * </p> - */ - @Override - protected double getSolverAbsoluteAccuracy() { - return solverAbsoluteAccuracy; - } - - /** - * {@inheritDoc} - * - * For lower limit {@code a}, upper limit {@code b} and mode {@code c}, the - * PDF is given by - * <ul> - * <li>{@code 2 * (x - a) / [(b - a) * (c - a)]} if {@code a <= x < c},</li> - * <li>{@code 2 / (b - a)} if {@code x = c},</li> - * <li>{@code 2 * (b - x) / [(b - a) * (b - c)]} if {@code c < x <= b},</li> - * <li>{@code 0} otherwise. - * </ul> - */ - public double density(double x) { - if (x < a) { - return 0; - } - if (a <= x && x < c) { - double divident = 2 * (x - a); - double divisor = (b - a) * (c - a); - return divident / divisor; - } - if (x == c) { - return 2 / (b - a); - } - if (c < x && x <= b) { - double divident = 2 * (b - x); - double divisor = (b - a) * (b - c); - return divident / divisor; - } - return 0; - } - - /** - * {@inheritDoc} - * - * For lower limit {@code a}, upper limit {@code b} and mode {@code c}, the - * CDF is given by - * <ul> - * <li>{@code 0} if {@code x < a},</li> - * <li>{@code (x - a)^2 / [(b - a) * (c - a)]} if {@code a <= x < c},</li> - * <li>{@code (c - a) / (b - a)} if {@code x = c},</li> - * <li>{@code 1 - (b - x)^2 / [(b - a) * (b - c)]} if {@code c < x <= b},</li> - * <li>{@code 1} if {@code x > b}.</li> - * </ul> - */ - public double cumulativeProbability(double x) { - if (x < a) { - return 0; - } - if (a <= x && x < c) { - double divident = (x - a) * (x - a); - double divisor = (b - a) * (c - a); - return divident / divisor; - } - if (x == c) { - return (c - a) / (b - a); - } - if (c < x && x <= b) { - double divident = (b - x) * (b - x); - double divisor = (b - a) * (b - c); - return 1 - (divident / divisor); - } - return 1; - } - - /** - * {@inheritDoc} - * - * For lower limit {@code a}, upper limit {@code b}, and mode {@code c}, - * the mean is {@code (a + b + c) / 3}. - */ - public double getNumericalMean() { - return (a + b + c) / 3; - } - - /** - * {@inheritDoc} - * - * For lower limit {@code a}, upper limit {@code b}, and mode {@code c}, - * the variance is {@code (a^2 + b^2 + c^2 - a * b - a * c - b * c) / 18}. - */ - public double getNumericalVariance() { - return (a * a + b * b + c * c - a * b - a * c - b * c) / 18; - } - - /** - * {@inheritDoc} - * - * The lower bound of the support is equal to the lower limit parameter - * {@code a} of the distribution. - * - * @return lower bound of the support - */ - public double getSupportLowerBound() { - return a; - } - - /** - * {@inheritDoc} - * - * The upper bound of the support is equal to the upper limit parameter - * {@code b} of the distribution. - * - * @return upper bound of the support - */ - public double getSupportUpperBound() { - return b; - } - - /** {@inheritDoc} */ - public boolean isSupportLowerBoundInclusive() { - return true; - } - - /** {@inheritDoc} */ - public boolean isSupportUpperBoundInclusive() { - return true; - } - - /** - * {@inheritDoc} - * - * The support of this distribution is connected. - * - * @return {@code true} - */ - public boolean isSupportConnected() { - return true; - } - - @Override - public double inverseCumulativeProbability(double p) - throws OutOfRangeException { - if (p < 0 || p > 1) { - throw new OutOfRangeException(p, 0, 1); - } - if (p == 0) { - return a; - } - if (p == 1) { - return b; - } - if (p < (c - a) / (b - a)) { - return a + FastMath.sqrt(p * (b - a) * (c - a)); - } - return b - FastMath.sqrt((1 - p) * (b - a) * (b - c)); - } -} http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/distribution/UniformIntegerDistribution.java ---------------------------------------------------------------------- diff --git a/src/main/java/org/apache/commons/math3/distribution/UniformIntegerDistribution.java b/src/main/java/org/apache/commons/math3/distribution/UniformIntegerDistribution.java deleted file mode 100644 index dd5ffa4..0000000 --- a/src/main/java/org/apache/commons/math3/distribution/UniformIntegerDistribution.java +++ /dev/null @@ -1,181 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -package org.apache.commons.math3.distribution; - -import org.apache.commons.math3.exception.NumberIsTooLargeException; -import org.apache.commons.math3.exception.util.LocalizedFormats; -import org.apache.commons.math3.random.RandomGenerator; -import org.apache.commons.math3.random.Well19937c; - -/** - * Implementation of the uniform integer distribution. - * - * @see <a href="http://en.wikipedia.org/wiki/Uniform_distribution_(discrete)" - * >Uniform distribution (discrete), at Wikipedia</a> - * - * @since 3.0 - */ -public class UniformIntegerDistribution extends AbstractIntegerDistribution { - /** Serializable version identifier. */ - private static final long serialVersionUID = 20120109L; - /** Lower bound (inclusive) of this distribution. */ - private final int lower; - /** Upper bound (inclusive) of this distribution. */ - private final int upper; - - /** - * Creates a new uniform integer distribution using the given lower and - * upper bounds (both inclusive). - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param lower Lower bound (inclusive) of this distribution. - * @param upper Upper bound (inclusive) of this distribution. - * @throws NumberIsTooLargeException if {@code lower >= upper}. - */ - public UniformIntegerDistribution(int lower, int upper) - throws NumberIsTooLargeException { - this(new Well19937c(), lower, upper); - } - - /** - * Creates a new uniform integer distribution using the given lower and - * upper bounds (both inclusive). - * - * @param rng Random number generator. - * @param lower Lower bound (inclusive) of this distribution. - * @param upper Upper bound (inclusive) of this distribution. - * @throws NumberIsTooLargeException if {@code lower > upper}. - * @since 3.1 - */ - public UniformIntegerDistribution(RandomGenerator rng, - int lower, - int upper) - throws NumberIsTooLargeException { - super(rng); - - if (lower > upper) { - throw new NumberIsTooLargeException( - LocalizedFormats.LOWER_BOUND_NOT_BELOW_UPPER_BOUND, - lower, upper, true); - } - this.lower = lower; - this.upper = upper; - } - - /** {@inheritDoc} */ - public double probability(int x) { - if (x < lower || x > upper) { - return 0; - } - return 1.0 / (upper - lower + 1); - } - - /** {@inheritDoc} */ - public double cumulativeProbability(int x) { - if (x < lower) { - return 0; - } - if (x > upper) { - return 1; - } - return (x - lower + 1.0) / (upper - lower + 1.0); - } - - /** - * {@inheritDoc} - * - * For lower bound {@code lower} and upper bound {@code upper}, the mean is - * {@code 0.5 * (lower + upper)}. - */ - public double getNumericalMean() { - return 0.5 * (lower + upper); - } - - /** - * {@inheritDoc} - * - * For lower bound {@code lower} and upper bound {@code upper}, and - * {@code n = upper - lower + 1}, the variance is {@code (n^2 - 1) / 12}. - */ - public double getNumericalVariance() { - double n = upper - lower + 1; - return (n * n - 1) / 12.0; - } - - /** - * {@inheritDoc} - * - * The lower bound of the support is equal to the lower bound parameter - * of the distribution. - * - * @return lower bound of the support - */ - public int getSupportLowerBound() { - return lower; - } - - /** - * {@inheritDoc} - * - * The upper bound of the support is equal to the upper bound parameter - * of the distribution. - * - * @return upper bound of the support - */ - public int getSupportUpperBound() { - return upper; - } - - /** - * {@inheritDoc} - * - * The support of this distribution is connected. - * - * @return {@code true} - */ - public boolean isSupportConnected() { - return true; - } - - /** {@inheritDoc} */ - @Override - public int sample() { - final int max = (upper - lower) + 1; - if (max <= 0) { - // The range is too wide to fit in a positive int (larger - // than 2^31); as it covers more than half the integer range, - // we use a simple rejection method. - while (true) { - final int r = random.nextInt(); - if (r >= lower && - r <= upper) { - return r; - } - } - } else { - // We can shift the range and directly generate a positive int. - return lower + random.nextInt(max); - } - } -} http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/distribution/UniformRealDistribution.java ---------------------------------------------------------------------- diff --git a/src/main/java/org/apache/commons/math3/distribution/UniformRealDistribution.java b/src/main/java/org/apache/commons/math3/distribution/UniformRealDistribution.java deleted file mode 100644 index 10c7712..0000000 --- a/src/main/java/org/apache/commons/math3/distribution/UniformRealDistribution.java +++ /dev/null @@ -1,243 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -package org.apache.commons.math3.distribution; - -import org.apache.commons.math3.exception.NumberIsTooLargeException; -import org.apache.commons.math3.exception.OutOfRangeException; -import org.apache.commons.math3.exception.util.LocalizedFormats; -import org.apache.commons.math3.random.RandomGenerator; -import org.apache.commons.math3.random.Well19937c; - -/** - * Implementation of the uniform real distribution. - * - * @see <a href="http://en.wikipedia.org/wiki/Uniform_distribution_(continuous)" - * >Uniform distribution (continuous), at Wikipedia</a> - * - * @since 3.0 - */ -public class UniformRealDistribution extends AbstractRealDistribution { - /** Default inverse cumulative probability accuracy. - * @deprecated as of 3.2 not used anymore, will be removed in 4.0 - */ - @Deprecated - public static final double DEFAULT_INVERSE_ABSOLUTE_ACCURACY = 1e-9; - /** Serializable version identifier. */ - private static final long serialVersionUID = 20120109L; - /** Lower bound of this distribution (inclusive). */ - private final double lower; - /** Upper bound of this distribution (exclusive). */ - private final double upper; - - /** - * Create a standard uniform real distribution with lower bound (inclusive) - * equal to zero and upper bound (exclusive) equal to one. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - */ - public UniformRealDistribution() { - this(0, 1); - } - - /** - * Create a uniform real distribution using the given lower and upper - * bounds. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param lower Lower bound of this distribution (inclusive). - * @param upper Upper bound of this distribution (exclusive). - * @throws NumberIsTooLargeException if {@code lower >= upper}. - */ - public UniformRealDistribution(double lower, double upper) - throws NumberIsTooLargeException { - this(new Well19937c(), lower, upper); - } - - /** - * Create a uniform distribution. - * - * @param lower Lower bound of this distribution (inclusive). - * @param upper Upper bound of this distribution (exclusive). - * @param inverseCumAccuracy Inverse cumulative probability accuracy. - * @throws NumberIsTooLargeException if {@code lower >= upper}. - * @deprecated as of 3.2, inverse CDF is now calculated analytically, use - * {@link #UniformRealDistribution(double, double)} instead. - */ - @Deprecated - public UniformRealDistribution(double lower, double upper, double inverseCumAccuracy) - throws NumberIsTooLargeException { - this(new Well19937c(), lower, upper); - } - - /** - * Creates a uniform distribution. - * - * @param rng Random number generator. - * @param lower Lower bound of this distribution (inclusive). - * @param upper Upper bound of this distribution (exclusive). - * @param inverseCumAccuracy Inverse cumulative probability accuracy. - * @throws NumberIsTooLargeException if {@code lower >= upper}. - * @since 3.1 - * @deprecated as of 3.2, inverse CDF is now calculated analytically, use - * {@link #UniformRealDistribution(RandomGenerator, double, double)} - * instead. - */ - @Deprecated - public UniformRealDistribution(RandomGenerator rng, - double lower, - double upper, - double inverseCumAccuracy){ - this(rng, lower, upper); - } - - /** - * Creates a uniform distribution. - * - * @param rng Random number generator. - * @param lower Lower bound of this distribution (inclusive). - * @param upper Upper bound of this distribution (exclusive). - * @throws NumberIsTooLargeException if {@code lower >= upper}. - * @since 3.1 - */ - public UniformRealDistribution(RandomGenerator rng, - double lower, - double upper) - throws NumberIsTooLargeException { - super(rng); - if (lower >= upper) { - throw new NumberIsTooLargeException( - LocalizedFormats.LOWER_BOUND_NOT_BELOW_UPPER_BOUND, - lower, upper, false); - } - - this.lower = lower; - this.upper = upper; - } - - /** {@inheritDoc} */ - public double density(double x) { - if (x < lower || x > upper) { - return 0.0; - } - return 1 / (upper - lower); - } - - /** {@inheritDoc} */ - public double cumulativeProbability(double x) { - if (x <= lower) { - return 0; - } - if (x >= upper) { - return 1; - } - return (x - lower) / (upper - lower); - } - - @Override - public double inverseCumulativeProbability(final double p) - throws OutOfRangeException { - if (p < 0.0 || p > 1.0) { - throw new OutOfRangeException(p, 0, 1); - } - return p * (upper - lower) + lower; - } - - /** - * {@inheritDoc} - * - * For lower bound {@code lower} and upper bound {@code upper}, the mean is - * {@code 0.5 * (lower + upper)}. - */ - public double getNumericalMean() { - return 0.5 * (lower + upper); - } - - /** - * {@inheritDoc} - * - * For lower bound {@code lower} and upper bound {@code upper}, the - * variance is {@code (upper - lower)^2 / 12}. - */ - public double getNumericalVariance() { - double ul = upper - lower; - return ul * ul / 12; - } - - /** - * {@inheritDoc} - * - * The lower bound of the support is equal to the lower bound parameter - * of the distribution. - * - * @return lower bound of the support - */ - public double getSupportLowerBound() { - return lower; - } - - /** - * {@inheritDoc} - * - * The upper bound of the support is equal to the upper bound parameter - * of the distribution. - * - * @return upper bound of the support - */ - public double getSupportUpperBound() { - return upper; - } - - /** {@inheritDoc} */ - public boolean isSupportLowerBoundInclusive() { - return true; - } - - /** {@inheritDoc} */ - public boolean isSupportUpperBoundInclusive() { - return true; - } - - /** - * {@inheritDoc} - * - * The support of this distribution is connected. - * - * @return {@code true} - */ - public boolean isSupportConnected() { - return true; - } - - /** {@inheritDoc} */ - @Override - public double sample() { - final double u = random.nextDouble(); - return u * upper + (1 - u) * lower; - } -} http://git-wip-us.apache.org/repos/asf/commons-math/blob/a7b4803f/src/main/java/org/apache/commons/math3/distribution/WeibullDistribution.java ---------------------------------------------------------------------- diff --git a/src/main/java/org/apache/commons/math3/distribution/WeibullDistribution.java b/src/main/java/org/apache/commons/math3/distribution/WeibullDistribution.java deleted file mode 100644 index f1fdabd..0000000 --- a/src/main/java/org/apache/commons/math3/distribution/WeibullDistribution.java +++ /dev/null @@ -1,353 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -package org.apache.commons.math3.distribution; - -import org.apache.commons.math3.exception.NotStrictlyPositiveException; -import org.apache.commons.math3.exception.OutOfRangeException; -import org.apache.commons.math3.exception.util.LocalizedFormats; -import org.apache.commons.math3.random.RandomGenerator; -import org.apache.commons.math3.random.Well19937c; -import org.apache.commons.math3.special.Gamma; -import org.apache.commons.math3.util.FastMath; - -/** - * Implementation of the Weibull distribution. This implementation uses the - * two parameter form of the distribution defined by - * <a href="http://mathworld.wolfram.com/WeibullDistribution.html";> - * Weibull Distribution</a>, equations (1) and (2). - * - * @see <a href="http://en.wikipedia.org/wiki/Weibull_distribution";>Weibull distribution (Wikipedia)</a> - * @see <a href="http://mathworld.wolfram.com/WeibullDistribution.html";>Weibull distribution (MathWorld)</a> - * @since 1.1 (changed to concrete class in 3.0) - */ -public class WeibullDistribution extends AbstractRealDistribution { - /** - * Default inverse cumulative probability accuracy. - * @since 2.1 - */ - public static final double DEFAULT_INVERSE_ABSOLUTE_ACCURACY = 1e-9; - /** Serializable version identifier. */ - private static final long serialVersionUID = 8589540077390120676L; - /** The shape parameter. */ - private final double shape; - /** The scale parameter. */ - private final double scale; - /** Inverse cumulative probability accuracy. */ - private final double solverAbsoluteAccuracy; - /** Cached numerical mean */ - private double numericalMean = Double.NaN; - /** Whether or not the numerical mean has been calculated */ - private boolean numericalMeanIsCalculated = false; - /** Cached numerical variance */ - private double numericalVariance = Double.NaN; - /** Whether or not the numerical variance has been calculated */ - private boolean numericalVarianceIsCalculated = false; - - /** - * Create a Weibull distribution with the given shape and scale and a - * location equal to zero. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param alpha Shape parameter. - * @param beta Scale parameter. - * @throws NotStrictlyPositiveException if {@code alpha <= 0} or - * {@code beta <= 0}. - */ - public WeibullDistribution(double alpha, double beta) - throws NotStrictlyPositiveException { - this(alpha, beta, DEFAULT_INVERSE_ABSOLUTE_ACCURACY); - } - - /** - * Create a Weibull distribution with the given shape, scale and inverse - * cumulative probability accuracy and a location equal to zero. - * <p> - * <b>Note:</b> this constructor will implicitly create an instance of - * {@link Well19937c} as random generator to be used for sampling only (see - * {@link #sample()} and {@link #sample(int)}). In case no sampling is - * needed for the created distribution, it is advised to pass {@code null} - * as random generator via the appropriate constructors to avoid the - * additional initialisation overhead. - * - * @param alpha Shape parameter. - * @param beta Scale parameter. - * @param inverseCumAccuracy Maximum absolute error in inverse - * cumulative probability estimates - * (defaults to {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY}). - * @throws NotStrictlyPositiveException if {@code alpha <= 0} or - * {@code beta <= 0}. - * @since 2.1 - */ - public WeibullDistribution(double alpha, double beta, - double inverseCumAccuracy) { - this(new Well19937c(), alpha, beta, inverseCumAccuracy); - } - - /** - * Creates a Weibull distribution. - * - * @param rng Random number generator. - * @param alpha Shape parameter. - * @param beta Scale parameter. - * @throws NotStrictlyPositiveException if {@code alpha <= 0} or {@code beta <= 0}. - * @since 3.3 - */ - public WeibullDistribution(RandomGenerator rng, double alpha, double beta) - throws NotStrictlyPositiveException { - this(rng, alpha, beta, DEFAULT_INVERSE_ABSOLUTE_ACCURACY); - } - - /** - * Creates a Weibull distribution. - * - * @param rng Random number generator. - * @param alpha Shape parameter. - * @param beta Scale parameter. - * @param inverseCumAccuracy Maximum absolute error in inverse - * cumulative probability estimates - * (defaults to {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY}). - * @throws NotStrictlyPositiveException if {@code alpha <= 0} or {@code beta <= 0}. - * @since 3.1 - */ - public WeibullDistribution(RandomGenerator rng, - double alpha, - double beta, - double inverseCumAccuracy) - throws NotStrictlyPositiveException { - super(rng); - - if (alpha <= 0) { - throw new NotStrictlyPositiveException(LocalizedFormats.SHAPE, - alpha); - } - if (beta <= 0) { - throw new NotStrictlyPositiveException(LocalizedFormats.SCALE, - beta); - } - scale = beta; - shape = alpha; - solverAbsoluteAccuracy = inverseCumAccuracy; - } - - /** - * Access the shape parameter, {@code alpha}. - * - * @return the shape parameter, {@code alpha}. - */ - public double getShape() { - return shape; - } - - /** - * Access the scale parameter, {@code beta}. - * - * @return the scale parameter, {@code beta}. - */ - public double getScale() { - return scale; - } - - /** {@inheritDoc} */ - public double density(double x) { - if (x < 0) { - return 0; - } - - final double xscale = x / scale; - final double xscalepow = FastMath.pow(xscale, shape - 1); - - /* - * FastMath.pow(x / scale, shape) = - * FastMath.pow(xscale, shape) = - * FastMath.pow(xscale, shape - 1) * xscale - */ - final double xscalepowshape = xscalepow * xscale; - - return (shape / scale) * xscalepow * FastMath.exp(-xscalepowshape); - } - - /** {@inheritDoc} */ - @Override - public double logDensity(double x) { - if (x < 0) { - return Double.NEGATIVE_INFINITY; - } - - final double xscale = x / scale; - final double logxscalepow = FastMath.log(xscale) * (shape - 1); - - /* - * FastMath.pow(x / scale, shape) = - * FastMath.pow(xscale, shape) = - * FastMath.pow(xscale, shape - 1) * xscale - */ - final double xscalepowshape = FastMath.exp(logxscalepow) * xscale; - - return FastMath.log(shape / scale) + logxscalepow - xscalepowshape; - } - - /** {@inheritDoc} */ - public double cumulativeProbability(double x) { - double ret; - if (x <= 0.0) { - ret = 0.0; - } else { - ret = 1.0 - FastMath.exp(-FastMath.pow(x / scale, shape)); - } - return ret; - } - - /** - * {@inheritDoc} - * - * Returns {@code 0} when {@code p == 0} and - * {@code Double.POSITIVE_INFINITY} when {@code p == 1}. - */ - @Override - public double inverseCumulativeProbability(double p) { - double ret; - if (p < 0.0 || p > 1.0) { - throw new OutOfRangeException(p, 0.0, 1.0); - } else if (p == 0) { - ret = 0.0; - } else if (p == 1) { - ret = Double.POSITIVE_INFINITY; - } else { - ret = scale * FastMath.pow(-FastMath.log1p(-p), 1.0 / shape); - } - return ret; - } - - /** - * Return the absolute accuracy setting of the solver used to estimate - * inverse cumulative probabilities. - * - * @return the solver absolute accuracy. - * @since 2.1 - */ - @Override - protected double getSolverAbsoluteAccuracy() { - return solverAbsoluteAccuracy; - } - - /** - * {@inheritDoc} - * - * The mean is {@code scale * Gamma(1 + (1 / shape))}, where {@code Gamma()} - * is the Gamma-function. - */ - public double getNumericalMean() { - if (!numericalMeanIsCalculated) { - numericalMean = calculateNumericalMean(); - numericalMeanIsCalculated = true; - } - return numericalMean; - } - - /** - * used by {@link #getNumericalMean()} - * - * @return the mean of this distribution - */ - protected double calculateNumericalMean() { - final double sh = getShape(); - final double sc = getScale(); - - return sc * FastMath.exp(Gamma.logGamma(1 + (1 / sh))); - } - - /** - * {@inheritDoc} - * - * The variance is {@code scale^2 * Gamma(1 + (2 / shape)) - mean^2} - * where {@code Gamma()} is the Gamma-function. - */ - public double getNumericalVariance() { - if (!numericalVarianceIsCalculated) { - numericalVariance = calculateNumericalVariance(); - numericalVarianceIsCalculated = true; - } - return numericalVariance; - } - - /** - * used by {@link #getNumericalVariance()} - * - * @return the variance of this distribution - */ - protected double calculateNumericalVariance() { - final double sh = getShape(); - final double sc = getScale(); - final double mn = getNumericalMean(); - - return (sc * sc) * FastMath.exp(Gamma.logGamma(1 + (2 / sh))) - - (mn * mn); - } - - /** - * {@inheritDoc} - * - * The lower bound of the support is always 0 no matter the parameters. - * - * @return lower bound of the support (always 0) - */ - public double getSupportLowerBound() { - return 0; - } - - /** - * {@inheritDoc} - * - * The upper bound of the support is always positive infinity - * no matter the parameters. - * - * @return upper bound of the support (always - * {@code Double.POSITIVE_INFINITY}) - */ - public double getSupportUpperBound() { - return Double.POSITIVE_INFINITY; - } - - /** {@inheritDoc} */ - public boolean isSupportLowerBoundInclusive() { - return true; - } - - /** {@inheritDoc} */ - public boolean isSupportUpperBoundInclusive() { - return false; - } - - /** - * {@inheritDoc} - * - * The support of this distribution is connected. - * - * @return {@code true} - */ - public boolean isSupportConnected() { - return true; - } -} -
