1 package nl.tudelft.simulation.jstats.distributions;
2
3 import org.djutils.exceptions.Throw;
4
5 import nl.tudelft.simulation.jstats.math.ProbMath;
6 import nl.tudelft.simulation.jstats.streams.StreamInterface;
7
8 /**
9 * The Negative Binomial distribution. It is also known as the Pascal distribution or Pólya distribution. It gives the
10 * probability of x failures where there are s-1 successes in a total of x+s-1 Bernoulli trials, and trial (x+s) is a success.
11 * The chance for success is p for each trial. For more information on this distribution see
12 * <a href="https://mathworld.wolfram.com/NegativeBinomialDistribution.html">
13 * https://mathworld.wolfram.com/NegativeBinomialDistribution.html </a>
14 * <p>
15 * Copyright (c) 2002-2024 Delft University of Technology, Jaffalaan 5, 2628 BX Delft, the Netherlands. All rights reserved. See
16 * for project information <a href="https://simulation.tudelft.nl/" target="_blank"> https://simulation.tudelft.nl</a>. The DSOL
17 * project is distributed under a three-clause BSD-style license, which can be found at
18 * <a href="https://https://simulation.tudelft.nl/dsol/docs/latest/license.html" target="_blank">
19 * https://https://simulation.tudelft.nl/dsol/docs/latest/license.html</a>.
20 * </p>
21 * @author <a href="https://www.linkedin.com/in/peterhmjacobs">Peter Jacobs </a>
22 * @author <a href="https://www.tudelft.nl/averbraeck">Alexander Verbraeck</a>
23 */
24 public class DistNegBinomial extends DistDiscrete
25 {
26 /** */
27 private static final long serialVersionUID = 1L;
28
29 /** s is the number of successes in the sequence of (x+n) trials, where trial (x+n) is a success. */
30 private int s;
31
32 /** p is the probability of success for each individual trial in the negative binomial distribution. */
33 private double p;
34
35 /** lnp is a helper variable equal to ln(1-p) to avoid repetitive calculation. */
36 private double lnp;
37
38 /**
39 * constructs a new negative binomial distribution.
40 * @param stream StreamInterface; the random number stream
41 * @param s int; the number of successes in the sequence of (x+n) trials, where trial (x+n) is a success
42 * @param p double; the probability of success for each individual trial in the negative binomial distribution
43 * @throws IllegalArgumentException when s <= 0 or p <= 0 or p >= 1
44 */
45 public DistNegBinomial(final StreamInterface stream, final int s, final double p)
46 {
47 super(stream);
48 Throw.when(s <= 0 || p <= 0.0 || p >= 1.0, IllegalArgumentException.class,
49 "Error NegBinomial - s<=0 or p<=0.0 or p>=1.0");
50 this.s = s;
51 this.p = p;
52 this.lnp = Math.log(1.0 - this.p);
53 }
54
55 /** {@inheritDoc} */
56 @Override
57 public long draw()
58 {
59 long x = 0;
60 for (int i = 0; i < this.s; i++)
61 {
62 double u = this.stream.nextDouble();
63 x = x + (long) (Math.floor(Math.log(u) / this.lnp));
64 }
65 return x;
66 }
67
68 /** {@inheritDoc} */
69 @Override
70 public double probability(final long observation)
71 {
72 if (observation >= 0)
73 {
74 return ProbMath.combinations(this.s + observation - 1, observation) * Math.pow(this.p, this.s)
75 * Math.pow(1 - this.p, observation);
76 }
77 return 0.0;
78 }
79
80 /**
81 * Return the number of successes in the sequence of (x+n) trials, where trial (x+n) is a success.
82 * @return int; the number of successes in the sequence of (x+n) trials, where trial (x+n) is a success
83 */
84 public int getS()
85 {
86 return this.s;
87 }
88
89 /**
90 * Return the probability of success for each individual trial in the negative binomial distribution.
91 * @return double; the probability of success for each individual trial in the negative binomial distribution
92 */
93 public double getP()
94 {
95 return this.p;
96 }
97
98 /** {@inheritDoc} */
99 @Override
100 public String toString()
101 {
102 return "NegBinomial(" + this.s + "," + this.p + ")";
103 }
104 }