1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16 package org.apache.commons.math.distribution;
17
18 import java.io.Serializable;
19
20 import org.apache.commons.math.MathException;
21 import org.apache.commons.math.special.Gamma;
22
23 /**
24 * The default implementation of {@link GammaDistribution}.
25 *
26 * @version $Revision: 355770 $ $Date: 2005-12-10 12:48:57 -0700 (Sat, 10 Dec 2005) $
27 */
28 public class GammaDistributionImpl extends AbstractContinuousDistribution
29 implements GammaDistribution, Serializable {
30
31 /** Serializable version identifier */
32 private static final long serialVersionUID = -3239549463135430361L;
33
34 /** The shape parameter. */
35 private double alpha;
36
37 /** The scale parameter. */
38 private double beta;
39
40 /**
41 * Create a new gamma distribution with the given alpha and beta values.
42 * @param alpha the shape parameter.
43 * @param beta the scale parameter.
44 */
45 public GammaDistributionImpl(double alpha, double beta) {
46 super();
47 setAlpha(alpha);
48 setBeta(beta);
49 }
50
51 /**
52 * For this disbution, X, this method returns P(X < x).
53 *
54 * The implementation of this method is based on:
55 * <ul>
56 * <li>
57 * <a href="http://mathworld.wolfram.com/Chi-SquaredDistribution.html">
58 * Chi-Squared Distribution</a>, equation (9).</li>
59 * <li>Casella, G., & Berger, R. (1990). <i>Statistical Inference</i>.
60 * Belmont, CA: Duxbury Press.</li>
61 * </ul>
62 *
63 * @param x the value at which the CDF is evaluated.
64 * @return CDF for this distribution.
65 * @throws MathException if the cumulative probability can not be
66 * computed due to convergence or other numerical errors.
67 */
68 public double cumulativeProbability(double x) throws MathException{
69 double ret;
70
71 if (x <= 0.0) {
72 ret = 0.0;
73 } else {
74 ret = Gamma.regularizedGammaP(getAlpha(), x / getBeta());
75 }
76
77 return ret;
78 }
79
80 /**
81 * For this distribution, X, this method returns the critical point x, such
82 * that P(X < x) = <code>p</code>.
83 * <p>
84 * Returns 0 for p=0 and <code>Double.POSITIVE_INFINITY</code> for p=1.
85 *
86 * @param p the desired probability
87 * @return x, such that P(X < x) = <code>p</code>
88 * @throws MathException if the inverse cumulative probability can not be
89 * computed due to convergence or other numerical errors.
90 * @throws IllegalArgumentException if <code>p</code> is not a valid
91 * probability.
92 */
93 public double inverseCumulativeProbability(final double p)
94 throws MathException {
95 if (p == 0) {
96 return 0d;
97 }
98 if (p == 1) {
99 return Double.POSITIVE_INFINITY;
100 }
101 return super.inverseCumulativeProbability(p);
102 }
103
104 /**
105 * Modify the shape parameter, alpha.
106 * @param alpha the new shape parameter.
107 * @throws IllegalArgumentException if <code>alpha</code> is not positive.
108 */
109 public void setAlpha(double alpha) {
110 if (alpha <= 0.0) {
111 throw new IllegalArgumentException("alpha must be positive");
112 }
113 this.alpha = alpha;
114 }
115
116 /**
117 * Access the shape parameter, alpha
118 * @return alpha.
119 */
120 public double getAlpha() {
121 return alpha;
122 }
123
124 /**
125 * Modify the scale parameter, beta.
126 * @param beta the new scale parameter.
127 * @throws IllegalArgumentException if <code>beta</code> is not positive.
128 */
129 public void setBeta(double beta) {
130 if (beta <= 0.0) {
131 throw new IllegalArgumentException("beta must be positive");
132 }
133 this.beta = beta;
134 }
135
136 /**
137 * Access the scale parameter, beta
138 * @return beta.
139 */
140 public double getBeta() {
141 return beta;
142 }
143
144 /**
145 * Access the domain value lower bound, based on <code>p</code>, used to
146 * bracket a CDF root. This method is used by
147 * {@link #inverseCumulativeProbability(double)} to find critical values.
148 *
149 * @param p the desired probability for the critical value
150 * @return domain value lower bound, i.e.
151 * P(X < <i>lower bound</i>) < <code>p</code>
152 */
153 protected double getDomainLowerBound(double p) {
154
155 return Double.MIN_VALUE;
156 }
157
158 /**
159 * Access the domain value upper bound, based on <code>p</code>, used to
160 * bracket a CDF root. This method is used by
161 * {@link #inverseCumulativeProbability(double)} to find critical values.
162 *
163 * @param p the desired probability for the critical value
164 * @return domain value upper bound, i.e.
165 * P(X < <i>upper bound</i>) > <code>p</code>
166 */
167 protected double getDomainUpperBound(double p) {
168
169
170
171
172 double ret;
173
174 if (p < .5) {
175
176 ret = getAlpha() * getBeta();
177 } else {
178
179 ret = Double.MAX_VALUE;
180 }
181
182 return ret;
183 }
184
185 /**
186 * Access the initial domain value, based on <code>p</code>, used to
187 * bracket a CDF root. This method is used by
188 * {@link #inverseCumulativeProbability(double)} to find critical values.
189 *
190 * @param p the desired probability for the critical value
191 * @return initial domain value
192 */
193 protected double getInitialDomain(double p) {
194
195
196
197 double ret;
198
199 if (p < .5) {
200
201 ret = getAlpha() * getBeta() * .5;
202 } else {
203
204 ret = getAlpha() * getBeta();
205 }
206
207 return ret;
208 }
209 }