1 /*
2 * Copyright 2004 The Apache Software Foundation.
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 * http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
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 import org.apache.commons.math.util.MathUtils;
23
24 /**
25 * Implementation for the {@link PoissonDistribution}.
26 *
27 * @version $Revision: 355770 $ $Date: 2005-12-10 12:48:57 -0700 (Sat, 10 Dec 2005) $
28 */
29 public class PoissonDistributionImpl extends AbstractIntegerDistribution
30 implements PoissonDistribution, Serializable {
31
32 /** Serializable version identifier */
33 private static final long serialVersionUID = -3349935121172596109L;
34
35 /**
36 * Holds the Poisson mean for the distribution.
37 */
38 private double mean;
39
40 /**
41 * Create a new Poisson distribution with the given the mean.
42 * The mean value must be positive; otherwise an
43 * <code>IllegalArgument</code> is thrown.
44 *
45 * @param p the Poisson mean
46 * @throws IllegalArgumentException if p ≤ 0
47 */
48 public PoissonDistributionImpl(double p) {
49 super();
50 setMean(p);
51 }
52
53 /**
54 * Get the Poisson mean for the distribution.
55 *
56 * @return the Poisson mean for the distribution.
57 */
58 public double getMean() {
59 return this.mean;
60 }
61
62 /**
63 * Set the Poisson mean for the distribution.
64 * The mean value must be positive; otherwise an
65 * <code>IllegalArgument</code> is thrown.
66 *
67 * @param p the Poisson mean value
68 * @throws IllegalArgumentException if p ≤ 0
69 */
70 public void setMean(double p) {
71 if (p <= 0) {
72 throw new IllegalArgumentException(
73 "The Poisson mean must be positive");
74 }
75 this.mean = p;
76 }
77
78 /**
79 * The probability mass function P(X = x) for a Poisson distribution.
80 *
81 * @param x the value at which the probability density function is evaluated.
82 * @return the value of the probability mass function at x
83 */
84 public double probability(int x) {
85 if (x < 0 || x == Integer.MAX_VALUE) {
86 return 0;
87 }
88 return Math.pow(getMean(), x) /
89 MathUtils.factorialDouble(x) * Math.exp(-mean);
90 }
91
92 /**
93 * The probability distribution function P(X <= x) for a Poisson distribution.
94 *
95 * @param x the value at which the PDF is evaluated.
96 * @return Poisson distribution function evaluated at x
97 * @throws MathException if the cumulative probability can not be
98 * computed due to convergence or other numerical errors.
99 */
100 public double cumulativeProbability(int x) throws MathException {
101 if (x < 0) {
102 return 0;
103 }
104 if (x == Integer.MAX_VALUE) {
105 return 1;
106 }
107 return Gamma.regularizedGammaQ((double)x + 1, mean,
108 1E-12, Integer.MAX_VALUE);
109 }
110
111 /**
112 * Calculates the Poisson distribution function using a normal
113 * approximation. The <code>N(mean, sqrt(mean))</code>
114 * distribution is used to approximate the Poisson distribution.
115 * <p>
116 * The computation uses "half-correction" -- evaluating the normal
117 * distribution function at <code>x + 0.5</code>
118 *
119 * @param x the upper bound, inclusive
120 * @return the distribution function value calculated using a normal approximation
121 * @throws MathException if an error occurs computing the normal approximation
122 */
123 public double normalApproximateProbability(int x) throws MathException {
124 NormalDistribution normal = DistributionFactory.newInstance()
125 .createNormalDistribution(getMean(),
126 Math.sqrt(getMean()));
127
128 // calculate the probability using half-correction
129 return normal.cumulativeProbability(x + 0.5);
130 }
131
132 /**
133 * Access the domain value lower bound, based on <code>p</code>, used to
134 * bracket a CDF root. This method is used by
135 * {@link #inverseCumulativeProbability(double)} to find critical values.
136 *
137 * @param p the desired probability for the critical value
138 * @return domain lower bound
139 */
140 protected int getDomainLowerBound(double p) {
141 return 0;
142 }
143
144 /**
145 * Access the domain value upper 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 upper bound
151 */
152 protected int getDomainUpperBound(double p) {
153 return Integer.MAX_VALUE;
154 }
155
156 }