1 | /*
|
---|
2 | * Licensed to the Apache Software Foundation (ASF) under one or more
|
---|
3 | * contributor license agreements. See the NOTICE file distributed with
|
---|
4 | * this work for additional information regarding copyright ownership.
|
---|
5 | * The ASF licenses this file to You under the Apache License, Version 2.0
|
---|
6 | * (the "License"); you may not use this file except in compliance with
|
---|
7 | * the License. You may obtain a copy of the License at
|
---|
8 | *
|
---|
9 | * http://www.apache.org/licenses/LICENSE-2.0
|
---|
10 | *
|
---|
11 | * Unless required by applicable law or agreed to in writing, software
|
---|
12 | * distributed under the License is distributed on an "AS IS" BASIS,
|
---|
13 | * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
|
---|
14 | * See the License for the specific language governing permissions and
|
---|
15 | * limitations under the License.
|
---|
16 | */
|
---|
17 | package agents.org.apache.commons.math.distribution;
|
---|
18 |
|
---|
19 | import java.io.Serializable;
|
---|
20 |
|
---|
21 | import agents.org.apache.commons.math.MathException;
|
---|
22 | import agents.org.apache.commons.math.MathRuntimeException;
|
---|
23 | import agents.org.apache.commons.math.exception.util.LocalizedFormats;
|
---|
24 | import agents.org.apache.commons.math.special.Beta;
|
---|
25 | import agents.org.apache.commons.math.util.FastMath;
|
---|
26 |
|
---|
27 | /**
|
---|
28 | * The default implementation of {@link BinomialDistribution}.
|
---|
29 | *
|
---|
30 | * @version $Revision: 1054524 $ $Date: 2011-01-03 05:59:18 +0100 (lun. 03 janv. 2011) $
|
---|
31 | */
|
---|
32 | public class BinomialDistributionImpl extends AbstractIntegerDistribution
|
---|
33 | implements BinomialDistribution, Serializable {
|
---|
34 |
|
---|
35 | /** Serializable version identifier */
|
---|
36 | private static final long serialVersionUID = 6751309484392813623L;
|
---|
37 |
|
---|
38 | /** The number of trials. */
|
---|
39 | private int numberOfTrials;
|
---|
40 |
|
---|
41 | /** The probability of success. */
|
---|
42 | private double probabilityOfSuccess;
|
---|
43 |
|
---|
44 | /**
|
---|
45 | * Create a binomial distribution with the given number of trials and
|
---|
46 | * probability of success.
|
---|
47 | *
|
---|
48 | * @param trials the number of trials.
|
---|
49 | * @param p the probability of success.
|
---|
50 | */
|
---|
51 | public BinomialDistributionImpl(int trials, double p) {
|
---|
52 | super();
|
---|
53 | setNumberOfTrialsInternal(trials);
|
---|
54 | setProbabilityOfSuccessInternal(p);
|
---|
55 | }
|
---|
56 |
|
---|
57 | /**
|
---|
58 | * Access the number of trials for this distribution.
|
---|
59 | *
|
---|
60 | * @return the number of trials.
|
---|
61 | */
|
---|
62 | public int getNumberOfTrials() {
|
---|
63 | return numberOfTrials;
|
---|
64 | }
|
---|
65 |
|
---|
66 | /**
|
---|
67 | * Access the probability of success for this distribution.
|
---|
68 | *
|
---|
69 | * @return the probability of success.
|
---|
70 | */
|
---|
71 | public double getProbabilityOfSuccess() {
|
---|
72 | return probabilityOfSuccess;
|
---|
73 | }
|
---|
74 |
|
---|
75 | /**
|
---|
76 | * Change the number of trials for this distribution.
|
---|
77 | *
|
---|
78 | * @param trials the new number of trials.
|
---|
79 | * @throws IllegalArgumentException if <code>trials</code> is not a valid
|
---|
80 | * number of trials.
|
---|
81 | * @deprecated as of 2.1 (class will become immutable in 3.0)
|
---|
82 | */
|
---|
83 | @Deprecated
|
---|
84 | public void setNumberOfTrials(int trials) {
|
---|
85 | setNumberOfTrialsInternal(trials);
|
---|
86 | }
|
---|
87 |
|
---|
88 | /**
|
---|
89 | * Change the number of trials for this distribution.
|
---|
90 | *
|
---|
91 | * @param trials the new number of trials.
|
---|
92 | * @throws IllegalArgumentException if <code>trials</code> is not a valid
|
---|
93 | * number of trials.
|
---|
94 | */
|
---|
95 | private void setNumberOfTrialsInternal(int trials) {
|
---|
96 | if (trials < 0) {
|
---|
97 | throw MathRuntimeException.createIllegalArgumentException(
|
---|
98 | LocalizedFormats.NEGATIVE_NUMBER_OF_TRIALS, trials);
|
---|
99 | }
|
---|
100 | numberOfTrials = trials;
|
---|
101 | }
|
---|
102 |
|
---|
103 | /**
|
---|
104 | * Change the probability of success for this distribution.
|
---|
105 | *
|
---|
106 | * @param p the new probability of success.
|
---|
107 | * @throws IllegalArgumentException if <code>p</code> is not a valid
|
---|
108 | * probability.
|
---|
109 | * @deprecated as of 2.1 (class will become immutable in 3.0)
|
---|
110 | */
|
---|
111 | @Deprecated
|
---|
112 | public void setProbabilityOfSuccess(double p) {
|
---|
113 | setProbabilityOfSuccessInternal(p);
|
---|
114 | }
|
---|
115 |
|
---|
116 | /**
|
---|
117 | * Change the probability of success for this distribution.
|
---|
118 | *
|
---|
119 | * @param p the new probability of success.
|
---|
120 | * @throws IllegalArgumentException if <code>p</code> is not a valid
|
---|
121 | * probability.
|
---|
122 | */
|
---|
123 | private void setProbabilityOfSuccessInternal(double p) {
|
---|
124 | if (p < 0.0 || p > 1.0) {
|
---|
125 | throw MathRuntimeException.createIllegalArgumentException(
|
---|
126 | LocalizedFormats.OUT_OF_RANGE_SIMPLE, p, 0.0, 1.0);
|
---|
127 | }
|
---|
128 | probabilityOfSuccess = p;
|
---|
129 | }
|
---|
130 |
|
---|
131 | /**
|
---|
132 | * Access the domain value lower bound, based on <code>p</code>, used to
|
---|
133 | * bracket a PDF root.
|
---|
134 | *
|
---|
135 | * @param p the desired probability for the critical value
|
---|
136 | * @return domain value lower bound, i.e. P(X < <i>lower bound</i>) <
|
---|
137 | * <code>p</code>
|
---|
138 | */
|
---|
139 | @Override
|
---|
140 | protected int getDomainLowerBound(double p) {
|
---|
141 | return -1;
|
---|
142 | }
|
---|
143 |
|
---|
144 | /**
|
---|
145 | * Access the domain value upper bound, based on <code>p</code>, used to
|
---|
146 | * bracket a PDF root.
|
---|
147 | *
|
---|
148 | * @param p the desired probability for the critical value
|
---|
149 | * @return domain value upper bound, i.e. P(X < <i>upper bound</i>) >
|
---|
150 | * <code>p</code>
|
---|
151 | */
|
---|
152 | @Override
|
---|
153 | protected int getDomainUpperBound(double p) {
|
---|
154 | return numberOfTrials;
|
---|
155 | }
|
---|
156 |
|
---|
157 | /**
|
---|
158 | * For this distribution, X, this method returns P(X ≤ x).
|
---|
159 | *
|
---|
160 | * @param x the value at which the PDF is evaluated.
|
---|
161 | * @return PDF for this distribution.
|
---|
162 | * @throws MathException if the cumulative probability can not be computed
|
---|
163 | * due to convergence or other numerical errors.
|
---|
164 | */
|
---|
165 | @Override
|
---|
166 | public double cumulativeProbability(int x) throws MathException {
|
---|
167 | double ret;
|
---|
168 | if (x < 0) {
|
---|
169 | ret = 0.0;
|
---|
170 | } else if (x >= numberOfTrials) {
|
---|
171 | ret = 1.0;
|
---|
172 | } else {
|
---|
173 | ret = 1.0 - Beta.regularizedBeta(getProbabilityOfSuccess(),
|
---|
174 | x + 1.0, numberOfTrials - x);
|
---|
175 | }
|
---|
176 | return ret;
|
---|
177 | }
|
---|
178 |
|
---|
179 | /**
|
---|
180 | * For this distribution, X, this method returns P(X = x).
|
---|
181 | *
|
---|
182 | * @param x the value at which the PMF is evaluated.
|
---|
183 | * @return PMF for this distribution.
|
---|
184 | */
|
---|
185 | public double probability(int x) {
|
---|
186 | double ret;
|
---|
187 | if (x < 0 || x > numberOfTrials) {
|
---|
188 | ret = 0.0;
|
---|
189 | } else {
|
---|
190 | ret = FastMath.exp(SaddlePointExpansion.logBinomialProbability(x,
|
---|
191 | numberOfTrials, probabilityOfSuccess,
|
---|
192 | 1.0 - probabilityOfSuccess));
|
---|
193 | }
|
---|
194 | return ret;
|
---|
195 | }
|
---|
196 |
|
---|
197 | /**
|
---|
198 | * For this distribution, X, this method returns the largest x, such that
|
---|
199 | * P(X ≤ x) ≤ <code>p</code>.
|
---|
200 | * <p>
|
---|
201 | * Returns <code>-1</code> for p=0 and <code>Integer.MAX_VALUE</code> for
|
---|
202 | * p=1.
|
---|
203 | * </p>
|
---|
204 | *
|
---|
205 | * @param p the desired probability
|
---|
206 | * @return the largest x such that P(X ≤ x) <= p
|
---|
207 | * @throws MathException if the inverse cumulative probability can not be
|
---|
208 | * computed due to convergence or other numerical errors.
|
---|
209 | * @throws IllegalArgumentException if p < 0 or p > 1
|
---|
210 | */
|
---|
211 | @Override
|
---|
212 | public int inverseCumulativeProbability(final double p)
|
---|
213 | throws MathException {
|
---|
214 | // handle extreme values explicitly
|
---|
215 | if (p == 0) {
|
---|
216 | return -1;
|
---|
217 | }
|
---|
218 | if (p == 1) {
|
---|
219 | return Integer.MAX_VALUE;
|
---|
220 | }
|
---|
221 |
|
---|
222 | // use default bisection impl
|
---|
223 | return super.inverseCumulativeProbability(p);
|
---|
224 | }
|
---|
225 |
|
---|
226 | /**
|
---|
227 | * Returns the lower bound of the support for the distribution.
|
---|
228 | *
|
---|
229 | * The lower bound of the support is always 0 no matter the number of trials
|
---|
230 | * and probability parameter.
|
---|
231 | *
|
---|
232 | * @return lower bound of the support (always 0)
|
---|
233 | * @since 2.2
|
---|
234 | */
|
---|
235 | public int getSupportLowerBound() {
|
---|
236 | return 0;
|
---|
237 | }
|
---|
238 |
|
---|
239 | /**
|
---|
240 | * Returns the upper bound of the support for the distribution.
|
---|
241 | *
|
---|
242 | * The upper bound of the support is the number of trials.
|
---|
243 | *
|
---|
244 | * @return upper bound of the support (equal to number of trials)
|
---|
245 | * @since 2.2
|
---|
246 | */
|
---|
247 | public int getSupportUpperBound() {
|
---|
248 | return getNumberOfTrials();
|
---|
249 | }
|
---|
250 |
|
---|
251 | /**
|
---|
252 | * Returns the mean.
|
---|
253 | *
|
---|
254 | * For <code>n</code> number of trials and
|
---|
255 | * probability parameter <code>p</code>, the mean is
|
---|
256 | * <code>n * p</code>
|
---|
257 | *
|
---|
258 | * @return the mean
|
---|
259 | * @since 2.2
|
---|
260 | */
|
---|
261 | public double getNumericalMean() {
|
---|
262 | return (double)getNumberOfTrials() * getProbabilityOfSuccess();
|
---|
263 | }
|
---|
264 |
|
---|
265 | /**
|
---|
266 | * Returns the variance.
|
---|
267 | *
|
---|
268 | * For <code>n</code> number of trials and
|
---|
269 | * probability parameter <code>p</code>, the variance is
|
---|
270 | * <code>n * p * (1 - p)</code>
|
---|
271 | *
|
---|
272 | * @return the variance
|
---|
273 | * @since 2.2
|
---|
274 | */
|
---|
275 | public double getNumericalVariance() {
|
---|
276 | final double p = getProbabilityOfSuccess();
|
---|
277 | return (double)getNumberOfTrials() * p * (1 - p);
|
---|
278 | }
|
---|
279 | }
|
---|