For this example, we will opt to #define two macros to control the error
and discrete handling policies. For this simple example, we want to avoid
throwing an exception (the default policy) and just return infinity.
We want to treat the distribution as if it was continuous, so we choose
a discrete_quantile policy of real, rather than the default policy integer_round_outwards.

The Geometric distribution describes the probability (p)
of a number of failures to get the first success in k
Bernoulli trials. (A Bernoulli
trial is one with only two possible outcomes, success of failure,
and p is the probability of success).

If the dice is thrown repeatedly until the first
time a three appears. The probablility distribution
of the number of times it is thrown not
getting a three (not-a-threes
number of failures to get a three) is a geometric
distribution with the success_fraction = 1/6 = 0.1666 ̇.

We therefore start by constructing a geometric distribution with the
one parameter success_fraction, the probability of success.

geometricg6(success_fraction);// type double by default.

To confirm, we can echo the success_fraction parameter of the distribution.

cout<<"success fraction of a six-sided dice is "<<g6.success_fraction()<<endl;

So the probability of getting a three at the first throw (zero failures)
is

cout<<pdf(g6,0)<<endl;// 0.1667cout<<cdf(g6,0)<<endl;// 0.1667

Note that the cdf and pdf are identical because the is only one throw.
If we want the probability of getting the first three
on the 2nd throw:

cout<<pdf(g6,1)<<endl;// 0.1389

If we want the probability of getting the first three
on the 1st or 2nd throw (allowing one failure):

cout<<"pdf(g6, 0) + pdf(g6, 1) = "<<pdf(g6,0)+pdf(g6,1)<<endl;

Or more conveniently, and more generally, we can use the Cumulative Distribution
Function CDF.

cout<<"cdf(g6, 1) = "<<cdf(g6,1)<<endl;// 0.3056

If we allow many more (12) throws, the probability of getting our three
gets very high:

The geometric distribution is related to the negative binomial negative_binomial_distribution(RealTyper,RealTypep);
with parameter r = 1. So we could get the same result
using the negative binomial, but using the geometric the results will
be faster, and may be more accurate.

Note too that Boost.Math geometric distribution is implemented as a continuous
function. Unlike other implementations (for example R) it uses
the number of failures as a real parameter,
not as an integer. If you want this integer behaviour, you may need to
enforce this by rounding the parameter you pass, probably rounding down,
to the nearest integer. For example, R returns the success fraction probability
for all values of failures from 0 to 0.999999 thus:

A company knows from warranty claims that 2% of their products will be
faulty, so the 'success_fraction' of finding a fault is 0.02. It wants
to interview a purchaser of faulty products to assess their 'user experience'.

To estimate how many customers they will probably need to contact in
order to find one who has suffered from the fault, we first construct
a geometric distribution with probability 0.02, and then chose a confidence,
say 80%, 95%, or 99% to finding a customer with a fault. Finally, we
probably want to round up the result to the integer above using the
ceil function. (We could
also use a policy, but that is hardly worthwhile for this simple application.)

(This also assumes that each customer only buys one product: if customers
bought more than one item, the probability of finding a customer with
a fault obviously improves.)

cout.precision(5);geometricg(0.02);// On average, 2 in 100 products are faulty.doublec=0.95;// 95% confidence.cout<<" quantile(g, "<<c<<") = "<<quantile(g,c)<<endl;cout<<"To be "<<c*100<<"% confident of finding we customer with a fault, need to survey "<<ceil(quantile(g,c))<<" customers."<<endl;// 148c=0.99;// Very confident.cout<<"To be "<<c*100<<"% confident of finding we customer with a fault, need to survey "<<ceil(quantile(g,c))<<" customers."<<endl;// 227c=0.80;// Only reasonably confident.cout<<"To be "<<c*100<<"% confident of finding we customer with a fault, need to survey "<<ceil(quantile(g,c))<<" customers."<<endl;// 79

According to Wikipedia, average pro basket ball players get free
throws in the baskets 70 to 80 % of the time, but some get as
high as 95%, and others as low as 50%. Suppose we want to compare the
probabilities of failing to get a score only on the first or on the fifth
shot? To start we will consider the average shooter, say 75%. So we construct
a geometric distribution with success_fraction parameter 75/100 = 0.75.

cout.precision(2);geometricgav(0.75);// Shooter averages 7.5 out of 10 in the basket.

What is probability of getting 1st try in the basket, that is with no
failures?

cout<<"Probability of score on 1st try = "<<pdf(gav,0)<<endl;// 0.75

This is, of course, the success_fraction probability 75%. What is the
probability that the shooter only scores on the fifth shot? So there
are 5-1 = 4 failures before the first success.

Now compare this with the poor and the best players success fraction.
We need to constructing new distributions with the different success
fractions, and then get the corresponding probability density functions
values:

Of course one man's failure is an other man's success. So a fault can
be defined as a 'success'.

If a fault occurs once after 100 flights, then one might naively say
that the risk of fault is obviously 1 in 100 = 1/100, a probability of
0.01.

This is the best estimate we can make, but while it is the truth, it
is not the whole truth, for it hides the big uncertainty when estimating
from a single event. "One swallow doesn't make a summer." To
show the magnitude of the uncertainty, the geometric (or the negative
binomial) distribution can be used.

If we chose the popular 95% confidence in the limits, corresponding to
an alpha of 0.05, because we are calculating a two-sided interval, we
must divide alpha by two.

So while we estimate the probability is 0.01, it might lie between 0.0003
and 0.04. Even if we relax our confidence to alpha = 90%, the bounds
only contract to 0.0005 and 0.03. And if we require a high confidence,
they widen to 0.00005 to 0.05.

In real life, there will usually be more than one event (fault or success),
when the negative binomial, which has the neccessary extra parameter,
will be needed.

As noted above, using a catch block is always a good idea, even if you
hope not to use it!

}catch(conststd::exception&e){// Since we have set an overflow policy of ignore_error,// an overflow exception should never be thrown.std::cout<<"\nMessage from thrown exception was:\n "<<e.what()<<std::endl;

For example, without a ignore domain error policy, if we asked for

pdf(g,-1)

for example, we would get an unhelpful abort, but with a catch:

Message from thrown exception was:
Error in function boost::math::pdf(const exponential_distribution<double>&, double):
Number of failures argument is -1, but must be >= 0 !