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.random.RandomDataImpl;
|
---|
25 | import agents.org.apache.commons.math.util.FastMath;
|
---|
26 |
|
---|
27 |
|
---|
28 | /**
|
---|
29 | * Base class for integer-valued discrete distributions. Default
|
---|
30 | * implementations are provided for some of the methods that do not vary
|
---|
31 | * from distribution to distribution.
|
---|
32 | *
|
---|
33 | * @version $Revision: 1067494 $ $Date: 2011-02-05 20:49:07 +0100 (sam. 05 févr. 2011) $
|
---|
34 | */
|
---|
35 | public abstract class AbstractIntegerDistribution extends AbstractDistribution
|
---|
36 | implements IntegerDistribution, Serializable {
|
---|
37 |
|
---|
38 | /** Serializable version identifier */
|
---|
39 | private static final long serialVersionUID = -1146319659338487221L;
|
---|
40 |
|
---|
41 | /**
|
---|
42 | * RandomData instance used to generate samples from the distribution
|
---|
43 | * @since 2.2
|
---|
44 | */
|
---|
45 | protected final RandomDataImpl randomData = new RandomDataImpl();
|
---|
46 |
|
---|
47 | /**
|
---|
48 | * Default constructor.
|
---|
49 | */
|
---|
50 | protected AbstractIntegerDistribution() {
|
---|
51 | super();
|
---|
52 | }
|
---|
53 |
|
---|
54 | /**
|
---|
55 | * For a random variable X whose values are distributed according
|
---|
56 | * to this distribution, this method returns P(X ≤ x). In other words,
|
---|
57 | * this method represents the (cumulative) distribution function, or
|
---|
58 | * CDF, for this distribution.
|
---|
59 | * <p>
|
---|
60 | * If <code>x</code> does not represent an integer value, the CDF is
|
---|
61 | * evaluated at the greatest integer less than x.
|
---|
62 | *
|
---|
63 | * @param x the value at which the distribution function is evaluated.
|
---|
64 | * @return cumulative probability that a random variable with this
|
---|
65 | * distribution takes a value less than or equal to <code>x</code>
|
---|
66 | * @throws MathException if the cumulative probability can not be
|
---|
67 | * computed due to convergence or other numerical errors.
|
---|
68 | */
|
---|
69 | public double cumulativeProbability(double x) throws MathException {
|
---|
70 | return cumulativeProbability((int) FastMath.floor(x));
|
---|
71 | }
|
---|
72 |
|
---|
73 | /**
|
---|
74 | * For a random variable X whose values are distributed according
|
---|
75 | * to this distribution, this method returns P(x0 ≤ X ≤ x1).
|
---|
76 | *
|
---|
77 | * @param x0 the (inclusive) lower bound
|
---|
78 | * @param x1 the (inclusive) upper bound
|
---|
79 | * @return the probability that a random variable with this distribution
|
---|
80 | * will take a value between <code>x0</code> and <code>x1</code>,
|
---|
81 | * including the endpoints.
|
---|
82 | * @throws MathException if the cumulative probability can not be
|
---|
83 | * computed due to convergence or other numerical errors.
|
---|
84 | * @throws IllegalArgumentException if <code>x0 > x1</code>
|
---|
85 | */
|
---|
86 | @Override
|
---|
87 | public double cumulativeProbability(double x0, double x1)
|
---|
88 | throws MathException {
|
---|
89 | if (x0 > x1) {
|
---|
90 | throw MathRuntimeException.createIllegalArgumentException(
|
---|
91 | LocalizedFormats.LOWER_ENDPOINT_ABOVE_UPPER_ENDPOINT, x0, x1);
|
---|
92 | }
|
---|
93 | if (FastMath.floor(x0) < x0) {
|
---|
94 | return cumulativeProbability(((int) FastMath.floor(x0)) + 1,
|
---|
95 | (int) FastMath.floor(x1)); // don't want to count mass below x0
|
---|
96 | } else { // x0 is mathematical integer, so use as is
|
---|
97 | return cumulativeProbability((int) FastMath.floor(x0),
|
---|
98 | (int) FastMath.floor(x1));
|
---|
99 | }
|
---|
100 | }
|
---|
101 |
|
---|
102 | /**
|
---|
103 | * For a random variable X whose values are distributed according
|
---|
104 | * to this distribution, this method returns P(X ≤ x). In other words,
|
---|
105 | * this method represents the probability distribution function, or PDF,
|
---|
106 | * for this distribution.
|
---|
107 | *
|
---|
108 | * @param x the value at which the PDF is evaluated.
|
---|
109 | * @return PDF for this distribution.
|
---|
110 | * @throws MathException if the cumulative probability can not be
|
---|
111 | * computed due to convergence or other numerical errors.
|
---|
112 | */
|
---|
113 | public abstract double cumulativeProbability(int x) throws MathException;
|
---|
114 |
|
---|
115 | /**
|
---|
116 | * For a random variable X whose values are distributed according
|
---|
117 | * to this distribution, this method returns P(X = x). In other words, this
|
---|
118 | * method represents the probability mass function, or PMF, for the distribution.
|
---|
119 | * <p>
|
---|
120 | * If <code>x</code> does not represent an integer value, 0 is returned.
|
---|
121 | *
|
---|
122 | * @param x the value at which the probability density function is evaluated
|
---|
123 | * @return the value of the probability density function at x
|
---|
124 | */
|
---|
125 | public double probability(double x) {
|
---|
126 | double fl = FastMath.floor(x);
|
---|
127 | if (fl == x) {
|
---|
128 | return this.probability((int) x);
|
---|
129 | } else {
|
---|
130 | return 0;
|
---|
131 | }
|
---|
132 | }
|
---|
133 |
|
---|
134 | /**
|
---|
135 | * For a random variable X whose values are distributed according
|
---|
136 | * to this distribution, this method returns P(x0 ≤ X ≤ x1).
|
---|
137 | *
|
---|
138 | * @param x0 the inclusive, lower bound
|
---|
139 | * @param x1 the inclusive, upper bound
|
---|
140 | * @return the cumulative probability.
|
---|
141 | * @throws MathException if the cumulative probability can not be
|
---|
142 | * computed due to convergence or other numerical errors.
|
---|
143 | * @throws IllegalArgumentException if x0 > x1
|
---|
144 | */
|
---|
145 | public double cumulativeProbability(int x0, int x1) throws MathException {
|
---|
146 | if (x0 > x1) {
|
---|
147 | throw MathRuntimeException.createIllegalArgumentException(
|
---|
148 | LocalizedFormats.LOWER_ENDPOINT_ABOVE_UPPER_ENDPOINT, x0, x1);
|
---|
149 | }
|
---|
150 | return cumulativeProbability(x1) - cumulativeProbability(x0 - 1);
|
---|
151 | }
|
---|
152 |
|
---|
153 | /**
|
---|
154 | * For a random variable X whose values are distributed according
|
---|
155 | * to this distribution, this method returns the largest x, such
|
---|
156 | * that P(X ≤ x) ≤ <code>p</code>.
|
---|
157 | *
|
---|
158 | * @param p the desired probability
|
---|
159 | * @return the largest x such that P(X ≤ x) <= p
|
---|
160 | * @throws MathException if the inverse cumulative probability can not be
|
---|
161 | * computed due to convergence or other numerical errors.
|
---|
162 | * @throws IllegalArgumentException if p < 0 or p > 1
|
---|
163 | */
|
---|
164 | public int inverseCumulativeProbability(final double p) throws MathException{
|
---|
165 | if (p < 0.0 || p > 1.0) {
|
---|
166 | throw MathRuntimeException.createIllegalArgumentException(
|
---|
167 | LocalizedFormats.OUT_OF_RANGE_SIMPLE, p, 0.0, 1.0);
|
---|
168 | }
|
---|
169 |
|
---|
170 | // by default, do simple bisection.
|
---|
171 | // subclasses can override if there is a better method.
|
---|
172 | int x0 = getDomainLowerBound(p);
|
---|
173 | int x1 = getDomainUpperBound(p);
|
---|
174 | double pm;
|
---|
175 | while (x0 < x1) {
|
---|
176 | int xm = x0 + (x1 - x0) / 2;
|
---|
177 | pm = checkedCumulativeProbability(xm);
|
---|
178 | if (pm > p) {
|
---|
179 | // update x1
|
---|
180 | if (xm == x1) {
|
---|
181 | // this can happen with integer division
|
---|
182 | // simply decrement x1
|
---|
183 | --x1;
|
---|
184 | } else {
|
---|
185 | // update x1 normally
|
---|
186 | x1 = xm;
|
---|
187 | }
|
---|
188 | } else {
|
---|
189 | // update x0
|
---|
190 | if (xm == x0) {
|
---|
191 | // this can happen with integer division
|
---|
192 | // simply increment x0
|
---|
193 | ++x0;
|
---|
194 | } else {
|
---|
195 | // update x0 normally
|
---|
196 | x0 = xm;
|
---|
197 | }
|
---|
198 | }
|
---|
199 | }
|
---|
200 |
|
---|
201 | // insure x0 is the correct critical point
|
---|
202 | pm = checkedCumulativeProbability(x0);
|
---|
203 | while (pm > p) {
|
---|
204 | --x0;
|
---|
205 | pm = checkedCumulativeProbability(x0);
|
---|
206 | }
|
---|
207 |
|
---|
208 | return x0;
|
---|
209 | }
|
---|
210 |
|
---|
211 | /**
|
---|
212 | * Reseeds the random generator used to generate samples.
|
---|
213 | *
|
---|
214 | * @param seed the new seed
|
---|
215 | * @since 2.2
|
---|
216 | */
|
---|
217 | public void reseedRandomGenerator(long seed) {
|
---|
218 | randomData.reSeed(seed);
|
---|
219 | }
|
---|
220 |
|
---|
221 | /**
|
---|
222 | * Generates a random value sampled from this distribution. The default
|
---|
223 | * implementation uses the
|
---|
224 | * <a href="http://en.wikipedia.org/wiki/Inverse_transform_sampling"> inversion method.</a>
|
---|
225 | *
|
---|
226 | * @return random value
|
---|
227 | * @since 2.2
|
---|
228 | * @throws MathException if an error occurs generating the random value
|
---|
229 | */
|
---|
230 | public int sample() throws MathException {
|
---|
231 | return randomData.nextInversionDeviate(this);
|
---|
232 | }
|
---|
233 |
|
---|
234 | /**
|
---|
235 | * Generates a random sample from the distribution. The default implementation
|
---|
236 | * generates the sample by calling {@link #sample()} in a loop.
|
---|
237 | *
|
---|
238 | * @param sampleSize number of random values to generate
|
---|
239 | * @since 2.2
|
---|
240 | * @return an array representing the random sample
|
---|
241 | * @throws MathException if an error occurs generating the sample
|
---|
242 | * @throws IllegalArgumentException if sampleSize is not positive
|
---|
243 | */
|
---|
244 | public int[] sample(int sampleSize) throws MathException {
|
---|
245 | if (sampleSize <= 0) {
|
---|
246 | MathRuntimeException.createIllegalArgumentException(LocalizedFormats.NOT_POSITIVE_SAMPLE_SIZE, sampleSize);
|
---|
247 | }
|
---|
248 | int[] out = new int[sampleSize];
|
---|
249 | for (int i = 0; i < sampleSize; i++) {
|
---|
250 | out[i] = sample();
|
---|
251 | }
|
---|
252 | return out;
|
---|
253 | }
|
---|
254 |
|
---|
255 | /**
|
---|
256 | * Computes the cumulative probability function and checks for NaN values returned.
|
---|
257 | * Throws MathException if the value is NaN. Rethrows any MathException encountered
|
---|
258 | * evaluating the cumulative probability function. Throws
|
---|
259 | * MathException if the cumulative probability function returns NaN.
|
---|
260 | *
|
---|
261 | * @param argument input value
|
---|
262 | * @return cumulative probability
|
---|
263 | * @throws MathException if the cumulative probability is NaN
|
---|
264 | */
|
---|
265 | private double checkedCumulativeProbability(int argument) throws MathException {
|
---|
266 | double result = Double.NaN;
|
---|
267 | result = cumulativeProbability(argument);
|
---|
268 | if (Double.isNaN(result)) {
|
---|
269 | throw new MathException(LocalizedFormats.DISCRETE_CUMULATIVE_PROBABILITY_RETURNED_NAN, argument);
|
---|
270 | }
|
---|
271 | return result;
|
---|
272 | }
|
---|
273 |
|
---|
274 | /**
|
---|
275 | * Access the domain value lower bound, based on <code>p</code>, used to
|
---|
276 | * bracket a PDF root. This method is used by
|
---|
277 | * {@link #inverseCumulativeProbability(double)} to find critical values.
|
---|
278 | *
|
---|
279 | * @param p the desired probability for the critical value
|
---|
280 | * @return domain value lower bound, i.e.
|
---|
281 | * P(X < <i>lower bound</i>) < <code>p</code>
|
---|
282 | */
|
---|
283 | protected abstract int getDomainLowerBound(double p);
|
---|
284 |
|
---|
285 | /**
|
---|
286 | * Access the domain value upper bound, based on <code>p</code>, used to
|
---|
287 | * bracket a PDF root. This method is used by
|
---|
288 | * {@link #inverseCumulativeProbability(double)} to find critical values.
|
---|
289 | *
|
---|
290 | * @param p the desired probability for the critical value
|
---|
291 | * @return domain value upper bound, i.e.
|
---|
292 | * P(X < <i>upper bound</i>) > <code>p</code>
|
---|
293 | */
|
---|
294 | protected abstract int getDomainUpperBound(double p);
|
---|
295 |
|
---|
296 | /**
|
---|
297 | * Use this method to get information about whether the lower bound
|
---|
298 | * of the support is inclusive or not. For discrete support,
|
---|
299 | * only true here is meaningful.
|
---|
300 | *
|
---|
301 | * @return true (always but at Integer.MIN_VALUE because of the nature of discrete support)
|
---|
302 | * @since 2.2
|
---|
303 | */
|
---|
304 | public boolean isSupportLowerBoundInclusive() {
|
---|
305 | return true;
|
---|
306 | }
|
---|
307 |
|
---|
308 | /**
|
---|
309 | * Use this method to get information about whether the upper bound
|
---|
310 | * of the support is inclusive or not. For discrete support,
|
---|
311 | * only true here is meaningful.
|
---|
312 | *
|
---|
313 | * @return true (always but at Integer.MAX_VALUE because of the nature of discrete support)
|
---|
314 | * @since 2.2
|
---|
315 | */
|
---|
316 | public boolean isSupportUpperBoundInclusive() {
|
---|
317 | return true;
|
---|
318 | }
|
---|
319 | }
|
---|